弹塑性有限元方法

第三章弹塑性有限元方法的实施
§3.1增量平衡方程和切线刚度矩阵
1、分段线性化的求解思想
塑性变形的特点决定了塑性本构关系的非线性和多值性，上面由塑性增量理论给 出了塑性应力一应变关系｛da ｝ = ［q,］｛dg ｝
说明当前应力状态不仅与当前应变有关，而且和达到这一变形状态的路径（加载历史）有 关。

这里包含了屈服准则、强化条件和加卸载准则。

由此，对物理非线性问题，通常采用分段线性化的纯增量法和逐次迭代的方法求解。

即将加载过程分成若干个增量步，选择其中任意一个增量步建立它的增量平衡方程并求 解，对整个过程的求解有普遍意义。

2、增量平衡方程和切线刚度矩阵
设t 时刻（加载至M 步终），结构（单元）在当前载荷（广义体力｛几｝和表面力｛£｝） 的作用下处于平衡状态，此时物体内一点的应力、应变状态为9｝、匹｝。

在此基础上，施 加一个载荷增量｛A/；.｝和｛«｝，即从时刻，则在体内必然引起一个位移增量｛△“｝ 和相应的｛Ao-｝. ｛△$｝,只要佩｝和｛M ｝足够小，就有｛"｝ = ［〈］｛△£｝。

倘若初始状态｛b ｝己知，加载过程己知，贝可以确定（即Jd 硝可以确定，然后 可在硬化曲线上得到材所对应的硬化系数）于是上面的方程成为线性的。

在+ M 这一 增量过程中，应用于虚功原理可得到如下虚功方程：
根据小变形几何关系△"二和 2=B\q,再由虚位移§（△?）的任意性，并设
P + AP= + y 0U+『川（£. +纣»柠，展开后，其中单元在/时刻载荷等效节点 力：P = j N T f v dV + J N T f s dS ；卜内增量载荷的等效力AP = \ N 丁恋dV + ） N 丁嗽dS ° 匕
Se 乙
其中
』［（（7+2）丁 一伉dU
-f(f s + Af s )r SAudS = O (1)
(2)
这样，由方程（1）可得平衡方程：J[B] {（y+\（y}dV = {P + \P}
岭
即：刁+$ = J B‘ bdV + J B l\(ydV-(P + AP) = O
(3)
因为/时刻（第i步终）结构处于平衡状态F,=P-\B T（yclV = O
匕
这样(2)式变为二AP= j B T A(jdV即：Z\F = j B1 AadV-AP = 0(4)
将{db} = 和△“代入上式得增量平衡方程：
J B l D cp BdV \q-NP = \F(5)对增量位移求导：
〃（血）一f B『D Bdv- K；妙）\ cp
(6)于是（5）式成为K：= AP(7)
K；为单元切向刚度矩阵。

集合所有单刚后得到结构总的增量平衡方程
K T ^q = \P(8)方程（8）是线性的，可以直接求解。

3.2硬化系数的数值表示
根据单一曲线定理，对于一般稳定性硬化材料，在其简单加载过程中，歹和g之间存在着一一对应的确定的函数关系歹=0）（科，这一关系可用单向拉伸实验來确定。

例如，对于Mises各向同性硬化材料
H, = da/ds p（8）在有限元分析中，作为初始参数应把这一曲线输入（用函数或数字的形式），在加载过程中弹塑性矩阵不断地修改，根据当前的应力或应变来确定。

目前，硬化曲线的输入格式有两种：
1）解析表达式
根据单一曲线定理，由单向拉伸试验曲线直接得出硬化曲线的解析式。

例如：（a） Mises各向同性线性硬化材料
单向拉伸曲线有：当cr < cr y cr = Es
当(J>a s cr = 6 + (w_ s s )E t (9)
(10)
其中： B = Am 斗t = E 〃G / 6)"i 2）根据离散的单拉实验数据，采用样条插值计算
（参看清华大学孟凡中教材：弹塑性有限变形理论和有限元方法）
3.3过渡单元弹塑性矩阵的确定
附：对于一般材料的硬化曲线的求法（求H'）
如单拉曲线 则硬化曲线
根据 6〜刍 ==》厅〜s p
===》®"=JdH=Jd 磅 其中单拉时等效应变为 匕岸届=彳（1+必
因为 务=£, s = s=-uw,平均应变为 q =扌（耳+6+£3）= *（1-2U ）S 所以/. ^=£ +
£0 ，当 U = *时 5=无
（b ） Mises 各向同性幕硬化材料
单向拉伸曲线有： 当cr ＜（y s a = Es
(11)
由屈服点条件：= A&
(12) 据（8）式得 得 A = Es^m
(18) 1. 三种变形状态
弹塑性变形体中，在一个载荷增量步内可能有三种变形状态：
1） 弹性区：加载前后均处于弹性状态，故采用弹性阵不变。

2） 塑性区：加载前后均处于塑性状态，其弹塑性矩阵［0』由塑性增量理论确定（与 当
前应力水平和塑性变形增量的总量有关）
3） 过渡区：加载前处于弹性状态，加载后进入塑性状态，所以，在这一过程中采用弹 性
矩阵［D ］或最终的［0」都不合适，必须寻找一个合适的弹塑性矩阵［万缈］。

2. 加权平均的弹塑性矩阵［万缈］
1）过渡单元在加载后的应力计算（以单拉状态为例）
在时间步内施加一个增量载荷后，讨论某单元的应力应变状态。

设某单元加（卸）载前的应力状态告，相应的应变各（4点）处于弹性状态（弹
性区间0€）。

加载后，按弹性计算得到应变增量到达B 点。

显然B 点不是实际 的应力状态，因为己经超过了 C 点，进入了塑性变形阶段，假设实际应到达D 点。

该增量步的弹塑性矩阵［万缈］是未知的，它的大小应该和该增量步内弹塑性应变所 占比例有关，只能经过迭代试算得出。

因为:
(14)
(15)
显然,0 </n < 1,且川=0时是全弹性,m=l 时是全塑性。

实际应力增量为
△q = △b ： + A CT /1 = + △of" = 0(1-171)^ +
=[(1- m)D + 〃©,]△§ = 2a
推广到一般的应力状态{A CT }, {△"为
△歹=[(1 一〃?)[£)] +加 D t .r △g = [Dp ［万订=（1一加）［切+ 〃?［2订
(13) 且设 (16) (17)
［万缈］一加权平均弹塑性矩阵；加一比例系数。

2）比例系数加的迭代公式
已知A点的％和6，同时由到达A点的路径确定g：和o■:
(18)
(19)
△刃=\£ = £q +\£-£'x
由定义:
3. 过渡单元加的确定
1） 确定是过渡单元。

即在第（M ）个增量步终（求解结束时）某单元是弹性的应力、
应变状态久和吕，且o-0<<（或进入第j 个增量步（△『内载荷增量△£）, 按弹性计算到达B 点，其应力,应变2。

空,可以确定该单元在第i 个增量 步内是过渡单元。

2） 关于加的迭代过程：
按弹性矩阵［Q ］计算该单元的切线刚度矩阵心，然后和其他单元集合成总刚K —代 入结构的增量平衡方程并求解得总位移向量,从{△§}：中提取该单元的 酬,并求出{△卅,及{对。

代入（19）式，计算出加⑴,再将加⑴ 代入（18）得订⑴，［0』中的［2订与当前应力和应变状态有关。

当前应力为:
{b} = {b}_ + {间⑴,{£} = {£ + {△・：）
（a ）按第1次迭代的计算值［0订⑴代入该单元计算切线刚度阵，并与其他单元集 合组装求解总的增量方程得{△§},）及相应的{△$},）, {A<T }；2）
及
匕}?}二{b}_ + {"}汽此时耳和®没有改变，再代入（19）式计算加叽将〃严和
{b严（若是硬化材料，还要根据当时塑性应变总量确定才的值）代入计算［2］
铁
（b）依次类推，求出［0」⑴；［0订宀.......直至前后两次的加值十分接近（到达给定的允许误差范围）停止迭代。

（C）将迭代终止时的［D.订作为该单元的弹塑性矩阵，求单刚集合，解方程，求出{△?};及相应的{△"：, {△b};将其累加到上一步的终值上作为下一步的初值。

总位移{q} = {q}i + {△“"
k}={d-i + {M-
{b}={bh+ {△□■}：
并记下每个单元的仅}和£：{刃={b}_ + {△・}：,以此作为（汁1）步的初始状态，
继续加载。

3）讨论
上面采用的是最简单的纯增量法，并取其中一个增量步a•步）内对〃7值的迭代，最终确定加权平均的弹塑性矩阵［^订
（1）釆用加权平均的弹塑性矩阵［nJ,在同一增量步内，对过渡单元的加值往往要迭代若干次，每次迭代都要重新计算单元的切线刚度阵，并重新组装总刚和解方程。

显然求解过程比较复杂繁琐，由此增加了许多工作量，但从提高精度和加速收敛两方面是大有好处的。

实践证明：即使加载步长比较大，势必在这一步内新进入屈服的单元（过渡元）比较
多，然而采用对〃7迭代计算出比较准确的［耳」后，仍能获得比较满意的结果。

如果不釆用对过渡单元迭代的办法，则为了保证解的精度，必须控制每个增量步的大
小，以保证每一步内新增加的塑性单元较少，否则将越來越偏离正确解，使求解失真，甚至发散。

（2）迭代收敛准则：加值的迭代从理论上讲要求两次值非常接近时方可结束。

而大量实践说明，一般迭代2~3次就可以，所以往往用迭代次数来控制即可，这是对一个过渡单元而言，而从整个结构来看，还要求在前后两次迭代中不再有新的塑性单元产生來决定是否可进入下一增量步的计算。

§ 3.4采用纯增量法作弹塑性有限元分析的步骤
以下仅限于简单加载过程（无反复加卸载过程）和Mises各向同性强化材料
1.开始，输入初始参数（几何；材料性质，卍，R；边界条件；外载荷P）
2.将外载荷一次加上作线弹性分析｛q｝T｛£｝T｛b｝Tbm3x （Ml.条件）
如果不存在塑性区则为弹性问题T直接输出结果结束！
否则bmax > bf 作弹塑性分析
3.计算弹性极限｛P｝。

设"了心/疔，则｛P｝=^
$ i ）e a
并可输出弹性极限载荷｛P｝e下的结果｛必、｛我、｛a｝e o
4.对剩余载荷｛P｝r = ｛P｝~｛P｝e作弹塑性分析
如果釆用等增量步格式，则将｛P｝,.等分为N个增量步，即每一增量步载荷为：△P =
*Z。

下面5.中是对N个增量步循环。

N
5.在i步上施加一个增量载荷己知当前状态下（人1步终），各单元的（o［高斯点）了，
2, 6。

判断三种类型的单元：1）弹性2）塑性3）过渡单元。

对本增量步内所有过渡单元经过2〜3次迭代得到合适的［0』，计算各单元的匕，并集合所有单元，形成总刚心，求解［K』｛M｝ = AP得｛△§｝,得到第i步的解。

⑷,=｛心+ ｛切和创=（4-1 + ｛塚；= 3｝円+ ｛呵
同时记录下各单元的当前状态。

< &
6.直至全部载荷施加完毕，输出结果，结束
讨论:
1）以上步骤只考虑简单加载情况，变形与加载路径无关，因此没有采用加卸载准则。

2）在△『（△P）增量步内方程是线性的，因为只要当前状态6&己知，［0』是可确定对于每一步［&•］｛△§｝ = ｛AP｝是线性的，但对于整个加载过程［K⑷］｛q｝ = ｛P｝仍是非线性的，这就是分段线性化的思想。