空间问题与板壳单元

1、 在用有限元法求解轴对称问题时，采用的单元一般为整圆环，如图所示。它们和子午面rz面相交的截面，可以是直边三角形或矩形，也可以是任意四边形、曲边三角形、曲边四边形等。各个单元之间以圆环形的铰相互连接，而每一个铰与子午面rz面的交点就称为节点。如图上的i、j、m等。所有单元将在子午面rz面上形成有限元网格，与在平面问题中形成的网格一样。因为在轴对称问题中采用的单元是一个整圆环，所以在计算单元的体积时要注意到这一点。下面以三角形截面单元为例，说明如何求解轴对称问题。 如图所示子午面rz面上的一个三角形单元，设单元上任一点的径向位移（沿r向位移）分量为 ，轴向位移（沿z向位移）为量为 。由于这两个位

2、移分量仅是r和z的函数，故令uw123456urzwrz 与平面问题一样，可将位移用形函数及节点位移表示为iijjmmiijjmmuN uN uN uwN wN wN w即eeeeijmu= NNNv uNIII式中： 二阶单位矩阵； 、 、 形函数矩阵，如下IiNjNmN 1()2iiiiNabrc z （i、j、m 轮换）1211iijjmmrzrzrz jjimmrzarz11jimzbz 11jimrcr式中： 由弹性力学知，轴对称问题中除平面内的应变分量 外，还有环向应变 。rzzr、 、因此其几何方程为rezrzururwzwurz 将位移代入上式，得eeeeeeeijmBBBB

3、引入iiiiabrc zhr则00001020iiiieiiiiiiiNrbNhrcNzcbNNzrB 应变矩阵 中的元素不全是常量，因此单元内的应变也不是常量，这是因为轴对称问题中采用的单元是圆环，径向的位移 必引起环向应变 ，而此应变的大小与点的位置有关。另外，由于 中含有 项，使单元的应变、应力及单元刚度矩阵的计算比平面问题复杂得多。1 reBueB根据弹性力学可知 ，对于轴对称问题，有ee D Terzrz 1011101111+1 210111 20002 1ED = 根据虚功原理或用最小位能原理，可以和平面问题一样推得其单元刚度矩阵的表达式为 ，在轴对称问题中，由于单元是一圆环，上

4、述积分式中的微分体积 可取为微分圆环的体积，即 ，故单元刚度矩阵为eeeTdVKBDBdV2dVrdrdz2eeTerdrdzKBDB 与平面问题一样，单元刚度矩阵 是一个 阶的方阵，矩阵 可分成三块，故 也可分成 个子矩阵，每个子矩阵为 阶的方阵，其表达式为3 32 2eBeK2eeTeststrdrdzKBDB( , ,)s ti j m 因为矩阵 与坐标有关，且坐标r处于分母上，因此积分不像平面问题中那么简单，常采用三种办法进行计算：1、显式积分；2、数值积分；3、简单的近似积分。一般采用第3种简单的积分，它不仅在程序上简单，而且还回避了节点在极轴上时带来的奇异问题。实践证明，在精度方面

5、它并不比精确的积分公式法差。 具体做法是取节点坐标平均值，即单元中心坐标eB33ijmijmrrrrzzzz并取iiiiabrc zhr在式刚度矩阵中以 代替 ，以 代替 可得ihihrir00102iieiiiibhccbB 这样就使得单元刚度矩阵中的被积函数化为常数，然后积分即可求得 ，具体表达式这里不再给出。estK 与平面问题一样，无论使用虚功原理或最小位能原理可以得到相同的载荷移置公式，其形式与平面问题相似。 1. 集中载荷 轴对称问题中的集中载荷实质上是沿着圆周线作用，均匀分布的一圈力。在子午面单元上任一点 处作用集中力 ，其中 、 为单位弧长上分布力的合力，其对应的等效节点载荷为

6、,CCrzTcrczcFFF =rcPrzP2eeTcccrFNF 2. 面力 设单元的分布面力为 ，则其对应的等效节点载荷为Trzqqq =2eeTqSrdSFNq 3. 体积力 设单元是分布体力为 ，其中 、 为单位体积的体力分量，其对应的等效节点载荷为TrzGGG =rGzG2eeTGArdrdzFN G 例如，在体力为自重的情况下，有 、 其中 为容重，于是有0rG zG 0000200020000202026TijmeGijmATijmATijmjmimijNNNrdrdzNNN-NNNrdrdzArrrrrrrrr F 1. 位移模式 如图所示的四面体单元，单元节点的编码为i，j，

7、m，n。每个节点具有三个位移分量，单元节点的位移列阵可表示为iTjeiiijjjmmmnnnmnuvwuvwuvwuvw 单元的位移模式采用线性多项式的形式，如下123456789101112uxyzvxyzwxyz 将单元上四个节点坐标代入上式，求解相应参数，则单元内任一点的位移可用节点位移和形函数表示，如下iijjmmnniijjmmnniijjmmnnuN uN uN uN uvN vN vN vN vwN wN wN wN w用矩阵表示为000000000000000000000000iijmnjeijmnmijmnneeijmnuNNNNvNNNNwNNNNNNNN uIIIIN 式

8、中：I I三阶单位阵； 、 、 、 四面体单元的形函数，如下iNjNmNnN1()61()61()61()6iiiiijjjjjmmmmmnnnnnNab xc yd zVNab xc yd zVNab xc yd zVNab xc yd zVjjjimmmnnnxyzaxyzxyz111jjimmnnyzbyzyz111jjimmnnxzcxzxz111jjimmnnxydxyxy1111iiijjjmmmnnnxyzxyzVxyzxyz式中：（i、j、m、n轮换） 式中：V四面体的体积。 2. 单元应变 根据几何方程，得eeeeijmnBBBBB 00000016000iiiiiiiiii

9、bcdVcbdcdbB（i、j、m、n轮换） 由于单元中的应变是常量，所以四面体单元是常应变单元。3. 单元应力由物理方程得单元的应力为eeeeeijmnDDBSSSSS式中： 空间问题的弹性矩阵，如下D10001110001110001111 211 2000002 11 2000002 11 2000002 1ED =式中： 四面体单元的应力矩阵，表示为S11111132222220600iiiiiiiiiiiiiiiibAcAdAbcAdAbAcdAA cA bVA dA cA dA bS （i、j、m、n轮换）式中： ， ， 。 11A21 22(1)A3(1)(1)(1 2 )EA

10、由于单元中的应力是常量，所以四面体单元是常应力单元。 4. 单元刚度矩阵 由虚功原理，按照平面问题的类似推导，可得空间问题的单元刚度矩阵为eeTeeTedxdydzVK =BDB= BDB单元刚度矩阵可表示为分块形式，如下iiijiminjijjjmjnemimjmmmnninjnmnnKKKKKKKKK =KKKKKKKK式中任一子块由下式计算147258369(1)36 (1+ )(12 )eTeststKKKEV =KKKVKKKK= BDB, ,s ti j m n122123124125261271281292()()()rsrsrsrsrsrsrsrsrsrsrsrsrsrsrsr

11、srsrsrsrsrsKb bA c cd dKAc bA b cKAd bA b dKAb cA c bKc cA d db bKAd cA c dKAb dA d dKAc dA d cKd dA b bc c式中：5. 等效节点载荷（1）集中载荷若集中载荷为 ，则等效节点载荷可表示为TccxcyczFFFFeeTccFNFeeTqdAFNq（3）体力若分布体力为 ，则等效节点载荷可表示为TxyzGGGG =eeTGdxdydzFNG（2）面力 对于单元的某一边界上的分布面力 ，则等效节点载荷可表示为qTxyzqqq 单 元 自 重 的 等 效 节 点 载 荷 可 由 上 式 计 算 ，

12、因为 ，所以等效节点载荷表示为0,xyzGGG 11110 00 00 00 0444411110 00 00 00 04444eeTGTTdxdydzVVVVVWWWWFNG 由上式可知，单元等效节点力可看作将单元自重W均匀移置到每个节点上。 1. 四面体单元 由上节可知，常应变四面体单元中的各点应力为常量，然而实际工程结构中的各点应力是随着坐标变化的，为了反映真实情况，可以假定高次位移模式，相应的单元称为高阶单元。例如，在四节点四面体单元的基础上增加节点数，形成四面体单元族，如图所示。图中的4节点、10节点和20节点四面体单元，与左侧的四面体比较，可以看出节点的数目与多项式的项数完全一致，

13、因此相应的位移模式可取为一次、二次和三次完全多项式。不难验证10节点二次四面体单元的应力、应变随着坐标线性变化，20节点三次四面体的应力、应变是坐标的二次函数，对于弹性力学问题其计算精度大大高于4节点常应变四面体单元。高次四面体单元 2. 六面体单元 在空间问题中经常使用六面体单元，常见的六面体单元有8节点、20 节点和32节点的立方体单元，如图5.5所示。8节点六面体单元精度明显高于同样节点数的四面体单元，20节点的六面体单元精度更高一些，但六面体形状难以适应工程结构的复杂外形，实际中使用不多。 空间六面体单元 板壳结构是工程结构中常见的一种结构，如图所示，中面为平面的称为板，中面为曲面的称

14、为壳。本节以四节点矩形弯曲薄板单元为例，说明板壳单元的基本理论。板壳结构( , )ww x ywuzx wvzy 1.应变在板的小挠度理论中，板单元的的位移表示为 板单元的应变表示为222222xeyxywuxxvwzzyyuvwyxx y 式中：222222Twwwxyx y 2.应力应力应变关系10110002(1)xxyyxyxyE 则应力为xyxyzDD式中： 弹性矩阵。D3.内力分量由弹性理论得22222312xhhhhyxyfMMz dzzdzMhMDDD式中： 薄板弯曲问题的弹性矩阵，如下fD32101012(1)(1)002fEhD 则薄板应力可由内力表示，如下312zhM 四

15、节点矩形薄板单元如图所示，节点位移为iiixiiyiiwwwywx （i、j、m、n轮换） 矩形薄板有4个节点，整个单元的位移列阵为TeTTTTijmn矩形薄板单元的位移模式可取为22312345672233389101112waa xa ya xa xya ya xa x ya xya ya x ya xyeew N 将薄板单元的坐标和节点位移代入上式，求解得到 ，并代回上式，则上式表示为112aa式中：ijmnNNNNN式中：iixiyiNNNN式中：222211(1)(2)811(1)(1)811(1)(1)8iiiiixiiiiyiiiiNNbNa （i、j、m、n轮换）式中：,iii

16、ixyxyaabb将位移代入应变表达式，可得到单元应变eeeeijmnBBBBB 式中：22222222222211122iieiiiiiaxzzybabx y NNNNBNN （i、j、m、n轮换）化简后得到单元应变矩阵，如下22223(1)01 3(1)3(1)1(1 3)04(334)(321)(321)iiiiieiiiiiiiiiiiibbazaaabbbaB则单元应力可由物理方程求得，如下eeee= DDS 单元刚度矩阵用分块形式表示为 将应变矩阵 、弹性矩阵 代入单元刚度一般表达式 中，得单元刚度矩阵为eBDeeTedVKB DB112112heeTehabd d dz KBDB

17、iiijiminjijjjmjnemimjmmmnninjnmnnKKKKKKKKKKKKKKKKK111213212223313233staaaaaaaaaK, ,s ti j m n式中子矩阵为33的矩阵，如下式中：222200001122222200122222001322220212231514455323515532351553235155iijiijjjbbbaaHaaabaaaHbbbbbaHaaaaaaHbbb 02200002222322003122322200003322(1)(35)5(3)(3)15()()323515515()()2(1)(35)5(3)(3)iijijjjiijijaaHbbaH abbbaHaaaaH abbaHaa 式中：30026012(1)ijijDEhHDab 至于单元节点等效载荷的求法以及总体刚阵的组装均与平面问题是类似的，这里不再给出。把四节点矩形平面应力单元和该节讨论的四节点矩形弯曲薄板单元结合起来，可以构造同时考虑中面的伸展和板弯曲的四节点矩形板单元，进而分析壳体单元。 如图（a）表示一个边长为1m的正方形薄板，四边固支，中心受垂直板面的