ABAQUS 反复加载卸载应力应变计算说明1、塑性模型采用 Armstrong and Frederick model (AF 模型)[1](1) 屈服准则(Mises屈服)3 ~ ~ ~ ~屈服函数为 f = 2(~ _<~) ： (~ _(~) - Y 2式中：S为应力偏量，为总背应力，Y为屈服极限(2) 流动准则(the associated flow rule)~p =迟d<~式中：~P为塑性应变对时间的微分，血为待定量[2], &为应力张量( 3) 硬化准则<1> 除)=3 h(i)~p - 匚(i)(~(i) p式中：pi)为背应力分量对时间的微分；h(i)，: (i)为材料常数为已知量，P为等效塑性应 变对时间的微分<2> Y = Y + 兰 r (i) (1 一 e-Q ⑴ p )0i=1式中：Y为对应于等效塑性应变p的屈服极限，Y为初始屈服极限为已知量，r(i)为材料常0数为已知量2、反复加载卸载应力应变计算过程说明假设在受载前，物体的初始应力、应变以及背应力均为零2.1 加载过程计算<1>外力不足以使得物体中的任何一点的Mises应力值大于屈服极限此时：塑性应变、背应力均保持为零，屈服极限保持不变。
应力由& = De : ~e计算， 总应变值等于弹性应变<2>外力使得物体中的任何一点的Mises应力值大于屈服极限为了说明ABAQUS是如何确定应变增量A~ ，有必要对ABAQUS求解材料非线性问n+1题进行简单介绍[3]ABAQUS首先将载荷分为若干个微小增量，如图1所示当结构收到一个微小增量AP 时，ABAQUS用与初始结构位移相对应的初始刚度矩阵%和载荷增量AP计算出结构在这一增量后的位移修正c、修正后的位移值U和相应的新的刚度矩阵KABAQUS用新的 a a a刚度矩阵计算结构的内力I，载荷P和I的差值为迭代的残余力R如果R在模型内的a a a a每一个自由度上的值都为零或小于一个给定的容差，如图1所示的a点，则结构处于平衡状 态即ABAQUS计算到的内外力是平衡的若假设整个物体只有一个单元，则位移修正c就为应变增量A~ ，刚度矩阵就是a n + 1UMAT 程序中的雅可比矩阵 根据塑性变形时的雅可比矩阵计算公式：Dp = De — 4G2L-1 : /其中De为弹性状态下的雅可比矩阵，可知Dp较De小UMAT计d算出的应力便是结构的内力Ia在 UMAT 程序中通过计算：3 ~ ~ ~ ~Ftri = (~tri —(~ )： (~tri —(~ ) — Y2n+1 2 n+1 n n+1 n式中：C~ tri =C~ + De : Af~n+1 n n+1并判断若Ft”是否大于零来决定节点是否达到屈服。
若Ftri小于零则节点未达到屈服,n +1 n + 1若 Ftri 大于等于零节点达到屈服，需要进行塑性迭代与重新计算雅可比矩阵n + 1图 1 增量法迭代原理ABAQUS 处理塑性加载过程问题可以描述为：首先读入上一状态的应力、应变、背应力及屈服强度，使用下标n表示（对应时间为t ）再根据上一步的刚度矩阵K和载荷增 nn量AP计算出给定的应变增量A~ ，然后在UMAT程序中求出t （ t二t +At ）时 n +1 n +1 n +1 n n +1刻满足所用塑性模型给出的屈服准则、流动准则、硬化准则下的应力、应变、背应力及屈服 强度UMAT程序利用已知的分，分⑴，~，~p，P，Y和给定的A~首先求出可满足 n n n n n n n+1所用塑性模型给出的屈服准则、流动准则、硬化准则的Ap ，再使用公式依次求出其余量n + 1等效塑性应变： P = P + APn +1 n n+1屈服极限：Y 二 Y +n +1 0i=1n+1塑性流动：nn+1(~s tri — 艺e (i)((i))n +1 n +1 nI i 1.(S tri —艺e (i) ((i)): (S 艺eI n+1 n+1 n n+1 n+1 n' i=1 i=1式中：e(i)n+11 — (i) Apn + 1〜 厅塑性应变增量：A〜p = Ap〜n+1 * 2 n+1 n+1~ ~ 2 ~背应力：((i) =e (i) (( (i) + h (i) A P )n +1 n + 1 n 3 n + 1+ ICJ trimC2应力：C~ =(~(i)+ ：Y nn+1 n+1 \l 3 n+1 n+1式中：C tri = ](C tri +C tri +C tri )m 3 11 22 33再计算出： E~ ， E~ Pnn加载过程按照此方式进行计算同样已知： C~ ， (~(i)， E~ ， E~p ，n n n n n n2.2 当外力卸载时P , Y并且给定：A~ 。
这里考虑卸载过程的第n+1一步迭代由于外力开始卸载，计算得到的(~tri —( )：(~tri —( ) < (S —( )：(S —()n +1 n n +1 n n n n n由于上一步迭代材料处于屈服状态，所以在UMAT程序运行结束时有下面关系成立：3 ~ ~ ~ ~ 3 ~ ~ ~ ~F = — (s —(~ ):(s —( ) — Y2 = 0 所以此时的 F =— (~tri —(~ )：(stri —(~ ) — Y 2 < 02 n n n n n n + 1 2 n + 1 n n + 1 n n根据程序判断不进行塑性计算，只计算弹性雅可比矩阵即按照弹性方式进行卸载各个参数有关系： ((i) =((i)， E p =E p ， p = p ， Y =Yn +1 n n +1 n n +1 n n +1 n弹性应变：Ee =Ee +Aen+1 n n +1总应变： E = E e + E Pn n +1 n +1按照此方式，直到卸载完成这里不讨论反向加载2.3 再次加载时由于（~⑴，~p , Y在卸载过程中不变化，所以在第二次加载过程中，读入的& , &⑴，n n n n n~~p，p，Y中（~⑴，~p 一定是不为零的，并且Y > Y。
对于&是否为零，与残余n n n n n n n 0 n应变是否引起了残余应力有关对于一个单元残余塑性应变是无法引起残余应力的；对于多 个单元构成的物体，当卸载完成后，一些单元的残余塑性应变使得它们与周围的单元之间产生作用力从而产生残余应力此时读入的是不为零的这是与第一次加载初始应力状态 n有区别的根据UMAT弹性试应力的计算公式&衍=<~ + De : A~ ，若物体中含有残余n+1 n n+1应力则在计算过程中直接代入程序计算~ 3 ~ ~ ~ ~由给定的应变增量As 计算出下一时刻的Ftrl = （stri -（~ ）：（s加） — Y 2若 n + 1 n+1 2 n+1 n n+1 nFtri小于零则节点未达到屈服按照2.1中的弹性加载过程进行计算;若Ftri大于等于零按照n +1 n + 12.1 中的塑性加载过程进行计算分析 若要此模型出现图2所示的应力应变曲线必须进行双向加载在单向反复加载时，若所受最大值是不变的则，在后续的加载过程中均为弹性变形图 2 应力应变曲线问题：我这样的分析合理吗？[1] 2002-implementation of cyclic plasticity models based on a general form of kinematic hardening[2] 有限单元法 王勛成 清华大学[3] 2005-ABAQUS用户材料子程序开发及应用（硕士论文-华中科技大学）。