有限元基本理论

凌俊华

对于许多力学和物理问题，人们可以给出他们的数学模型，即应遵循的控制方程和相应的定解条件，由于方程的复杂性和非线性性质，和求解域的复杂性，直接求解方程得到精确解是不可能的，因此使用数值分析方法求解近似解称为途径。
已经发展的微分方程的数值分析方法分为两类，一类是有限差分法，步骤：划网格、在结点上用差分方程来近似微分方程，网格越密，精度愈高。直接求解基本方程和相应定解条件的近似解。一类是有限元它是从与其等效的积分形式出发。等效积分的形式是加权余量法、他适用于普遍的方程形式。如果原问题的方程的泛函存在，他的等效积分的伽辽金法可以归结为某个泛函的变分。也就是求解泛函的驻值问题（里兹法）。
然而有限元法区别于传统的加权余量法和求解泛函驻值的变分法，该法不是在整个求解域上假设近似函数，而是在各个单元上分片假设近似函数。整体离散，分片近似的思想。
前一部分内容主要是微分方程和相应的定解条件的求解方法。
1、        对于任意的函数，可以转化为等效积分形式，并且进行分部积分可以得到另一种积分形式。等效积分的弱形式。适当提高对任意函数的连续性要求，以降低对微分方程场函数u的连续性要求所建立的等效积分形式称为微分方程的等效积分的弱形式。定解问题的微分方程称为强形式。强形式要求函数U必须可微，且可微的次数至少等于微分方程的次数，原始微分方程往往对解提出了过分平滑的要求，而弱形式较原始方程有更逼近的解。
2、        未知场函数U可以使用近似函数来表示。U=Na，A为待定参数，N为试探函数（基函数、形函数）的已知函数。近似解通常选择使之满足强制边界条件和连续性要求。显然近似解不能精确地满足微分方程和全部的边界条件。它们产生残差R及R拔。通过选择待定系数A，强迫余量在某种平均意义上等于0。余量的加权积分为0就得到一组求解方程，用以求解近似解的待定系数A。从而得到原问题的解答。
采用使余量的加权积分为0来求得微分方程的近似解的方法称为加权余量法。Weighted Residual method。权函数的选择方法。二阶精度高。伽辽金法得到的求解方程的系数矩阵是对称的，并且当微分方程存在相应的泛函时，伽辽金法和变分法往往具有一致的结果。
3、        微分算子和变分算子，泛函的构造，原问题微分方程和边界条件的等效积分的伽辽金提法等效于它的变分原理。即原问题的微分方程和边界条件等效于泛函的变分等于0，亦即泛函取驻值。泛函的极值性对判断解的近似性质有意义，利用它可以对解的上下界作出估计。由变分原理求近似解的标准过程-里兹方法。具体步骤是：未知函数的近似解仍由一族带有待定参数的试探函数来近似表示。。即U=na.将U代入泛函得到用近似函数表达的泛函表达式。泛函的变分为0相当于将泛函对所包含的待定参数进行全微分，并令所得的方程等于0。得到一组与待定系数的个数相等的方程组。这种求解方法叫做里兹法。如果泛函中u和它的导数的最高方次为二次，则称泛函为二次泛函。大量的工程和物理问题中的泛函。对于二次泛函，泛函的变分为0，可以将方程简化为一组线性方程。方程的系数矩阵为对称矩阵。
4、        由于里兹法以变分原理为基础，其收敛性有严格的理论基础，得到的求解方程的系数矩阵是对称的，而且在场函数事先满足强制边界条件的情况下，解通常具有明显的上、下界等性质。因此长期以来在物理和力学的微分方程的近似解法中占有重要的位置。实际应用中的困难1、求解域复杂的情况下，选取满足边界条件的试探函数很难。2、为提高精度，需增加试探函数的项数，增加了求解的复杂性。而同样建立于变分原理基础上的有限元法，虽然本质上和里兹法是类似的，由于近似函数是定义在子域上，因此可以克服上述两方面的困难。
5、        弹性力学基本方程的矩阵形式。弹性体在载荷作用下，任一点的应力状态和应变状态。三大方程：平衡方程、几何方程、物理方程（弹性矩阵和柔度矩阵）。弹性体的力的边界条件和位移边界条件。弹性力学方程式：满足微分方程和边界条件。弹性体的应变能和余能。弹性体的应变能密度和余能密度，在线弹性力学中是相等的，上下的面积。
6、        平衡方程和几何方程的等效积分弱形式-虚功原理
变形体的虚功原理描述如下：体系外力的虚功和体系内力的虚功之和为0.
虚功原理是虚位移原理和虚应力原理的总称。平衡方程和力的边界条件，几何方程和位移边界条件的等效积分的弱形式。
虚功原理在推导过程中均为涉及物理方程，因此它适用于非线弹性和弹塑性等不同的力学问题。利用最小位能原理可得到位移近似解的弹性变形能是精确变形能的下界，即近似的位移场在总体上偏小。计算模型偏硬。而利用最小余能原理得到的应力近似解的弹性余能是精确解余能的上界，即近似的应力解在整体上偏大，结构的计算模型偏于柔软。

7、        构造位移函数就是构造近似函数，位移函数等于插值函数矩阵（形函数矩阵）乘以单元结点位移列阵。确定了单元位移后，可以利用几何方程和物理方程求得单元的应变和应力，当单元的结点位移确定后，由应变矩阵得到的单元应变是常量。也就是说在载荷作用下单元中各点具有同样的应变值。因此，3结点单元也称为常应变单元。与应变矩阵相同，应力矩阵也是常量矩阵。即3结点单元中各点的应力和应变是相同的。单元的选择就是位移函数的选择，就是应力和应变以及结点位移在单元内的分布选择。
基于弹性力学最小位能原理形成的有限元求解方程的一般原理，在具体计算时涉及单元刚度矩阵的形成和单元等效结点载荷列阵的形成，以及集合单元刚度矩阵和单元等效结点载荷列阵形成结构刚度矩阵和结构等效载荷列阵，还有给定位移边界条件的引入等问题。
由单元到结构的集成。
8、        单元刚度矩阵的力学意义和性质，利用最小位能原理建立一个单元的求解方程.KA=P+F，单元等效结点载荷、相邻单元对该单元作用力。统称为点力。每个结点在各个自由度方向上各有一个平衡方程。方程的左端是通过单元结点位移表示的单元结点内力,方程右端是单元结点力。单元刚度矩阵第一列的物理意义是，当a1=1时，其它节点位移都为0时，需要在单元各节点位移方向上施加的结点力大小。因此单元刚度矩阵中任一元素Kij表示，当单元第j个节点位移单位位移而其他位移为0时，需要在第i个节点位移方向上施加的节点力大小。用自由度方向更简洁。单元刚度越大，使节点产生单位位移的节点力就越大。因此，单元刚度矩阵中的每个元素反映了单元刚度的大小，称为刚度系数。
1、        对称性，各种形式的单元都普遍具有对称的性质。
2、        奇异性：当a1=1,其它结点位移为0是，考虑单元两个方向力的平衡和绕任一点力矩的平衡。刚度系数之间存在一定关系,即3节点三角形单元6*6阶的刚度矩阵，只有3行或3列是独立的，因此刚度矩阵是奇异的，即系数矩阵的行列式为0.虽然可以在给定位移条件下，计算出单元的节点力。当给出节点力，却不能唯一地确定单元的节点位移。这是因为单元还可以有任意的刚体位移。
3、        主元恒正，主对角线上的元素恒为正，要使节点位移a=1，施加在a方向的节点力必须与位移a同向，这是结构稳定性的要求。
9、        集成。将单元刚度矩阵K扩大到与结构刚度矩阵同阶，以便进行矩相加。将单元刚度矩阵中各子块按照单元结点的实际编码安放到扩大的矩阵中，单元对结构刚度矩阵K的哪些刚度系数有贡献。单元等效节点载荷和结构载荷列阵的关系也类似。结构刚度矩阵的特点：
物理意义：结构第j个自由度的位移为单位值，而其他自由度的位移为0时，需在第i个自由度方向上施加的节点力大小。与单元不同的地方是，结构是单元的集合体，每个单元对结构起一定的作用。
由于单元刚度矩阵是对称的，奇异的，主元恒正。他们集成的结构刚度矩阵除了具有同样的性质外，还具有稀疏性和非零元素呈带状分布。
任意给定结构的节点位移所得到的结构节点力总体上是满足力和力矩平衡的。反之，给定任意满足结构力和力矩平衡的节点载荷P，由于K的奇异性却不能得到结构的位移解。这是因为结构仍可能发生任意的刚体位移，为消除K的奇异性，结构至少需要给出限制刚体位移的约束条件。
由于每个节点周围只有几个单元，因此虽然总体单元数和节点数很多，结构刚度矩阵的阶数很高，但刚度系数中的非零系数却很少，这就是刚度矩阵的大型和稀疏性，只要节点编号是合理的，这些稀疏的非零元素将集中在以对角线为中心的一条带状区域内。
10、        引入位移边界条件
最小位能原理是具有附加条件的变分原理，他要求场函数u满足几何方程和位移边界条件。现在离散模型的近似场函数在单元内满足几何方程。但选择的场函数的试探函数却没有提出在边界上满足位移边界条件的要求。
1、        直接代入法
优缺点：将已知节点位移消去，得到一组修正方程，用以求解其它待定节点的位移。若结构总体位移为n个，已知节点位移m个，则得到一组n-m个待定节点位移的修正方程组。这种方法重新组合了方程，是原方程的阶数降低，节点位移顺序性已被破坏，给程序编制带来麻烦。
2、        对角元素改1法（对于0位移），对角元素乘大数法（对于非0位移）。
对于固支边界条件，在K中将与零节点位移相对应的行列中，将主对角线元素改为1，其它元素改为0，在载荷列阵中将与零节点位移对应的元素改为0.这种方法不改变原来方程的阶数和节点未知量的顺序编号。对于节点位移为给定值aj时，第j个方程做如下修改，对角元素乘以大数a（e10）,并将载荷表示为kjjaja,这种方法简单。
在求得结构的节点位移后，可以利用单元的几何方程和物理方程求解单元的应变和应力。由于导数运算的结果，精度低于位移。表现在单元交界面上应力不连续，力的边界条件也不能精确满足，特别对于3节点三角形单元，由于它是常应力单元，相邻单元的应力会出现明显的跳跃现象。
如何选择合适的单元进行计算，涉及求解问题的类型、对计算精度的要求及经济性等多方面的因素。
(广义坐标在求近似函数时为未知量，在求应力应变时为已知量，包含常数项和一次项，刚体位移和常应变、与单元内坐标有关的应变).
广义坐标有限元位移模式的选择和插值函数的构造。
1、        广义坐标是由节点场变量确定的，因此它的个数应与节点的自由度数相等，eg，三角形3节点单元具有6个自由度，因此两个方向各取3项多项式。
2、        选择多项式时，常数项和坐标的一次项必须完备，对于三角形3节点单元，s11和s22都为常数，也反应了当单元趋于一点时，单元应变趋于常应变。常数项引入了与坐标无关的位移量值-刚体位移-否则将产生非0应变。
3、        多项式的选取应由低阶到高阶，对于每边有两个端点的应保证一次完全多项式，每边有三个端点的应保证二次完全多项式，完备性要求应变是与x,y有关。
当边界上位移不连续时，则有限元解就不可能收敛于精确解,单元的边界上位移要协调。当泛函中出现位移函数的m阶导数时，则单元位移函数至少是m阶完全多形式，并且在单元和内部具有Cm-1级连续性。
等参元计算中数值积分阶次的选择。单元的阶次和积分的阶次。二阶单元比一阶单在单元的边上多了一个节点。单元边界上位移不在是线性变化，而是二次变化，形状选择与结构的构型有关，三角形比较适合不规则形状，而四边形则比较适合规则形状。单元的阶次的选择与求解域内应力的变化特点有关，应力梯度大的区域，单元阶次应较高，否则即使网格很密也难达到理想的结果。积分阶次是精度的问题，完全积分精度高，减缩积分精度低。使用最小位能原理建立的位移有限元法，其解答具有下限性质。即有限元的计算模型具有实际结构偏大的整体刚度。选择减缩积分使有限元计算模型的刚度偏小，因此可能有利于提高计算精度。
保证结构的总体刚度矩阵K是非奇异的。完全积分没有必要，减缩积分可能出现0能量模式,而且二阶减缩积分比一阶减缩积分抵抗0能量的模式更强。

清泉石上流

挺好的。

sunshine_hust_

楼主理解的深刻，可不可以推荐些书？

bluekooozcb

赞一个，说的一针见血

山高我为峰srf

太深刻了，看来我得补补我的数学知识了。。。

郝Smlie

楼主说的非常好，赞一个!

fei123456

不懂。。。飘过   

lxcb123

太深奥了, 啊啊啊啊啊

kun_sleepy

好强

MorganMo

有点像是笔记，很多地方比较简洁， 让我较快地回忆了相应的知识点，个人感觉更像是文字版的思维导图，相当精彩。

偶然一笔

不行啊！看不懂啊我！

hbhu

感谢分享！ 有限元看起来稍学一下就能操作，但真正理解及能够解决问题还是有一定门槛的。

qiaokelijiejie

多谢楼主分享