1、第3章 平面问题的有限元法,3.1 结构的离散化,3.2 三角形常应变单元的位移模式和形函数,3.5 整体分析,3.6 等效节点载荷计算,3.8 有限元分析的实例,3.3 单元刚度矩阵,3.4 单元位移函数的选择原则,3.7 约束条件的处理,将连续体变换为离散化结构 将连续体划分为有限多个、有限大小的单元， 并使这些单元仅在节点处连结起来，构成所谓“离散化结构”。,3.1 结构的离散化,离散化要注意:,1.单元形状的选择:,平面问题的单元，按其几何 特性可分为两类：,以三节点三角形为基础； 以任意四边形为基础。,较高精度的三角形等参数单元；运用非常广泛的四边形等参数单元。,这两类都可以增加节点
2、也构成一系列单元：,首选三角形单元和等参数单元。,2.对称性的利用,利用结构和载荷的对称性：如结构和载荷都对于某轴对称，可以取一半来分析；若对于x轴和y轴都对称，可以取四分之一来分析。,3.单元的划分原则,通常集中载荷的作用点、分布载荷强度的突变点、分布载荷与自由边界的分界点，支承点都应取为节点,单元的形状和尺寸可以根据要求进行调整。对于重要或应力变化急剧的部位，单元应划分得小些；对于次要和应力变化缓慢的部位，单元可划分得大些；中间地带以大小逐渐变化的单元来过渡。,单元的划分原则,单元数量要根据计算精度和计算机的容量来决定。在保证精度的前提下，尽可能减少单元数量。,不要把不同厚度或不同材料的区
3、域划分在一个单元里。,单元的划分原则,根据误差分析，应力及位移的误差都和单元的最小内角正弦成反比，所以单元的边长力求接近相等。即单元的三（四）条边长尽量不要悬殊太大。,4.节点的编号,应尽量使同一单元的节点编号相差小些，以减少整体刚度矩阵的半带宽，节约计算机存储。,上图，节点顺短边编号为好。,3.2 三角形常应变单元的位移模式和形函数,首先以平面单元中最基本的三节点三角形单元为例，介绍有限元法。,单元分析的步骤可表示如下：,节点位移,内部各点位移,应变,应力,节点力,单元分析,分为四步求出相邻各量之间的转换关系,综合起来,得出由节点位移求节点力的转换关系:,单元刚度矩阵,位移模式,1.位移模式
4、,单元的若干个节点有基本未知量，即,位移模式: 单元内任一点的位移表达式，假定为坐标的简单函数。 反映单元的位移分布形态，是单元内的插值函数。 在节点处等于该节点位移。,位移模式可表示为：,N为形态矩阵（形函数矩阵）,平面问题每个节点位移分量有两个，所以整个单元有6个节点位移分量，即6个自由度。,单元节点位移列阵:,三角形单元有6个自由度，内部任一点的位移是由6个节点 位移分量完全确定的，位移模式中应含有6个待定系数， 所以位移模式可取为：,位移函数一般用多项式来构造。,位移模式: 单元内任一点的位移表达式，假定为坐标的简单函数。 反映单元的位移分布形态。,在弹性体内，位移变化非常复杂。有限元
5、法将整个弹性体分割成许多小单元，在每个单元内采用简单的函数来近似表达单元的真实位移，将各单元连接起来，便可近似表达整个弹性体的真实位移函数。,这种化整为零、化繁为简的方法，正是有限元法的精华。,假设节点i，j，m的坐标分别为(xi,yi),(xj,yj),(xm,ym),2.形函数,联立求解左边3个方程，得：,其中A为三角形单元的面积,注意：为了使得出的面积值不为负值，节点i,j,m的次序必须是逆时针。至于将那个节点作为起始点i则没有关系。,同理，求解右边的三个方程，得到a4，a5，a6，解得：,i，j，m轮换,整理后得：,其中Ni，Nj，Nm是坐标的线性函数，反应了单元的位移形态，称为形（状
6、）函数。,写成矩阵形式,式中：I 二阶单位阵，N 形函数矩阵,3.三角形面积坐标,定义：在三角形内任一点P，向三个角点（节点）连线，将原三角形分割成三个子三角形，设子三角形的面积分别是：Ai，Aj，Am，则：,即面积坐标定义为子三角形与原三角形面积之比； 记为：P（Li，Lj，Lm）。,面积坐标的性质：,1.,Li，Lj，Lm中只有两个是独立的。,2.三角形三个角点处,3.三条边上,i-j:Lm=0,j-m:Li=0,m-i:Lj=0,形心处：,推论：三角形内一条平行于三角形任一边的直线上的各点，具有相同的与该边对应的坐标值。,面积坐标与直角坐标的转换：,（i,j,m),(i,j,m),因此：
7、,即三角形面积坐标就是三角形相应的形函数。,所以，位移模式也可以用面积坐标表示为：,(i,j,m),将面积坐标的表达式：,写成矩阵形式：,求逆得：,第1行展开为面积坐标性质1，第2行和第3行展开即为局部的面积坐标和整体直角坐标的关系：,例 题,下图为一平面应力的直角三角形单元，直角边长均为a,厚度为t，弹性模量为E,泊松比=0.3,求形函数。,1.单元应变,应变矩阵为常量，单元内应变是常数,3.3 单元刚度矩阵,2.单元应力,S称为应力转换矩阵,应用平面应力问题的弹性矩阵：,应变矩阵为常量，单元内应力也是常数，相邻单元的应变与应力将产生突变，但位移是连续的。,能量转换与守恒定律，是自然界基本的
8、运动规律之一。,实功原理：处于平衡状态的可变形固体，在受外力作用而变形时外力对其相应的位移所做的功（实功），等于积蓄在物体中的应变能（实应变能）。,能量法的优点：与坐标系的选择无关，因而应用极为广泛。,能量法与数学工具变分法的结合，导出虚位移（虚功）原理，使得用数学分析的方法解决力学问题的理论得到发展而更趋完善。,3.虚位移（功）原理,单元节点力列阵：,单元节点虚位移列阵：,节点力在虚位移所做的功：,简写为：,4.单元刚度矩阵,单元虚应变：,单元内应力在虚应变上所做的功（虚应变能）：,其中：t为单元厚度,单元应力：,单元刚度矩阵ke取决于单元的大小、方向和弹性常数，而与单元的位置无关，即不随单
9、元或坐标轴的平行移动而改变。,对于三角形常应变单元：,单元刚度矩阵为对称矩阵。,例 题,下图为一平面应力的直角三角形单元，直角边长均为a,厚度为t，弹性模量为E,泊松比=0.3,求单元刚度矩阵。,理论力学中质点、质点系（刚体）的虚位移原理；,材料力学中杆件的虚位移原理。,弹性力学中的虚位移（虚功）原理：,在外力作用下处于平衡状态的变形体，当给与该物体微小位移时，外力总虚功在数值上等于变形体的总虚应变能。,虚：微小的、任意的、可能的，变分的思路,实功是力在自己产生位移上所做的功，虚功是力在别的（人为的）因素产生的位移上做的功。所谓”虚“并不是虚无，而是可能、虚设的意思。,“虚”的表达：,虚位移（
10、虚功）原理：,3.4 单元位移函数的选择原则,三角形常应变单元简单，精度较差，要提高精度： 1.增加单元数目和节点数目； 2.采用更高精度的单元。,FEM中的一系列工作，都是以位移模式为基础的。所以当单元趋于很小时，即x, y0时，为了使FEM之解逼近于真解，即为了保证FEM收敛性,位移模式应满足下列条件：,1. 位移模式必须能反映单元的刚体位移。,单元位移包含两部分：本单元的形变引起的位移；其他单元的形变引起的位移，即刚体位移。在位移函数中，常数项即提供刚体位移。,2. 位移模式必须能反映单元的常量应变。,单元应变包含两部分：变量应变和常量应变。位移函数的一次项提供常量应变。,当单元0时，单
11、元中的位移和应变都趋近于基本量刚体位移和常量位移。,3. 位移模式应尽可能反映位移的连续性,使相邻单元之间的位移保持连续，即受力后，相邻单元在公共边界上，即既不互相脱离，也不互相嵌入。 使相邻单元在公共节点处具有相同的位移。,使单元内部的位移保持连续。位移函数取坐标的单值连续函数。,满足条件1、2的单元，称为完备单元；满足条件3的单元，称为协调单元。,常采用“帕斯卡三角形”来选取位移模式代数多项式的形式。,按照帕斯卡三角形选择位移模式的原则：,1.多项式的阶次及项数，由单元的节点数目和自由度数目来决定。保证多项式中的待定系数同单元的自由度数目相一致，以避免在确定待定系数时增加困难。,2.当高次
12、多项式只选取一部分项时，应遵循“对称性”原则，即取其最高次中的位置对称的相应项，以保证在各坐标轴方向上具有相同的精度。,3.应满足完备性和协调性要求。,3节点三角形单元：,6节点三角形单元：,4节点四边形单元：,3.5 整体分析,结构的整体分析是将离散后的所有单元通过节点连接成原结构，进行分析。分析过程是将所有单元平衡方程组集成整体平衡方程，引进边界条件后求解整体节点位移向量。,整体平衡方程：,F=K,K为整体刚度矩阵,设弹性体被划分为N个三角形单元和n个节点，则结构就有2n个自由度。,K2n2n,整体刚度矩阵的组装：,例：求下面结构的整体刚度矩阵,解：1）结构离散，单元和节点编码 用三角形单
13、元把该结构分成4个单元，6个节点,节点两种编码：一是节点总码；二是节点局部码，每个三角形单元的三个节点按逆时针方向的顺序各自编码为i，j,m。,单元1：节点号码1，2，3,单元2：节点号码2，5，3,单元3：节点号码5，6，3,单元4：节点号码2，4，5,2）分别写出各个单元的分块刚度矩阵：,单元1：节点号码1，2，3,单元2：节点号码2，5，3,单元3：节点号码5，6，3,单元4：节点号码2，4，5,3）组装整体刚度矩阵,利用单元分块矩阵中，各子块的节点和单元信息，直接把单元刚度的各元素送入总体刚度矩阵的相应行列上，并同总体刚度矩阵该元素的已有值相加。“对号入座”,组装一般规则：,1)当Kr
14、s中r=s时，该点被哪几个单元所共有，则整体刚度矩阵中的子矩阵Krs就是这几个单元的刚度矩阵中的子矩阵Krse的相加。,2)当Krs中rs时，若rs边是组合体的内边，则整体刚度矩阵中的子矩阵Krs就是共用该边的两相邻单元刚度矩阵中的子矩阵Krse的相加。,3)当Krs中r和s不同属于任何单元时，整体刚度矩阵中的子矩阵Krs=0 。,整体刚度矩阵的性质：,1）整体刚度矩阵是对称矩阵。,2）整体刚度矩阵每一个元素的物理意义：,3）整体刚度矩阵的主对角线上的元素总是正的。,4）整体刚度矩阵是一个奇异阵。只有排除刚体位移后，K才是正定的，其逆矩阵才存在。,在 F=K中，令节点1在x方向的位移u1=1，
15、而其余节点位移均为0，则：,5）整体刚度矩阵是一个稀疏阵。,离散后结构的任一节点，只和与它相连的元素发生联系，所以K存在大量的零元素，而非零元素往往分布在主对角线的附近。,带形矩阵,半带宽：在半个斜带形区域内，每行具有的元素个数，用d表示。,半带宽d=（相邻节点码的最大差值+1）2,半带存储：利用带形矩阵的特点和矩阵的对称性，计算机中可以只存储上半带的元素。,在同一网格中，如果采用不同的编码方式，则相应的半带宽也可能不同。,应采取合理的节点编码方式（使相邻节点码尽可能小），以便得到最小的半带宽，从而节约计算机存储容量。,不同的编码方式，相邻节点的最大差值分别为4，6，8，半带宽分别为10，14
16、，18。,3.6 等效节点载荷计算,根据有限元法的思想，所有有关的量都要转换为节点的量。,结构所受的载荷也必须转换为等效的节点载荷。,整体刚度方程中的载荷列阵F，是由弹性体全部单元等效节点力集合而成，而单元的等效节点力，是由作用在单元上的集中力、表面力和体积力分别移植到节点上，再逐点加以合成求得。,1.单元自重：,下面用上述公式计算几种常用载荷作用下的等效节点力。,三角形单元i,j,m的厚度为t，重度为，面积为A，则体积力：,节点力为：,由形函数的性质得：,则：,受自重载荷作用下的等效节点力为单元重量的1/3。,2.均布面力：,三角形单元i,j,m的ij边上作用有均匀的分布力，集度为：,单元节
17、点力为：,由形函数性质：,把作用于ij边上的均布面力按静力等效平均分配到该边两端的节点上。,3.线性分布面力：,三角形单元i,j,m的ij边上作用有三角形分布表面力 设j点表面力为0，i点集度为：,4.集中力：,集中力G作用与ij边上作用,总载荷的2/3分配给i点，1/3分配给j点。,整体刚度矩阵的奇异性，可以通过引入边界约束条件来排除弹性体的刚体位移，以达到求解的目的。引用边界条件后，待求节点未知量的数目和方程的数目可相应的减少。,3.7 约束条件的处理,引入节点位移最常用的方法有以下两种：,计算机常用的方法是，以某种方法引入已知的节点位移（包括零约束位移），而保持非常原有的数目不变，只是修
18、正K和F中的某些元素，以避免计算机存储做大的变动。,设已知u1=1，u2=3，则,若已知节点i在y方向位移vi，则令K中的元素K（2i）（2i）为1，第2i行和第2i列的其余元素都为零。F中的第2i个元素则用位移vi的已知值代入，F中的其他各行元素都减去节点位移的已知值与原来K中这行的相应元素的乘积。,若已知节点i在x方向位移ui，则令K中的元素K（2i-1）（2i-1）为1，第2i-1行和第2i-1列的其余元素都为零。F中的第2i-1个元素则用位移ui的已知值代入，F中的其他各行元素都减去节点位移的已知值与原来K中这行的相应元素的乘积。,1. 化1置0法,2. 乘大数法,将K中与已知节点位移
19、相关的主对角线元素乘上一个计算机可接受的充分大的数，同时将F中的对应元素换上已知节点位移与对角线元素及同一个大数的乘积。,设已知u1=1，u2=3，则,3.8 有限元分析的实例,有限元法的解题过程,2.结构的离散化。包括单元划分、节点和单元编号、节点坐标计算。,3.等效节点力的计算。 按单元逐个进行分析，计算体积力、表面力和集中力的等效节点力，进行叠加，得到每个单元的等效节点力载荷。,对每个节点，所有环绕该节点的单元节点力求和，得到整个结构的节点力载荷列阵。,1.力学模型的确定。根据工程实际情况确定问题的力学模型，并按一定比例绘制结构图，确定尺寸、载荷和约束情况等。,4.计算各单元的刚度矩阵。
20、由各单元的常数bi、bj、bm、ci、cj、cm及单元面积、弹性常数，计算各单元的刚度矩阵。,5.组装整体刚度矩阵。将各个单元的刚度矩阵组集在一起，形成整体刚度矩阵。,6.建立整体平衡方程，引入约束条件，求解节点位移。,7.求解单元应力。根据节点位移e表示的节点应力的公式，计算单元应力。,有限元分析实例,两端固定的矩形深梁,跨度为2a,梁高为a,厚度为t,已知E,=0,承受均布压力q,用有限元法求解此平面应力问题。,解：由对称性，取梁的一半进行分析。载荷和约束情况如图。,1）结构离散：划分为2个三角形单元，4个节点。,单元1：节点号码2，3，1,单元2：节点号码3，2，4,2）单元刚度矩阵的计算,单元1的刚度矩阵：,单元2的刚度矩阵：,2 3 1,3 2 4,23 1,32 4,3）整体刚度矩阵的组装,4）处理载荷和位移边界条件，求解节点位移,Fy1=0， Fy3=qa/2,所以：,位移边界条件：,采用“化1置0法”，再划去相应的行列：,解得：,5）应力计算,应力矩阵,得单元节点位移向量,计算应力,作 业,1. 下图为一固定端梁受集中力P作用，试用图所示单元划分，求出节点位移（按平面应力问题计算，取E为常量，=1/6，t=1)。,2. 三角形单元ijm的jm边作用有右图所示线性分布面载荷，求节点载荷向量。,
