材料成形计算机模拟 第二章 线性弹性力学问题的有限元法

第二章线性弹性力学问题的有限元法本章内容提要本章介绍有限元的基础知识,仅考虑线弹性力学问题。针对弹性力学平面问题，通过最简单的三角形常应变单元,阐明有限元用于弹性体应力分析的基本原理和方法,然后简要地介绍轴对称问题和三维问题的有限元法。简要介绍等参单元以及数值积分方法.本章最后讨论了有限元的一般化，说明有限元是一种求解微分方程边值问题和初值问题的普遍适用的离散化方法。§2-1弹性力学的基本方程一、线弹性静力分析问题的假设二、基本方程三、利用虚功原理建立有限元方程一、线弹性静力分析问题的假设1位移梯度是小量，因此应变与位移之间的关系是线性的；2物体始终保持为弹性状态，应力与应变之间的关系是线性的；3边界条件中不包含接触边界条件。二、基本方程1、应力平衡方程根据静力学平衡条件导出应力分量之间的关系式0jiijbx在静载荷作用下的物体，是在变形态构形上达到平衡的，此时物体中任一点和任一部分都应处于平衡状态。推导过程请参考《金属塑性成形过程模拟》李尚健主编机械工业出版社P21为单位体积的材料上作用的体积力ib2、几何方程（应变与位移的关系式)根据变形体的连续性和均匀性，导出应变与位移之间的关系式1()2jijjiuuxxi3、本构方程根据实验与假设导出应变与应力分量之间的关系式，也叫物理方程eijijklklC为弹性张量eijklC本构关系反映物质宏观性质的数学模型，又称本构方程。最熟知的本构关系有胡克定律、牛顿粘性定律、理想气体状态方程、热传导方程等。对于弹性体来说，材料的本构关系，即应力应变关系，是通过弹性张量来建立的，也就通常所说的广义虎克定律。三、利用虚功原理建立有限元方程1、虚功的概念变形体在产生虚位移过程中真实力(外力，应力)所作的功，称为虚功。虚位移是一种加在变形体上假想的微小的位移。虚位移必须满足变形体的约束条件并保持其连续性。而虚位移的类型(移动，转动)和方向是任意的。所谓微小的就是在产生虚位移过程中，各力的作用线保持不变。2、虚功原理如果变形体处于平衡状态，对其给一任意微小虚位移，则外力所作的总虚功必须等于内力(应力合力)所作的总虚功，这就是变形体的虚功原理。3、虚功原理的应用在弹性力学中，利用虚功原理可以对某些比较复杂的问题求得近似的解析解，而在有限元法中可对某些更复杂的问题求得近似的数值解。通常用虚功方程代替平衡方程与几何方程。虚功原理是力学中的一个普通原理，它不仅可用于线弹性问题，而且可以用于非线性弹性及弹塑性等材料的非线性问题。4、有限元方程的建立利用虚功原理来建立有限元方程。设物体处于平衡状态，给物体内各点以任意的虚位移，它仅是坐标x的函数，同时在Su上满足。根据虚功原理，外力在虚位移上所作的虚功等于因虚位移引起的虚应变能，即pijijiiiiVVSdVbudVpudS§2-2弹性力学平面问题的有限元列式一、平面问题3节点三角形单元的插值函数二、单元刚度方程三、单元等效节点载荷列阵四、整体刚度矩阵和总节点载荷列阵的集成五、引入边界条件六、实施步骤与注意事项一、平面问题3节点三角形单元的插值函数在有限单元法中，选择节点位移作为基本未知量时称为位移法；选择节点力作为基本未知量时称为力法；取一部分节点力和一部分节点位移作为基本未知量时称为混合法。位移法易于实现计算自动化，所以，在有限单元法中位移法应用范围最广。当采用位移法时，物体或结构物离散化之后，就可把单元总的一些物理量如位移，应变和应力等由节点位移来表示。这时可以对单元中位移的分布采用一些能逼近原函数的近似函数予以描述。通常，有限元法我们就将位移表示为坐标变量的简单函数，这种函数称为位移模式或位移函数。（一）单元的位移模式及插值函数在一个单元中，每个节点有两个位移分量，节点i的位移用列向量记为iiiuuv(,,)ijm(,,)ijm表示其他节点的相应分量。每个单元共有3个节点，6个位移自由度，可表示如下单元节点位移列阵()()eTTTTTijmiijjmmuuuuuvuvuv上标“T”表示转置。有限元法中单元位移模式或位移函数一般采用多项式作为近似函数。为使求得面积不致成为负值，节点局部编码1、2、3（对应的整体编码为i，j，m）的次序必须是逆时针排列。3节点三角形单元位移模式iijjmmiijjmmuNuNuNuvNvNvNv1()2iiiiNabxcyA(,,)ijm其中jNiNmN是坐标函数，它们反映单元的位移状态，因而称为形函数。3节点三角形单元位移模式写成矩阵形式为ieeijmjijmmuuuINININuNNNuNuvuN为形函数矩阵；I为二阶单位矩阵；u为单元位移函数列阵；ue为单元节点位移列阵。一、平面问题3节点三角形单元的插值函数（二）形函数的性质和面积坐标1、三角形单元的面积坐标点P在三角形中的相对位置可表示为(,,)ijmPLLL面积坐标与三节点三角形的形函数完全相同。其中iiALA(,,)ijm式中，A是原三角形的面积；称为面积坐标。,,ijmLLL2、形函数的性质（1）在节点上形函数的值满足1(,)0当当iijijjiNxyji,,)（ijm,iixxyyiuu上式是插值函数必须满足的条件，这样才能保证当时（即在节点i）有2、形函数的性质（2）在单元中任一点有当发生刚体位移时，该式保证了单元中的位移处处都等于该刚体位移值，这是求得正确解答形函数必须满足的要求。（3）3节点三角形单元的形函数是线性的，对于线性单元，相邻单元公共节点位移相等就保证了相邻单元公共边上位移的连续性。1ijmNNN3、面积坐标与直角坐标的转换关系112iiiijjjjmmmmLabcLabcxAabcyL1111ijjjjmmmmLxxxxLyyyyL用直角坐标表示的多项式（插值函数等）也可以表示为面积坐标的同阶多项式，反之亦然。一、平面问题3节点三角形单元的插值函数（三）应变矩阵和应力矩阵ε＝Lu＝LNue=L(NiNjNm)ue＝(BiBjBm)ue=Bue式中L称为平面问题的微分算子。B称为应变矩阵。其分块子矩阵是0102iiiiibBcAcb(,,)ijm式中(,,)、iibcijm都是由节点坐标确定的常数，因此单元的应变为常数，3节点三角形单元称为常应变单元。eeeeCCBuSueC是弹性矩阵，S称为应力矩阵。二、单元刚度方程采用虚功原理建立有限元刚度方程。eeepTeeTTTcVVSBCBdVuNbdVNpdSNPeTeeVBCBdVKeTebVNbdVPepTesSNpdSPTeccNPPeeeebscPPPP令eKebPesPecPeP称为单元刚度矩阵；为体积力引起的等效节点载荷；为面力引起的等效节点载荷；为集中力引起的等效节点载荷；称为单元等效节点载荷列阵。单元刚度矩阵具有如下特性：1）对称性；2）奇异性；3）主元恒为正，即0iik三、单元等效节点载荷列阵在有限元分析中，要针对节点自由度建立求解方程，因此单元所受到的面积力和体积力都要由作用于其节点的等效载荷来表示。体积力引起的等效节点载荷面积力引起的等效节点载荷集中力引起的等效节点载荷eeTbAPNbtdxdy（eA为三角形单元的面积）epeTsLPNptdl（l表示边界线段）eTccPNP对于平面问题：几种常见载荷的计算方法：1、均匀等厚单元的自重2、均布测压q四、整体刚度矩阵和总节点载荷列阵的集成（一）集成方法假设弹性体被离散成ne个三角形单元和n个节点，则整个弹性体的总虚功应是ne个单元虚功的总和。将扩阶后的全部ne个单元的单元刚度方程叠加起来，则弹性体的整体刚度方程可写成KU=P。K是所有单元的刚度矩阵之和，称为整体刚度矩阵或总刚度矩阵。P是由作用于弹性体上的载荷（包括体力、面力和集中力）分别移置到有关单元的节点上，并逐个节点加以合成得到的总等效节点力，按照节点编号由小到大排列组成的，所以它称为总节力列阵或载荷列阵。在实际编程计算的过程中不采用上述先将单元刚度矩阵和单元等效节点载荷列阵扩阶，然后再相加的方法，在计算得到Ke和Pe的各元素后，只需按照单元中各节点的整体自由度编码，“对号入座”地叠加到整体刚度矩阵和整体节点载荷列阵地相应位置上即可。这个过程又被称为组装或集成。（二）整体刚度矩阵的特点1）与单元刚度矩阵类似，整体刚度矩阵K中的一列元素kij的物理意义是要使弹性体的节点自由度i发生单位位移，而其他节点自由度都保持零位移的状态下，需要施加的与节点自由度j对应的节点力。2）整体刚度矩阵K的主元素总是正的。3）整体刚度矩阵K是对称矩阵。4）整体刚度矩阵K是一个稀疏矩阵，如果遵守一定的节点编号规则，可使非零元素都集中于主对角线附近而呈带状。5）与单元刚度矩阵类似，整体刚度矩阵K也是一个奇异矩阵。在引入位移边界条件排除刚体位移后，它是正定矩阵。五、引入边界条件应用虚功原理时，要求虚位移满足几何方程和位移边界条件。对于求位移场的问题，至少要约束系统的刚性位移，以消除整体刚度矩阵的奇异性。1、对角元素改1法（给定位移是零位移时）2、对角元素乘大数法（对任何给定位移零值或非零值都适用）引入位移边界条件消除整体刚度矩阵的奇异性以后，即可由整体刚度方程解出整体节点位移向量UU=K-1P求得节点位移U以后，将各单元的节点位移ue分别代入式（3-25）和式（3-27），即可求得各单元的应变和应力。这一步骤可称为单元分析Ⅱ。六、实施步骤与注意事项（一）有限元分析的步骤从应用的角度看，可以将整个分析过程分为建模、分网、加载、求解和后处理几个阶段。（二）注意事项1、对称性的利用在建模、分网和加载时，要充分利用物体几何形状和所受载荷的对称和反对称特性，确定是取整个物体还是取部分物体作为计算模型。2、节点的选择和单元的划分节点的选择：通常集中载荷的作用点、分布载荷强度的突变点、分布载荷与自由边界的分界点、支撑点、几何形状突变点等都应取作为节点。单元的划分：厚度不同或材料不同的区域不要划在同一个单元。在保证精度的前提下，力求采用较少的单元。在划单元时，对应力变化急剧的区域要分得细一些，应力变化平缓的区域可以分得粗一些。单元的三条边的长度不要悬殊太大，尽可能接近相等。3、节点的编号在节点编号中，应注意使同一单元的相邻节点的号码差尽可能小，以便缩小钢锯矩阵的带宽，节约计算机存储。4、应力计算结果的整理由位移法得到的位移解在全域是连续的，应变和应力解在单元内部是连续的而在单元间是不连续的，即在单元边界上发生突跳。因此需要对计算得到的应力进行处理，以改善所得到的结果。最简单的处理应力结果的方法是取相邻单元或围绕节点各单元应力的平均值。七、收敛准则和多项式位移模式阶次的选择（一）收敛准则（二）多项式位移模式阶次的选择