UMAT/VUMAT入门攻略-Davvife留个Simwe的纪念帖

tangguixiang

谢谢lz  正想深入学习下

敦诚

好吧，支持美女版主，提个问题。vumat中一般会涉及大变形，这个时候应变总量无法通过简单地线性叠加来实现了，请问您一般怎么解决？

Davvife

好吧，支持美女版主，提个问题。vumat中一般会涉及大变形，这个时候应变总量无法通过简单地线性叠加来实现了，请问您一般怎么解决？
敦诚 发表于 2011-7-26 19:01

谢谢敦大的支持啊！我真的还遇到过这个问题啊！嘻嘻嘻
这是一个好问题，很多人遇到过，但是始终没有人详细地解释过，这是为什么呢？嘻嘻。因为这个问题涉及到一个比较重要的概念应力客观率，和本构积分的客观性。小女子不才，在编写梁的vumat子程序中涉及到过这个问题，当时也是一头雾水，后来翻阅了很多文献和书籍总算是整理出来一些思路，这里斗胆地说一下。
首先是应力率的客观性，大家知道流体力学里面一般很少设计什么客观性的，应为流体力学是采用cauchy应力作为变形体的一个控制指标，而cauchy应力就是客观地，也就是不随着坐标轴变化。但是应力率就不是这样了，它不是客观地，举个例子来说，一格杆件，在初始x方向受到100n的拉力，假设中间部分的应力大小为1Mpa，然后我们旋转这个杆件，让它围着z轴发生发生90°的旋转，这样的话原来的x轴变为了y轴，y轴变为了x轴，那么也就是说原来x方向的应力从1mpa变为了0，而y轴的应力从0变为了1Mpa，可是如果我们站在这个杆件上面，就是发现其实这个杆件的变形没有增加，那么应力为什么发生变化了呢？呵呵呵
就是因为应力率不为0，而应力率又不是由于变形产生的，而是由于旋转产生的，所以概括地说来就是旋转导致了应力发生变化，很多学者早就对于这一问题进行了研究，提出了应力客观率的概念，最常用的为Jaumann率和Naghdi率还有Truesdell率（有兴趣可以看看《连续体和结构的非线性有限元》这本书，描述的很详细）。这些客观率可以很好地提出旋转部分导致的应力的变化。
现在说说Vumat，它本身输入的应变增量是已经剔除了旋转部分的应变增量了，但是由于一个杆件时时可能旋转，而总的应变却是全部应变的和（包括了旋转和变形），那么这个时候随意地将应变相加可能就与自己想要的情况，不同，所以一般我会采用将得到的总应力反算应变总量的办法，由于一开始说的应力是客观地，那么使用它显然更保险一些啊！
这些都是小女子的一些经验，可能有存在问题的地方，希望各位大牛指正啊！

jian4dan

楼主你好，我觉得让人快且好的学会umat .等，应该把理论推导公式，和源代码都发表出来，如果理论发表难，最好也不出处写出来。这样大家才好看程序，不然发一堆原代码，我想没几个人看，尤其是对初学者更看不懂。所以提出帮人帮到底，呵呵。是不是有点贪心，我觉得比较合理 9# Davvife

Davvife

呵呵，我觉得您确实挺贪心的！
不过，我能理解。不过我还得请大家共同参与，讨论和交流，来完成这个事情。因为我也是笨丫头一个，对理论部分也不是全然了解，慢慢啃呢
众人拾柴火焰高，我还是希望大家助人也助己，共同出谋划策也好，提供材料也好，一起来做。
如果单凭我一己之力，我怕难当大任，让大家失望啊....。

26# jian4dan

zgw99

http://forum.simwe.com/thread-995894-1-1.html
这个帖子中，对于推导率相关的晶体塑性求解Δγ时的，β的求解对角元为什么没有2倍。不知道LZ能不能解释下，谢谢了

whh1235068

2# 北鹰南飞

这个问题应该用UEL来解决吧，觉得UEL中的“应变梯度理论”也可以解决你的问题吧。
上传一个文件。

jinzhuce

支持一下LZ，希望你早日完成大业！！

kmliangl

求学习教程，资料，。。。

heraculia

支持楼主，顺便和楼主请教个问题，我现在在试着调用vumat来做损伤
材料的强化部分采用线性硬化方式，这个是与http://forum.simwe.com/viewthrea ... umat%2B%CB%F0%C9%CB中用的方法基本相同。而后的损伤也是借助abaqus的ductile damage 部分的原理来写。线性损伤值damage=UPL/UPL0  UPL=特征长度xPEEQ增量 UPL0为设定的固定值   或者damage=（1-exp（-3*（UPL/UPL0）））/（1-exp（-3））等形式。  实际的应力=不考虑损伤的应力X（1-damage）
现在遇到几个问题
1、实际应力的计算可能不准确，因为不考虑损伤的应力无法准确计算，将stressnew转存为statenew后 与 下一分析步的 stressold 不相等  
2、在单向拉伸过程中，我现在已经基本上能实现整个过程，从硬化到损伤到断裂判别，但是换到悬臂弯曲截面模拟时，在第一个单元产生破坏后计算 由于 相邻单元的过大变形而终止，不知如何处理。
相关准则图示和我的结构分析图 随后奉上。

语言表述能力较差，请楼主见谅。期待楼主的指导，谢谢。

heraculia

这个图 是 材料的损伤基本思路 是帮助文档中的图

该图是我悬臂弯曲测试的 情况，右上角的三角位置就是扭曲单元，计算被迫终止。

这个是单向拉伸试验的情况

由于本身不是搞力学的，所以理论功底较差，如果有什么错误希望楼主指出，另外如有需要我可奉上子程序。

lungo999

借这个好贴问各位高手一个问题，不知道大家有没有遇到同样的问题。本人正在写的UMAT里面有3个独立的内变量，虽然是J2理论，但是内变量的演化有点复杂，现阶段已经写了一个初步的代码，并用一个3D单元进行了测试，在单调加载和荷载子步较小的（0.001）循环加载时，结果都蛮好的，但是在荷载子步较大的（0.01）循环加载时，应力-应变图出现抖动，而在多个单元的测试时，发现在较大塑性应变时（0.4）较难收敛，不知道大家是否可以给些建议和意见，非常感谢。

附件为对单个单元加载时的应力应变图，荷载子步分别为=0.001和0.01，在0.01时应力反向加载时出现抖动。



2011.8.22更新问题：

最近还在DEBUGGING。发现对于复杂的材料模型，一致切线模量矩阵是非常重要的，错误的切线模量矩阵不但不会收敛，也可能会造成上述的应力抖动现象。本人现在发现了一致切线模量矩阵的一个bug，现在已经可以算荷载子步为0.01时的单个单元了，且应力-应变曲线平滑无抖动。不过在多个单元测试时发现在大塑性应变下还是有应力抖动，现在发现问题在于牛顿迭代，继续DEBUGGING.

敦诚

32# heraculia

我解释一下，损伤存在的时候你不能说应力是否准确，因为这个东西是靠实验判断的，而不是解析解可以解得的。还有我们在实验中也不是直接测量得到的应力，而是得到应变，我们能得到的一般是位移载荷曲线曲线中有了下降段，我们就认为这个是材料所反映特性，至于损伤位置应力实际为多少，其实谁也不知道。
第二个问题。
弯曲梁的结果一般不要那么做，如果一个结构没有存在应力集中区域，或者说变形局部化带的话（比如有缺口），那么首先破坏的位置就是端部（位移载荷和力载荷施加的部位），或者接触部位。所以你的弯曲模型一个角失效以后，如果那个位置还有力边界的话，就会导致周围部位急速应力集中，立刻出发巨大变形（就好比这个地方有个裂缝存在一样），所以大量单元大变形，计算停止。
至于你的拉伸算例就是因为形状导致中间部位应力集中，所以坏也沿着一条线来实现，而且这个位置没有力的作用，不会像弯曲梁的形式那样，一触即溃。

heraculia

35# 敦诚
谢谢 您的回复，拜读过不少您的帖子，收益颇多。
可能还是有一些我没有说明白的地方。
1、对于stress不准确的部分，是这样的。
比如我计算的应力，我先存为statenew变量，然后让stressnew=statenew，按理说这样应该对结果是没有什么影响的。但是这样计算的结果，和我直接将应力存为stressnew得到的结果不一样，随着总加载时间的增加，差别越来越大。（后期相差达到0.5倍） 我想这是不是与单元之间的协调关系有关？是否在我将stressnew传递到abaqus后台计算后到下一个分析步的stressold的过程中，stress被重新计算过。
2、关于我选择悬臂梁的问题，实属被逼无奈，这个是导师的要求。其实在abaqus自带的材料模型中，我模拟的单向拉伸试验已经能够很好的体现了材料的应力应变关系，（至少我是这么认为的）。
关于您提得几个问题
（1）导师想看到单元一个接一个被破坏所产生的施载反力下降趋势。本来在abaqus自带的材料模型中，也基本上可以实现，但是由于手头资料有限，只有一条单向拉伸曲线，无法准确定义材料压缩的相关参数，这样造成悬臂梁弯曲时，梁的根部底端也会产生压缩破坏，因此子程序中我采用了一个静水压力是否大于零的判定，如果静水压力小于零，即判定压缩状态，仍然材料硬化斜率计算，不产生破坏。这也是我采用子程序来进行计算的主要原因。
（2）为了避免施载点被破坏的可能，我在施加载荷位置绑定了刚体，将位置载荷加到刚体上，实际上也对结构产生了约束，但是也是没有办法的事情，导师非要我做这个东西啊。
（3）在不考虑压缩问题的时候我用abaqus自带的ductile damage损伤模型，做悬臂梁仿真时，在单元产生破坏后也是仍然可以正常计算到全部断裂的，因此我觉得可能还是我的子程序模型有一些问题。

本身自己不是搞力学的，导师非要做一个看似不现实的模型，我也没有办法，只能硬着头皮上。顺便充实一下自己，学习一些力学和有限元方面的知识。

最后，感谢各位版主的辛勤努力，让我在这学到了很多东西。

将我总结和使用的vumat模型奉上，希望各位大侠能帮忙看看有什么问题，这个就是我计算悬臂梁终止的模型，不知道在哪方面能够有所修改。写的有点乱，因为现在还在调试阶段，没有将程序语言整理好，另外里面我加了一些注解。

敦诚

我对于您的问题是这样认为的
1.首先是stressnew的问题，对于小变形来说，基本是stress算到多少，后处理现实的就是多少，如果存在旋转的话，那么可能两者不一致，因为应力=应力+应力增量，这个是简单的公司、、公式，可是在vumat中你计算的只是这个应力增量的一部分，还有一部分是共旋所引起的，这个你可以参考一下应力客观率的介绍，比如Jaumann率，这个时候应力增量有三部分组成第一部份为变形导致的应力增加，后面两部分为旋转导致的，这两个部分应力增量abaqus会在你计算完stressnew以后给你增加上。
2.一般对于使用vumat的还都到载荷作用的区域，我会选择将这一位置设置为弹性体，这是个比较常用的方式，如果你的计算感兴趣区域不再载荷作用位置。这个解释起来很简单，圣维南原理解释得很清楚，还有这个一般不会影响你的计算结果，如果你想问为什么，这个其实可以开个帖子来说，这里不详细说明了。
3.这个问题就是个人观点了，首先对于下降段出现的问题，本身就是控制函数发生了变化，从椭圆问题退化到双曲问题，我是只这个函数在空间上的控制方程的变化，这个时候一般会有震颤（因为存在间断解），但是存在震颤并不是人家abaqus计算的不对，而是本来就应该有震颤，好比流体激波行为，差分出来的解本来就是有震颤的。但是这个震颤可以消除就是通过数值耗散的方法，但是数值耗散其实是一中数值误差，他将原本正确的震荡解变为了平滑的解。但是这样曲线看起来好看多了。至于为什么abaqu的只带模型好使，就是因为人家增加了数值耗散。
这个问题我在这个帖子中回了，你可以看看http://forum.simwe.com/viewthrea ... F0%C9%CB&page=2。

lyk_302

在看umat帮助时候有这么一段话：
Large volume changes with geometric nonlinearity：
....
For rate-form constitutive laws, the exact consistent Jacobian is given by


能不能以线弹性
为例子给出Jacobian的最后形式。
我已经搜索别的帖子，没有发现相关的问题的答案。

谢谢！

lungo999

关于大家所说的UMAT的资料问题，个人建议如下：找几本经典的固体力学本构方程相关的书籍，花1-2个月把所有公式推一遍，然后照着ABAQUS把线性，各向同性强化，随动强化的例题自己照着写一遍，搞懂所有的变量，大概就对UMAT有了一定的了解了。然后再开始自己的模型，上手会顺利很多。《Introduction to computional plasticity》对简单的本构模型进行了推导，但是在复杂模型时不怎么管用，所以建议用来理解概念，《连续体和结构的非线性有限元》里面的公式具有一般性，即使对于复杂模型也可以用，所以非常推荐，虽然只有中文版的可以看。另外Simo等人的书和论文可以好好研究下，对于提高迭代的效率很有帮助。

wenyuan721

顶，学习中。。。

wenyuan721

大侠们，可以给一点儿算钢筋混凝土开裂的建议么！。。。

curing2011

学习下，很感谢！！！