第9章 非线性问题的有限单元法

第9章非线性问题的有限单元法9.1 非线性问题概述前面章节讨论的都是线性问题，但在很多实际问题中，线弹性力学中的基本方程已不能满足，需要用非线性有限单元法。

非线性问题的基本特征是变化的结构刚度，它可以分为三大类：材料非线性、几何非线性、状态非线性。

1. 材料非线性(塑性, 超弹性, 蠕变)材料非线性指的是材料的物理定律是非线性的。

它又可分为非线性弹性问题和非线性弹塑性问题两大类。

例如在结构的形状有不连续变化（如缺口、裂纹等）的部位存在应力集中，当外载荷到达一定数值时该部位首先进入塑性，这时在该部位线弹性的应力应变关系不再适用，虽然结构的其他大部分区域仍保持弹性。

2. 几何非线性(大应变, 大挠度, 应力刚化)几何非线性是有结构变形的大位移引起的。

例如钓鱼杆，在轻微的垂向载荷作用下，会产生很大的变形。

随着垂向载荷的增加，杆不断的弯曲，以至于动力臂明显减少，结构刚度增加。

3. 状态非线性(接触, 单元死活)状态非线性是一种与状态相关的非线性行为。

例如，只承受张力的电缆的松弛与张紧；轴承与轴承套的接触与脱开；冻土的冻结与融化。

这些系统的刚度随着它们状态的变化而发生显著变化。

9.2 非线性有限元问题的求解方法对于线性方程组，由于刚度方程是常数矩阵，可以直接求解，但对于非线性方程组，由于刚度方程是某个未知量的函数则不能直接求解。

以下将简要介绍借助于重复求解线性方程组以得到非线性方程组解答的一些常用方法。

1.迭代法迭代法与直接法不同，它不是求方程组的直接解，而是用某一近似值代人，逐步迭代，使近似值逐渐逼近，当达到允许的规定误差时，就取这些近似值为方程组的解。

与直接法相比，迭代法的计算程序较简单，但迭代法耗用的机时较直接法长。

它不必存贮带宽以内的零元素，因此存贮量大大减少，且计算中舍入误差的积累也较小。

以平面问题为例，迭代法的存贮量一般只需直接法的14左右。

在求解非线性方程组时，一般采用迭代法。

2. 牛顿—拉斐逊方法ANSYS程序的方程求解器计算一系列的联立线性方程来预测工程系统的响应。

然而，非线性结构的行为不能直接用这样一系列的线性方程表示。

需要一系列的带校正的线性近似来求解非线性问题。

一种近似的非线性救求解是将载荷分成一系列的载荷增量，即逐步递增载荷和平衡迭代。

它可以在几个载荷步内或者在一个载荷步的几个子步内施加载荷增量。

在每一个增量的求解完成后，继续进行下一个载荷增量之前程序调整刚度矩阵以反映结构刚度的非线性变化。

但是纯粹的增量近似会随着每一个载荷增量积累的误差，导种结果最终失去平衡，如图9-1（a）所示所示。

.（a）纯粹增量式解（b)全牛顿－拉普森迭代求解（2个载荷增量）图9-1 纯粹增量近似与牛顿－拉普森近似的关系。

ANSYS程序通过使用牛顿－拉普森平衡迭代（Newton-Raphson—简称NR）克服了这种困难，它迫使在每一个载荷增量的末端解达到平衡收敛（在某个容限范围内）。

图9-1（b)描述了在单自由度非线性分析中牛顿－拉普森平衡迭代的使用。

在每次求解前，NR方法估算出残差矢量，这个矢量是回复力（对应于单元应力的载荷）和所加载荷的差值。

程序然后使用非平衡载荷进行线性求解，且核查收敛性。

如果不满足收敛准则，重新估算非平衡载荷，修改刚度矩阵，获得新解。

持续这种迭代过程直到问题收敛。

ANSYS程序提供了一系列命令来增强问题的收敛性，如自适应下降（Adaptive descent）、线性搜索（Line search）、自动载荷步（Automatic load stepping）及二分法（Bisection）等，可被激活来加强问题的收敛性，如果不能得到收敛，那么程序或者继续计算下一个载荷前或者终止（根据用户的指示）。

3. 弧长法对某些物理意义上不稳定系统的非线性静态分析，如果你仅仅使用NR方法，正切刚度矩阵可能变为降秩短阵，导致严重的收敛问题。

这样的情况包括独立实体从固定表面分离的静态接触分析，结构或者完全崩溃或者“突然变成”另一个稳定形状的非线性弯曲问题。

对这样的情况，你可以激活另外一种迭代方法：弧长方法（Arc-length method），来帮助稳定求解。

弧长方法导致NR平衡迭代沿一段弧收敛，从而即使当正切刚度矩阵的倾斜为零或负值时，也能阻止发散。

这种迭代方法以图形表示在图9-2中。

图9-2 传统的NR方法与弧长法的比较4. 非线性求解的组织级别非线性求解被分成三个操作级别：载荷步（Load steps）、子步（Substeps or time steps）、平衡迭代（Equilibrium iterations）。

“顶层”级别由在一定“时间”范围内明确定义的载荷步组成，并假定载荷在载荷步内是线性地变化的。

在每一个载荷步内，为了逐步加载可以控制程序来执行多次求解（子步或时间步）。

在每一个子步内，程序将进行一系列的平衡迭代以获得收敛的解。

图9-3说明了一段用于非线性分析的典型的载荷历史。

图9-3 载荷步、子步、及“时间”5. 收敛容限当确定收敛准则时，ANSYS程序会列出一系列的选择：将收敛检查建立在力、力矩、位移、转动或这些项目的任意组合上。

另外，每一个项目可以有不同的收敛容限值。

对多自由度问题，同样也有收敛准则的选择问题。

以力为基础的收敛提供了收敛的绝对量度，而以位移为基础的收敛仅提供了表观收敛的相对量度。

因此，如果需要总是使用以力为基础（或以力矩为基础的）收敛容限。

如果需要可以增加以位移为基础的（或以转动为基础的）收敛检查，但是通常不单独使用它们。

图9-4说明了一种单独使用位移收敛检查导致出错情况。

在第二次迭代后计算出的位移很小可能被认为是收敛的解，尽管问题仍旧远离真正的解。

要防止这样的错误，应当使用力收敛检查。

图9-4 完全依赖位移收敛检查6. 保守行为与非保守行为：过程依赖性如果通过外载输入系统的总能量当载荷移去时复原，那么这个系统是保守的。

如果能量被系统消耗（如由于塑性应变或滑动摩擦），则是非保守的。

如图9-5所示的非守恒系统。

图9-5 非守恒（过程相关的）过程一个保守系统的分析是与过程无关的：通常可以任何顺序和以任何数目的增量加载而不影响最终结果。

相反地，一个非保守系统的分析是过程相关的；必须紧紧跟随系统的实际加载历史，以获得精确的结果。

如果对于给定的载荷范围，可以有多于一个的解是有效的（如在突然转变分析中）这样的分析也可能是过程相关的。

过程相关问题通常要求缓慢加载（也就是，使用许多子步）到最终的载荷值。

7. 子步较多的子步能得到较好的精度，但它是以增加运行时间为代价的，ANSYS提供了两种方法来控制子步数。

❶子步数或时间步长，即用户可以指定实际的子步数，也可以通过指定时间步长来控制子步数。

❷自动时间步长。

ANSYS软件会根据结构的特性和系统的响应，自动检查时间步长。

8. 子步数如果结构在它的整个加载历史过程中显示出高度的非线性特点，而且用户也对结构的行为非常了解，以至于可以保证能得到收敛的解，这样用户就可以自己来确定时间步长的大小，并且对所有的载荷步都可以使用同一时间步。

9. 自动时间步如果结构的行为将从线性变化到非线性，或者用户想要在系统响应的非线性部分变化时间步长，这时可以激活自动时间步长，以便系统根据需要自动调节步长，可以获得精度和代价之间的良好平衡。

但是，如果用户不能保证分析问题能否成功收敛，这时可以使用自动时间分步来激活ANSYS软件的二分法。

二分法提供了一种对收敛失败自动矫正的方法，无论何时只要平衡迭代收敛失败，二分法将把时间步长分成两半，然后从最后收敛的子步自动重启动，如果已二分的时间步长再次收敛失败，二分法将再次分割时间步长然后再重启动，如此循环下去，直到获得收敛或者到达用户指定的最小时间步长。

10. 载荷和位移方向当结构经历大变形时应该考虑到载荷将发生了什么变化。

在许多情况中，无论结构如何变形施加在系统中的载荷保持恒定的方向。

而在另一些情况中，力将改变方向，随着单元方向的改变而变化。

ANSYS程序对这两种情况都可以建模，依赖于所施加的载荷类型。

加速度和集中力将不管单元方向的改变而保持它们最初的方向，表面载荷作用在变形单元表面的法向，且可被用来模拟“跟随”力。

图9-6说明了恒力和跟随力。

注意：在大变形分析中不修正节点坐标系方向。

因此计算出的位移在最初的方向上输出。

图9-6 变形前后载荷方向9.3 材料非线性问题的有限单元法许多与材料有关的参数可以使结构刚度在分析期间改变。

塑性、非线性弹性、超弹性材料、混凝土材料的非线性应力—应变关系。

可以使结构刚度在不同载荷水平下（以及在不同温度下）改变。

蠕变、粘塑性和粘弹性可以引起与时间、温度和应力相关的非线性。

膨胀可以引起作为温度、时间、中子流水平（或其他类似量）函数的应变。

ANSYS程序应可以考虑多种材料非线性特性：1.率不相关塑性指材料中产生的不可恢复的即时应变。

2.率相关塑性也可称之为粘塑性，材料的塑性应变大小将是加载速度与时间的函数。

3.材料的蠕变行为也是率相关的。

产生随时间变化的不可恢复的应变，但蠕变的时间尺度要比率相关塑性大的多。

4.非线性弹性允许材料的非线性应力—应变关系，但应变是可以恢复的。

5.超弹性材料应力—应变关系由一个应变能密度势函数定义，用于模拟橡胶，泡菜类材料，变形是可以恢复的。

6.粘弹性是一种率相关的材+料特性，这种材料应变中包含了弹性应变和粘性应变。

7.混凝土材料具有模拟断裂和压碎的能力。

8.膨胀是指材料在中子流作用下的体积扩大效应。

9.3.1 塑性分析1. 塑性理论简介许多常用的工程材料，在应力水平低于比例极限时，应力—应变关系为线性的。

超过这一极限后，应力—应变关系变为非线性，但却不一定是非弹性的。

以不可恢复的应变为特征的塑性，则在应力超过屈服点开始出现。

由于屈服极限与比例极限相差很小，ANSYS程序在塑性分析中，假设这两点相同，如图9-7所示塑性是一种非保守的（不可逆的），与路径无关的现象。

换句话说，载荷施加的顺序，以及什么时候发生塑性响应，影响最终求解结果。

如果用户预计在分析中会出现塑性响应，则应把载荷处理成一系列的小增量载荷步或时间步，以使模型尽可能符合载荷响应路径。

最大塑性应变是在输出（Jobname.OUT）文件的子步信息中打印的。

在一个子步中，如果执行了大量的平衡迭代，或得到大于15%的塑性应变增量，则塑性将激活自动时间步选项[AUTOTS]{GUI：Main menu | Solution | Sol’n Control：Basic Tab或Main Menu | Solution | Unabridged Menu | Time/Frequenc | Time and Substps }。

如果取了太大的时间步，则程序将二分时间步，并重新求解。