第二章 弹性力学平面问题有限元法1

1、 唐文亭唐文亭 1338927295813389272958 2011.9 2011.9 内容结构内容结构第一章第一章 概述概述第六章第六章 轴对称问题的有限单元法轴对称问题的有限单元法第七章第七章 软件介绍及计算示例软件介绍及计算示例第五章第五章 空间问题的有限单元法空间问题的有限单元法第四章第四章 等参元等参元第三章第三章 杆系结构静力分析的有限单元法杆系结构静力分析的有限单元法F第二章第二章 弹性力学平面问题有限元法弹性力学平面问题有限元法有限元法的基本思想及优越性有限元法的基本思想及优越性在应用有限元法时，我们首先将一个连续的弹性体看作由在应用有限元法时，我们首先将一个连续的弹性体看作

2、由许多尺寸有限的小单元许多尺寸有限的小单元-有限元组成。有限元组成。 这就是所谓这就是所谓区域划分区域划分，在数学上称为，在数学上称为“离散化离散化”。2. 2. 根据计算对象的简化模型，单元的形状取成平面三角形根据计算对象的简化模型，单元的形状取成平面三角形或四边形，四面体或六面体等。单元与单元之间，通过若或四边形，四面体或六面体等。单元与单元之间，通过若干个称为干个称为“节点节点”的点铰接相连，由此组合成整体。的点铰接相连，由此组合成整体。3. 3. 以一个个小单元为计算单位，首先进行单元分析，然后以一个个小单元为计算单位，首先进行单元分析，然后把它们组装起来，进行整体分析，最后求出结构的

3、近似解。把它们组装起来，进行整体分析，最后求出结构的近似解。这种把复杂结构看成有限个单元组成的整体，就是有限元法的这种把复杂结构看成有限个单元组成的整体，就是有限元法的基本思想基本思想。因此，有限元法是以变分原理和分片插值为基础的。得到的是近似解。真实系统真实系统有限元模型有限元模型10节点：节点：空间中的坐标位置，具有一空间中的坐标位置，具有一定自由度和存在相互物理作用定自由度和存在相互物理作用单元单元: : 一组节点自由度间相互作用的一组节点自由度间相互作用的数值矩阵描述数值矩阵描述( (称为刚度或系数矩阵称为刚度或系数矩阵) )单元有线、面或实体以及二维或三维单元有线、面或实体以及二维或

4、三维的单元等种类。的单元等种类。载荷载荷12节点自由度是随连接该节点节点自由度是随连接该节点单元类型单元类型变化的变化的JIIJJKLILKIPOMNKJIL三维杆单元三维杆单元 (铰接铰接)UX, UY, UZ三维梁单元三维梁单元二维或轴对称实体单元二维或轴对称实体单元UX, UY三维四边形壳单元三维四边形壳单元UX, UY, UZ,三维实体热单元三维实体热单元TEMPJPOMNKJIL三维实体结构单元三维实体结构单元ROTX, ROTY, ROTZROTX, ROTY, ROTZUX, UY, UZ,UX, UY, UZ应力计算方程求解整体分析区域 剖分 单元分析 严格地说，任何弹性体都是

5、处于三维受力状态，因而都是空严格地说，任何弹性体都是处于三维受力状态，因而都是空间 问 题 ， 但 是 在 一 定 条 件 下 ， 许 多 空 间 问间 问 题 ， 但 是 在 一 定 条 件 下 ， 许 多 空 间 问题都可以简化成平面问题。题都可以简化成平面问题。 平面问题可以分为两类：平面问题可以分为两类：平面应力问题平面应力问题和和平面应变问题平面应变问题。 平面问题应力状态2.1 2.1 平面应力问题和平面应变问题平面应力问题和平面应变问题 平面应力问题 如图所示的深梁结构，其厚度方向的尺寸远比其它两个方向的如图所示的深梁结构，其厚度方向的尺寸远比其它两个方向的尺寸小得多，可视为一薄

6、板。它只承受作用在其平面内的载荷，且尺寸小得多，可视为一薄板。它只承受作用在其平面内的载荷，且沿厚度方向不变，计算时以中性面为研究对象。其沿厚度方向不变，计算时以中性面为研究对象。其力学特点力学特点是：是：, 0, 0, 0zyyzzxxzz0z平面应力问题的应力应变转换矩阵即弹性矩阵为：平面应力问题的应力应变转换矩阵即弹性矩阵为：。 2100010112ED平面应力问题平面应力问题 平面应变问题平面应变问题平面应变问题 图示为一圆形涵洞的横截面。其长度方图示为一圆形涵洞的横截面。其长度方向上的尺寸远比其它两个方向上的尺寸大得向上的尺寸远比其它两个方向上的尺寸大得多，同样，载荷作用在多，同样，

7、载荷作用在xy坐标面内，且沿坐标面内，且沿z轴轴方向均匀分布。其力学特点是：方向均匀分布。其力学特点是：0,0,0zxzyz但一般情况下: 0z平面应变问题的弹性矩阵只需将式平面应变问题的弹性矩阵只需将式(4-1)中的中的 E 换成换成 21E换成换成 ， 1即可。即可。 )1 (22100011011)21)(1 ()1 (uED。无论是平面应力问题还是平面应变问题的应力无论是平面应力问题还是平面应变问题的应力 与 应变应变 之间的关系均为： 0D Txyyx Txyyx其中其中： 0为初应变。式中(a) 三结点三角形单元 (b) 四结点正方形单元 (c) 四结点矩形单元 (d) 四结点四边

8、形单元平面问题单元的主要类型图2-1剖分要一直进行到弹性区域的边界上剖分要一直进行到弹性区域的边界上 当边界是直线段时，就取其为三角形单元的一条边；当边界是曲线时，则在每小段上用相应的直线近似地代替曲线而作为三角形单元的一边，如图2-1单元的大小和数目要根据精度的要求和计算机容量来确定单元的大小和数目要根据精度的要求和计算机容量来确定。图2-2a图2-2b 具体进行剖分时，一般应注意以下几点具体进行剖分时，一般应注意以下几点：图2-3a图2-3b图2-4 (a) 均匀受力板力学模型 (b) 力学模型离散化平面问题有限单元法的计算力学模型( ,)iii x y(,)jjj xy(,)mmm xy

9、iuivmujumvjv图2-5 jjjuv 对于平面问题，三个节点的位移分别为： mmmuv iiiuv mvmuivjvjuiu(,)iii xy(,)jjj xy(,)mmm xy iiiejjjmmmuvuvuv(2 1) a123456uxyvxy126, b123123123iiijjjmmmuxyuxyuxy456456456iiijjjmmmvxyvxyvxy c123121212iijjmmiijjmmiijjmma uauaub ub ubuc uc ucu d456121212iijjmmiijjmmiijjmma vavavb vb vbvc vcvcv e,ijmmj

10、ijmijmjm iimjmijmimijjimijmijax yx ybyycxxax yxybyycxxaxyx ybyycxx(22)12jmijmijiimmjx yx yx yx yx yx y111211()21()2iijjmmijjiijmxyxyxyb cb caaa(23)111222111222i ijjm mi ijjm mi ijjm mi ij jm mi ij jm mi ij jm muauaua ububub uxcucuc uyvavava vbv bvb vxcvcvc vy (24)111222111222ijiiijjjmmmiiijmjjmmijmm

11、uabxc yab xc yab xc yvabuuuvvxc yab xc yab xc y v f121212eiiiiejjjjemmmmNab xc yNabxcyNabxcy(25)eeeiijjmmeeeiijjmmuN uN uN uvN vN vN v(26) 000000iieeejijmeeejijmmmuvuNNNufvNNNvuv 222eeeeijmI NI NI N21001I 222eeeeijmNI NI NI N111222111222iijjmmiijjmmiijjmmi ijjm mi ijjm mi ijjm muaua ua ubub ub uxcuc

12、 uc uyvava va vbvb vb vxcvc vc vy (24)(2-9a)由由(2-4)(2-4)式，即可得到用节点位移表示单元任一点的应变表达式：式，即可得到用节点位移表示单元任一点的应变表达式：1()21()21() () ()2xi ijjm myi ijjm mxyi ii ijjjjm mm mubub ub uxvcvc vc vyuvcubvc ub vc ub vyx00010002iijjmmuvbbbmxijucccymijvcbcbcbxymmiijjuveeg轾犏犏犏犏犏犏犏犏臌轾犏犏犏轾犏犏犏犏=犏犏D犏犏犏犏臌犏犏犏臌(2-9b)00010002iij

13、jmmuvbbbmxijucccymijvcbcbcbxymmiijjuveeg轾犏犏犏犏犏犏犏犏臌轾犏犏犏轾犏犏犏犏=犏犏D犏犏犏犏臌犏犏犏臌1()21()21() () ()2xi ijjm myi ijjm mxyi ii ijjjjm mm mubub ub uxvcvc vc vyuvcubvc ub vc ub vyx0001 0002bbbmijBcccmijcbcbcbmmiijj轾犏犏=犏D 犏犏犏臌(2-11)(2-10)称为单元的称为单元的应变矩阵应变矩阵。3 3）节点位移）节点位移 e e和单元应力分量和单元应力分量 的关系的关系 由弹性力学已知，平面应力情况下，由弹性

14、力学已知，平面应力情况下， 应力与应变之间的关系应力与应变之间的关系可表达为：可表达为：22()()(0)11xyxEExymmgmsee=+-22()()( 0)11xyEEyxymmgmsee-=+(0)(0)2( 1)xyxyExytgmee+=+1010210210 xyxyxyxyEemmssstmem g轾犏犏犏犏犏犏犏犏犏犏轾犏犏犏犏犏犏犏犏轾犏犏犏犏犏犏犏犏臌犏臌=-=1010210210 xyxyxyxyEemmssstmem g轾犏犏犏犏犏犏犏犏犏犏轾犏犏犏犏犏犏犏犏轾犏犏犏犏犏犏犏犏臌犏臌=-=22()()(0)11xyxEExymmgmsee=+-22()()( 0)1

15、1xyEEyxymmgmsee-=+(0)(0)2( 1)xyxyExytgmee+=+10 10211002EDmmmm轾犏犏犏犏犏犏犏犏犏犏犏犏臌=-(2-13a)101( 1) 101( 1)( 1 2 )1 2002( 1)EDmmmmmmmmm轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌-=-+-( 2 13 )b- iiiUFV轾犏=犏犏臌jjjUFV轾犏=犏犏臌mmmUFV轾犏=犏犏臌0 xyijmiUmUjUiVjVmViiijejjmmmUVFUFFVFUV轾犏犏犏轾犏犏犏犏=犏犏犏犏犏犏臌犏犏犏臌(215)-111213141516212223242526313233343536414

16、24344454651525354555661iiijjmmiiijjmmjiijjmmjiijjmmmiijjmmmiUk uk vk uk vk uk vVk uk vk uk vk uk vUk uk vk uk vk uk vVk uk vk uk vk uk vUk uk vk uk vk uk vVk uk6263646566ijjmmvk uk vk uk v111213141516212223242526313233343536414243444546555152535456616263646566kkkkkkkkkkkkkkkkkkkkkkkkUiViUjVjUmkkkkkk

17、kkkkkkVm轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏轾犏犏犏犏犏犏犏犏犏犏犏犏犏=uiviujvjumvm轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌( 2 16 )a-111213141516212223242526313233343536414243444546555152535456616263646566ekkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkk轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌轾犏臌= 称为单元刚度矩阵单元刚度矩阵。11 112 123 124 135 136 1UkVkUkVkUkVk轾轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏

18、犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌臌=12311k41k31k21k61k51kijjikk 1 11 21 31 41 51 62 22 32 42 52 63 33 43 53 64 44 54 65 55 66 6ekkkkkkkkkkkkkkkkkkkkkk轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=对称 (a) 实际力系 (b) 虚设位移弹性体虚功原理的应用*uiviujvjumvmd轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌*=* *xyxyuxvyvuxyeeeg*轾犏犏轾犏犏犏犏= 犏犏犏犏犏犏抖犏犏臌+犏抖犏臌 *TeiiiijjjjmmmmuUvVu

19、UvVuUvVF( )a( )b*Txxyyxyxytdxdytdxdy ( )c *TTeFtdxdy( )d *()TeeTFBD Btdxdy *TTTeeFBD Btdxdy( ) e *TTTeeFBDBtdxdy( )f TeeFBDB tdxdy eeeFk TekBDB tdxdy( )g( )h(2 18 )adxdy 6 623 63 36 3000001001000004 ( 1)10 0020eiiiijjkjjmmmmbccbbbbmijbcEtcccmijcbcbcbcbmmi ijjbccbmmmm轾犏犏犏轾犏犏犏犏犏犏犏犏犏犏轾犏犏犏犏臌犏犏犏犏犏犏犏犏犏犏臌犏

20、犏犏犏臌轾犏犏=犏D -犏-犏臌(2 18 )b0 11,1 011 0 1,0 0 00 0 0,0 1 1ijmijmjmijmimijmijbyycxxbyycxxbyycxx 0001 0002bbbmijBcccmijcbcbcbmmiijj轾犏犏犏犏犏犏犏犏犏犏犏犏臌=V101000010001110110轾犏犏犏犏犏犏犏犏犏臌-=-101. 60. 40 100. 41. 60211000. 6002EDmmmm轾犏轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌犏犏臌=-1. 60. 40101000 0. 41. 60010001000. 6110110SD B轾轾犏犏犏犏犏犏犏犏犏

21、犏犏犏犏犏犏犏犏犏臌臌-=-1. 60. 41. 6000. 40. 41. 60. 4001. 60. 60. 600. 60. 60轾犏犏犏犏犏犏犏犏犏臌-= -101011100000001010轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌-=1. 60. 41. 6000. 40. 41. 60. 4001. 60. 5 0. 40. 60. 600. 60. 60轾犏犏犏犏犏犏犏犏犏臌-创-0.440.200.320.120.120.080.440.080.120.120.320.32000.320.120.1200.1200.32轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏

22、犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌-=图( 2 9)-0. 1200. 120. 1200. 120. 320. 080. 320. 0800. 440. 200. 320. 120. 440. 080. 120. 3200. 12 ek轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌-=-答：,ijmmjijmijmjm iimjmijmimijjimijmijax yx ybyycxxax yxybyycxxaxyx ybyycxx kkkiiijimekkkkjijjjmkkkm mm im j轾轾轾轾犏犏犏犏犏犏犏犏臌臌犏犏臌犏犏轾轾轾犏犏犏犏犏犏犏犏犏犏犏犏臌臌臌犏犏轾犏轾轾犏犏犏犏犏犏臌

23、犏臌犏犏臌臌=(219)-,2 224 (1)E tkrsm=-V2 211221122b bc cb cc brsr sr srsc bb cc cb brsr sr srsmmmmmm轾犏犏犏犏犏犏犏臌-+-+(), , ;, ,ri jm si jm=(220 )a-( 1)4( 1)( 12 )Etkrsmmm-=+-V1 21 22( 1)12( 1)1 21 212( 1)2( 1)r sr sr sr sr sr sr sr sb bc cb cc bc bb cc cb bmmmmmmmmmmmm轾-犏+犏-犏-+犏-犏臌(220 )b- 通常，我们总将集中力的作用点取为节点，

24、不需要移通常，我们总将集中力的作用点取为节点，不需要移置。置。 因此，下面只讨论：因此，下面只讨论： 体力和面力的移置体力和面力的移置mbbj13bcbi图(211)-eiYejYemYeiY图(211)-mjbieiY113eiWY 图(211)-01eiWX 0eiX 图(212)- TeeeeeeeiijjmmRXYXYXY0101013TW2ejqlR 2eiqlR 图(213)-22eiejqlRqlRij2ejqlR 2ejqlR 图(213)-0emR 01 1 022TTelqlRqq 0 02TeixiyjxjylRqqqqiq 202 1 066TTeiiiqllRqq 2

25、20 06TeixiyjxjylRqqqq3eiilRq6ejilRqiqiqjq26eiijlRqq26ejjilRqq0emR 2206TTeeeeijmijjilRRRRqqqq TeeeeeeeiijjmmRXYXYXY或写成分量形式：或写成分量形式：(2) (2) (2) (2) 0 06Tixjxiyjyjxixjyiylqqqqqqqq 式中式中 分别是分别是 在在 x , y x , y 方向上的分量，其方向与方向上的分量，其方向与 x , y x , y 轴正向一致为正，反轴正向一致为正，反之为负。之为负。,ixiyjxjyqqqq,ijqqiqjq(2)6eiijlRqq(

26、2)6ejjilRqq2.3 2.3 整体分析整体分析1.1.基本方程基本方程 总刚度矩阵总刚度矩阵KK的形成的形成(1)(1)节点力的组合节点力的组合 以上我们分析了一个单元的情况，现在进而研究单元以上我们分析了一个单元的情况，现在进而研究单元的组合。为了说明问题，今选用一个包含的组合。为了说明问题，今选用一个包含9 9个节点个节点8 8个单元个单元的平面问题来分析。如图的平面问题来分析。如图2-162-16所示，除所示，除1 1、3 3、7 7、9 9四个节四个节点外，其余五个节点均联接着四个单元。点外，其余五个节点均联接着四个单元。 对于联结着对于联结着 n n 个单元的节点，一个节点的

27、位移当个单元的节点，一个节点的位移当然涉及到然涉及到 n n 个单元，与节点位移相应的节点力将是个单元，与节点位移相应的节点力将是 n n 个单元的综合效应。个单元的综合效应。 如节点如节点5 5的位移涉及到的位移涉及到(2)(2)、(3)(3)、(6)(6)、(7)(7)单元。与单元。与节点节点5 5相应的节点力将与上述四个单元有关。相应的节点力将与上述四个单元有关。 首先，我们按首先，我们按(2-16a)(2-16a)式逐个建立单元节点力和节点式逐个建立单元节点力和节点位移的关系。下面，选写其中位移的关系。下面，选写其中(2)(2)、(3)(3)、(6)(6)、(7)(7)四个四个单元。单

28、元。对于单元对于单元(2)(2)，节点编号，节点编号2 2、5 5 、4 493949399, 1334393, 10373824344494, 1047482510, 310, 410, 910, 1010, 710, 847374797, 10777848384898, 108787809985kkkkkkUkkkkkkUkkkkkkVVkkkkkkUkkkkkkVkkkkkk轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=225544uvuvuv轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌犏犏犏犏犏犏犏犏犏犏对于单元对于单元(3)(3)，节点编号，节点编号2 2，6 6，5 533343, 1

29、13, 12393, 10243444, 114, 12494, 10211, 311, 411, 1111, 1211, 911, 106612, 312, 412, 1112, 1212, 912, 10593949, 119, 12999, 10510kkkkkkUkkkkkkVkkkkkkUVkkkkkkUkkkkkkVk轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=226655, 310, 410, 1110, 1210, 910, 10uvuvuvkkkkk轾犏犏犏犏轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌犏犏犏犏犏臌对于单元对于单元(6),(6)

30、,节点编号节点编号4 4，5 5，8 87778797, 107, 157, 1648788898, 108, 158, 1649798999, 109, 159, 165510, 710, 810, 910, 1010, 1510, 16815, 715, 815, 915, 1015, 1515, 16816kkkkkkUkkkkkkVkkkkkkUVkkkkkkUkkkkkkVk轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=445588, 716, 816, 916, 1016, 1516, 16uvuvuvkkkkk轾犏犏犏犏轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏

31、犏犏犏犏犏犏臌犏犏犏犏犏臌对于单元对于单元(7),(7),节点编号节点编号5 5，6 6，8 8999, 109, 119, 129, 159, 16510, 910, 1010, 1110, 1210, 1510, 16511, 911, 1011, 1111, 1211, 1511, 166612, 912, 1012, 1112, 1212, 1512, 16815, 9158kkkkkkUkkkkkkVkkkkkkUVkkkkkkUkkV轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=55668, 1015, 1115, 1215, 1515, 16816, 916, 1016, 11

32、16, 1216, 1516, 16uvuvukkkkvkkkkkk轾犏犏犏犏轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌犏犏犏犏犏臌 同样可写出同样可写出(1)(1)、(4)(4)、(5)(5)、(8)(8)单元的节点力和节单元的节点力和节点位移的关系。点位移的关系。( (学习者可自己完成）。学习者可自己完成）。 对于第对于第5 5个节点，其个节点，其X X方向的节点力方向的节点力 即即: : U5=U5(2)+ U5(3)+ U5(6)+ U5(7) 4551eUU=999793949,10989,119993949,129,109,99,1597989,109

33、,169(2)(2)(2)(2)(2)(2)(2)5225544(3)(3)(3)(3)(3)(3)(3)5226655(6)(6)(6)(6)(6)(6)(6)5445588(7)5UkukvkukvkukvUkukvkukvkukvUkukvkukvkukvUk9,119,1599,109,129,16(7)(7)(7)(7)(7)(7)556688ukvkukvkukv( )b( )a将上式代入将上式代入(a)(a)939394949797989899999,9999,109,109,109,109,119,119,129,129,15(2)(3)(2)(3)522(2)(6)(2)(6

34、)44(2)(3)(6)(7)(2)(3)(6)(7)55(3)(7)(3)(7)66(6()()()()()()()()(Ukkukkvkkukkvkkkkukkkkvkkukkvk9,159,169,16)(7)(6)(7)889329429749849959,1059,1169,1269,1589,168)()kukkvK uKvKuKvKuKvKuKvKuKv式中式中: :939394949797989899999,9999,109,109,109,109,119,119,129,12(2)(3)93(2)(3)94(2)(6)97(2)(6)98(2)(3)(6)(7)99(2)(3

35、)(6)(7)9,10(3)(7)9,11(3)(9,12()()()()()()()(KkkKkkKkkKkkKkkkkKkkkkKkkKkk9,159,159,169,167)(6)(7)9,15(6)(7)9,16)()()KkkKkk所以式中的所以式中的eijijKk1111213141718121222324272823132333435363738393,103,113,1223344556677889900000000000000000000000000000UKKKKKKVKKKKKKUKKKKKKKKKKKKVUVUVUVUVUVUVUV4142434445464748494

36、,104,114,12535455565,115,12636465666,116,12717273747778797,107,137,147,157,16818283848788898,1000000000000000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK8,138,148,158,1693949798999,109,119,129,159,1610,310,410,710,810,910,1010,1110,1210,1510,1611,311,411,511,611,911,1011,

37、1111,1211,1511,1611,17000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK11,1812,312,412,512,612,912,1012,1112,1212,1512,1612,1712,1813,713,813,1313,1413,1513,1614,714,814,1314,1414,1514,1615,715,815,915,1015000000000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKK,1115,1215,1315,141

38、5,1515,1615,1715,1816,716,816,916,1016,1116,1216,1316,1416, 1516,1616,1716,1817,1117,1217,1517,1617,1717,1818,1118,1218,151000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKK112233445566778898,1618,1718,189uvuvuvuvuvuvuvuvuKKv . . . .111213141, 21. . . .212223242, 22. . . .313233343, 2341424

39、5212KKKKKnFKKKKKnFKKKKKFnKKFFFnFn轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌=-MM. . . .43444, 2. . . .515253545, 2. . . . .21, 121, 2. . . .2 , 12 , 22 , 32 , 42 , 2KKKnKKKKKnKKnnnKKKKKnnnnn n轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌-MM12.2nddd轾犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏犏臌上式可简写成： F F=K K式中K K称为总刚度矩阵。它是一个对

40、称正定阵(其对角线上各元素为正值，且有元素Kij=Kji)。总刚度矩阵由单元刚度矩阵累加而成，每个单元刚度矩阵对总刚度矩阵都有一定贡献。 K K=k ke (221)- 假定结构离散化后有 n 个节点，每个节点有两个方程。因此总刚度矩阵为 (2n x 2n)的矩阵。 可将单元刚度矩阵用补零的方法由6x6扩大到(2n x 2n) 的方阵(图中虚点上的元素均为0） 如图所示，则单元刚度矩阵中各元素在总刚度矩阵中的位置即可确定。 例如将第3单元刚度矩阵中的元素填入总刚度矩阵（亦即将该单元刚度矩阵用补零的方法扩大成总刚度矩阵）33343, 113, 12393, 1043444, 114, 12494

41、, 1011, 311, 411, 1111, 1211, 911, 1012, 312, 412, 1112, 1212, 912, 1093949, 119, 12999, 1010, 310, 410, 1110, 1210, 910, 10kkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkkk轾犏犏犏犏臌21,2121,221,2121,221,2121,22 ,212 ,22 ,212 ,22 ,212 ,221,2121,221,2121,221,2121,22 ,212 ,212 ,21iiiiijijimimiiiiijijimimjijijjjjjmjmji

42、jijjkkkkkkkkkkkkkkkkkkkkkk2 ,22 ,212 ,221,2121,221,2121,221,2121,22,212,22,212,22,212,2jjjmjmmimimjmjmmmmmimimjmjmmmmkkkkkkkkkkkkkkijijmm 在计算程序编制中，我们的做法是先都按局部编码，使用同一过程体算出各单元的刚度矩阵；然后转换成总体编码，最后将相同编码的元素合并成总刚度矩阵中的元素。转换过程示意如下：111213212223313233eeeeeeeeekkkkkkkkkeeeiiijimeeejijjjmeeemimjmmkkkkkkkkkiiijimj

43、ijjjmmimjmmikkkjkkkmkkk局部编码局部编码编码转换编码转换总体编码总体编码相相同同编编码码合合并并形成总刚度矩阵系数形成总刚度矩阵系数ijm(2)节点载荷组合 当进行单元组合时，除了考虑节点力的组合外，同时还应进行节点载荷的组合。 设结构上的载荷（如体力、面力、集中力等）均已移置到节点上，则单元节点载荷列阵为：YmeXmeYieXieXjeYje 图2-17xy eieeiieejjeejmememXYRXRRYRXY对于联结着两个以上单元的节点，把相同方向上的载荷迭加起来，显然有：eiieXXeiieYY这就是该节点载荷在x和y方向的两个分量。（ 表示对环绕节点i的单元求

44、和）若各节点上的Xi和Yi（i=1，2，3n）均已经求出，并按节点编码的顺序排列起来，就得到弹性体总的节点载荷列阵： 12321122.TTnTnnRRRRRXYXYXY(222)-(3)平衡方程 在求得了节点外力矩阵以后，我们就可以写出位移法中位移分量必须满足的平衡条件。在有限元法中，也就是节点位移必须满足的节点平衡条件。 根据公式（2-21）和（2-22），表示所有节点内力与外力平衡的数学表达式为： FKR通常写成： 等式右端为总节点载荷列阵，当弹性体上的外载荷确定时，它是已知的； 总刚度矩阵K由单元刚度矩阵ke集合而成。因此，式（2-23a）是一个以 为未知量，以K为系数的线性代数方程组

45、，这是有限元法的基本方程。 KR(223 )a-其矩阵展开式为：11121,21121222,2112 ,12 ,22 ,2.nnnnnnnnkkkuxkkkvykkkvy 顺便提一下，为了编程方便，对各矩阵元素的下标，均按其所在位置标定：11121,21121222,2222 ,12 ,22 ,222.nnnnnnnnkkkRkkkRkkkR即 为节点1，x方向的位移u1， 为节点1，y 方向位移v1 ,.余此类推。 R1为节点1，x方向的节点载荷x1， R2为节点1，y方向的节点载荷y1,.余此类推。1221112211.TTnnTTnnuvvRRRxyy式中122.关于总刚度矩阵的性质和

46、常用的存贮方法： 在有限元法中，结构总刚度矩阵的性质和常用的存贮方法： 在有限元法中，结构总刚度矩阵的阶数为节点总数和自由度数的乘积。 如平面问题，自由度为2，若有 n 个节点，则总刚度矩阵为:(2n x 2n) 阶。 若一个具有200个节点的小型平面问题，其总刚度矩阵的元素就有400 x400=160000个，为一般中小型电子计算机的内存容量所不允许的。因此，如何缩小总刚度矩阵所需的存贮单元，是有限元法程序编制中一个需要考虑的突出问题。我们可以根据总刚度矩阵的某些特性，寻求节省存贮量的途径。1)总刚度矩阵是对称阵。 利用对称性，我们只需要存贮总刚度矩阵的上三角形部分（或下三角形部分）中的元素

47、。2)当合理编排节点编码时，总刚度矩阵可呈带状的稀疏阵。其中有不少元素为零，而非零元素对称地分布于主对角线的两旁，形成一带状阵。下面介绍两种常见的压缩存贮方法：(1)半宽带存放： 仍以图2-16为例，分析与该图相应的总刚度矩阵，其元素的分布表示于图2-18中1111213141718121222324272823132333435363738393,103,113,1223344556677889900000000000000000000000000000UKKKKKKVKKKKKKUKKKKKKKKKKKKVUVUVUVUVUVUVUV4142434445464748494,104,114,

48、12535455565,115,12636465666,116,12717273747778797,107,137,147,157,16818283848788898,1000000000000000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK8,138,148,158,1693949798999,109,119,129,159,1610,310,410,710,810,910,1010,1110,1210,1510,1611,311,411,511,611,911,1011,1111,1211

49、,1511,1611,17000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK11,1812,312,412,512,612,912,1012,1112,1212,1512,1612,1712,1813,713,813,1313,1413,1513,1614,714,814,1314,1414,1514,1615,715,815,915,1015000000000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKK,1115,1215,1315,1415,1515,16

50、15,1715,1816,716,816,916,1016,1116,1216,1316,1416,1516,1616,1716,1817,1117,1217,1517,1617,1717,1818,1118,1218,151000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKK12345678910111213141516178,1618,1718,1818KK 当把对角线向右平移到第十个元素时，剩下的就全部是零元素了。 这是因为节点编排时，就一个单元而言，节点编码最大相差五个号码，而每个节点又有两个方向的位移所致。 我们把自对

51、角线开始向右侧平移至最后一根带有非零元素的斜线为止，两线之间一行内包含的元素个数（包括两线之间的零元素）称为最大半带宽。 而刚度矩阵中每一行的半带宽均取决于一个单元中节点编号差和节点的自由度（即有几个方向的位移）。 如以NBD表示最大半带宽，以d来表示各单元中最大节点差值。（对于三角形单元，也就是相邻节点编号的最大差值），则对具有两个自由度的平面问题，最大半带宽的计算公式为：NBD=（d+1）x2 对于以上所述的这样一个对称的，具有NBD半带宽的总刚度矩阵，由于它在半带宽以外有许多零元素，而对角线一侧的元素又和另一侧元素对称相同，因此为了节省计算机内存，可以采用长方形压缩存储法只将半带宽内的元

52、素存放起来，如图（2-18b）所示，这种存贮总刚度矩阵的方法也称为“半带宽存贮”或“等带宽存贮”。 图2-181112131417182223242728333435363738393,103,113,124445464748494,104,114,1255565,115,12666,116,127778797,100000000000000000000000000000000000000000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK07,137,147,157,1688898,108,138,148,158,16999,

53、109,119,129,159,1610,1010,1110,1210,1510,1611 ,1111 ,1211 ,1511 ,1611 ,1711 ,1812,1212,1512,1612,1712,1813,00000000000000000000KKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKKK1313,1413,1513,1614,1414,1514,1615,1515,1615,1715,1816,1616,1716,1817,1717,1818,180000KKKKKKKKKKKKKKKK KSK (a) (b) 如节点总数为n，而每一节点有两个方程，则方程总数

54、为：NEQ=2n 我们可以用一个两维数组SK1：NEQ;1：NBD即图（2-18b）来存放总刚度矩阵中必须保存的元素。总刚度矩阵K和长方形压缩存贮的数组SK之间有如下的对应关系：矩阵K数组SK对角线第一列r行r行r行s列元素r行(s-r+1)列元素由K形成SK时，元素的行号相同。新的列号等于原先的列号减去行号加1，即:新列号=原列号-行号+1 为了节约内存，我们力求减小半带宽NBD，这就是要求在编排节点号码时，应使相邻节点的节点差值d尽可能小。 SK右下角三角块中的元素不和K中的任何元素对应，那些元素应是零。(224)- 考虑到总刚度矩阵K中各行的半带宽并不相同，有时，由于结构几何形状等原因，

55、某些行的半带宽特别大，而其它行又较小，这种情况下如以NBD为半带宽（即等带宽）存贮，就可能把许多零元素也存贮起来，这对节省计算机存储量来说是不利的。(2)一维压缩存贮法一维压缩存贮法是将总刚度矩阵K的下三角形中每一行从第一个非零元素开始，逐行存放入一维数组K1：S中，S是元素的总个数）举例如下，设有一对称正定阵：534603007200008009008对称在一维数组K1：13中依次存放为： 5，3，4，6，0，3，7，2，8，9，0，0，8共13个数。 但是，仅使用这样一个一维数组并不能将元素在K中的位置确定下来。 为此，还须将主对角线上的元素在一维压缩存贮中的序号用另一个一维数组N1：2n

56、0存放起来（n0是节点总数）。如对于上述矩阵数组N1：6中存放的是：1，3，6，8，9，13 即指出，在一维数组K1：13中，第1，第3，第6.个元素是对角元，显然，这两个数组完全确定了各元素在K中的位置。2.4 2.4 方程求解方程求解1.引入位移约束条件 上一节，我们建立了有限元法的基本方程式： 有了这个方程还不能立即求解节点位移，因为到现在为止，我们还没有考虑到弹性体的几何边界条件，即边界位移的约束条件. KR 很明显，如弹性体的边界没有位移约束，则在外载荷的作用下，它将有产生刚体运动的可能性，反映在基本方程上，其系数矩阵K将是一个奇异阵（对应的行列式的值等于零）。逆矩阵不存在，方程将具

57、有不定解。 在进行应力分析时，为了使解具有唯一性（即排除刚体运动），必须根据弹性体具体的边界位移约束条件，对基本方程加以处理，方能求解。1）引入约束条件的原因（1）结构实际上可能存在若干约束条件。它们应该加以考虑，否则计算的结果将与实际不符。如一端固定，一端铰支的静不定梁，如图2-19所示。 图2-19节点1，2，3及15既不允许有x方向的位移也不允许有y方向的位移。(2)有些约束是由于考虑到结构和载荷的对称性（或反对称性）可取其中一部分作为计算对象而附加的。如图2-20所示一对角受压的方形薄板。 图2-20 由于结构和载荷的均对称，可以只计算薄板的四分之一。 图2-20b因为变形对称于对角线

58、，水平对角线nn上不可能有y方向的位移，所以附加的约束条件应是： （节点4，5，6为垂直的可动铰支座）；4560vvv 垂直对角线mm不可能有x方向的位移，所以附加的约束条件应是: （节点1，2，4为水平的可动铰支座），薄板的中心点O，任何方向的位移均为零，故节点4处为固定铰支座。1240uuu2）引入约束条件的处理： （1）零位移约束条件的处理 举例说明，对于在剖分后有n个节点的弹性体，假定第n个节点处有约束，其位移为零， 即： un=0 vn=0此外，还已知第（n-1）个节点沿y轴方向有约束，其y方向位移为零， 即： vn-1=0 （这对平面问题来说，是不产生刚体运动的最低限度的边界位移约

59、束条件）。则应对基本方程作这样的处理： 把刚度矩阵K的最后三行和最后三列划去，得到一个（2n-3）阶的方阵 ，称为总刚度矩阵的缩聚或降阶。研究表明，降阶以后的总刚度矩阵K是一个非奇异的对称正定矩阵。 相应地，分别将节点位移列阵和节点载荷列阵R的最后三行划去，得到（2n-3）维列阵 和 。 这样，经过零位移边界条件处理后的基本方程就成为: KRK R(225)- 这样得到的结果，不但满足平衡条件和相容条件，而且满足全部边界条件，按上面的假定及处理方法，方程的具体形式是：11121,231121222,232223,123,223,232323.nnnnnnnnkkkRkkkRkkkR式中Kij是

60、第i行，j列的元素。 如果零位移的节点编号不在最后而在中间，也可以用同样的方法划去想相应的行和列，而得到降阶后的式（2-25）： KR 但是在计算机的计算程序中，我们采用的方法是使总刚度矩阵K中要划去的行和列除对角线元素充成1以外，其余的元素均充成零。 如此得到的矩阵不降阶，仍为2n阶方阵，记为 ，并且也是一个非奇异的对称正定矩阵。 与此同时，将节点载荷列阵中相应的元素也充成零，如此得到的2n维列阵记为 ,因而基本方程就变为： R KR（2-26）K 式中 为包含所有节点位移的列阵。显然，在求解式（2-26）时，原来节点位移为零的仍保持为零。式（2-26）与式（2-25）实际上是等价的。 按前