非线性有限元分析

非线性有限元分析1 概述在科学技术领域，对于许多力学问题和物理问题，人们已经得到了它们所应遵循的基本方程（常微分方程或偏微分方程）和相应的定解条件（边界条件）。

但能够用解析方法求出精确解的只是少数方程性质比较简单，并且几何形状相当规则的问题。

对于大多数工程实际问题，由于方程的某些特征的非线性性质，或由于求解区域的几何形状比较复杂，则不能得到解析的答案。

这类问题的解决通常有两种途径。

一是引入简化假设，将方程和几何边界简化为能够处理的情况，从而得到问题在简化状态下的解答。

但是这种方法只是在有限的情况下是可行的，因为过多的简化可能导致误差很大甚至是错误的解答。

因此人们多年来一直在致力于寻找和发展另一种求解途径和方法——数值解法。

特别是五十多年来，随着电子计算机的飞速发展和广泛应用，数值分析方法已成为求解科学技术问题的主要工具。

已经发展的数值分析方法可以分为两大类。

一类以有限差分法为代表，主要特点是直接求解基本方程和相应定解条件的近似解。

其具体解法是将求解区域划分为网格，然后在网格的结点上用差分方程来近似微分方程，当采用较多结点时，近似解的精度可以得到改善。

但是当用于求解几何形状复杂的问题时，有限差分法的精度将降低，甚至发生困难。

另一类数值分析方法是首先建立和原问题基本方程及相应定解条件相等效的积分提法，然后再建立近似解法并求解。

如果原问题的方程具有某些特定的性质，则它的等效积分提法可以归结为某个泛函的变分，相应的近似解法实际上就是求解泛函的驻值问题。

诸如里兹法，配点法，最小二乘法，伽辽金法，力矩法等都属于这一类方法。

但此类方法也只能局限于几何形状规则的问题，原因在于它们都是在整个求解区域上假设近似函数，因此，对于几何形状复杂的问题，不可能建立合乎要求的近似函数。

1960年，R.W.CLOUGH发表了有限单元法的第一篇文献“The Finite Element Method in Plane Stress Analysis”，这同时也标志着有限单元法（FEM）的问世。

有限单元法的基本思想是将连续的求解区域离散为一组有限个，且按一定方式相互联接在一起的单元的组合体。

由于单元能按不同的联结方式进行组合，且单元本身又可以有不同形状，因此可以模型化几何形状复杂的求解域。

并且可以利用在每一个单元假设的近似函数来分片地表示全求解域上待求的未知场函数，从而使一个连续的无限自由度问题变成离散的有限自由度问题。

现已证明，有限单元法是基于变分原理的里兹法的另一种形式，从而使里兹法分析的所有理论基础都适用于有限单元法，确认了有限单元法是处理连续介质问题的一种普遍方法。

利用变分原理建立有限元方程和经典里兹法的主要区别是有限单元法假设的近似函数不是在全求解域而是在单元上规定的，而且事先不要求满足任何边界条件，因此可以用来处理很复杂的连续介质问题。

在短短四十余年的时间里，有限单元的分析方法已经迅速地发展为适合于使用各种类型计算机解决复杂工程问题的一种相当普及的方法。

如今，有限元广泛地应用于各个学科门类，已经成为工程师和科研人员用于解决实际工程问题，进行科学研究不可或缺的有力工具。

有限单元法的应用围已由弹性力学平面问题扩展到空间问题，板壳问题，由静力平衡问题扩展到稳定问题，动力问题和波动问题。

分析的对象从弹性材料扩展到塑性，粘弹性，粘塑性和复合材料等，从固体力学扩展到流体力学，传热学等连续介质力学领域。

在工程分析中的作用已从分析和校核扩展到优化设计并和计算机辅助设计技术相结合。

各种各样商业化的大型通用有限元软件层出不穷，不断推出新。

可以预见，随着现代力学，计算数学，计算机技术等学科的发展，有限单元法作为一个具有巩固理论基础和广泛应用围的数值分析工具，必将得到进一步的完善和发展。

2 非线性问题的类型和求解特点2．1 非线性问题的类型2. 1. 1 线性分析的含义在有限元分析中的线性假设包含下列含义：即结点位移为无限小量，材料为线弹性，加载时边界条件的性质保持不变。

于是，静力平衡方程可以表示为：[]{}{}R U K = (2.1)其中，[]K 为刚度矩阵，{}R 为荷载矢量。

由于[]K 和{}R 的元素为常数，故位移响应{}U 是荷载矢量{}R 的线性函数。

也就是说，如果{}R 变为{}R α，则{}U 变为{}U α，其中，α为常数。

这就是所谓的线性有限元分析。

如果上述假设中的任何一条不能得到满足，那么就属于非线性有限元分析。

2. 1. 2 非线性分析的必要性结构力学问题，从本质上讲都是非线性的，线性假设只是实际工程问题的一种简化。

当然，任何实际工程问题的求解都避免不了适当地简化，简化是否合理主要应根据求解效果和实际经验来判断。

对于目前工程实际中的很多问题，如地震作用下结构的弹塑性动力响应，高层建筑抗风，大跨度网壳结构动力稳定性，索膜结构找形荷载与裁减分析，大型桥梁风致振动等问题的研究，仅仅假设为线性问题是很不够的，常常需要进一步考虑为非线性问题。

因此，对各种工程结构的非线性分析就是必不可少且日趋重要了。

对于结构力学的非线性问题来说，有限单元法是最为有效的数值分析方法。

2. 1. 3 非线性问题的类型通常，把非线性问题分为两大类，即分为几何非线性和材料非线性。

但从建立基本方程和程序设计的方便出发，又可分为三种类型：1．材料非线性：非线性效应仅由应力应变关系的非线性引起，位移分量仍假设为无限小量，故仍可采用工程应力和工程应变来描述，即仅材料为非线性。

非线性的应力应变关系是结构非线性的常见原因，许多因素都可以影响材料的应力应变性质，包括加载历史（如在弹塑性响应状况下），环境状况（如温度），加载的时间总量（如在蠕变响应状况下）等。

2．几何非线性：如果结构经受大变形，则变化了的几何形状可能会引起结构的非线性响应，这又可以分为两种情形：第一种情形，大位移小应变。

只是物体经历了大的刚体平动和转动，固连于物体坐标系中的应变分量仍假设为无限小。

此时的应力应变关系则根据实际材料和实际问题可以是线性的也可以是非线性的。

第二种情形，大位移大应变。

也即最一般的的情况，此时结构的平动位移，转动位移和应变都不再是无限小量，本构关系也是非线性的。

3．状态非线性：除以上两种非线性问题之外，还有一种非线性问题，即由于系统刚度和边界条件的性质随物体的运动发生变化所引起的非线性响应。

例如，一根只能拉伸的钢索可能是松散的,也可能是绷紧的；轴承套可能是接触的,也可能是不接触的； 冻土可能是冻结的,也可能是融化的。

这些系统的刚度和边界条件由于系统状态的改变在不同的值之间突然变化。

状态改变也许和载荷直接有关，也可能由某种外部原因引起。

最为典型的就是接触问题，接触是状态非线性类型中一个特殊而重要的子集。

通常情况下，状态非线性问题可以在上述材料非线性和几何非线性类型中的每一种同时出现，从而使得问题的分析变得更为复杂。

2．2 非线性问题的求解特点2. 2. 1 非线性分析的基本问题非线性分析的基本问题是求出在当前荷载作用下的平衡状态。

如果作用的荷载被描述成时间的函数，则物体有限元离散系统的平衡方程可以表示为：{}{}0=-F R t t (2.2)其中，矢量{}R t由t 时刻外荷载的结点力分量所构成，而矢量{}F t 则表示t 时刻的单元应力所引起的结点力分量。

平衡方程(2.2)应针对t 时刻的几何位形建立，并应计入所有的非线性效应。

如果是动力分析，矢量{}R t中还应当包括惯性力和阻尼力。

在求解非线性问题时，(2.2)式应在全部加载历史中成立。

变量t 的引入并不意味着一定是动力问题。

在静力分析中，t 不具有真实“时间”的含义，它的不同取值只是表示相应于不同位形的不同的荷载水平。

但是，在动力分析或具有时间效应的静力分析中，变量t 就有了它本来的“时间”的含义。

2. 2. 2 非线性方程组的增量逐步解法对于许多工程结构，我们所关心的常常是在特定的荷载水平下，或相应的时间物体中的应力和变形。

实际问题根据其解法可以分为两大类型。

第一类问题无需计算中间变形过程，可直接求解在给定荷载下的平衡位形。

但是，如果问题的几何性质或材料性质与路径相关或与时间相关，即该问题依赖于变形历史，则中间变形过程的计算是不可缺少的，这就是第二类问题。

从本质上来说，非线性问题是第二类问题。

此时，往往采用增量分析的方法。

增量逐步解法的基本思想是：假定t 时刻的解为已知，要求t ＋Δt 时刻的解，其中，Δt 是适当选择的时间增量。

在t ＋Δt 时刻，式(2.2)写成为：{}{}0=-∆+∆+F R t t t t (2.3)这里，左上标表示为t ＋Δt 时刻的量。

由于t 时刻的解为已知，因此，可以写为： {}{}{}F F F t tt +=∆+ (2.4) 式中，{}F 表示t 到t ＋Δt 时间间隔，由于单元应力增量所引起的结点力增量矢量。

这一矢量可以近似表示为：{}[]{}U K F t ≈ (2.5)式中，[]K t为相应于t 时刻材料和几何条件的切线刚度矩阵。

{}U 为Δt 时间间隔中的结点位移增量，现在它还是未知的。

将式(2.4)和(2.5)代入式(2.3)中，得到：[]{}{}{}F R U K t t t t -=∆+ (2.6)上式中只有位移增量{}U 为未知，一旦解出，即可算得t ＋Δt 时刻的位移： {}{}{}U U U t t t +=∆+ (2.7)根据{}U t t ∆+，就容易算出t ＋Δt 时刻的应力及{}F t t ∆+，{}K t t ∆+，于是马上可以着手下一步的计算。

但要指出的是，式(2.5)是一个近似表达式，因此t ＋Δt 时刻的解也是近似的，如果急于求成的作下去，最终结果可能出现不可忽视的重大误差以致于达到荒谬的地步。

解决这一困难的办法是以花费计算时间为代价，即在t 到t ＋Δt 时步中进行足够次数的迭代，以保证最终的解获得足够的精度。

2. 2. 3 引入修正Newton －Raphson 迭代格式的增量逐步解法现在更多采用的方法是在每一个荷载增量步中，使用Newton －Raphson 迭代法或修正的Newton －Raphson 迭代法。

由于后者不需要每次迭代时都计算切线刚度矩阵，因此在实际中具有更广泛的应用。

现对该方法做简单的介绍。

在t 时刻到t ＋Δt 时刻的时步中，修正Newton －Raphson 法的迭代公式可以表示为：[]{}(){}{}()1-∆+∆+-=∆i t t t t i t F R U K (2.8){}(){}(){}()i i t t i t t U U U ∆+=-∆+∆+1 (2.9)其中，i 表示迭代步数，依次取1，2，3，…，其迭代所用的初始值正是t 时刻的解，即： {}(){}{}(){}F F U U t t t t t t ==∆+∆+00, (2.10)式(2.8)的右端项：{}{}()1-∆+∆+-i t t t t F R 称为第i 步迭代前的不平衡荷载。