非线性有限元方法及实例分析

《非线性有限元方法及实例分析》由会员分享，可在线阅读，更多相关《非线性有限元方法及实例分析（9页珍藏版）》请在金锄头文库上搜索。

1、http:/ 非线性有限元方法及实例分析非线性有限元方法及实例分析 梁军 河海大学水利水电工程学院，南京（210098） 摘摘 要要： 对在地下工程稳定性分析中常用的非线性方程组的求解方法进行研究， 讨论了非线 性计算的迭代收敛准则，并利用非线性有限元方法分析了一个钢棒单轴拉伸的实例。 关键词关键词：非线性有限元，方程组求解，实例分析 1 引引 言言 有限单元法已成为一种强有力的数值解法来解决工程中遇到的大量问题， 其应用范围从固体到流体，从静力到动力，从力学问题到非力学问题。有限元的线性分析已经设计工具被广泛采用。 但对于绝大多数水利工程中遇到的实际问题如地下洞室等， 将其作为非线性问题加以

2、考虑更符合实际情况。根据产生非线性的原因，非线性问题主要有 3 种类型1： 1材料非线性问题（简称材料非线性或物理非线性） 2几何非线性问题 3接触非线性问题（简称接触非线性或边界非线性） 2 非线性方程组的求解非线性方程组的求解 在非线性力学中，无论是哪一类非线性问题，经过有限元离散后，它们都归结为求解一个非线性代数方程组2： () ()()00021212211=nnnn（1.1） 其中n,21是未知量，n,21是n,21的非线性函数，引用矢量记号 Tn21（1.2） Tn21=（1.3） 上述方程组（1.1）可表示为 ( )0=（1.4） 可以将它改写为 ( )( )( )0=RKRF（

3、1.5） 其中( )K是一个的矩阵，其元素是矢量的函数，nnijkR为已知矢量。在位移有限元中，代表未知的结点位移，( )F是等效结点力，R为等效结点荷载，方程( )0=表示结点平衡方程。 在线弹性有限元中，线性方程组 -1- http:/ 0RK（1.6） 可以毫无困难地求解，但对线性方程组( )0=则不行。一般来说，难以求得其精确解，通常采用数值解法，把非线性问题转化为一系列线性问题。为了使这一系列线性解收敛于非线性解，曾经有过许多方法，但这些解法都有一定的局限性。某解法对某一类非线性问题有效，但对另一类问题可能不合适。因而，根据问题性质正确选用求解方法成为非线性有限元的一个极重要的问题。

4、 常见的求解非线性方程组的数值方法有迭代法、 增量法和混合法345。 2.1. 迭代法迭代法 在每次迭代过程中都施加全部荷载， 但逐步修改位移和应变， 使之满足非线性的应力应变关系。迭代法还分为直接迭代法、Newton-Raphson 方法、修正的 Newton-Raphson 方法和拟 Newton 法。 直接迭代法 对非线性方程组 ( )0RK（1.7） 设其初始的近似解为，由此确定近似的0K矩阵 ( )00KK （1.8） 可得改进的近似解 ()RK101（1.9） 重复这一过程，以第i次近似解求出第1+i次近似解的迭代公式为 ( )iiKK（1.10） ()RKii11+（1.11）

5、直到变得充分小，即近似收敛时，终止迭代。 iii=+1对于单变量问题，这一迭代过程是收敛的，单对于多自由度情况，由于未知量通过矩阵 K耦合，迭代过程可能不收敛。 在迭代过程中，得到的近似解一般不会满足( )0RK即 ( )( )0RKiii（1.12） ( )作为对平衡偏离的一种度量，称为失衡力。 Newton-Raphson 方法 Newton-Raphson 方法是求解非线性方程组 -2- http:/ ( )( )( )0=RKRF（1.13） 的一个著名方法，简称 Newton 法。 设( )为具有一阶导数的连续函数，是方程（1.13）的第i次近似解。若 i ( )( )0=RFiii

6、（1.14） 希望能找到一个更好的、方程（1.13）的近似解为 iii+=+1（1.14） 将（1.15）代入（1.13） ，并在附近按一阶 Tayor 级数展开，则i( )在 i处的线性近似公式为 ii ii +=+1（1.16） 引入记号( )i i Ti TKK =Newton 法的迭代公式可归纳为 ()() (ii Tii TiFRKK=11) （1.17） ii i TFK = =（1.18） iii+=+1（1.19） Newton 法的收敛性是比较好的，但对于某些非线性问题，如理想塑性和塑性软化问题，在迭代过程中可能事奇异或病态的，于是的求逆就会出现困难。为此，可引入一个阻尼因子

7、TKTK，使矩阵或者成为非奇异的，或者使它的病态减弱。 IKii T+修正的 Newton-Raphson 方法 上述两种方法求解非线性方程组时，在迭代过程中的每一步都需要重新计算。如将Newton 法迭代公式中的改用初始矩阵i TKi TK( )00TTKK =， 就成为修正的 Newton-Raphson 法。此时，仅第一步迭代需要完全求解线性方程组，并将三角分解后的存储起来，以后的每一步迭代都采用公式 0 TK()i TiK10=（1.20） 修正的 Newton 法的每一步迭代所用的计算时间较少，但迭代的收敛速度降低。 拟 Newton 法 拟 Newton 法的主要特点是每次迭代后用

8、一个简单的方法修正K，K的修正要满足以下的拟牛顿方程 ()()( )iiiiiK=+111（1.21） -3- http:/ 仿照位移的迭代公式来建立劲度拟矩阵的迭代公式： ()()()1111+=iiiKKK（1.22） 只要由和求出。 ii()1iK综上所述， 迭代法就是用总荷载作用下不平衡的线性解去逼近平衡的非线性解， 迭代的过程就是消除失衡力的过程。对于所采用的迭代方法，Newton 法收敛最快，修正的 Newton法最慢。但对于某个具体问题，往往需要进行数值实验，才能判断哪个方法较好。 2.2增量法增量法 增量法与迭代法的最大不同就在于用线性方法求解非线性方程组时， 对荷载增量进行线

9、性化处理。基本思想是将荷载分成许多小的荷载部分（增量） ，每次施加一个荷载增量。此时，假定方程是线性的，劲度矩阵K为常矩阵，对每级增量求出位移增量，对它累加，就可以得到总位移。增量法分为 Euler 法和修正的 Euler 法。 Euler 法 设R为总荷载，引入荷载因子参数，令RR=，非线性 方程组称为 ()( )( )0，RFRF（1.23） 迭代公式成为 ()mmmTmRRK= 111（1.24） 修正的 Euler 法 将由 Euler 法第级荷载增量求得的mm作为中间结果，记为m，它与前一级结果1m加权平均为 () 11mmm+=（1.25） 式中为加权系数，由m确定可得修正的 Eu

10、ler 法的基本公式 mTK,即 mmTmRK= 1 ,mmm+=1（1.26） 2.3混合法混合法 对于一个实际的非线性问题往往用一种非线性解法是得不到满意得结果的， 由此混合使用增量法和迭代法，则称之为混合法或逐步迭代法。一般在总体上采用 Euler 法，而在同一级荷载增量内，采用迭代法。 3 收敛准则收敛准则 在采用迭代法和混合法求解非线性方程组时， 必须给出迭代的收敛准则， 否则就无法终-4- http:/ 止计算。目前，在非线性有限元计算中常用的迭代收敛准则有位移准则、失衡力准则和能量准则。 3.1位移准则位移准则 i di（1.27） 式中符号表示范数，d为位移收敛容许误差，是事先

11、指定的一个很小的正数，通常取%5%1 . 0d3.2失衡力准则失衡力准则 ( )Rqi（1.28） q是失衡力容许误差。当材料表现处明显的软化时，或材料接近理想塑性时，失衡力的微小变化将引起位移增量的很大误差，此时不能采用失衡力误差。 3.3能量准则能量准则 能量准则同时考虑位移增量和失衡力， 是把每次迭代后的内能增量与初始内能增量相比较。内能增量指失衡力在位移增量上所做的功。能量准则可表示为 ()()()00TeiTi（1.29） e是能量收敛容许误差。 4 非线性分析实例非线性分析实例 4.1问题的提出问题的提出 一直径为 10mm，长为 100mm 的钢棒，将其沿轴向拉伸 20mm，求其

12、最后的变形情况和应力分布。弹性模量，泊松比3200EEX =3 . 0=NUXY。塑性时的应力应变关系如下表 1 表 1 塑性时的应力应变关系 数据点 1 2 3 4 5 6 7 应变 0.002 0.003 0.004 0.005 0.006 0.008 0.01 应力 400 416.3 429 438 445 456.8 465.3 数据点 8 9 10 11 12 13 14 应变 0.015 0.02 0.03 0.05 0.1 0.2 0.3 应力 480.8 491.8 507.4 527.5 555.7 586 603 -5- http:/ 由于模型和载荷都是轴对称的，因此建模

13、时可以只取钢棒轴向截面的 1/4，用平面轴对称单元来实现。 4.2结果分析结果分析 经过计算得出以下成果： 1. 应力应变曲线如图 1 所示。 2. 钢棒单轴拉伸等效应力分布如图 2。 3. 钢棒 X 方向应力分布等值线如图 3 所示。 4. 钢棒 Y 方向应力分布等值线如图 4 所示。 分析钢棒受力后的结果，可以看出： （1）对于钢棒的单轴拉伸，由于模型和载荷都是轴对称的，因此建模时可以只取钢棒轴向截面的 1/4，用平面轴对称单元来实现。而且为了产生颈缩现象，在建模时钢棒端部和中部截面作了 5的误差，使其诱导出颈缩现象。 （2）从钢棒单轴拉伸等效应力图上可以看出，最大主应力出现在受拉一侧的拉

14、伸受力处，其值为 603.4Mpa；最小主应力出现在钢棒的中部，其值为 38.37Mpa。 对于钢棒 X 方向应力分布可以看出，最大和最小应力的位置发生了变化，出现在钢棒中下部，靠近受拉一侧。 图 1 应力应变曲线 -6- http:/ 图 2 钢棒单轴拉伸等效应力分布 图 3 钢棒 X 方向应力分布等值线 -7- http:/ 图 4 钢棒 Y 方向应力分布等值线 5 小结小结 本文主要对在地下工程稳定性分析中常用的非线性方程组的求解方法进行研究，列出了一些基本的求解方法，讨论了非线性计算的迭代收敛准则。通过一个钢棒单轴拉伸的材料非线性分析实例来说明非线性分析的过程。 -8- http:/ 参考文献参考文献 1 任青文.非线性有限单元法（上，下）M.河海大学.2003 2 吕和祥，蒋和洋.非线性有限元M.北京：化学工业出版社.1992，10 3 韦博成，万方焕，朱宏图译.非线性回归分析及应用M.北京：中国统计出版社.1997 4 Shi G.H. & Goodman,R.E., Underground support design using block theory to determine key block bolts req