第2章有限元分析基础

第2章有限元法基础
第1节有限单法的形成
一、有限元法的形成
在工程技术领域内，经常会遇到两类典型的问题。

其中的第一类问题，可以归结为有限个已知单元体的组合。

例如，材料力学中的连续梁、建筑结构框架和桁架结构。

我们把这类问题，称为离散系统。

尽管离散系统是可解的，但是求解这类复杂的离散系统，要依靠计算机技术；第二类问题，通常可以建立它们应遵循的基本方程，即微分方程和相应的边界条件。

例如弹性力学问题、热传导问题和电磁场问题等。

由于建立基本方程所研究的对象通常是无限小的单元，这类问题称为连续系统。

尽管已经建立了连续系统的基本方程，由于边界条件的限制，通常只能得到少数简单问题的精确解答。

对于许多实际的工程问题，还无法给出精确的解答。

为解决这个困难，工程师们和数学家们提出了许多近似方法。

在寻找连续系统求解方法的过程中，工程师和数学家从两个不同的路线得到了相同的结果，即有限元法。

有限元法的形成可以追溯到二十世纪50年代，来源于固体力学中矩阵结构法的发展和工程师对结构相似性的直觉判断。

从固体力学的角度来看，桁架结构等标准离散系统与人为地分割成有限个分区后的连续系统在结构上存在相似性。

1956年M.J.Turner，R.W.Clough，H.C.Martin，L.J.Topp在纽约举行的航空学会年会上介绍了一种新的计算方法，将矩阵位移法推广到求解平面应力问题。

他们把结构划分成一个个三角形和矩形的“单元”，利用单元中近似位移函数，求得单元节点力与节点位移关系的单元刚度矩阵。

1954—1955年，J.H.Argyris在航空工程杂志上发表了一组能量原理和结构分析论文。

1960年，Clough在他的名为“The finite element in plane stress analysis”的论文中首次提出了有限元（finite element）这一术语。

数学家们则发展了微分方程的近似解法，包括有限差分方法、变分原理和加权余量法。

在1963年前后，经过J.F.Besseling，R.J.Melosh，R.E.Jones，R.H.Gallaher，T.H.Pian（卞学磺）等许多人的工作，认识到有限元法就是变分原理中Ritz近似法的一种变形，发展了用各种不同变分原理导出的有限元计算公式。

1965年O.C.Zienkiewicz和Y.K.Cheung（张佑启）发现只要能写成变分形式的所有场问题，都可以用与固体力学有限元法的相同步骤求解。

1969年B.A.Szabo和G.C.Lee指出可以用加权余量法特别是Galerkin法，导出标准的有限元过程来求解非结构问题。

我国的力学工作者为有限元方法的初期发展做出了许多贡献，其中比较著名的有：陈伯屏（结构矩阵方法），钱令希（余能原理），钱伟长（广义变分原理），胡海昌（广义变分原理），冯康（有限单元法理论）。

遗憾的是，从1966年开始的近十年期间，我国的研究工作
受到阻碍。

有限元法不仅能应用于结构分析，还能解决归结为场问题的工程问题，从二十世纪六十年代中期以来，有限元法得到了巨大的发展，为工程设计和优化提供了有力的工具。

二、有限元法的基本思路
有限元法的基本思路可以归结为：假想的把一连续体分割成数目有限的小体（有限元中称为单元），彼此间只在数目有限的指定点（有限元中称为节点）相互连结，组成一个单元的集合体以代替原来的连续体，再在节点上引进等效力以代替实际作用于单元上的外力。

有限元法的实质是：把有无限个自由度的连续体，理想化为只有有限个自由度的单元集合体，使问题简化为适合于数值解法的结构型问题。

选择一个简单的函数来近似地表示位移分量的分布规律，建立位移和节点力之间的关系，这个函数就是第1章讲的单元位移模式。

下面用在自重作用下的等截面直杆来说明有限元法的思路。

图2-1 受自重作用的等截面直杆 图2-2 离散后的直杆
受自重作用的等截面直杆如图2-1所示，杆的长度为L ，截面积为A ，弹性模量为E ，单位长度的重量为q ，杆的内力为N 。

试求：杆的位移分布，杆的应变和应力。

等截面直杆在自重作用下的有限元法解法。

1．离散化
如图2-2所示，将直杆划分成n 个有限段，有限段之间通过一个铰接点连接。

称两段之间的连接点为节点，称每个有限段为单元。

第i 个单元的长度为L i ，包含第i ，i+1个节点。

2．用单元节点位移表示单元内部位移
第i 个单元中的位移用所包含的节点位移来表示，
)()(1i i
i
i i x x L u u u x u --+
=+ （2- 1）
其中u i 为第i 节点的位移，x i 为第i 节点的坐标。

第i 个单元的应变为εi ，应力为σi ，内
力为N i ：
i
i
i i L u u dx du -=
=
+1ε （2- 2）
i
i i i i L u u E E )
(1-=
=+εσ
（2-3）
i
i i i i L u u EA A N )
(1-=
=+σ
（2- 4）
3．把外载荷集中到节点上 如图2-3所示，把第i 单元和第i+1单元重量的一半 q*(L i +L i+1)/2，集中到第i+1节点上。

4．建立节点的力平衡方程
对于第i+1节点，由力的平衡方程可得： 2
)
(11+++=
-i i i i L L q N N
（2- 5令1
+=i i i L L
λ，并将（2- 5）代入得：
221)11(2)1(i i
i i i i i L EA q u u u λλλ+=
-++-++ （2 - 6根据约束条件，01=u 。

对于第n+1个节点，2
n
n qL N =
EA
qL u u n n n 221
=+-+
（2-7）
建立所有节点的力平衡方程，可以得到由n+1个方程构成的方程组，可解出n+1个未知的节点位移。

由此，可由位移法通过几何方程求出各单元的应变，通过物理方程求出各单元的应力。

第2节 有限元法的计算步骤
有限元法的计算步骤归纳为以下四个基本步骤：模型选取、网格划分、单元分析和整体分析。

一、力学模型的选取
根据受力体的特点把它归结为平面问题、平面应变问题、平面应力问题、轴对称问题、空间问题、板、梁、杆或组合体等等问题，并注意利用受力体的对称或反对称性质。

二、单元的选取、结构的离散化（网格划分）
有限元法的基础是用有限个单元体的集合来代替原有的连续体。

因此首先要对弹性体进行必要的简化，再将弹性体划分为有限个单元组成的离散体。

单元之间通过单元节点相连接。

由单元、节点、节点连线构成的集合称为网格。

通常把平面问题划分成三角形或四边形单元的网格，三维实体划分成4面体或6面体单元的网格，如图2-4、图2-5所示。

图2-4 平面问题的三节点三角形单元划分
图2-5 四节点四面体单元及三维实体的四面体单元划分
三、单元分析
对于弹性力学问题，单元分析，就是建立各个单元的节点位移和节点力之间的关系式。

由于将单元的节点位移作为基本变量，进行单元分析首先要为单元内部的位移确定一个近似表达式，然后计算单元的应变、应力，再建立单元中节点力与节点位移的关系式。

以平面问题的三角形3节点单元为例。

如图2-6所示，单元e 有三个节点I 、J 、M ，每个节点有两个位移u 、v 和两个节点力U 、V 。

单元的所有节点位移、节点力，可以表示为节点位移向量（vector ）：
节点位移{
}⎪⎪⎪⎭
⎪⎪⎪⎬⎫
⎪⎪⎪⎩⎪⎪⎪⎨⎧=m m j j i i e
v u v u v u δ
节点力{}
⎪⎪⎪⎭
⎪⎪⎪⎬⎫⎪⎪⎪⎩⎪⎪⎪⎨⎧=m m j j i i e
V U V U V U F 图2-6 三节点三角形单元
单元的节点位移和节点力之间的关系用张量（tensor ）来表示
{}[]{}e e e K F δ=
（2-8）
上式就是表征单元的节点力和节点位移之间关系的刚度方程，[k]e 就是单元刚度矩阵，对平面问题单元刚度矩阵是与单元各节点坐标及材料常数E 和μ有限的常数矩阵。

四、整体分析
对由各个单元组成的整体进行分析，建立节点外载荷与节点位移的关系，以解出节点位移，这个过程为整体分析。

再以弹性力学的平面问题为例，假设弹性体被划分为N 个单元和n 个节点，对每个单元按前述方法进行分析计算，便可得到N 组形如（2-8）式的方程。

将这些方程集合起来，就可得到表征整个弹性体的平衡关系式。

为此，我们先引入整个弹性体的节点位移列阵 {δ}2n ×1，它是由各节点位移按节点号码以从小到大的顺序排列组成，即
I
I
j
j
{}[]T i
i i v u =δ{}[]
T
N
e e i N e e i T
i i i V U Y X R ⎥⎦
⎤⎢⎣⎡==∑∑==11其中子矩阵
(i =1,2, …, n ) 是节点i 的位移分量。

继而再引入整个弹性体的载荷列阵{R }2n ×1 ，它是移置到节点上的等效节点载荷依节点号码从小到大的顺序排列组成，即
其中子矩阵
(i =1,2, …, n )
是节点i 上的等效节点载荷。

各单元的节点力列阵经过这样的扩充之后就可以进行相加，把全部单元的节点力列阵叠加在一起，便可得到弹性体的载荷列阵，即
这是由于相邻单元公共边内力引起的等效节点力，在叠加过程中必然会全部相互抵消，所以只剩下载荷所引起的等效节点力。

同样，将（2-8）式的六阶方阵[k ]加以扩充，使之成为2n 阶的方阵。

考虑到[k ]扩充以后，除了对应的i , j , m 双行和双列上的九个子矩阵之外，其余元素均为零，故（2-8）式中的单元位移列阵{δ}e 2n ×1 便可用整体的位移列阵{δ}2n ×1 来替代。

这样，（2-8）式可改写为
把上式对N 个单元进行求和叠加，得
上式左边就是弹性体所有单元刚度矩阵的总和，称为弹性体的整体刚度矩阵（或简称为
总刚），记为[K ]。

由此，便可得到关于节点位移的所有2n 个线性方程，即
[K ]{δ}={R } (2-10)
在2-10式中，整体刚度矩阵[K ]中的各个元素均可由单元刚度矩阵求出，{R}载荷列向量为各节点等效节点力均可求出，由此通过迭代法或消元法解2n 个线性方程可求出各节点的位移。

然后，由位移法可求出各节点的应变和应力。

第３节 载荷向节点的移置方法
在进行有限元分析时，受力体要受到体力、集中载荷、分布载荷等各种载荷的作用，这些载荷必须向节点移置。

载荷移置的依据是：圣维南原理。

{}[]
9)
-(2 2112T
T
n T
T n δδδδ =⨯{}[]
T
T n T T n R R R R 2112=⨯{}{}[]
T
T
n
T
T
N
e e R R R R R 211
==∑=[]{}{}e n n n n R k 1
21222⨯⨯⨯=δ[]{
}{}∑∑===⎪⎭
⎫
⎝⎛N e e N e R k 11δ
图2-7 重力向节点的转移
一、体力向节点的移置
先求应当移置到节点i 的垂直载荷y i 。

根据圣维南原理，假设节点i 沿y 方向作单位位移，而其余两节点不动，由于三角形单元采用线形位移模式，所以当m ，j 两节点的位移为零时，b 点的位移必然为零。

显然，mb=bj ，bc=1/3*i 。

所以当i 点的位移为1时，c 点的位移为1/3。

按静力等效原则原载荷W 的虚功应当等于y i 的虚功，因此得到y i =1/3*W 。

同理，y j =1/3*W ，y m =1/3*W 这样，就把不作用在节点上的重力载荷转移置到三个节点上去了。

二、集中载荷向节点的移置
如图2-8所示，单元ijm 的ij 边长为L ，在ij 边上沿x 轴方向有一集中载荷p ，其作用点距i ，j 分别为Li 及Lj ，用与上述相同的方法，可将载荷p 移置到i ，j 节点上去，x i =L j /L*p ，x j = L i /L*p 。

三、分布载荷向节点的移置
如图2-9所示，单元ijm 的ij 边长为L ，在ij 边上沿x 轴方向有一分布载荷，j 点为零，i 点为q ，把分布载荷移置到i ，j 节点上去。

先把分布载荷等效为集中载荷，大小为1/2*qL ，其作用点距i ，j 分别为1/3*L 及2/3*L ，用与上述相同的方法，可将分布载荷q 移置到i ，j 节点上去，x i =1/3*q ，x j =1/6*q 。

图2-8 集中载荷的移置
图2-9 分布载荷的移置
i
i
X
X
注意：如果集中载荷不与坐标轴同向，可把该载荷向x 、y 、z 轴分解，然后向各节点移置；如果分布载荷不与坐标轴同向，可把该分布载荷先等效为集中载荷，然后按集中载荷向节点移置的方法向节点移置。

第４节 有限元分析应注意的问题
在对受力体进行有限元分析时，要遵循一定的原则，应注意以下问题。

1．根据分析对象的特点选择合适的单元
常用的平面问题单元有三节点三角形、四节点矩形、六节点三角形、八节点矩形；常用的空间单元有四节点四面体单元、八节点六面体单元等。

2．单元节点编号的原则
单元编号原则是编号从1开始，按每次递增1顺序编号，不能有间断；节点编号原则为编号从1开始，按每次递增1顺序编号，不能有间断。

3．单元划分的疏密对计算精度和计算成本的影响
单元划分的密，则计算精度高，但计算成本高（对计算机硬件的要求高，计算时间长）；反之，单元划分的稀疏，则计算精度低，但计算成本低（对计算机硬件的要求低，计算时间短）。

因此，在划分单元时要根据计算精度的要求适度控制单元划分的疏密。

4．划分网格应注意的问题
边界曲折、应力集中、应力变化大的部位，单元划分应细化，否则，可划分的稀疏一些。

单元由细到疏应当逐步过渡。

如图2-10所示。

a) 复杂形体的网格划分 b) 边界曲折形体的单元划分
5．同一单元边长的要求
对三角形单元，三条边长应尽量接近，不能出现钝角，以免计算结果出现大的误差；对矩形单元，长度和宽度不宜相差过大；对空间单元单元的各条边边长也应尽量接近，否则，边长相差越大，则产生的误差也越大。

6．单元之间的要求
图2-10 网格的划分
任意一个单元的角点必须同时也是相邻单元的角点，而不能是相邻单元边上的内点，划分单元必须遵守此原则。

7．受力体中不同材料、不同厚度的处理
如果计算对象具有不同的厚度或不同的弹性系数，则厚度或弹性系数突变之处应是单元的边线。

8．应力集中和应力突变的处理
应力分布有突变之处，或者受有应力集中载荷处布置节点，其附近的单元也应划分得细些。

9．对称性的利用
应充分利用受力体结构的对称性（几何形状和支承条件对某轴对称，同时截面和材料性质也对称此轴）。

第５节有限元法的应用
有限元法不仅能应用于结构分析，还能解决归结为场问题的工程问题，从二十世纪六十年代中期以来，有限元法得到了巨大的发展，特别是近20年来有限元理论的发展和有限元分析软件的应用为工程设计和优化提供了有力的工具。

一、算法与有限元软件
从二十世纪60年代中期以来，进行了大量的理论研究，不但拓展了有限元法的应用领域，还开发了许多通用或专用的有限元分析软件。

理论研究的一个重要领域是计算方法的研究，主要有：
1．大型线性方程组的解法；
2．非线性问题的解法；
3．动力问题计算方法。

目前应用较多的通用有限元软件如表2-1所列：
表2-1 常用通用有限元软件
在结构分析中，由于ANSYS软件的功能强大，被广泛应用。

ANSYS是世界上著名的大型通用有限元计算软件，它包括热、电、磁、流体和结构等诸多模块，具有强大的求解器和前、后处理功能，为我们解决复杂、庞大的工程项目和致力于高水平的科研攻关提供了一个优良的工作环境，更使我们从繁琐、单调的常规有限元编程中解脱出来。

ANSYS本身不仅具有较为完善的分析功能，同时也为用户自己进行二次开发提供了友好的开发环境。

ANSYS程序自身有着较为强大三维建模能力，仅靠ANSYS的GUI(图形界面)就可建立各种复杂的几何模型；此外，ANSYS还提供较为灵活的图形接口及数据接口。

因而，利用这些功能，可以实现不同分析软件之间的模型转换。

ANSYS/Multiphysics是ANSYS产品的“旗舰”，它包括所有工程学科的所有性能，ANSYS 产品家族如图2-11所示。

ANSYS/Multiphysics有三个主要的组成产品：
● ANSYS/Mechanical-ANSYS/机械-结构及热；
● ANSYS/Emag-ANSYS电磁学；
● ANSYS/FLOTRAN-ANSYS计算流体动力学。

其它产品：
● ANSYS/LS-DYNA -高度非线性结构问题；
● Design Space –CAD环境下，适合快速分析容易使用的设计和分析工具；
● ANSYS/ProFEA –Pro/ENGINEER的ANSYS分析接口。

图2-11 ANSYS产品家族
二、ANSYS软件进行有限元分析的基本方法
有限元分析的基本方法
1．建立实际工程问题的计算模型
（1）利用几何、载荷的对称性简化模型；
（2）建立等效模型。

2．选择适当的分析工具
侧重考虑以下几个方面：
（1）物理场耦合问题；
（2）大变形；
（3）网格重划分。

3．前处理（Preprocessing）
（1）单元属性定义（单元类型、实常数、材料属性）；
（2）建立几何模型（Geometric Modeling，自下而上，或基本单元组合）；
（3）有限单元划分（Meshing）与网格控制。

4．求解（Solution）
（1）施加约束（Constraint）和载荷（Load）；
（2）求解方法选择；
（3）计算参数设定；
（4）求解solve。

5）后处理（Postprocessing）
后处理的目的在于分析计算模型是否合理，提出结论。

可进行以下工作：
（1）用可视化方法查看分析结果（等值线、等值面、色块图）分析计算结果，包括位移、应力、应变、温度等；
（2）最大最小值分析；
（3）检验结果。