第6章 非线性有限元法(几何非线性)

1、第六章第六章 非线性有限元法（几何非线性）非线性有限元法（几何非线性） 1、变形体的运动描述 x3 x1 x2 P0 t0=0 tn+1=tn+tn tn Pn Pn+1 A0 An+1 An 变形体上的质点的运动状态变形体上的质点的运动状态 可以随不同的坐标选取以下几可以随不同的坐标选取以下几 种描述方法：种描述方法： 1 1、全拉格朗日列式法全拉格朗日列式法（T.LT.L列式列式 法法Total Lagrangian Formulation)Total Lagrangian Formulation)： 选取选取t t0 0=0=0时刻未变形物体的构时刻未变形物体的构 形形A A0 0作为参

2、照构形进行分析。作为参照构形进行分析。 2 2、修正拉格朗日列式法修正拉格朗日列式法（U.LU.L列式法列式法Updated Lagrangian FormulationUpdated Lagrangian Formulation）： 选取选取t tn n时刻的物体构形时刻的物体构形A An n作为参照构形。由于作为参照构形。由于A An n随计算而变化，因随计算而变化，因 此其构形和坐标值也是变化的，即与此其构形和坐标值也是变化的，即与t t有关。有关。t tn n为非线性增量求解时增为非线性增量求解时增 量步的开始时刻。量步的开始时刻。 3 3、欧拉描述法欧拉描述法(Eulerian Fo

3、rmulation)(Eulerian Formulation)： 独立变量是质点当前时刻的位置独立变量是质点当前时刻的位置x xn+1 n+1与时间 与时间t tn+1 n+1。 。 几何非线性的有限元方程一几何非线性的有限元方程一 般采用般采用T.LT.L或或U.LU.L列式法建立！列式法建立！ 、变形梯度张量 x3 x1 x2 PP P P 初始初始/ /未变形未变形 变形后变形后 位移位移u u x x x x 1 1、首先采用、首先采用LagrangianLagrangian方法，方法， 将一个物体的加载过程划分为将一个物体的加载过程划分为 一系列平衡状态。一系列平衡状态。 iii

4、uxx位移方程位移方程 初始状态与变形后状态之间坐初始状态与变形后状态之间坐 标关系为：标关系为： 2 2、然后，考虑材料方向矢量，这个矢量、然后，考虑材料方向矢量，这个矢量 描述物体内一段无限小的单元。描述物体内一段无限小的单元。 jijj j i i xdFxd x x dx x3 x1 x2 i xP ii dxxP ii xdxQ i xP 式中，式中，F Fij ij称为变形梯度张量。称为变形梯度张量。 j i ij x x F 初始状态与变形后状态之间材料方向矢量初始状态与变形后状态之间材料方向矢量 的关系：的关系： 、变形梯度张量 j i j i j i ij x u x x x

5、 x F jjiiijij FFFFI 2 1 2 由位移方程，得：由位移方程，得： j i ijij x u F ijijii FFI 1 由二阶张量特性，变形梯度张量由二阶张量特性，变形梯度张量 的三个不变量为：的三个不变量为： VJdVdFdV ij det JFI ij det 3 由于由于F Fij ij表示从初始状态到变表示从初始状态到变 形后状态的一个映射，其逆映射形后状态的一个映射，其逆映射 Fij-1一定存在，即：一定存在，即： AdFJNdAn ijij 1 或写为：或写为： 体积映射体积映射: : 面积映射：面积映射： 变形前面积变形前面积dA dA Ni ( (初始面积

6、法向矢量初始面积法向矢量) ) 变形后面积变形后面积dAdA ni(变形后面积法向矢量变形后面积法向矢量) ) 映射映射F Fij ij 逆映射逆映射F F-1 -1ij ij F Fij ij是一个二阶张量。是一个二阶张量。 j i ij j i ij x u x x F 1 、应变与变形测度 由于变形梯度张量由于变形梯度张量F Fij ij中包含了刚体运动，因此不能直接用于定中包含了刚体运动，因此不能直接用于定 义应变测度。而材料方向矢量则不包含刚体运动，因此它的平方值义应变测度。而材料方向矢量则不包含刚体运动，因此它的平方值 可以作为衡量从某一状态到变形后状态的一个测度，定义为：可以作为

7、衡量从某一状态到变形后状态的一个测度，定义为： ii xdxdsd 2 初始状态初始状态: : 一个应变测度应该能反映出材料一段一个应变测度应该能反映出材料一段 长度发生的改变。因此，应变张量可以由长度发生的改变。因此，应变张量可以由 下式定义：下式定义： iiii xdxddxdxsdds 22 x3 x1 x2 i xP ii dxxP ii xdxQ i xP iidx dxds 2 变形后状态：变形后状态： 提醒：提醒：由于由于GreenGreen应变张量表达式中的变形梯度张量对应于初始状应变张量表达式中的变形梯度张量对应于初始状 态，因此该应变张量也应在初始状态下计算。态，因此该应变

8、张量也应在初始状态下计算。 、应变与变形测度 、AlmanshiAlmanshi应变张量应变张量1 1、Green Green 应变张量应变张量 GreenGreen应变张量采用应变张量采用LagrangianLagrangian运运 动描述方法，即按初始状态下的动描述方法，即按初始状态下的 构形定义应变张量。构形定义应变张量。 iiijiiijkjki iijkjkii iiii xdxdexdxdFF xdxdxdFFxd xdxddxdxsdds 2 22 ijkjkiij FFe 2 1 式中，式中，e eij ij称为称为GreenGreen应变张量应变张量或或 Green-Lagr

9、angianGreen-Lagrangian应变张量应变张量。 AlmanshiAlmanshi应变张量采用应变张量采用EularEular运动运动 描述方法，即按当前状态下的构描述方法，即按当前状态下的构 形定义应变张量。形定义应变张量。 iiijiikjkiij jkjkiiii iiii dxdxEdxdxFF dxFFdxdxdx xdxddxdxsdds 2 11 11 22 11 2 1 kjkiijij FFE 式中，式中，E Eij ij称为称为AlmanshiAlmanshi应变张量应变张量 或或Almanshi Almanshi EularEular应变张量应变张量。 由于

10、大变形问题有由于大变形问题有 限元方程主要采用限元方程主要采用 T.LT.L列式法列式法或或U.LU.L列式列式 法法建立，因此应在初建立，因此应在初 始状态下定义应变张始状态下定义应变张 量，即采用量，即采用GreenGreen应应 变张量。变张量。 可以证明可以证明GreenGreen应变张量和应变张量和AlmanshiAlmanshi应变张量都是二阶对称张量。应变张量都是二阶对称张量。 、应变与变形测度 2 2、Green Green Lagrangian Lagrangian应变张量应变张量e eij ij与小应变张量与小应变张量ij ij的关系的关系 将变形梯度张量表达式代入到将变形

11、梯度张量表达式代入到 GreenGreen应变张量公式中，得：应变张量公式中，得： ijij j k i k i j j i ij j k i k i j j i ij ij j k kj i k kiij x u x u x u x u x u x u x u x u x u x u e 2 1 2 1 2 1 2 1 i j j i ij x u x u 2 1 式中：式中： 为小变形应变张量；为小变形应变张量； kjkiij FFC j k i k ij x u x u 2 1 2 2、GreenGreen变形张量也可写为：变形张量也可写为： 为非线性二次项为非线性二次项 1 1、Gre

12、enGreen应变张量应变张量 为小应变张量与一个非线性二为小应变张量与一个非线性二 次项之和，这意味所有大变形次项之和，这意味所有大变形 分析都是非线性的。分析都是非线性的。 ijijij Ce 2 1 ijijij e 式中，式中，C Cij ij是是CauchyCauchy变形张量变形张量 由于由于CauchyCauchy变形张量是正定对称变形张量是正定对称 阵，因此该张量有三个实特征值；阵，因此该张量有三个实特征值； 这些特征值的平方根记为材料的这些特征值的平方根记为材料的 主轴拉伸。主轴拉伸。 、大变形的应力测度 1 1、柯西应力张量、柯西应力张量(Cauchys stress (C

13、auchys stress tensor)tensor) 取三维空间笛卡儿坐标系，在取三维空间笛卡儿坐标系，在t t时刻时刻 的现时构形中截取一个四面体素，斜面的现时构形中截取一个四面体素，斜面 的法线为的法线为n n，另外三个面元与所取坐标，另外三个面元与所取坐标 面平行。由四面体素的平衡条件得出其面平行。由四面体素的平衡条件得出其 上的应力为：上的应力为： jiji nn n i i n 3 x 2 x 1 x 这里这里ij=ji便是便是柯西应力张量柯西应力张量，它是二阶对称张量。，它是二阶对称张量。 、柯西、柯西(Cauchy)(Cauchy)应力张量是一种采用欧拉描述法应力张量是一种采

14、用欧拉描述法( (是以质点的瞬时是以质点的瞬时 坐标坐标x xk k和时间和时间t t作为自变量描述作为自变量描述) )定义在定义在t t时刻的现时构形上的应力张时刻的现时构形上的应力张 量量ij ij，又称，又称欧拉应力张量欧拉应力张量。 、在大变形、在大变形( (有限变形有限变形) )情况下，由于变形前的初始构形和变形后情况下，由于变形前的初始构形和变形后 的现时构形差别较大，柯西的现时构形差别较大，柯西(Cauchy)(Cauchy)应力张量难于适应。应力张量难于适应。 柯西应力是定义在现柯西应力是定义在现 时构形（变形后状态时构形（变形后状态 下）的单位面积上的下）的单位面积上的 力，

15、是与变形相关的力，是与变形相关的 真实应力。真实应力。 3、大变形的应力测度 2 2、一阶、一阶Piola-KirchoffPiola-Kirchoff应力张量应力张量 一阶一阶Piola-KirchoffPiola-Kirchoff应力张量的定义是建立在总应力张量的定义是建立在总 力相等的基础上。即：在参考状态下该应力张量力相等的基础上。即：在参考状态下该应力张量 能给出与变形后状态下柯西应力张量相同的力。能给出与变形后状态下柯西应力张量相同的力。 变形后状态下：变形后状态下：dAndP jiji 称为一阶称为一阶Piola-KirchoffPiola-Kirchoff应力张量应力张量或或名

16、义应力名义应力 参考后状态下：参考后状态下： AdNTdP jiji 变形前面积变形前面积dA dA Ni ( (参考面积法向矢量参考面积法向矢量) ) 变形后面积变形后面积dAdA ni(变形后面积法向矢量变形后面积法向矢量) ) AdFJNdAn ijij 1 将面积映射关系：将面积映射关系： 代入上式，得：代入上式，得： 1 jkikij JFT AdNTdAn jijjij AdNTAdFJN jijkjkij 1 同样，柯西应力张量也可以由一同样，柯西应力张量也可以由一 阶阶Piola-KirchoffPiola-Kirchoff应力张量表示：应力张量表示：ikjkij TFJ 1

17、从该式可以看出，一阶从该式可以看出，一阶Piola-KirchoffPiola-Kirchoff应力应力 张量张量提供了以参考状态表示实际力的形式。但提供了以参考状态表示实际力的形式。但 是，直接应用一阶是，直接应用一阶Piola-KirchoffPiola-Kirchoff应力张量应力张量可能可能 存在以下两个困难：存在以下两个困难： 1 1、从能量角度上，、从能量角度上，T Tij ij不适合与不适合与GreenGreen应变张量应变张量 共同使用。因为共同使用。因为T Tij ij乘以乘以GreenGreen应变张量不会产应变张量不会产 生与生与CauchyCauchy应力张量与小应变张

18、量相同的能量应力张量与小应变张量相同的能量 密度。密度。 2 2、T Tij ij不对称，因而较难应用到有限元分析中。不对称，因而较难应用到有限元分析中。 、大变形的应力测度 3 3、二阶、二阶Piola-KirchoffPiola-Kirchoff应力张量应力张量 如不采用变形后状态如不采用变形后状态dPdP推导应力张量，而推导应力张量，而 是将作用在变形后状态下的是将作用在变形后状态下的dPdP映射到未变形映射到未变形 状态上（映射是采用逆变形梯度张量），即：状态上（映射是采用逆变形梯度张量），即： jiji dPFPd 1 AdNSPd jiji 这样可以定义另一个应力张量这样可以定义另

19、一个应力张量S S，它给出了，它给出了 未变形状态下作用在未变形面积上的总力：未变形状态下作用在未变形面积上的总力： 现在，变换柯西应力张量，使：现在，变换柯西应力张量，使： dAnFdPFPd kjkijjiji 11 AdNJFFPd rrkjkiji 11 AdFJNdAn ijij 1 将面积映射关系将面积映射关系 代入上式：代入上式： ( 1 )( 1 ) ( 2 )( 2 ) ( 3 )( 3 ) ( 4 )( 4 ) 对比对比(2)(2)、(4)(4)式可得：式可得： 11 rkjkijij JFFS rkijijjk FJSF 1 S Sij ij称为称为二阶二阶Piola-K

20、irchoffPiola-Kirchoff应力张量应力张量 或或伪应力伪应力 同样，由上式可得：同样，由上式可得： 二阶二阶Piola-KirchoffPiola-Kirchoff应力张量应力张量S Sij ij 的性质：的性质： S Sij ij是对称阵；是对称阵； S Sij ij在能量角度下与在能量角度下与GreenGreen应变张应变张 量协调，即：量协调，即： 该表达式的优点在于等式右边是该表达式的优点在于等式右边是 在参考状态下计算的。在参考状态下计算的。 1. 1. S Sij ij与与T Tij ij有以下关系：有以下关系： ijijijij eS rjijri TSF 二阶二

21、阶Piola-KirchoffPiola-Kirchoff应力张量应力张量 的物理意义是明确的：真实的物理意义是明确的：真实 的力元可以看成是由的力元可以看成是由S Sij ij定义定义 的力元经与变形相同的方式的力元经与变形相同的方式 被被“拉长和转动拉长和转动”后得到的。后得到的。 、大变形的应力测度 4 4、三个应力张量的比较、三个应力张量的比较 张量张量 作用力作用力 作用面积作用面积 柯西应力张量柯西应力张量ij ij 变形后状态下的力 变形后状态下的力 变形后状态下的面积变形后状态下的面积 一阶一阶P-KP-K应力张量变形后状态下的力未变形状态下的面积应力张量变形后状态下的力未变形

22、状态下的面积 二阶二阶P-KP-K应力张量未变形状态下的力未变形状态下的面积应力张量未变形状态下的力未变形状态下的面积 因此，虽然二阶因此，虽然二阶P-KP-K应力张量有其应用上的优点，但其本身的物理应力张量有其应用上的优点，但其本身的物理 意义很难理解。它主要是起到求解大变形问题的桥梁作用，通过意义很难理解。它主要是起到求解大变形问题的桥梁作用，通过 它计算出柯西应力张量。它计算出柯西应力张量。 、几何非线性有限元方程的建立 如前所述，几何非线性的有限元方程一般采用如前所述，几何非线性的有限元方程一般采用T.LT.L或或U.LU.L列式法建立：列式法建立： 1 1、全拉格朗日列式法全拉格朗日

23、列式法（T.LT.L列式法列式法) )： 选取选取t t0 0=0=0时刻未变形物体的构形时刻未变形物体的构形A A0 0作为参照构形进行分析。作为参照构形进行分析。 2 2、修正拉格朗日列式法修正拉格朗日列式法（U.LU.L列式法）列式法）： 选取选取t tn n时刻的物体构形时刻的物体构形A An n作为参照构形。由于作为参照构形。由于A An n随计算而变化，因随计算而变化，因 此其构形和坐标值也是变化的，即与此其构形和坐标值也是变化的，即与t t有关。有关。t tn n为非线性增量求解时增为非线性增量求解时增 量步的开始时刻。即增量分析。量步的开始时刻。即增量分析。 x3 x1 x2

24、P0 t0=0 tn+1=tn+tn tn Pn Pn+1 A0 An+1 An 图示物体同时作用有体积力图示物体同时作用有体积力fib和和 面力面力fiS，在时刻，在时刻t tn+1 n+1=t =tn n+ +t tn n的平的平 衡方程可以按虚功原理建立：衡方程可以按虚功原理建立： S * i S i V * i b i V ijij dSufdVufdV * 提醒：提醒：该方程此时不可解，因为应力该方程此时不可解，因为应力 和应变在变形后状态下表示未知。和应变在变形后状态下表示未知。 、几何非线性有限元方程的建立 2 2、在外力作用点和方向都不改、在外力作用点和方向都不改 变的条件下，

25、也可以将体积力变的条件下，也可以将体积力fib 和面力和面力fiS定义到初始状态下：定义到初始状态下： dSfSdfdSfSdf b i b i S i S i S * i S i V * i b i V ijij dSufdVufdV * S * i S i V * i b i V ijij SdufVdufVdeS * 提醒：提醒：上式给出的虚功方程是从上式给出的虚功方程是从变形后状态下变形后状态下的虚功方程转换的虚功方程转换 而来，因此是准确的，但是已经完全定义在初始状态下了。而来，因此是准确的，但是已经完全定义在初始状态下了。 为了求解，需将以上变形后状态下表示的虚功方程转换到为了求解

26、，需将以上变形后状态下表示的虚功方程转换到 初始状态下表达。初始状态下表达。 V ijij V ijij VdeSdV * 1 1、采用二阶、采用二阶PiolaPiola应力张量和应力张量和 GreenGreen应变张量将虚应变能转换应变张量将虚应变能转换 到初始状态下表示：到初始状态下表示： 将以上关系代入到虚功方程中：将以上关系代入到虚功方程中： 得：得： ( a )( a ) 、几何非线性有限元方程的建立 ij R 2 1 表示该张量对应的时刻：表示该张量对应的时刻：1 1代表初始代表初始 状态时刻，状态时刻，2 2为变形后状态时刻；为变形后状态时刻； 如该标识缺省，则表示从初始状态如该

27、标识缺省，则表示从初始状态 变化到变形后状态该张量的增量。变化到变形后状态该张量的增量。 代表定义该张量所对应的构形：代表定义该张量所对应的构形： 1 1为初始状态构形，为初始状态构形，2 2为变形后为变形后 状态构形；如该标识缺省，则状态构形；如该标识缺省，则 为初始状态构形。为初始状态构形。 在利用增量法（在利用增量法（修正拉格朗日列式法）修正拉格朗日列式法）求解时，为了分析的方求解时，为了分析的方 便，在张量符号的左侧引入上下标，分别该张量对应时刻以及定义便，在张量符号的左侧引入上下标，分别该张量对应时刻以及定义 该张量的构形：该张量的构形： S * i S i V * i b i V

28、ijij SdufVdufVdeS 2 1 2 1 *2 1 2 1 当引入以上表示后，当引入以上表示后， 按按t t1 1+ +t t时刻构形建立的虚功方程可以写为：时刻构形建立的虚功方程可以写为： 或写为：或写为： QVdeS V ijij 2 1 *2 1 2 1 式中，式中， 表示外力所做的虚功。表示外力所做的虚功。 S * i S i V * i b i SdufVdufQ 2 1 2 1 2 1 、几何非线性有限元方程的建立 引入此前引入此前GreenGreen应变张量表达式，可得：应变张量表达式，可得： ijijijijijij ee i j j i ij x u x u 2 1

29、 线线性性项项为为： ijijij SSS 1 1 1 2 1 j k i k j k i k ij x u x u x u x u 2 1 非线性项为：非线性项为： 再将变形后状态下再将变形后状态下KirchoffKirchoff应力张量表示为未变应力张量表示为未变 形状态的形状态的KirchoffKirchoff应力张量加上一个应力增量：应力张量加上一个应力增量： ( a )( a ) ( b )( b ) ij S 1 1 ijij S 11 1 注意，式注意，式 ( b )( b )中中 为作用在未变形构形上并以未变形状态下表示为作用在未变形构形上并以未变形状态下表示 的的Kircho

30、ffKirchoff应力张量，实际上就是柯西应力张量：应力张量，实际上就是柯西应力张量： 。 虚功方程：虚功方程： QVdeS V ijij 2 1 *2 1 2 1 ijijij SS 1 12 1 ( c )( c ) 、几何非线性有限元方程的建立 S * i S i V * i b i SdufVduf 1 1 1 1 VdQ V * ijij 1 11 1 为为t tn n时刻初始构形上时刻初始构形上 外力所做虚功。外力所做虚功。 QVdeSVdQ VdSVdSVdVdVdeS V ijij V ijij V ijij V ijij V ijij V ijij V ijijijij 2

31、 1 * 11 * 1 11 1 * 11 * 11 * 1 1* 1 1* 1 * 11 1 将以上将以上 ( a )( a )、( c )( c )两式代入到虚功方程中，可得：两式代入到虚功方程中，可得： QQVdVdeS V ijij V ijij 1 1 2 1 * 1 1* 11 即变形后状态下的虚功方程为：即变形后状态下的虚功方程为： 式中：式中： 为为t tn n+ +t t时刻初始构形上外力所做的虚功。时刻初始构形上外力所做的虚功。 S * i S i V * i b i SdufVdufQ 2 1 2 1 2 1 这里，虚功方程中由于包含了非线性二这里，虚功方程中由于包含了非线性二 次项，因此方程是非线性方程。这个方次项，因此方程是非线性方程。这个方 程还不能直接求解。为了求解这个方程，程还不能直接求解。为了求解这个方程， 需要将方程线性化。需要将方程线性化。 6、非线性平衡增量方