第8章 材料非线性问题的有限元法

T
(7.20)
' x x cp ' xy xy ' y y cp ' yz yz ' z z cp ' zx zx
(7.21)
式中，ζcp＝(ζx+ζy+ζz)／3称为平均应力。等效应力可以用应 力偏量表示为
非线性问题经有限元法离散后，得到如下形式的一组代数方程
Ψ ( ) K ({ }){ } { f } 0
即刚度方程[K]是节点位移向量{δ}的函数。 材料非线性问题是由材料非线性应力应变关系引起的，通常表 现为非线性弹性问题和弹塑性问题，此外还有与时间有关的应 力应变关系。 非线性弹性问题和弹塑性问题的塑性阶段呈现非线性物理性质。 加载过程时，这两类问题的非线性性质是一样的。不同之处在 于两点：一是弹塑性材料有一个从弹性到塑性的折点，二是卸 载过程两者有完全不同的路径。 在常应力状态下，变形随时间变化的特性成为粘性，变形随时 间变化的现象称为徐变（蠕变）。这类问题包括粘弹性问题、 粘弹塑性问题、徐变问题。
1 米赛斯屈服准则 米赛斯屈服准则认为；材科在复杂应力状态下的形状改变 能达到了单向拉伸屈服时的形状改变能，材料开始屈服。于是 得到米赛斯屈服条件是
[( 1 2 ) 2 ( 2 3 ) 2 ( 3 1 ) 2 ] / 2 T
等效应力为
(7.19)
式中，ζ1、ζ2、ζ3为主应力，ζT是单向拉伸时的屈服极限。
第7章 非线性有限元
——材料非线性问题
1. 2. 3. 4. 5. 材料非线性问题的求解方法 塑性应力应变关系 弹塑性矩阵的表达式 弹塑性问题的求解方法 弹塑性问题的实例计算
7.1 材料非线性问题的求解方法
前面各章中，我们所讨论的问题都是线弹性力学问题。在线弹 性力学中，位移与应变的关系（几何方程）是线性的，应变与 应力的关系（本构方程）也是线性的。
引进假想的线性弹性应力
{ }el [ D]{ }
那末
{ 0 } { }－{ }el
(7.12)
把(7.11)式代入(7.1)式，有
T T [ B ] [ D ] [ B ] d V { } { R } [ B ] {ζ 0}dV
T [ K ] [ B ] 设 0 [ D][ B]dV ，即由线弹性 矩阵所定义的结构整体刚度矩阵。于是 上式可以写成
A1B1＝( KT )（ 1 Δ） 2＝R F (1 )
把上式和(e)式比较可以看出，δ2就是位移的第二次近似。 如此不断重复，则得迭代公式如下
( KT ) （ ＝R F ( n )； n1 n （Δ） n Δ） n 1 n 1
因此，图7.2(a) 就是求解 (7.3) 式 的图解表示。由于KT表示结构的切线 刚度，因而牛顿—拉斐逊方法也称为 切线刚度法。 同样，修正的牛顿－拉斐逊方法 可以用图7.2 (b) 表示。由于每次迭代 不改变它的刚度值，所以也称为等刚 度法。
3 2 2 2 2 2 ＝ '2 ' ' 2 ( ' ' ' x y z xy yz zx ) 2
若记：
{ '} { ' x ' y ' z
dY dF ＝ KT d d
式中， KT 是曲线 F ＝ Kδ 的斜率，代表切 线刚度。第二步，从 B1 点作曲线 F ＝ Kδ 的切线交直线F＝R于A2点，取A2点的横 坐标是δ2。从图中看出
A1B1 /( 2 1 )＝( KT )1
由于 得
A1B1＝R F (1 )； 2 1 ( Δ )2
3 初应变法 某些问题，例如蠕变问题，它的应变值是由应力值决定。 因此，物理方程可写成
{ } f ({ })
类似于初应力法，设用具有初应变的线弹性物理方程
(7.15)
{ } [ D （ ] { } { 0 }）
(7.16)
所代替。式中{ε0}是初应变列阵；调整{ε}0的值，使它在给定 的应力值{ζ}用(7.15)式或(7.16)式可以得到同样的应变值{ε}。于 是有
K ({δ}){δ} {R} 0
(7.3)
1 牛顿-拉斐逊（Newton-Raphson）法 任何具有一阶导数的连续函数Y(x)，在xn点作一阶泰勒级数展开， 它在xn点的线性近似公式是
dY Y ( x) Y ( xn ) ( x xn ) dx n
＝ [( 1 2 ) 2 ( 2 3 ) 2 ( 3 1 ) 2 ] / 2
2 2 2 2 ＝ ( x y ) 2 ( y z ) 2 ( z x ) 2 6( xy yz zx ) 2
则米赛斯屈服条件是 应力偏量 为：
{ 0 }＝{ } [ D]1{ } f ({ }) [ D]1{ }
el 1 令 { } ＝[ D] { } ，于是
{ 0 }＝{ } { }el
(7.17)
把(7.16)式代入(7.1)式，有 [ K 0 ]{ } {R} [ B]T [ D]{ 0 }dV
因此，非线性方程Y(x)=0在xn附近的近似方程是线性方程
dY Y ( x) Y ( xn ) ( x xn )＝0 dx n
它的解是
dY Δx Y ( xn ) /( ) n； dx
xn 1 xn＋Δx n 1
这就是牛顿－拉斐逊方法的迭代公式。
牛顿－拉斐逊方法的迭代过程如图7.1(a)所示，它要求在每 次迭代时计算 Y ' ( x) dY ( x) / dx ，因此计算工作量巨大。 修正的牛顿－拉斐逊方法 迭代公式是
1 表征材料应力应变关系的本构方程是线性的； 2 描述应变与位移关系的几何方程是线性的； 3 以变形前的状态建立的平衡方程仍适用于变形后的体系． 但是，工程中的许多问题的位移与应变、应变与应力的关系不 满足上述线性关系，呈非线性状态。通常把不满足条件 1 的称 为材料非线性，把不满足条件2，3的称为几何非线性。
在每次迭代时Y’(x)值是不变的，迭代过程如图7.1(b)所示。
dY Δx Y ( xn ) /( ) 0； dx
xn 1 x＋Δx n 1
牛顿－拉斐逊方法求解平衡方程的迭代过程 结构的平衡方程式为
K ({δ}){δ} {R} 0
(7.4)
简单起见，考虑单自由度系统。设Y(δ)=K(δ)δ－R＝0，因为 K是δ的函数，即K＝K(δ)。令F(δ)＝Kδ，于是
7.2 塑性应力应变关系
1 材料的塑性性质 简单拉伸及薄壁筒扭转实验所得到的应力应变曲线是我们 研究材料塑性性质的基本资料，图7.7(a)是低碳钢的拉伸曲线。 实验知，应力增加到屈服极限时，应力应变曲线上出现屈服阶 段。过了屈服阶段以后，大多数材料要使继续增加变形，必须
使应力进一步增加，即dζ/dε>0; 这所谓强化现象或称加工硬化。 加果屈服阶段很长，我们可以简化为如图7.7 (b) 所示的曲线— —dζ/dε=0；称为理想塑性或完全塑性。
考虑由于增量d{δ}引起{Y}的变化。因为 {R}与{δ}无关，则 d{Y ({ })} [ B]T d{ζ}dV
( [ B]T [ DT ({ })][ B]dV )d{ } [ KT ]d{ } 切线刚度矩阵 T
利用牛顿—拉斐逊方法，得到迭代公式 [ KT ]n { }n 1 {Y }n
T
[ K T ] [ B] [ DT ({ })][ B]dV
(7.8)
[ B] {ζ}n dV {R}； { }n 1 { }n { }n 1
(7.9)
3 初应力法 设材料的物理方程取为 { } f ({ }) (7.10) 即由给定的应变值确定相应的应力值。上式可用具有初应力的 线弹性物理方程
{ }n [ D]({ }n { 0 }n1 )
接着利用应力应变曲线，即是由(7.15)式得出对应于{ζ}n的 应变值。这样，就能算出初应变值
el 1 { 0 }n＝{ } { } { } [ D ] { }n n n n
作为下次迭代之用。把上式代入(7.8)式，求解代数方程就可以求 得位移的第n+1次近似值{δ}n+1。重复迭代直至收敛。
对于材料非线性问题进行有限元分析，由于考虑的是小变形， 平衡方程和几何关系依然成立，即
T [ B ] {ζ}dV {R} {ε} [ B]{δ}
(7.1)
但是物理方程是非线性的，可以写成如下的一般形式
f ({ζ},{ε}) 0
(7.2)
必须注意，由于小变形的关系，应力形式的平衡方程仍然是线 性的，但是以结点位移列阵{δ} 表示的平衡方程则不再是线性的 了。因为应力{ζ}和应变{ε} 之间是非线性的，从而应力{ζ}与位 移{δ}之间也是非线性的；于是(7.1)式可以写成
{ } [ D]{ } { 0 }
(7.11)
所代替。式中{ζ0}是初应力列阵；[D]是线性弹性矩阵，就是非 线性材料在{δ}＝0时的切线弹性矩阵。 调整初应力值{ζ0} ，使它在给定的应变值{ε}用 (7.10) 式或 (7.11)式可以得到相同的应力值{ζ}。有
{ 0 } { }－[ D]{ }＝f ({ }) [ D]{ }
Y ( ) F{δ} {R} 0
用牛顿－拉斐逊方法求非线性方程Y(δ)＝0的根，(7.1)式的 迭代公式可以写成
dY ( ) n Δ n 1 R F ( n ) ； n 1 n＋Δ n 1 d
(e)
我们来叙述牛顿—拉斐逊方法求解(7.3)式的图解表示。 在图7.2中，给出了牛顿—拉斐逊迭代方法。曲线F＝Kδ和 直线F＝R的交点A的横坐标是(7.3)式的精确解。迭代开始时按线 性理论求解位移δ1 作为第一次近似值，即图7.2a 中A1 点的横坐 标。如果载荷R不因变形而改变它的大小和方向，则有
[ K 0 ]{ } {R} [ B]T {ζ 0 }dV
把上式写成如下的迭代公式
(7.13)来自百度文库
[ K 0 ]{ }n 1 {R} {R}n {R}n [ B] ({ζ}n {ζ} )dV
T el n
(7.14)
在整个迭代过程中刚度矩阵 [K0] 保持不变，因此也称为等 刚度法。
可以把上式写成如下迭代公式
[ K 0 ]{ }n 1 {R} {R}n {R}n [ B]T [ D]{ 0 }dV
(7.18)
式中，{R}n称为矫正载荷，它近 似地等于抵消初应变所需要的载 荷。
迭代过程是调整初应变的过程，现把它叙述于下。
如果己知位移的第n次近似值{δ}n ，可以利用(7.2)式算出应 变值{ε}n。由{ε}n以及前一次迭代所得到的初应变{ε}n－1，就可 以利用(7.16)式求出应力值
把(7.3)式写成迭代公式
(7.5)
[ K ]n1{ }n {R}
迭代步骤 如图7.3。
(7.6)
(2) 切线刚度法 设材料的应力应变关系表示为增量形式
d{ } [ DT ({ })]d{ } [ DT ({ })] 是切线弹 就可以利用切线刚度法， 性矩阵。把公式（7.1）改写成 {Y ({ })} [ B]T {ζ}dV {R} 0 (7.7)
2 变刚度法 (1) 割线刚度法 如果材料的应力应变关系能够表示成如 下形式
{ } [ D{ }]{ }
于是由(1.2)式，上式可以写成
{ } [ D({ })][ B]{ } [ D({ })][ B]{ }
把上式代入(7.1)式，并利用(7.3)式，得
[ K ({R})] [ B]T [ D({ })][ B]dV