第2章_弹性力学基础及有限元法的基本原理1

1、第第2 2章章 有限元法的基本原理有限元法的基本原理平面问题有限元法平面问题有限元法（1 1）力学应用）力学应用引言（2 2）力学回顾）力学回顾力学回顾力学回顾力学回顾力学回顾第一节第一节 弹性力学基础知识弹性力学基础知识 1.1 1.1 弹性力学中的物理量弹性力学中的物理量1 1、载荷、载荷定义：作用在弹性体上的力（力矩），又称外力。定义：作用在弹性体上的力（力矩），又称外力。 载荷可分为：载荷可分为：体力、面力、集中力体力、面力、集中力 1 1）体力）体力 定义：分布于整个弹性体体积内的外力。如：定义：分布于整个弹性体体积内的外力。如：重力重力 可以分解为三个坐标系上的分量，用向量表示为：

2、可以分解为三个坐标系上的分量，用向量表示为： TvvxvyvzPppp2 2）面力）面力定义：作用于弹性体表面上的外力，如：定义：作用于弹性体表面上的外力，如：流体压力流体压力可以分解为三个坐标系上的分量，用向量表示为：可以分解为三个坐标系上的分量，用向量表示为：3 3）集中力）集中力定义：集中在某一点上的外力，如：定义：集中在某一点上的外力，如：牵引力牵引力可以分解为三个坐标系上的分量，用向量表示为：可以分解为三个坐标系上的分量，用向量表示为： TssxsyszPppp TccxcyczPppp2 2、应力（注意下标）、应力（注意下标）定义：弹性体内某一点作用于定义：弹性体内某一点作用于某个

3、截面某个截面单位面积单位面积上的内力，上的内力，反映了反映了内力在截面上的分布密度内力在截面上的分布密度。微元体表面上的应力：一个正应力微元体表面上的应力：一个正应力 （拉压）（拉压） 两个切应力两个切应力 （剪切）（剪切）切应力互等定律：切应力互等定律：符号标注规则：第一个角标是所在面垂直于哪个轴；符号标注规则：第一个角标是所在面垂直于哪个轴； 第二个是指向哪个轴第二个是指向哪个轴弹性体内某一点的应力状态由六个应力所决定弹性体内某一点的应力状态由六个应力所决定应力向量可以表示为：应力向量可以表示为：, , xyyxyzzyzxxz Txyzxyyzzx3 3、应变（对应于应力）、应变（对应于

4、应力）定义：微元体体发生变形后，定义：微元体体发生变形后，单位长度的变形量单位长度的变形量。对应于应力，应变向量可以表示为：对应于应力，应变向量可以表示为：4 4、位移、位移定义：弹性体定义：弹性体内质点位置的变化内质点位置的变化。位移向量可以表示为：。位移向量可以表示为： Txyzxyyzzx Tduvw 1.2 弹性力学的基本方程弹性力学的基本方程主要是描述应力、应变、位移及外力间的相互关系主要是描述应力、应变、位移及外力间的相互关系1、平衡方程、平衡方程 （应力间的关系应力间的关系）000yxxzxvxxyyzyvyyzxzzvzpxyzpxyzpxyz 0TvLPL微分算子2、几何方程

5、（应变与位移的关系）、几何方程（应变与位移的关系） 000000000 xyzxyyzzxxuxvyyuwzzvuvyxwyxvwzyzywuxzzx Ld3、物理方程（应力与应变之间的关系）、物理方程（应力与应变之间的关系）111111xxyzyyzxzzxyxyxyyzyzzxzxEEEGGG 2(1)EGEG其中: 为为为剪切弹性模量且:杨氏弹性模量泊松比 100011100011100011(1)1 2(1)(1 2 )000002(1)1 2000002(1)1 2000002(1)DED因此物理方程可以简写为： 弹性矩阵E杨氏弹性模量1D问题的弹性模量泊松比是指材料在单向受拉或受压

6、时，横向正应变与轴向正应变的绝对值的比值，也叫横向变形系数，它是反映材料横向变形的弹性常数。可以这样记忆：空气的泊松比为0，45#钢0.3，水的泊松比为0.5，中间的可以推出。若在弹性范围内加载，横向应变x与纵向应变y之间存在下列关系：x=- y式中为材料的一个弹性常数，称为泊松比。泊松比是量纲为一的量。 未知数未知数 应力应力 6个个+应变应变 6个个+位移位移 3个个=15个个 方程个数方程个数 平衡方程平衡方程 3个个+几何方程几何方程6个个+物理方程物理方程6个个=15个个原则上可以根据原则上可以根据15个方程求出个方程求出15个未知物理量个未知物理量但实际求解时先求出一部分再通过方程

7、求解剩下的。但实际求解时先求出一部分再通过方程求解剩下的。 目前有限元法主要采用的是目前有限元法主要采用的是位移法位移法，以三个位移，以三个位移分量为分量为基本未知量基本未知量。位移位移- -应变应变- -应力，应力和外力平衡应力，应力和外力平衡 1.3 虚位移原理虚位移原理（1）虚功与虚应变能）虚功与虚应变能 弹性体在外力作用下变形，外力对弹性体做功，弹性体在外力作用下变形，外力对弹性体做功，所做的功以所做的功以应变能应变能的形式储存于弹性体中。的形式储存于弹性体中。弹性体弹性体单位体积单位体积的的应变能应变能为：为：虚位移虚位移定义：在约束条件允许的范围内弹性体可能发生的定义：在约束条件允

8、许的范围内弹性体可能发生的 任意任意微小位移微小位移。 虚位移与时间及外载荷无关虚位移与时间及外载荷无关实际位移是在外载荷作用下可能的虚位移实际位移是在外载荷作用下可能的虚位移 12TU弹性体在平衡状态下发生虚位移弹性体在平衡状态下发生虚位移1 1）外力所做的虚功为：）外力所做的虚功为：2 2）应力在虚应变上所做的虚功，也就是存储在弹性）应力在虚应变上所做的虚功，也就是存储在弹性 体内的虚应变能为：体内的虚应变能为： TWfRWfR其中：为虚功，为虚位移，为外力。 TVUdV（2）虚位移原理）虚位移原理表述：如果在在虚位移发生之前弹性体是平衡的，表述：如果在在虚位移发生之前弹性体是平衡的， 那

9、么在虚位移发生时那么在虚位移发生时外力在虚位移上所做的外力在虚位移上所做的 功就等于弹性体的虚应变能功就等于弹性体的虚应变能。即：。即：当外力的形式是多样的时，外力的虚功等于：当外力的形式是多样的时，外力的虚功等于：WU TTTcvsvsWfPfP dVfP dS 1.4 1.4 平面问题定义平面问题定义 严格地讲，任何结构都是空间的。对于某些特殊情严格地讲，任何结构都是空间的。对于某些特殊情况，空间问题可以转化为平面问题。况，空间问题可以转化为平面问题。（1 1）平面应力问题）平面应力问题满足条件：满足条件：1 1）几何条件）几何条件 厚度尺寸远远小于截面尺寸；厚度尺寸远远小于截面尺寸；2

10、2）载荷条件）载荷条件 载荷平行于板平面且载荷平行于板平面且沿厚度方向均匀沿厚度方向均匀分布分布，而，而板平面板平面不受任何外力作用。不受任何外力作用。此时，应力应变分量变为：几何方程 物理方程00 1zxzyzzxzyzxy以及 TxyxyTxyxy 00 xyxyxuvxyxx 2101011002DED简化为8个方程（2）平面应变问题满足条件1）几何条件 结构呈等截面的细长形；2）载荷条件 载荷垂直于厚度方向(平行横截面)且沿厚度均匀分布，两个端面不受力。此时，00 zxyzzyzzxzxy 以及应力应变分量变为： 物理方程 TxyxyTxyxy 10111011 211 2002 1D

11、ED简化为8个方程 1.5 1.5 弹性力学的参量及方程汇总弹性力学的参量及方程汇总1 1）参量）参量位移：位移：应力：应力：应变：应变：2 2）方程）方程平衡方程：平衡方程：几何方程：几何方程：物理方程：物理方程： Tduvw Txyzxyyzzx Txyzxyyzzx 0TvLP Ld D TvvxvyvzPppp其中：其中： 100011100011100011(1)1 2(1)(1 2 )000002(1)1 2000002(1)1 2000002(1)ED 000000000 xyzLyxzyzx111111xxyzyyzxzzxyxyxyyzyzzxzxEEEGGG 12）没有正应

12、力没有正应变）没有正应变没有正应力3）没有应变没有位移4）没有位移没有应变第第2 2章作业一：章作业一：手画微分体应力分量图。手画微分体应力分量图。 回顾回顾 弹性力学的参量及方程汇总弹性力学的参量及方程汇总1 1）参量）参量位移：位移：应力：应力：应变：应变：2 2）方程）方程平衡方程：平衡方程：几何方程：几何方程：物理方程：物理方程：3 3）虚位移、虚功原理）虚位移、虚功原理4 4）平面应力问题）平面应力问题 ；平面应力问题；平面应力问题 Tduvw Txyzxyyzzx Txyzxyyzzx 0TvLP Ld D TvvxvyvzPpppWU0 zzxzy0 zyzzx第二节第二节 平面

13、问题的有限元法平面问题的有限元法 00,zzTdu v平面应力平面问题平面应变共同点：三个方向的位移只有两个是独立的即： 2.1 2.1 平面问题的有限元分析步骤（平面应力问题）平面问题的有限元分析步骤（平面应力问题）1 1、结构离散、结构离散离散：将一个连续的弹性体分割为一定形状和数量离散：将一个连续的弹性体分割为一定形状和数量 的的单元的组合单元的组合 单元也称为网格单元也称为网格 连续体连续体 有限个单元的组合体有限个单元的组合体可用于离散的单元：可用于离散的单元： 三角形单元；三角形单元； 矩形单元；矩形单元； 不规则四边形单元。不规则四边形单元。DOF节点的自由度：节点所具有的节点的

14、自由度：节点所具有的位移分量的数量位移分量的数量。一个单元所有节点的自由度总和称为单元自由度。一个单元所有节点的自由度总和称为单元自由度。（1）单元参数只能通过节点传递到相邻单元）单元参数只能通过节点传递到相邻单元（2）单元和节点必须统一编号）单元和节点必须统一编号2.2 2.2 单元分析（位移、应力、应变）单元分析（位移、应力、应变）任务：形成单元刚度矩阵，建立单元特性方程任务：形成单元刚度矩阵，建立单元特性方程因此必须建立坐标系，如下图：因此必须建立坐标系，如下图：1 1）位移函数）位移函数分片插值分片插值 假设一种函数来表示单元位移分布假设一种函数来表示单元位移分布一般选取多项式（简单而

15、且易求导）一般选取多项式（简单而且易求导）对于三角形单节点单元（对于三角形单节点单元（DOF=6DOF=6）2212345622123456( , )( , )u x yxyxxyyv x yxyxxyy123456( , ) ( , )u x yxyv x yxy（2-1）126123456123456123456,6, iijjmmiiiiiijjjjjjmmmmmmu vu vuvuxyvxyuxyvxyuxyvxy 是 个待定系数(广义坐标)，可以表示为节点坐标和位移的函数.设三个节点的位移分别为、和，将其与节点坐标一起代入（2-1）可得： (2-2)根据克莱姆法则（ Cramer，瑞

16、士数学家），可求得：: , , , , , , 11 221 2ijmmjijmijmjmiimjmijmimjmmjmijmijiijjmmiijjmmiijjmmax yx ybyycxxax yx ybyycxxax yx ybyycxxaua ua ubub ub uAAcuc uc uA123设求解方程（2-2）得：1 211 22iijjmmiijjmmiijjmmava va vAbvb vb vcvc vc vAA45611: 12112 12 iijjmmiiiijjjjmmmmiiiijjjjmmmmxyAxyxyuab xc y uAab xc y uab xc y uv

17、ab xc y vAab xc y vab xc y v其中为三角形单元的面积则位移函数可以表示为：121 (23)212 (24)iiiijjjjmmmmiijjmmiijjmmNab xc yANab xc yANab xc yAuN uN uN uvN vN vN v引入形函数：则单元位移可以表示为： 000000ijmijmTeiijmmejNNNdNNNNNquvuvuvq用矩阵表示为：其中：为形函数矩阵为单元节点位移列阵注：（1）形函数只与节点坐标有关而与节点位移无关（2）单元的位移函数就可以表示为：单元坐标的函数与节点位移列阵的乘积。10, iiiiiuvuNvNNi 形函数的物

18、理意义：当，而另两个节点位移为 的时候 因此，是当节点 在某坐标上发生单位位移而其他节点的位移为0时，单元内的位移分布说明1 iijjmmiijjmmuN uN uN uvN vN vN v1( ,)1(,)(,)0(,)1( ,)(,)0(,)1( ,)(,)0( , )( , )( , )1iiiiijjimmjjjjiijmmmmmmiiijjijmNiN x yN xyN xyNxyNx yNxyNxyNx yN xyN x yNx yNx y 形函数的性质：）在 节点上值为1，而在其他节点处为0,2）单元的任意一点处，三个形函数之和为1(刚体位移)3)单元每一条边的形函数只( , )

19、1,( , ),( , )0iiijmjijixxxxN x yNx yNx yxxxx 与该边的节点位置有关 而与其它节点位置无关，如在i,j边上:说明2位移函数位移函数应该满足以下几个条件应该满足以下几个条件（1 1）包括常数项（保证刚体位移）包括常数项（保证刚体位移）（2 2）包括一次项（保证常应变）包括一次项（保证常应变）（3 3）保证位移的连续性（性质）保证位移的连续性（性质3 3保证）保证）（4 4）各项几何同性（）各项几何同性（x, yx, y应该是可以互换的应该是可以互换的）( , )( , )0( , )( , )mniijjiijjNx yNx yu x yN uN uv

20、x yN vN v iijjmmiijjmmuN uN uN uvN vN vN v123456( , ) ( , )u x yxyv x yxy满足上述三个条件的目的是满足有限元的收敛性满足上述三个条件的目的是满足有限元的收敛性（1 1）和（）和（2 2）是收敛的必要条件）是收敛的必要条件 完备性条件完备性条件（3 3）是收敛的充分条件）是收敛的充分条件 协调条件协调条件注意：非协调单元的解不一定不收敛注意：非协调单元的解不一定不收敛2)单元应力和应变单元应力和应变将位移表达式（将位移表达式（2-3）和（）和（2-4）代入几何方程得：）代入几何方程得： 263523560000010002x

21、yxyiiijmejijmjiijjmmmmxuvyyxuvbbbucccBqvAcb cb cbuv 将、代入上式，可得： 000000ijmijmijmiijjmmBLNbbbcccBBBcb cb cb其中：应变矩阵则应变向量可以表示为：则应变向量可以表示为： DD、BB均为常数矩阵，因此三角形三节点单元为常均为常数矩阵，因此三角形三节点单元为常应力单元。应力单元。 2 (, ,)2 11122eeijmllllllllDDBqSqSDBSSSbcESbcDBli j mAcb3 3）单元刚度矩阵）单元刚度矩阵设作用在单元节点上的单元节点力列阵为：设作用在单元节点上的单元节点力列阵为：而

22、节点发生的虚位移为：而节点发生的虚位移为：则节点力在虚位移上做的虚功为：则节点力在虚位移上做的虚功为： TTeijmixiyjxjymxmyFFFFFFFFFFTeiijjmmquvuvuv iixiiyjjxjjymmxmmyeTeWu Fv Fu Fv Fu Fv FqF TTvTeTeTeTTeTeTeTTeTeeTeTeTeeTeUdVtdxdytBqqBDBqUqBtdxdyqBDBqtdxdyqBDBqtdxdyqBDBqtAWUqFqBDBqtA单元的虚应变能：其中： 为板的厚度根据根据虚位移原理：考虑到虚位移的任意性，等式两边消 eTTeeeeFkqqBDBqtA去得：（2-5

23、） 2( , ,)2 21122 114 122 iiijimejijjjmmimjmmrsTrsrsrsrsrsrsrsrsrsrsDBkkkkkkkkkkkr si j mkBDB tAb bc cb cc bEtAc bb cc cb b将、的表达式代入（2-5）可得：其中为分块矩阵将（将（2-52-5）展开得：）展开得：单元刚度矩阵的物理意义：单元刚度矩阵的物理意义： 在一个节点处产生单位位移而其他点为零时，在一个节点处产生单位位移而其他点为零时，在该节点上需要的在该节点上需要的外力的大小外力的大小。单元某个元素的影响单元某个元素的影响 (26)iiiiijjimmjjiijjjjmm

24、mmiimjjmmmFkqkqkqFkqkqkqFkqkqkq 单元刚度矩阵的特性单元刚度矩阵的特性（1 1）对称性）对称性： （弹性力学互等定理）（弹性力学互等定理）（2 2）奇异性：）奇异性： （刚体位移）（刚体位移） eeTkk 0eeeekkqF不等于零，无力下单元无位移，与弹性力学理论相悖。不等于零，无力下单元无位移，与弹性力学理论相悖。2.3 2.3 总刚度矩阵的集成总刚度矩阵的集成 通过单元特性方程通过单元特性方程 并不能求出单并不能求出单元节点位移元节点位移 。因为。因为 包含单元间的作用包含单元间的作用力。因此，必须将每个单元的特性方程相加消除力。因此，必须将每个单元的特性方

25、程相加消除内力的影响内力的影响。这就是总刚度矩阵集成的目的。这就是总刚度矩阵集成的目的。1 1）总刚集成原理）总刚集成原理 在整个结构中，一个节点为几个单元共有。在整个结构中，一个节点为几个单元共有。在第在第i i个节点处的平衡方程为：个节点处的平衡方程为： eeekqF eq eF , , eeiissiees i j mFkqR 1, ,112121, ,inneissiies i j minnnsies i j mRikqRKqRqq qqRR RRKKk 为作用在节点 上的外载荷对所有节点都应用平衡方程，则有：或者：其中：是所有节点位移分量组成的列阵是作用在所有节点上的载荷组成的列阵为

26、总刚度矩阵 且例：总刚的形成过程例：总刚的形成过程 1 2 3 4 5 6 1 2 3 4 5 62 2）总刚度矩阵集成过程）总刚度矩阵集成过程（1 1）扩阶过程（可由转换矩阵完成）扩阶过程（可由转换矩阵完成）（2 2）叠加过程：）叠加过程： 1eneiKK 6 2 1 0000000000002 2TeeKGkGGnijmnIGIII 其中：为阶转换矩阵其中： 是的单位矩阵总刚矩阵的特点：总刚矩阵的特点：（1 1）对称性）对称性 节省存储容量节省存储容量（2 2）稀疏性）稀疏性 可能存在大量零元素可能存在大量零元素（3 3）带状性）带状性 半带宽与节点的编号有关半带宽与节点的编号有关（4 4

27、）奇异性）奇异性 保证刚体位移保证刚体位移2.42.4 载荷移置载荷移置 移置可能在局部产生误差，但不会影响整个结移置可能在局部产生误差，但不会影响整个结构的力学特性。构的力学特性。1 1）集中力的移置（虚功等效）集中力的移置（虚功等效） ccTccxcyeiicejjcPmmcPPppRNPRRNPRNP对于集中力移置后的等效节点载荷为：2 2）面力的移置）面力的移置 sssTssxsyTesPeisiejjsPmPmsPppRNP tdlNP tdlRRRNP tdlRNP tdl对于面力移置后的等效节点载荷为：3 3）体力的移置）体力的移置 vvvTvvxvyTevPeiviejjvPm

28、PmvPppRNP tdxdyNP tdxdyRRRNP tdxdyRNP tdxdy对于体力移置后的等效节点载荷为：2.5 2.5 约束处理约束处理1 1）边界位移为零）边界位移为零2 2）边界位移为已知量）边界位移为已知量2.6 2.6 求解线性方程组求解线性方程组 2.7 2.7 计算其它物理量计算其它物理量2.82.8 计算结果处理计算结果处理2.92.9 结果显示、打印、分析结果显示、打印、分析 KqR%-function k=Triangle2D3Node_Stiffness(E,NU,t,xi,yi,xj,yj,xm,ym,ID)%该函数计算单元的刚度矩阵%输入弹性模量E，泊松比

29、NU，厚度t%输入三个节点i、j、m的坐标xi,yi,xj,yj,xm,ym%输入平面问题性质指示参数ID(1为平面应力，2为平面应变)%输出单元刚度矩阵k(6X6)%-A = (xi*(yj-ym) + xj*(ym-yi) + xm*(yi-yj)/2;betai = yj-ym;betaj = ym-yi;betam = yi-yj;gammai = xm-xj;gammaj = xi-xm;gammam = xj-xi;B = betai 0 betaj 0 betam 0 ; 0 gammai 0 gammaj 0 gammam ; gammai betai gammaj betaj gammam betam/(2*A);if ID = 1 D = (E/(1-NU*NU)*1 NU 0 ; NU 1 0 ; 0 0 (1-NU)/2;elseif ID = 2 D = (E/(1+NU)/(1-2*NU)*1-NU NU 0 ; NU 1-NU 0 ; 0 0 (1-2*NU)/2;endk= t*A*B*D*B;（1） Triangle2D3Node_Stiffness第三节第三节 平面平面3 3节点三角形单元及节点三角形单元及MATLABMATLAB编程编程%-function z = Triangle2D3Node_Assembly(KK,k,i,j,