有限单元法概述及其基本概念的分析

2012-08-01孙洪铁

山西建筑 2012年14期

孙洪铁

(中国天辰工程有限公司，天津 300400)

1 有限单元法概述

有限单元法是一种获得工程问题近似解的数值分析方法，它通过离散的有限个单元集合体来代替真实的结构，通过对单元的分析，进而得到整个结构的近似解，在工程设计中，这种近似解已能满足工程的需要，有足够的精确度。有限单元法是最初结构力学中杆系结构的矩阵位移法的一种推广，并应用于二维和三维问题。然而，与杆系结构不同，二维和三维实体没有明显的联结点，这就需要建立许多人为的节点，将连续体离散化为许多任意形状的单元，用此方法，连续体便可用有限个自由度的系统来近似表示，并求得其近似解。

2 有限元分析过程

从求解未知量角度来看，有限元法可分为应力场，位移场及混合场三种场变量，本文以应用较为普遍的位移场为例来阐述有限元的分析过程，其主要分析步骤为:连续体的离散化，选择位移函数，有限元模型的建立，集合离散化单元形成系统方程组，求解节点位移及通过位移计算单元的应变和应力。

2.1 连续体的离散化

离散化主要目标是将物体分成充分小的单元，使得简单的位移模型即能足够近似的表示精确解，在这个过程中，将人为通过网格划分出有限个单元，单元与单元之间以节点相连。以二维问题为例，节点的设置要遵循以下几点要求:集中载荷处，分布载荷突变处，几何形状不连续处，材料性质突变处，几何形状的凹角处等，单元的形状主要有三角形，矩形等，三角形要求不要出现钝角;矩形尽量不要采用长宽比较大的矩形，否则都会影响近似解的精度。

2.2 选择位移函数

有限元法的基本原理是分块近似，即把所要研究的区域，分割成有限个子区域，然后假设一些比较简单的函数来表示每个子区域中的解，这些假定的函数被称为位移函数，位移场或位移模式。

一般多项式被作为位移函数，如二维位移函数一般为:

不仅因为其在微分与积分计算中比较容易，而且任一阶的多项式可以近似表示真实解，我们知道无限多次多项式与精确解相对应，我们在不同的阶次将其截断，就得到不同程度的近似解(如图1a)～图1c)所示);其次有限元法作为一种数值计算，所选择的位移函数一定要保证其收敛或趋向于问题的精确解，必须要满足以下三个条件:1)位移函数在单元内必须连续，并且相邻单元的位移必须协调。单元内的连续性，通过选择具备连续性的多项式即能满足，相邻单元之间的位移协调，要求单元之间不能开裂、重叠，这要求单元边界的位移仅与该边界上各节点的位移有关。2)位移函数必须包含单元的刚体位移，即单元一定要产生但不引起应力的那部分位移，如平面应力单元会在平面内任意方向产生均匀移动和转动而不引起应变。在多项式位移函数中，常数项提供单元的刚体位移。3)位移函数必须包含单元的常应变状态，从物理上我们可以理解常应变状态的必要性，设想表示结构集合体的单元越来越多，则在极限情况下，每个单元趋近于一个非常小的尺寸，单元内的应变接近为常量，所以位移函数中必须包含常应变，才能使计算收敛于正确解，多项式一次项提供了单元的常应变。以上1)条称为有限单元的协调性;2)，3)条称做有限单元的完备性，如单元既完备又协调，则收敛是单调的，即以某种模来度量的分析结果的精度会随着单元数目的增加而不断提高，如果只具备完备但不协调，则分析结果在极限时也可能收敛于“精确”结果，但一般收敛是不单调的，但在实际工程中，非协调元有时比与它密切相关的协调单元要好，原因在于采用的近似解的性质，我们都知道，利用假定的位移函数得到的近似结构比实际结构要更刚一些，但由于允许单元分离，重叠使这种近似结构变柔了，这两种影响相互抵消，常常得到好的结果，比如非协调元在板单元中的运用。选择位移函数的另一个因素，即该函数与局部坐标系的方位无关，即在任何一组相对于单元来说的方向固定的荷载作用下，单元的反应(指在和单元一起移动的坐标系中的单元应变能或应变)不依赖于单元本身及它的荷载在全局坐标xy中的方向，单元的这种性质叫做单元的几何各向同性或几何不变性，对于一般的多项式位移函数，按帕斯卡三角形(如图1d)所示)对称选取即可。

图1 位移函数

2.3 有限元模型的建立

作为位移场中有限单元的建立，我们所要求的是连接有限个单元体的节点的位移与节点力之间的转化关系，即确定单元的刚度矩阵。有限元矩阵的建立一般分广义坐标有限模型和等参有限元模型，前者在2.2节中位移函数式(1)我们已经看到，其中α1，α2，α3，…，αm被称为广义坐标，通过广义坐标将单元内任一点的位移与单元节点位移相联系，并且通过一个非奇异矩阵来求得以节点位移及节点坐标表示的单元内任意一点的位移。而等参有限元模型建立的基本思想是:直接通过使用插值函数(或称形函数)来得到单元内任意一点的位移和单元节点位移之间的关系。在实际分析中，使用等参有限元有时更为有效，可以避免在广义坐标计算中出现奇异矩阵的可能性，也减少了矩阵运算的次数;同时二次或高阶等参单元既可以模拟直边也可以模拟曲边，在模拟曲线边界的结构时非常有用。然而广义坐标有限元模型的建立，能更好的让我们加深对有限元法的理解，因此本文以广义坐标法为例来介绍有限单元模型的建立，即求解单元刚度矩阵的过程。从推导方法来看，分为三类:1)直接法，该方法易于理解，适用于简单的问题。2)变分法，把有限元归结为求泛函的极值问题，比如固体力学中的最小势能原理的应用，它使有限元法建立在更加坚实的数学基础上，扩大了有限元法的应用范围。3)加权余数法，不需要利用泛函的概念，而直接从基本的微分方程出发，求出近似解，对于不存在泛函的工程领域，提供了有效的解决方法。

以下我们将通过比较直观易于理解的直接法来了解一下有限元模型的建立过程，以弹性力学平面应力问题为例(见图2):首先假定线性多项式为位移函数，使广义坐标α个数等于单元节点位移个数:

可以写成:

可简写为:

此为单元内部点位移与广义坐标之间的转化关系，为了求出节点位移与单元内部点位移的关系，需要求出节点位移{δ}(e)与{α}的转换关系。

三角形单元的节点位移可以表示为:

{δ}(e)= ［{δ1}，{δ2}，{δ3}］T=［υ1v1υ2v2υ3v3］T。

三角形单元的节点的力向量可以表示为:

{F}(e)= ［{F1}，{F2}，{F3}］T=［u1υ1u2υ2u3υ3］T。

我们要找到{δ}(e)与{F}(e)的转换关系，即{F}(e)=［k］(e)×{δ}(e)，其中，［k］(e)为单元的刚度矩阵。以节点的水平分量为例代入式(2)，可以写成:

节点1:υ1=α1+α2x1+α3y1;节点2:υ2=α1+α2x2+α3y2;节点3:υ3=α1+α2x3+α3y3。

由上面三个方程，易得 α1，α2，α3;同理，通过节点竖直分力，可得α4，α5，α6，最终可以得到一个广义坐标与节点位移的关系式:

其中，［A］为一个只与三角形单元节点坐标有关的矩阵。将式(5)代入式(4)，得:

至此单元内部任意点位移已由节点坐标及节点位移表示。接下来可以求解单元的应变与位移的关系，在弹性力学平面问题中由:

求得:

转换矩阵［B］是一个仅仅与三角形单元的几何性质有关的常量，被称为几何矩阵，由上式我们确定单元位移场的应变状态。

为引入单元材料性质的影响，需要通过应力与应变的本构方程，求解单元应力:

转换矩阵［D］被称为弹性矩阵，对平面应力问题:

接下来，由单元应力推算节点力，需要用到平衡方程，在有限元法中通常用虚功原理来代替平衡方程。在平面应力问题中虚功原理表达式为:

如图3所示中，令实际受力状态在虚设位移状态上作虚功，由{ε*}=［B］{δ*}(e)得:

代入式(10)中，由于{δ*}为任意的，整理后得:

因在常应变三角形单元中［B］，{σ}均为常量，则:

结合式(7)及式(8)得:

其中，［k］(e)为单元的刚度矩阵。至此，我们已经完成了单元分析，为单元的集成提供了必要的条件。

图2 平面应力三角形单元节点位移及节点力

图3 平面应力三角形单元节点力虚设位移

2.4 有限单元的集成及形成系统方程组

单元刚度集成的方法基于各有限单元体连接节点处的协调性要求，即有限单元在相互连接的节点处，必须具有相同的位移，这里的位移为广义位移，它可以是位移、转角。出于这种考虑，我们采用直接刚度法建立整体刚度矩阵［K］，因为单元刚度矩阵［k］(e)与整体刚度矩阵［K］的阶数并不相同，需要将单元刚度矩阵扩大成与整体刚度矩阵同阶数，扩大后的矩阵叫做单元的贡献矩阵，它们表示每个单元单独变形时对整体刚度矩阵提供的贡献。然后根据局部编码和整体编码的对应关系，进行单元刚度矩阵的叠加，节点荷载也采用同样的方式进行叠加。在有限元法中，通常在整体刚度矩阵［K］形成后，引入已知的边界条件，而后形成可求解的有限元系统方程组:{F}=［K］{δ}，其中，［K］为整体刚度矩阵;{F}为结构的节点力;{δ}为结构节点位移。

2.5 结构位移及单元应力的求解

有限元系统方程组一旦建立，我们可以根据适用于计算机处理的方程组的解法，如高斯消元法等，求解出结构各有限单元节点处的位移，然后按式(7)及式(8)即可算出有限单元任意点处应力及应变分量，然而我们必须指出一点，单元的应力并不能保证单个单元的平衡条件，即在单元的节点处内应力和所加载荷是不一定平衡的，我们通常用应力平均值来代表单元的应力，最常用的平均法是采用单元形心处的应力。