吐血悬赏--关于平面应力子程序的问题

tianyaxiaoke · 发表于 2011-5-13 10:33:48

最近一直研究在平面应力状态下如何使用VUMAT子程序建立材料本构模型，先从最简单的Mises屈服方程做起，从已有的全应力状态子程序改写成平面应力状态，参考了很多资料，程序算出的结果与CAE自带的Mises准则结果比较，差别很大。先请教：
1. 常规壳单元（比如S4R单元）是否一定要用平面应力的UMAT或者VUMAT才行？这个是前提，个人感觉问题答案：是。
2. 在问题1答案肯定的基础上，平面应力条件stres33=0如何写入VUMAT？我认为直接在程序里把stressNew(k,3)=0,这样结果可能不对。
3. 我的思路是：法向应变首先按弹性方式计算出，并且代入。等效应变和应力更新按vumat中常规方法估算，法向应力更新时直接赋值为0。法向应变增量修正为法向弹性应变与法向塑性应变之和。这样做对吗？请高手指教。如能帮着解决问题，可将本人全部仿真币奉献！

shawn2008 · 发表于 2011-5-13 10:33:49

第一个是正确的，对于S4R单元，VUMAT中直接分量是3个，而UMAT中直接分量是2个。

对楼主的问题，只需write输出ndi和ndir便可知道答案。很多时候，自己去找答案，往往比等待别人的回答要深刻和快的多。

shawn2008 · 发表于 2011-5-13 12:25:01

1、有时为了将开发的子程序应用到更多的单元，可以先写一个通用的三维应力状态的本构模型子程序，然后在UMAT或VUMAT中根据传入的ndir和nshr，将abaqus主程序传入的应力、应变张量转换到三维状态，然后调用通用的本构模型子程序，得到结果后，再变换回原来的形式，返回给UMAT或VUMAT。

2、的确是直接将应力增量赋值为0.0

3、没理解你的问题，不过法向应变需要更新是因为我们要看单元厚度的变化。

tianyaxiaoke · 发表于 2011-5-13 15:51:12

2# shawn2008

tianyaxiaoke · 发表于 2011-5-13 16:44:26

本帖最后由 tianyaxiaoke 于 2011-7-4 17:45 编辑

十分感谢shawn的热心解答。
1.上面第1个问题，我是想问，采用常规壳单元时一定要用平面应力子程序吗？
您的意思是可以调用三维全应力子程序，再将结算结果转化为平面应力？
2. 第二问题，我在dynaform帮助文件中看到，在UMAT中要通过迭代计算出合适的法向应变增量，使法向应力更新为零。但在vumat中，等效应变增量是近似给出的，没有通过迭代运算，以此值来计算各塑性应变增量，进而算出法向应变增量。应力更新时法向应力显然不为零，所以想问，在Vumat中是直接将法向增量赋值为零对吗？还是也要迭代算出合适的法向应变增量？
3.第三个问题，我写子程序时也参考过您以前的相关帖子。但是子程序计算结果与cae计算结果差别很大。
现把子程序和inp文件（abaqus6.9版本）贴上来，烦请指点。

wyz1983121 · 发表于 2011-5-13 16:49:54

有点晕啊
我也想学习学习

tianyaxiaoke · 发表于 2011-5-13 17:05:26

第三个问题对我来说是最着急最重要的，恳请各位大侠帮助解答。

shawn2008 · 发表于 2011-5-13 19:59:20

本帖最后由 shawn2008 于 2011-5-13 20:01 编辑

4# tianyaxiaoke

1、不同类型的单元，其假设条件和单元位移模式是有差异的。就S4R单元来说，UMAT或VUMAT子程序中传入的变量和要求返回的变量，它们的形式是与平面应力状态相对应的。不过你可以在UMAT或VUMAT子程序中调用自己开发的三维全应力本构子程序。在调用你自己的三维子程序的前后，需要对变量形式进行转换。

比如，UMAT中对平面应力单元，传入的应力STRESS(NSTEN)只有三个分量，你需要将其扩展为三维状态的应力，

STRESS_3D(1) = STRESS(1)
STRESS_3D(2) = STRESS(2)
STRESS_3D(3) = 0.0
STRESS_3D(4) = 0.0
STRESS_3D(5) = 0.0
STRESS_3D(6) = STRESS(3)

最后，返回给UMAT的时候也要将STRESS_3D()转换成STRESS()

2、我也看晕了，不太明白你说的。不过，法向应变增量是按照 de33 = de33^e + de33^p 计算的，当你更新完平面内的应力之后，可以得到 de11^e，de22^e，de11^p和de22^p。利用广义虎克定律和法向应力为零的假设，可以通过de11^e和de22^e计算得到de33^e；利用塑性变形不可压缩假设，de33^p = -de11^p - de22^p。最后就可以计算得到法向方向的应变增量了，用于描述板料厚度的变化。

tianyaxiaoke · 发表于 2011-5-13 20:20:05

PS 7# shawn2008 感谢您的悉心解答。我大概明白了些，能帮我看看我上传的平面应力子程序中是否正确吗？我是按您上面说的方法编写的，但结果与cae计算的相比，应力应变都差异很大

shawn2008 · 发表于 2011-5-13 20:40:11

我只是在等程序的调试结果中偶尔看一眼帖子，恐怕没时间读代码。

cumtmbqq1 · 发表于 2011-5-13 21:11:42

个人觉得，你的第一个问题的前提不成立
首先，对于壳单元，有限元理论的节点自由度肯定是3个，而平面应力问题只有2个，而且对于壳单元，所针对的都是空间问题，不能简化成为平面应力问题。

tianyaxiaoke · 发表于 2011-5-13 21:44:51

10# cumtmbqq1
那请问，壳单元要调用全应力子程序吗？我尝试把一个用在实体单元上很好的子程序用在壳单元上，但差别与cae计算的也非常大

tianyaxiaoke · 发表于 2011-5-13 21:49:47

10# cumtmbqq1
这是我在参考资料中找到的，好像就是说壳单元要用平面应力子程序把。欢迎回帖讨论

tianyaxiaoke · 发表于 2011-5-14 09:30:59

abaqus中壳单元NDI=3，还是NDI=2？换句话说，就是下面两种方式哪种正确？
1.
StressNew(k,1)=s11；StressNew(k,2)=s22；StressNew(k,3)=s33；StressNew(k,4)=s12；
StrainInc(k,1)=ᇫe11；StrainInc(k,2)=ᇫe22；StrainInc(k,3)=ᇫe33；StrainInc(k,4)=ᇫe12
2.
StressNew(k,1)=s11；StressNew(k,2)=s22；StressNew(k,3)=s12；
StrainInc(k,1)=ᇫe11；StrainInc(k,2)=ᇫe22；StrainInc(k,3=ᇫe12

tianyaxiaoke · 发表于 2011-6-17 16:16:13

呵呵，感觉好寂静呀

tianyaxiaoke · 发表于 2011-6-30 08:30:18

15# shawn2008
验证啦，确实如您所说的，感谢指教！

zc78 · 发表于 2011-7-17 10:31:08

学习了，人多力量大.

tianyaxiaoke · 发表于 2011-9-9 17:48:56

时间过去这么久啦，估计帖子早就沉啦，尽管结果不是很理想，还是感谢上面的热心人！

三言两语 · 发表于 2012-11-14 21:41:58

shawn2008 发表于 2011-5-13 19:59
4# tianyaxiaoke

1、不同类型的单元，其假设条件和单元位移模式是有差异的。就S4R单元来说，UMAT或VUMAT ...

谢谢！

asisherry · 发表于 2013-1-14 09:43:21

tianyaxiaoke 发表于 2011-5-14 09:30
abaqus中壳单元NDI=3，还是NDI=2？换句话说，就是下面两种方式哪种正确？
1.
StressNew(k,1)=s11；StressNe ...

请问，您说的 1,2两种情况，那种是对的呢？abaqus中壳单元NDI=3，还是NDI=2？谢谢赐教

账号		自动登录	找回密码
密码			注册

账号		自动登录	找回密码
密码			立即注册

[子程序] 吐血悬赏--关于平面应力子程序的问题

最佳答案

相关帖子

评分

评分

本帖子中包含更多资源