6板壳有限元

板壳有限元分析概述弹性薄板矩形（R12）单元弹性薄板三角形（T9）单元薄板拟协调元与广义协调元8结点Hencky(Mindlin)板单元平面壳体单元与曲面壳元实际结构分析的若干问题薄板弯曲问题与平面应力问题不同，薄板弯曲问题是具有图示几何特征的结构在横向荷载作用下的分析。xyz平面应力：0yzxzz弹性薄板基本知识弹性薄板基本概念xyzuvw中面所谓薄板是指板厚h比板最小尺寸b在如下范围的平板1111~~1008085hb板面位移如图所示。当挠度w小于板厚h时，有如下克希霍夫(G.kirchhoff)假定成立：a)板中面是中性面，没有变形。b)中面法线变形后仍为挠曲面法线，长度不变。c)忽略应力z和应变z。平分厚度的平面称中面。弹性薄板基本知识由克希霍夫假定c)忽略应变z可推得w与z无关，;;(,)wwuzvzwwxyxy由此结果可得22222;;2xyxyx向曲率y向曲率扭率上述三个量完全确定板的变形，因此称它们组成的矩阵为形变矩阵，记作[]，也即T22222;;2位移只与挠度w有关由b)可知zx和yz等于零,再加上中面无变形,最终可得弹性薄板基本知识由此可得薄板应变矩阵为[]=z[]。薄板内力和总势能1)设平面应力弹性矩阵为[D]’，则薄板应力矩阵为[]=-z[D]’[]。xyzdydzxdydzxy/2/2-zdydzhhxy/2/2-zdydzhhx/2/2-zdxdzhhyx/2/2-zdxdzhhy2)薄板内力微元体如图所示。内力结果在下页弹性薄板基本知识xyzdydzxdydzxy/2/2-zdydzhhxy/2/2-zdydzhhx/2/2-zdxdzhhyx/2/2-zdxdzhhy/2/2-zdzdyhhx'xM/2/2-zdzdxhhy'yM/2/2-zdzdyhhxy'xyM/2/2-zdzdxhhyx'yxM扭矩弯矩由图可得弹性薄板基本知识由此可得薄板单位长度内力为Mx、My、Mxy=Myx(dx=dy=1)，依此顺序排列的列阵称内力矩阵，记作[M]。/2TT'/2TT11ddd2211χχdχd22hVAhAAUVDzADAMA将应力应变关系代入并对z进行积分，可得[M]=[D][]式中[D]=（h3/12）[D]’称作薄板的弹性矩阵。3)薄板的应变能弹性薄板基本知识4)薄板的总势能设薄板受z方向分布荷载q(x,y)作用，则线弹性薄板的总势能为TT1χχd-(,)(,)d21χd-(,)(,)d2AAAADAqxywxyAMAqxywxyA上式就是下面作有限元分析的理论依据。作业：试写出各向同性弹性体的[D]矩阵。试写出正交各向异性弹性体的[D]矩阵。对于各向同性体，弹性矩阵：弹性薄板基本知识1)各向同性板抗弯刚度：2)正交各向异性板抗弯刚度：式中：弹性薄板矩形（R12）单元薄板单元位移模式xyzw3y3x3Q1My1Mx112431)结点位移和结点力矩阵TixiyiidwTixiyiiFQMMTTTTT1234edddddTTTTT1234eFFFFF图中还给出了各结点位移和结点力的示意图。;xywwyx设局部编号1、2、3、4，x、y方向长度分别为2a、2b的矩形板单元如图所示。弹性薄板矩形（R12）单元薄板的形函数可以用广义坐标法，也可以用试凑法得到。由于单元自由度为12，因此可有12个广义坐标，位移模式可设为如下不完全四次多项式利用12个结点位移条件，由广义坐标法可建立形函数，显然十分麻烦。22312345672233389101112waaxayaxaxyayaxaxyaxyayaxyaxy？2)形函数的确定xyzw3y3x3Q1My1Mx11243龙驭球院士指出利用对称性可减少确定广义坐标的工作量，但是仍然不很方便。为此,介绍试凑法，首先引入自然坐标=x/a,=y/b。四次项的选取为了保证坐标的对称性，且曲率与扭率同阶次。弹性薄板矩形（R12）单元2-1)试凑形函数N1利用所有点N1的导数为零条件，P.122~123经式(2)~(12)的推导,可得))(1)(1(221edcbaN考虑到挠度是非完全四次式，为使自动满足它点为零N1(j)=0，可设)-2)(1)(1(221dN再由本点处位移的条件,可得d=-1/8,由此)/8-2)(1)(1(221N由形函数性质，对N1有:N1(1)=1；N1(j)=0，j=2,3,4N1对x，y的偏导数在各结点处均为零。xyzw3y3x3Q1My1Mx11243弹性薄板矩形（R12）单元课堂练习：试凑形函数Nx1？22122222(1)(1)()xNabcde考虑到挠度是非完全四次式，设由形函数性质，对Nx1有:Nx1(i)=0，i=1,2,3,4Nx1对x，y（本点除外）的偏导数在各结点均为零。1111214,10(1)10(2)0(3),(4)xixxixiNNNybNx222220bcdae最后利用本点1，确定a2=b/8，代回21(1)(1)(1-)/8xNb弹性薄板矩形（R12）单元记0=i；0=i仿N1可得:)/8-1)(1)(1(200ixibN对于转角xi相关的形函数，同样思路推导可得)/82()1)(1(220000iN)/8-1)(1)(1(200iyiaN对于转角yi相关的形函数，可推导得，2-2)其他形函数Nj、Nix、Niyxyzw3y3x3Q1My1Mx11243弹性薄板矩形（R12）单元3)薄板的挠度场有了每一结点的形函数,记eiiidNdNw41则薄板的挠度场可由结点位移表示为yixiiiNNNN41NNNxyzw3y3x3Q1My1Mx112434)单元间位移的协调性可以证明,上述w在边线上任意一点的挠度和转角都是三次多项式。弹性薄板矩形（R12）单元因此，边线的挠度和转角可由两端点的挠度和沿边线导数对应的转角唯一地确定。eiiidNdNw41由此可见，由形函数所建立的挠度场但是，边界法向转角只有两端两个法向转角位移条件，xyzw3y3x3Q1My1Mx11243当然无法唯一地确定，所以相邻单元法向转角位移不协调。是非完全协调的。(教材上有更严密的数学证明)1）相邻单元公共边挠度：位移的非完全协调性证明：结点位移协调：代入可以验证位移满足协调条件：2）相邻单元公共边切向转角：弹性薄板矩形（R12）单元3）相邻单元公共边法向转角：位移的非完全协调性证明：该转角的确定包含了单元全部结点位移参数，由于非公共边上结点位移的协调关系不能保证，因此一般综上所述，本节构造的位移场不能完全满足收敛的协调性准则，具体为挠度及切向转角跨单元协调，法向转角跨单元不协调，因此该单元不是完全协调元。请大家思考上述结论的力学实质？弹性薄板矩形（R12）单元弹性薄板矩形（R12）单元5)非完全协调元的收敛性eiiidNdNw41对于薄板等位移场非完全协调的位移模式,如何才能保证收敛呢?Irons给出了小片检验准则：用待检验的单元组成一小片，在无荷载、单元结点位移满足“常应变”状态位移条件时，如果各结点能够保持平衡且获得“常应力”受力状态，则这种位移模式对应的单元在如此网格下一定收敛。当程序无法计算已知支座位移问题时弹性薄板矩形（R12）单元小片检验的具体做法：•取某一单元小片，并在小片的边界上给出对应于完全二次多项式的边界条件。•按此位移模式进行（在无荷载作用下的）单元、整体分析，并在上述位移边界条件下求解。•若所求得的结点位移构造的小片上的挠度为一完全二次多项式，则单元的位移模式通过分片检验。•取某一单元小片，对小片的每一结点给以对应于完全二次多项式的结点位移。•每一单元按Fe=keδe求单元结点力，式中ke为对应所考察位移模式的单元刚度矩阵。•检验小片内部结点处是否均满足结点的自平衡条件，若条件均成立，则单元位移模式能通过小片检验。P.125给出了R12检查例子。请大家自己看一看。弹性薄板矩形（R12）单元薄板的单元列式BBBBB1234将位移模式代入可得形变矩阵为;;T222222在本章前面已导出：;;eeNdBdxyxyT222222xyzw3y3x3Q1My1Mx11243T2222222112;;iiababBN弹性薄板矩形（R12）单元2)薄板单元刚度方程式中eeeekdFFE由总势能的一阶变分为零可得eAkBDBATdAeFNqxyATE(,)d1)总势能用结点位移表示(-(,))eAAedΠBDBAqxyNAdTTdd2弹性薄板三角形（T9）单元同平面问题一样，对于复杂边界条件，三角形单元更理想。位移模式的选取1、直角坐标系下的位移模式方案zyxw1θy1θx1123三角形单元有9个结点位移，由帕斯卡三角形可知，要达到完整三次多项式要10项，而舍去任意一三次项将无法保证坐标不变性。教材上介绍了多种处理的方案并指出了他们的缺点。2、Zienkiewicz(监克维奇)提出了面积坐标位移模式1.直角坐标下单元位移模式：Tocher方案：当▲单元有两边分别平行于x轴和y轴时，上述位移模式中的待定系数将无法确定，因此离散时，网格划分有局限性。Adini方案：弹性薄板三角形（T9）单元Bell方案：增加单元内部位移参数——三角形形心挠度。整体分析前需要消去内部自由度（静力凝聚），Zienkiewicz指出这种单元不能保证收敛。1.直角坐标下单元位移模式：弹性薄板三角形（T9）单元Adini方案：舍去了二次项xy，致使常扭率无法保证，单元过刚、位移偏小，因此分析结果只有一阶精度。2wxy1.面积坐标下单元位移模式：以上三种方案均存在一定的问题，Zienkiewicz采用面积坐标解决了直角坐标下遇到的困难。利用9个已知的结点位移参数确定上述广义坐标参数，代入原式即得形函数。弹性薄板三角形（T9）单元脚标轮换规则,。1231基于面积坐标性质，上述位移模式还可改写为：弹性薄板三角形（T9）单元1.面积坐标下单元位移模式：作业：试证明其结果和监克维奇给出的结果完全相同。第一步：以作为结点自由度，求与之相对应的形函数；两步法确定形函数：将式中的换成，然后整理即可得到对应于直角坐标结点位移参数表达的形函数。第二步：利用关系式弹性薄板三角形（T9）单元1.面积坐标下单元位移模式：单元位移场坐标转换关系式123TTTT123T(1,2,3)eeiixiyieiixiyiwNNNiwNδNNNδNδδδδδ两步法比普通广义坐标法直接求待定系数简单。弹性薄板三角形（T9）单元两步法确定形函数：1.面积坐标下单元位移模式：弹性薄板三角形（T9）单元导数间的关系为其中三角形面积有了上述这些准备，由虚位移或势能