0.5有限元原理基础知识

24 0．5 有限元原理基础知识 弹性力学问题的本质是求解偏微分方程的边值问题。弹性力学的基本方程是以偏微分方程组来 示的。在应力时，一般总是从构件的连续性出发，依据无限小单元的物理～数学模型建立微 分方程式，然后求解微分方程式获得问题的解。从理论上讲，弹性力学能解决一切弹性体的应力和 应变问题。但在工程实际中，一般构件的形状、受力状态、边界条件都比较复杂，所以除少数的典 型问题外，对大多数工程实际问题，由于偏微分方程边值问题的复杂性，往往无法用弹性力学的基 本方程直接进行解析求解，只能采取各种近似方法或者渐近方法通过数值计算来求得其近似解。变 分原理就是将弹性力学的基本方程－—偏微分方程的边值问题转换为代数方程求解的一种方法。在 弹性力学中，主要有四个变分原理：虚位移原理、最小势能原理、虚应力原理和最小余能原理。 其中，最小势能原理如下陈述：在所有的可能位移场中，真实位移场的总势能取最小值。所以 这一原理称为最小势能原理。 U WΠ = − { } [ ]{ } { } { } { } { }∫∫∫ −−= S T V T V T dSpudVFudVD εε 2 1 数学描述即总势能的一阶变分为零，而且二阶变分是正定的（大于零）。 00 2 >Π=Π δδ , 即 { } [ ]{ } { } { } { } { }1 0 2 T T T V V S D dV u F dV u p dSδ δ ε ε   Π = − − =    ∫ ∫ ∫ 必须强调指出的是，真实位移与其他的可能位移之间的差别在于是否满足静力平衡条件，所以 说最小势能原理是用变分形式表达的平衡条件。 通过总势能的一阶变分为零，可以推导出平衡微分方程和面力边界条件，这和虚功原理是相同 的，即最小势能原理也等价于平衡微分方程和面力边界条件。 最小势能原理的主要用途并非推导平衡微分方程和面力边界条件，它是弹性力学问题近似解法 的基础。如果要使得某个原理要应用于实际问题，必须有对应的求解方法。本节介绍基于最小势能 原理的两种近似解法: 瑞利-里兹（Rayleigh-Ritz）法和伽辽金（Гапёркин）法。 根据最小势能原理，如果能够列出所有的几何可能位移，那么使总势能Π取最小值的那一组位 移就是真实位移。问题是列出所有几何可能的位移是非常困难的，甚至是不可能的。 因此，对于实际问题的计算，只能凭借经验和直觉缩小寻找范围，在这个范围内的一族几何可 能的位移中，找到一组位移使得总势能Π最小。 虽然这一组位移一般的说并不是真实的，但是可以肯定，它是在这个缩小的给定范围内部，与 真实位移最为接近的一组位移，由此解答可以作为近似解。 从上述思想出发，在一般情况下，可以将位移分量选择为如下的形式 0 1 0 1 0 1 n m m m n m m m n m m m u u A u v v B v w w C w = = =  = +   = +   = +  ∑ ∑ ∑ 25 其中， , ,m m mA B C 均为任意的常数； 0 0 0, ,u v w 以及 , ,m m mu v w 都是坐标的已知函数，并且在位移边 界 uS 上，有 这样构造的位移试函数，不论系数 , ,m m mA B C 取何值，总是满足位移边界条件的。而且对于连 续函数，必然满足几何方程。因此满足几何可能位移的条件。 现在的问题是将要如何选择待定系数 , ,m m mA B C ，使得总势能Π在位移表达式表示的这一族位移中 取最小值。 为此，将位移表达式代入几何方程求得应变分量，然后代入总势能Π的表达式，注意到应变能 密度函数是应变分量的齐二次函数，因此总势能Π表达式的第一个积分成为待定系数 , ,m m mA B C 的 齐二次函数，而第二和第三个积分为 , ,m m mA B C 的一次函数。于是，总势能Π原本是自变函数的泛 函，现在成为待定系数 , ,m m mA B C 的二次函数。 这样就把求解泛函的极值问题，转化成为求解函数的极值问题。总势能 Π取极值的条件为 0 mA ∂Π = ∂ ， 0 mB ∂Π = ∂ ， 0 mC ∂Π = ∂ 总势能 Π 0 0x m x m m V V S U dV F u dV p u dS A σ ∂ − − = ∂ ∫∫∫ ∫∫∫ ∫∫ 取极值的条件又可以写作 0 0y m y m m V V S U dV F v dV p v dS B σ ∂ − − = ∂ ∫∫∫ ∫∫∫ ∫∫ 0 0z m z m m V V S U dV F w dV p w dS C σ ∂ − − = ∂ ∫∫∫ ∫∫∫ ∫∫ 上述公式是一组以 , ,m m mA B C （m=1，2，3…）为未知数的线性非齐次代数方程组，求解方程可 得待定系数，回代就可以得到近似位移解答。这一方法称为瑞利—里茨法。 例： 两端简支的等截面梁，受均匀分布载荷q作用如图所示。试求解梁的挠度 ( )w x 。 26 解：首先使用瑞利—里茨法求解。 为了满足梁的位移边界条件，即简支梁两端的约束条件: 在 x=0和 l 处，w=0，取位移试函数， 即挠曲线方程为 ∑= m m l xmsinCw π 问题的总势能为 22 2 0 02 l lEI d w dx qwdx dx   Π = −    ∫ ∫ 即 4 4 2 3 1,3,5 2 4 m m m m CEI qlm C l m π π = Π = −∑ ∑  根据 0 mC ∂Π = ∂ , 所以 所以： 4 5 5 4 m qlC EI mπ = （m为奇数） 0mC = （m为偶数） 回代到位移公式，可得 4 5 5 1,3,5, 4 1 sin m ql m xw EI m l π π = = ∑  。 挠曲线表达式是无穷级数，它给出了本问题的精确解答。 这个级数收敛很快，只要取少数几项 就可以得到足够的精度。 有限元原理是目前工程上应用最为广泛的结构数值分析方法，它的理论基础仍然是弹性力学的 变分原理。那么，为什么变分原理在工程上的应用有限，而有限元原理却应用广泛。有限元原理与 一般的变分原理求解方法有什么不同呢。问题在于变分原理用于弹性体分析时，不论是瑞利-里茨法 还是伽辽金法，采用整体建立位移试函数或者应力试函数的方法。由于试函数要满足一定的条件， 导致对于实际工程问题求解仍然困难重重。 有限元方法选取的试函数不是整体的，而是在弹性体内分区（单元）完成的，因此试函数形式 简单统一。当然，这使得转换的代数方程阶数比较高。但是，面对强大的计算机处理能力，线性方 程组的求解不再有任何困难。因此，有限元原理成为目前工程结构分析的重要工具。 近年来，随着现代科学技术的发展，特别是计算机技术的迅速发展和广泛应用，使得有限元方 法首先在弹性力学和结构力学领域发展起来。以有限元方法为代表的计算力学的发展，迅速改变了 弹性力学理论和方法在工程应用领域的处境。以计算机的强大计算能力为后盾开发的大型通用有限 元程序，可以求解数十万自由度的线性代数方程组，目前已经成为工程技术人员手中强大的结构分 析工具。在此基础之上，CAD, CAE等技术的应用使得计算机不仅成为数值分析的工具，而且成为设 计分析的工具。 本节将从位移变分方程引出有限元方法的基本概念。对于最小势能原理，物体的总势能为 { } { } { } { }0 T T V V S U dV u F dV u p dS σ Π = − −∫∫∫ ∫∫∫ ∫∫ 27 { } [ ]{ } { } { }1 2 T T V D u F dVε ε = −  ∫∫∫ { } { } T S u p dS σ −∫∫ 如果将物体分解为若干个有限尺寸的单元，则物体总势能为所有单元体总势能的和。有 1 m e e= Π = Π∑ ( )1, 2,3, ,e m=  其中 e为单元序号，m为单元总数。而任意一个单元体的总势能为 { } [ ]{ } { } { }12 e TTe e V D u F dVε ε Π = −  ∫∫∫ { } { }( )e Te S u p dS σ − ∫∫ 这里Ve,(S)e分别表示第e个单元的体积和面力边界。显然如果选取的位移试函数是连续的，物 体的总势能可以用单元总势能的和表示。有限元分析中，物体的位移是由单元位移确定的，因此位 移连续需要分单元选取的位移试函数保证单元的边界位移与所有相邻单元位移相同。有限元原理就 是采用分单元选取位移函数，而选取的位移函数在弹性体内又是连续的这一基本思想的变分方法。 所以有限元方法的理论基础是最小势能原理。当然，这一思想同样可以应用于最小余能原理，乃至 广义变分原理等。 可以看出，微分形式的求解方法与积分形式的求解方法有着本质上的不同，正是引入了试函数， 使得求解的难度大大降低。由于工程问题非常复杂，要求所采用的方法具有较好的性、较低的 难度、较低的函数连续性要求、较明确的物理概念、较好的通用性。而基于最小势能原理的求解方 法具有较明显的综合优势，因此，可以在该原理的基础上发展出能广泛适用于工程中任意复杂问题 的求解方法，但必须处理的技术难点是： ·复杂物体的几何描述 ·试函数(许可位移场)的确定与选取(规范化形式) ·全场试函数的表达 因此，在具备大规模计算能力的前提下，将复杂的几何物体等效离散为一系列的形状的几 何体，再在标准的几何体上研究规范化的试函数表达及其全场试函数的构建，然后利用最小势能原 理建立起力学问题的线性方程组，这就是有限元方法的基本思路。在计算机技术高度发展的今天， 有限元方法在工程中得到最广泛的应用，也发展了上千种具有完善功能的有限元分析商品化软件。 28 0．6有限元方法的基本思想 实际上，可以认为有限单元法的概念是源于结构理论。对于一个杆系结构，通常是由许多结 构单元(杆件)所组成。结构中的这些单元仅在有限个节点上彼此相连。对每个单元而言，诸如力与 位移之间的关系这样的结构特性都是用节点上所确认的自由度来惟一地予以规定；而整体杆系结构 的特性则可通过组集这些单元特性来加以描述。 图 0.6.1所示的是一个由二根杆件组成的铰接桁架，杆件的截面积为 A，弹性模量为 E，长度分 别为 1l 和 2l 。该桁架在各铰接点处受有外力 321321 Y,Y,Y,X,X,X 。因杆件在节点处是铰接，不承受 弯矩，只能承受轴向力，所以，每个节点的力和位移各有两个分量，即每个节点均具有两个自由度， 而每个单元则有四个自由度。为此，必须用四个方程来描述每个单元的力与位移之间的关系。对于 单元①，有 1 1 1 1 1 1 11 1 12 1 13 2 14 2U k u k v k u k v′ ′ ′ ′= + + + 1 1 1 1 1 1 21 1 22 1 23 2 24 2V k u k v k u k v′ ′ ′ ′= + + + 1 1 1 1 1 2 31 1 32 1 33 2 34 2U k u k v k u k v′ ′ ′ ′= + + + （0.6.1） 1 1 1 1 1 1 41 1 42 1 43 2 44 2V k u k v k u k v′ ′ ′ ′= + + + 式中, 12 1 2 1 1 1 1 V,U,V,U 分别为节点 1和节点 2施于单元①的节点力沿坐标方向的分量； 1 2 1 2 1 1 1 1 v,u,v,u 分 别为节点 1和节点 2的位移沿相应坐标方向的分量。此处，上标是单元编号，下标为节点编号。 图 0.6.1铰接桁架 可将式(0.6.1)写成矩阵形式 1 111 12 13 14 1 1 1 1 21 22 23 241 1 1 1 2 231 32 33 34 1 1 2 2 41 42 43 44 k k k kU u k k k kV v U uk k k k V vk k k k  ′ ′ ′ ′         ′ ′ ′ ′    =     ′ ′ ′ ′            ′ ′ ′ ′   (0.6.2) 29 或简记为: { } [ ] { }111 δkR = (0.6.3) 式中,{ }1δ { }Tv,u,v,u 12121111= 被称为单元① 的节点位移向量；{ } =1R { }TVUVU 12121111 称为单元① 的节点力向量； [ ]1k 为单元①的刚度矩阵，其中 11 12 44, , ,k k k′ ′ ′ 为刚度系数。 刚度系数的物理意义如下： 若令只在节点 1处有沿 x方向的单位位移，其他位移不存在，则有 111 =u 及 012 1 2 1 1 === vuv 则得 1 1 11U k′= ， 1 1 21V k′= ， 1 2 31U k′= ， 1 2 41V k′= 这表明，当节点 l 沿 x方向产生一单位位移( 111 =u )，而其余节点、其余自由度方向上的位移为零 时，刚度系数就等于施加于该单元①各自由度方向上的力。这些力组成一个平衡力系，它们表示单 元①抵抗位移 11u 的刚度。这些力的值很容易根据材料力学知识求得：当位移 1 1 1 =u ，其余节点位移 都等于零时（见下图），单元的长度将缩短 1 cosl θ∆ = − 。 图 0.6.2 杆单元 在此状态下，根据材料力学知识，在节点 1 处需要向单元①施加的轴向压力为 1 1 1 cosEA EAl l l θ  ∆ =    ，这就是节点 1作用于单元①上的力，它在 x和 y方向的分量是 =11U 2 11 1 cosEAk l θ′ = 30 =11V 21 1 cos sinEAk l θ θ′ = 对于节点 2作用于单元①的两个自由度方向上的力，则与之大小相等、方向相反。即 =12U 2 31 1 cosEAk l θ′ = − =12V 41 1 cos sinEAk l θ θ′ = − 令 111 =v 及 012 1 2 1 1 === vuu 则得 1 1 12U k′= ， 1 1 22V k′= ， 1 2 32U k′= ， 1 2 42V k′= 再继续作类似的分析（见附页），便可得到其余的刚度系数．有 12 1sin cosk EA lθ θ′ = 2 22 1sink EA lθ′ = 由于作用力和反作用力的关系 32 1sin cosk EA lθ θ′ = − 2 42 1sink EA lθ′ = − 对节点 2做类似的分析，可得到其它 8个刚度系数，整理得： [ ] 2 2 2 2 1 2 2 1 2 2 cos cos sin cos cos sin cos sin sin cos sin sin cos cos sin cos cos sin cos sin sin cos sin sin EAk l θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ  − −  − − =  − −   − −  (0.6.4) 同理，可以求出作用于单元②的节点力和位移之间的关系----刚度矩阵 =               2 3 2 3 2 2 2 2 V U V U [ ]2k 2 2 2 2 2 3 2 3 u v u v               31 [ ]2 2 0 0 0 0 0 1 0 1 0 0 0 0 0 1 0 1 EAk l    − =    −  11 12 13 14 21 22 23 24 31 32 33 34 41 42 43 44 k k k k k k k k k k k k k k k k  ′′ ′′ ′′ ′′    ′′ ′′ ′′ ′′ =   ′′ ′′ ′′ ′′    ′′ ′′ ′′ ′′   (0.6.5) 对两个单元的力----位移关系加以集合得到本整体结构的力----位移关系。为了获得整体结构 的力与位移之间的关系，需要引入整体结构的节点位移分量 11 v,u , 22 v,u , 33 v,u 和单元节点位移分 量 11 ii v,u , 22 ii v,u , 33 ii v,u (上标 i为单元编号)之间的协调关系，即 1 11 uu = 1 11 vv = 1 2 2 2 2u u u= = 1 2 2 2 2v v v= = 2 3 3u u= 2 3 3v v= 另外，根据力的平衡条件，作用在节点上的外力应该等于与该节点相连的各单元所受到的节点力之 和。因此可得到结构的力与位移的关系 ( )1 2 21 1 1 1 2 2 1 cos cos sin cos cos sinEAX U u v u v l θ θ θ θ θ θ= = + − − ( )1 2 21 1 1 1 2 2 1 cos sin sin cos sin sinEAY V u v u v l θ θ θ θ θ θ= = + − − ( )1 2 2 22 2 2 1 1 2 2 1 cos cos sin cos cos sinEAX U U u v u v l θ θ θ θ θ θ= + = − − + + ( )1 2 2 22 2 2 1 1 2 2 1 cos sin sin cos sin sinEAY V V u v u v l θ θ θ θ θ θ= + = − − + + (0.6.6) ( )2 3 2 EA v v l + − 2 3 3 0X U= = 2 3 3Y V= = ( )2 3 2 EA v v l − − 这些方程就是整体结构的力----位移关系，写成矩阵形式．有 32 1 1 2 2 3 3 X Y X Y X Y          =           2 2 1 1 1 1 2 2 1 1 1 1 2 2 1 1 1 1 2 2 1 1 1 1 2 2 2 2 cos cos sin cos cos sin 0 0 cos sin sin cos sin sin 0 0 cos cos sin cos cos sin 0 0 cos sin sin cos sin sin 1 0 1 0 0 0 0 0 0 0 0 0 1 0 1 l l l l l l l l l l l l EA l l l l l l l l θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ  − −  − −   − −   − − + −      −   1 1 2 2 3 3 u v u v u v                     (0.6.7) 或简记为 { } [ ]{ }δKR = (0.6.8) 上式表明了节点力与节点位移之间的关系, 这就是有限单元法所要建立的基本方程组。上式中{ }R 是由作用在节点上的外载荷所组成的向量。称为载荷向量(或载荷列阵)；{ }δ 是由基本未知量—— 节点位移所组成的向量；矩阵 [ ]K 表示了作用在结构各单元上的总体效果,所以称为结构的整体刚度 矩阵。 [ ] 2 2 1 1 1 1 2 2 1 1 1 1 2 2 1 1 1 1 2 2 1 1 1 1 2 2 2 2 cos cos sin cos cos sin 0 0 cos sin sin cos sin sin 0 0 cos cos sin cos cos sin 0 0 cos sin sin cos sin sin 1 0 1 0 0 0 0 0 0 0 0 0 1 0 1 l l l l l l l l l l l l K EA l l l l l l l l θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ θ  − −  − −   − − =   − − + −      −   (0.6.9) 在本问题中，我们须加入边界条件，才能求出问题的解。根据 2点的力边界条件 PY,X −== 22 0 以及 1 点和 3 点的位移边界条件 03311 ==== vuvu ，确定未知的位移 22 v,u 和未知的力 3311 Y,X,Y,X 。 由式(0.6.8){ } [ ]{ }δKR = 得 [ ]                     − =                     3 3 1 1 2 2 0 0 0 0 0 Y X P Y X v u K (0.6.10) 由矩阵(0.6.9)可知，其行列式为零，故总刚矩阵 [ ]K 是一个奇异阵，上式需要合适的方法求解。 33 对应(0.6.10)式，可知 [ ]K 中的元素只有 44433433 k,k,k,k 与 22 v,u 有关，既       − =             Pv u kk kk 0 2 2 4443 3433 (0.6.11) 或写成       − =             + Pv u llsinlsincos lsincoslcos EA 0 1 2 2 21 2 1 11 2 θθθ θθθ (0.6.12) 即可以求出       − =       1 2 2 2 θtg EA Pl v u (0.6.13) 将(0.6.13)代回(0.6.10)求未知的 3311 Y,X,Y,X ，有               =                   3 3 1 1 2 2 64 54 24 14 63 53 23 13 Y X Y X v u k k k k k k k k 即               =                     − − − − − 3 3 1 1 2 2 2 1 2 1 1 1 2 1 0 0 0 Y X Y X v u l lsin lsincos lsincos lcos θ θθ θθ θ 得               =               2 2 3 3 1 1 1 0 0 0 l EA Pl Y X Y X 不难看出，结构的整体刚度矩阵是由各单元的刚度矩阵叠加组成的。在上式的 [ ]K 矩阵中，左 上方的虚线部分恰好是单元①的刚度矩阵，而右下方的虚线部分正是单元②的刚度矩阵；两虚线框 的重叠部分中的元素，则是两个单元刚度矩阵在同一位置处的元素之和。 结构的整体刚度矩阵具有许多特性： 1. 它是一个对称矩阵； 2. 对角线上的主元素 iik 总是非负的，因为作用力的方向将与它引起的对应位移的方向相同； 3. 整体刚度矩阵是奇异的。根据行列式的性质可以知道，矩阵[K]的对应行列式的值等于零，所以 它是奇异的。故求解时，必须引入几何边界条件，消除刚体位移，方可求出未知位移。 总刚矩阵的物理意义： 结构中第 s个自由度的单位变形所引起的第 r个节点力。 34 此时，方程组(0.2.7)还不能立即用来求解结点位移，其物理原因是结构的几何约束尚未设置， 可能产生刚体位移。只有加上几何边界条件，对刚度矩阵加以修改，排除刚体位移后，才能解出全 部位移分量。 以上我们通过这一简单的例子，说明了刚度矩阵的物理意义和它的性质以及从单元刚度矩阵集 合成整体刚度矩阵的概念。这些性质带有普遍性，可以推广到连续体问题中去。建立整体结构的刚 度矩阵是运用有限单元法求解问题的核心内容，一旦获得了整体刚度矩阵，就等于列出了有限单元 法的基本方程。而建立整体刚度矩阵的问题，又可归结为求单元的刚度矩阵问题。 对于一个连续体的求解问题，有限单元法的实质就是将具有无限多个自由度的连续体，理想化 为只有有限个自由度的单元集合体，单元之间仅在节点处相连接，从而使问题简化为适合于数值求 解的结构型问题。这样，只要确定了单元的力学特性，就可以按结构分析的方法来进行求解。 有限单元法的分析过程 通过以上的简单论述，可以把有限单元法的分析过程归纳为以下几个方面。 1．结构的离散化 结构的离散化是进行有限单元法分析的第一步。数学上，把无限自由度处理成有限自由度的过程 叫做”离散化”。有限单元法中的结构离散化过程，简单地说，就是将分析的对象划分为有限个单元 体，并在单元上选定一定数量的点作为节点，各单元体之间仅在指定的节点处相连。有限单元法的 整个分析过程就是针对这种单元集合体来进行的。单元的划分，通常需要考虑分析对象的结构形状 和受载情况。对于前面所研究讨论的桁架，其单元的划分比较简单，因为分析对象本身就是由一系 列杆件相互连接而成，所以可直接取每根杆件作为一个单元。 但是，对于其他非杆件的机械结构物，如齿轮、轧机机架等，为了能有效地逼近实际的分析对象， 就必须认真考虑划分、选择何种类型单元以及划分的单元数目等等。对于一些比较复杂的结构， 有时还要采用几种不同类型的单元来进行离散化。许多大型有限元分析软件都备有多达几十种单元 类型的单元库，供分析计算人员选用。常用的主要有杆单元、平面单元、块单元、等参元、壳单元 等，以后将陆续介绍、讨论有关这方面的具体实施方法。 35 36 2．位移模式的选择 有限单元法是应用局部的近似解来求得整个问题的解的一种方法。根据分块近似的思想，可以选 择一个简单的函数来近似地构造每一单元内的近似解。本书中讲授的有限单元法是以节点位移为基 37 本未知量，所以为了能用节点位移表示单元体的位移、应变和应力，在分析求解时，必须对单元中 位移的分布作出一定的假设，即选择一个简单的函数来近似地表示单元位移分量随坐标变化的分布 规律，这种函数称为位移模式。 位移模式的选择是有限单元法分析中的关键。由于多项式的数学运算比较简单、易于处理，所 以通常是选用多项式作为位移模式。多项式的项数和阶数的选择，一般要考虑单元的自由度数，和 解答的收敛性要求等，它的阶次应包含常数项和线性项。以后要对此作详细的讨论。 根据所选定的位移模式，就可以导出用结点位移表示单元内任一点位移的关系式，其矩阵形式是 { } [ ]{ }e eu N δ= (0.6.14) 式中 { }eu 为单元内任一点的位移列阵，{ }eδ 为单元的节点位移列阵， [ ]N 称为形函数矩阵，它的 元素是位置坐标的函数。 在此，我们顺便指出：有限单元法比起经典的近似法具有明显的优越性。例如，在经典的里兹 法中，要求选取一个函数来近似地描述整个求解区域中的位移，并须满足边界条件；而在有限单元 法中则采用分块近似，只需对一种单元选择一个近似位移函数。此时，不必考虑位移边界条件，只 须考虑单元之间位移的连续性就可以了。这样做当然比起在整个区域中选取一个连续函数要简单得 多，特别是对于复杂的几何形状或者材料性质。对于作用载荷有突变的结构，采用分段函数，比起 采用连续性较强的整段函数来近似精确的位移函数更为适宜。 3．单元的力学特性分析 位移模式选定以后，就可以进行单元力学特性的分析。分析单元的力学特性主要包括以下三部 分内容： (1)通过几何方程由（0.6.14）建立单元应变与节点位移的关系式； { } [ ]{ }eB δε = （0.6.15） { }ε 是单元内任一点的应变列阵， [ ]B 称为单元应变矩阵。 (2)利用物理方程，由应变的表达式（0.6.15）导出单元应力与节点位移的关系式； { } [ ][ ]{ }eBD δσ = （0.6.16） { }σ 是单元内任一点的应力列阵， [ ]D 是单元材料有关的弹性矩阵。 (3)由虚功原理推出作用于单元上的节点力与节点位移之间的关系式，及单元的刚度方程。 { } [ ]{ }ee kR δ= （0.6.17） [ ]k 称为单元刚度矩阵，在以后将推得 [ ] [ ] [ ][ ]∫∫∫= dxdydzBDBk T （0.6.18） 上式的积分应遍及整个单元的体积。 在以上三项中，导出单元刚度矩阵是单元特性分析的核心内容。 4．等效节点力的计算 弹性体经过离散化后，假定力是通过结点从一个单元传递到另一个单元，但是作为实际的连续 体，力是从单元的公共边界传递到另一个单元的。因而，这种作用在单元边界上的表面力以及作用 38 在单元上的体积力、集中力等都需要等效移置到结点上去，也就是用等效的结点力来替代所有作用 在单元上的力。移置的方法是按照作用在单元上的力与等效结点力，在任何虚位移上的虚功都相等 的原则(下称虚功等效)进行的，具体的计算将在以后介绍。 5．集合所有单元的刚度方程，建立整个结构的平衡方程 建立整体结构的平衡方程也叫做结构的整体分析，这个集合过程包括有两方面的内容。一是由各 个单元的刚度矩阵集合成整个物体的整体刚度矩阵；二是将作用于各单元的等效结点力列阵集合成 整体结构的节点载荷向量，即总的载荷列阵。最常用的集合刚度矩阵的方法是直接刚度法，即按节 点编号对号入座，直接利用单元刚度矩阵中的刚度系数子阵进行叠加。一般来说，从单元到整体的 组集过程主要是依据两点：一是所有相邻的单元在公共节点处的位移相等；二是所有各节点必须满 足平衡条件。于是得到以整体刚度矩阵 [ ]K 、载荷列阵 [ ]R 以及整个物体的结点位移列阵{ }δ 表示的 整个结构的平衡方程 [ ]{ } { }RK =δ （0.6.19） 这些方程还应在考虑了几何边界条件并作适当的修改之后，才能够解出所有的未知结点位移。 6．求解未知的节点位移及单元应力 在上述组集整体刚度矩阵时，没有考虑整体结构的平衡条件，所以组集得到的整体刚度矩阵是一 个奇异矩阵，尚不能对平衡方程直接进行求解。只有在引入边界约束条件、对所建立的平衡方程加 以适当的修改之后，方可根据方程组的具体特点选择恰当的计算方法来求得节点位移，继而求出单 元应变和应力。应注意的是，引入边界条件修改平衡方程实质上就是消除整体结构的刚体位移。 因此，有限元原理就是利用在单元内部的构造位移试函数，根据单元和整个弹性体总势能最小 的原则确定位移场。不难看出这与位移变分方程的类似之处，只是位移变分方程需要整体确定位移 试函数，而有限元原理是在单元内部根据单元节点位移确定的。因此有限元原理适应几乎任何形状 的弹性体。 尽管有限元原理选取的位移函数属于简单函数，但是可以证明分析是收敛的，因此随着单元的 逐渐减小，计算模型趋于原结构，因此计算精度可以保证的。 下面给出有限元分析的主要过程： 1. 建立计算模型（结构离散化，划分单元，节点编号等）； 2. 计算单元刚度矩阵和单元的节点力列阵； 3. 形成整体刚度矩阵； 4. 计算整体节点力； 5. 引入位移边界条件； 6. 求解整体平衡方程； 7. 回代求解单元应力等。 0．5 有限元原理基础知识 0．6有限元方法的基本思想