第九章板壳结构有限元案例

第九章板壳结构有限元板壳结构基本知识板壳结构在工程上应用十分广泛。在设计分析中采用板壳单元进行结构分析，可以得到足够的精度和良好的效果。板壳结构基本知识厚度方向的尺寸小于长度和宽度方向尺寸的结构。其中，表面为平面的成为板，表面为曲面的称为壳。平板：分薄板和厚板。载荷作用在垂直于板面的方向。对于薄板小挠度问题，它的变形完全由横向变形确定；对于薄板大挠度问题，则属于几何非线性问题。对于厚板，应考虑横向剪切变形的影响。壳体：壳体的变形除了横向弯曲变形外，同时存在中面变形。因此可以认为壳体是平面应力问题和平板弯曲问题的组合。当然，对于厚壳结构，仍需要横向剪切变形的影响。1111~~1008085hb薄板结构有限元薄板基础理论知识薄平板，取其中性面为坐标面，z轴垂直于中性面。其中t为板厚。当板受有垂直于板中性面的外力时，板的中性面将发生弯扭变形，从而变成一个曲面。板变形的同时，在板的横截面上将存在内力——弯矩和扭矩。薄板基础理论知识对于薄板问题采用如下假设：（1）直法线假设：薄板中面法线变形后仍保持为法线且长度不变。（2）忽略板中面的法线应力分量，且不计其引起的应变。（3）薄板中面内的各点没有平行于中面的位移，即中面不变形。由第（2）条可知挠度w与z无关，由第（1）条可知zx和yz等于零，另外根据第（3）条中面无变形薄板基础理论知识薄板弯曲问题只需要考虑三个分量。根据几何方程，应变可表示为对于薄板问题，一般采用形变分量表示x向曲率y向曲率扭率薄板基础理论知识相应的内力可表示为：[D]为平面应力问题的弹性矩阵：薄板基础理论知识图示为板的一个微元体。为方便计，取x和y的方向的宽度均为1。在垂直于x轴的横截面上的正应力σx可合成为一个力偶，从而构成该横截面上的弯矩（单位宽度上的弯矩）Mx。同理，在垂直于y轴的横截面上的正应力σy合成弯矩My，剪应力τxy合成扭矩Mxy，剪应力τyx合成扭矩Myx,由于剪应力互等内力列向量为Mxy=Myx薄板基础理论知识内力可以根据应力进行计算得到使用记号于是平面应力问题中的弹性矩阵薄板基础理论知识进行反向回代，可以得到在板的上、下表面处，z=±0.5t，于是应力为薄板基础理论知识假设发生虚位移，则薄板内部会发生虚应变如果薄板在z方向承受分布荷载此时薄板内部产生应力与之平衡，则可以采用虚功原理应力做的虚功为外力做的虚功为在后面我们会利用虚功原理来建立有限元控制方程。薄板矩形单元设局部编号1、2、3、4，x、y方向长度分别为2a、2b的矩形板单元如图所示。每个结点的位移分量为则一个单元的位移向量和载荷向量为每个结点的载荷分量为薄板矩形单元下面开始尝试建立形函数。一个单元有12个位移分量，那么位移函数应该为利用12个结点位移条件，由广义坐标法可建立形函数，显然十分麻烦。因此形函数的建立采用拉格朗日插值函数形成，完成这项工作首先需要将其转化为一个2×2的正方形，对于矩形单元，这项操作并不困难。四次项的选取为了保证坐标的对称性，且曲率与扭率同阶次。薄板矩形单元下面开始尝试建立形函数。建立的形函数形式如下：每个分块的薄板矩形单元薄板矩形单元形函数的性质对N1有:N1(1)=1；N1(j)=0，j=2,3,4另外，N1对x，y的偏导数在各结点处均为零。？？？于是，位移函数可表达为：薄板矩形单元薄板矩形单元应变离散薄板矩形单元薄板矩形单元应力离散那么，相应的，内力矩薄板矩形单元薄板矩形单元的单元刚度矩阵，其形式也为通用的展开进行积分单元刚度矩阵由16个子矩阵组成，其表示如下薄板矩形单元具体的元素计算为：式中：薄板矩形单元结点载荷向量的计算：积分展开，得假设板单元受横向均布载荷p作用，则等效结点力为如果承受的分布荷载随位置（x,y）变化，积分工作量较大薄板矩形单元应用实例受中心集中力的四边支承板的计算结果(边长为1，厚度为0.01，弹模为1，泊松比为0.3)单元数（1/4板）四边固定板中心挠度wD/PL2边中点弯矩M/P2×20.00614-0.11784×40.00580-0.12336×60.00571-0.1245理论解0.00560-0.1257薄板三角形单元三角形单元能较好地适应斜边界，实际中广泛应用。单元的结点位移仍然为结点处的挠度wi和绕x，y轴的转角θxi、θyi，独立变量为wi。三角形单元位移模式应包含9个参数。如果在直角坐标系下建立位移模式，则完全三次多项式需要10个参数zyxw1θy1θx1123若以此为基础构造位移函数，则必须去掉一项。无法保证对称。薄板三角形单元三角形单元采用直角坐标系建立位移模式的尝试：Tocher方案单元有两边分别平行于x轴和y轴时，上述位移模式中的待定系数将无法确定，因此离散时，网格划分有局限性。Adini方案舍去了二次项xy，致使常扭率无法保证，单元过刚、位移偏小，因此分析结果只有一阶精度。Bell方案增加单元内部位移参数——三角形形心挠度。整体分析前需要消去内部自由度（静力凝聚），Zienkiewicz指出这种单元不能保证收敛。薄板三角形单元Zienkiewicz采用面积坐标解决了直角坐标下遇到的困难。面积坐标采用面积坐标表达的位移模式为：基于面积坐标性质，上述位移模式还可改写为：薄板三角形单元进行一系列的推导后，可得到形函数的具体计算式为：薄板三角形单元建立了位移模式之后，那么剩下的工作就并不复杂：位移模式→应变离散→应力离散→刚度矩阵→载荷向量→约束处理→求解。但位移模式是建立在面积坐标上的，相关的计算怎么进行？zyxw1θy1θx1123导数间的关系为那么薄板三角形单元应用实例四边简支板的中心挠度系数计算单元数（1/4板）板中心挠度wD/qL42×20.0042494×40.0041538×80.004098解析解0.004042薄板单元关于薄板单元，要提醒大家注意的：不管是三角形单元还是矩形单元，事实上其都是非完全协调元。对左图所示的相邻单元公共边挠度公共边切向转角公共边法向转角位移场不能完全满足收敛的协调性准则，具体为挠度及切向转角跨单元协调，法向转角跨单元不协调，因此该单元不是完全协调元。要解决这一问题，可以通过增加结点数或结点自由度来构造高阶协调元，目前国内较成熟的工作有大连理工大学唐立民教授等提出的拟协调元和龙驭球院士提出的广义协调元等，这方面仍然处于研究之中。厚板结构有限元厚板基础理论知识对于厚板，基本假设：a．板的挠度w微小；b．板中性面法线在变形后仍保持直线，但不再垂直变形后的中性面；c．垂直于中性面的应力可以忽略。中性面法线在变形后不再垂直变形后的中性面则意味着还存在横向剪切变形，此时即需要采用Hencky理论进行分析，这种情况下板任意一点有三个下面的变形事实上对于薄板，也是这样三个变形，只不过对于厚板厚板基础理论知识对于厚板，板的曲率和扭率为厚板与薄板的区别在于，还要考虑由于剪切而产生的应变会产生的转角差：厚板基础理论知识对于厚板，其微元体的平衡为采用与薄板类似的积分可得到厚板结构有限元对于厚板，板内任意一点(x,y,z)的位移为：其位移可以采用两种方式表达。方式一：对应的应变矩阵为厚板结构有限元方式二：对应的应变矩阵为事实上这只是写法的区别，没有实质影响厚板结构有限元由应变矩阵获得应力矩阵事实上这只是写法的区别，没有实质影响厚板结构有限元以最常用的8结点厚板单元给大家进行介绍首先需要将一个厚板单元进行等参变换，注意其是二维问题中面形状和厚度厚板结构有限元熟悉的二维8结点等参元形函数计算方法：结点位移采用第一种方式表示，则：22234158单元的位移可以采用形函数和结点位移表示为：其矩阵形式为：厚板结构有限元应变的表达厚板结构有限元应力的表达分块形式：厚板结构有限元应力的表达分块形式：厚板结构有限元单元刚度矩阵具体数据计算如下：厚板结构有限元结点荷载的等效应用举例设单元表面作用有均布荷载q(x,y)，等效结点荷载为承受均布荷载q的方板，四边简支。4×4网格，挠度=？h/L有限元厚板薄板0.010.044380.044390.044370.10.046280.046320.044370.20.052020.052170.044370.30.061600.061920.044370.40.075000.075570.04437厚板结构有限元一个问题：薄板单元是非协调元，所以总是让我们感觉不适，厚板单元认为两个转角与挠度独立，而且根据前面等参元的使用我们知道厚板单元一定是协调的，那么为什么不就采用厚板单元进行薄板的计算呢？用厚板单元进行薄板的计算在数学构造上并没有太大的问题，无非是现在转角是挠度的偏导数，所以不需要对一个矩形设置8个结点，4个结点已经足够表达。计算结果表明：当采用2×2高斯数值积分时，厚板单元也可用于薄板的分析，但是太薄时将产生剪切闭锁现象（进行刚度矩阵计算时与剪切变形相关的剪切刚度时会出现无穷大而导致单元刚度矩阵变成奇异矩阵）。壳结构有限元壳结构基础理论知识壳体的中性面是一个曲面，壳单元受力状态及应力状态见图。在作结构分析时，一般采用平面单元（板）或者曲面单元处理。平面单元是平面应力单元和平面弯曲单元的组合体，它依赖于平板理论。在几何上以平板代替壳体，结构模拟是一种近似。但是，这种单元简单，只要结构离散化分合理，完全可以满足工程上的要求。壳结构基础理论知识曲面单元能够更好地模拟真实结构，相应得到的计算结果会更有效。但是，曲面壳体的变形与平板变形有所区别。壳体的中性面变形不能忽略，在壳体中的内力包括弯曲内力和中性面内力。对于曲面单元，现常采用考虑横向剪切变形的超参数曲面壳单元。曲面壳元往往较难满足完备性和协调性要求，这里不作具体介绍。壳结构基础理论知识任何单曲或双曲薄壳，在单元较小时均可用薄板单元组成的单向或双向折板体系来近似，也就是采用平面壳单元进行分析。平面壳单元可以视为平面应力单元与板弯曲单元的组合体。平面应力单元（亦称膜单元）仅仅能够承受作用于平面内的载荷，不能够承受其它载荷。假设z方向上的位移w=0，每一结点仅存在沿x轴和y轴的位移板弯曲单元仅仅承受弯曲载荷，此类单元只有沿坐标z方向的位移，实际上就是利用前面已经接触过的薄板单元。平面壳单元有限元根据前面的假定，那么单元上任意一点(x,y,z)的位移为平面应力位移薄板弯曲位移注意：上面的位移表达式是基于局部坐标系建立的，不然则不成立。平面壳单元有限元根据上述位移关系，单元的应变矩阵为注意：上面的位移表达式是基于局部坐标系建立的，不然则不成立。平面壳单元有限元为进行以下的单元分析，定义单元结点i的位移列阵为i结点力列阵为单元结点位移、结点力矩阵为均为0在局部坐标系下其实位移列阵和力列阵的最后一项没有意义，但是考虑到最后还是要在整体坐标下进行计算，所以占一格。平面壳单元有限元有了刚才的位移列阵和力列阵，采用虚功原理，可以逐步完成有限元格式的建立过程。以经典的三结点平面壳单元为例局部坐标系下整体坐标系下平面壳单元有限元以经典的三结点平面壳单元为例三角形平面壳单元的3个结点i，j，m共有15个自由度，位移函数可采用如下形式平面壳单元有限元以经典的三结点平面壳单元为例将三个结点的位移代入进去，则可以反推出单元位移=形函数×结点位移的三个表达式（u,v,w）。根据位移函数的表达形式，不难看出其就是平面应力单元和薄板弯曲单元的结合。后续分析过程较复杂，因此在这里只做文字性叙述注意事项。单元位移表达式（u,v,w）建立后，下面的工作就是进行应变计算。但是注意up,vp并不是u,v平面壳单元有限元以经典的三结点平面壳单元为例得到应变后，接着就是计算应力，刚才的应变拆分成了两项。单元内应力分量可以简单地将相应的平面应力单元和薄板弯曲单元的应力分量进行叠加。平面应力单元所产生的应力薄板弯曲单元所产生的应力平面壳单元有限元以经典的三结点平面壳单元为例在上述的分析过程中我们将得到应变矩阵[B]和[D]，那么就可以和往常一样的建立起单元刚度矩阵[k]，但是此时建立起来的单刚矩阵是在局部坐标下的，所以需要转换成整体坐标。组装完毕后得到总刚和平衡方程：同样的也要先进行约束的处理才能进