第八章一般壳体问题的有限元法课件

Click to edit Master title style,Click to edit Master text styles,Second level,Third level,Fourth level,Fifth level,11/7/2009,#,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,单击此处编辑母版标题样式,1,*,壳体实质上是从平板演变而来的，它的中面是一个曲面。在分析壳中应力时，虽然平板的基本假定同样有效，但是壳体的变形有着很大程度的不同，它除了弯曲变形外还存在中面变形。因而，壳中内力包括有弯曲内力和中面内力。,应用有限单元法分析壳体结构时，广泛地采用了平面单元和曲面的单元这两类壳体单元。本章首先介绍平面单元，它是平面应力问题和平板弯曲问题的组合；这种单元虽然简单，但是相当有效。然后讨论一个考虑横向剪切影响的曲面单元，称为八结点40个自由度的一般壳单元，可以适用于厚壳和薄壳。,第八章 一般壳体问题的有限元法,壳体实质上是从平板演变而来的，它的中面是一个曲面。在分,将壳体曲面划分为有限个单元，它们都是曲面单元。但是在单元细分时，用平面单元组成的一个单向或双向折板来近似壳体的几何形状将会得到良好的结果。通常对于任意形状的壳体，采用三角形单元比较方便，如图8-1所示。如果在壳体上容易找到同一平面上的四个点，可以采用平面四边形单元。例如具有正交边界的柱面壳体，如图8-2所示。,图8-1 任意壳体作为平面三角形单元的集合 图8-2 圆柱壳作为平面矩形单元的集合,8-1 平面壳体单元,图8-1 任意壳体作为平面三角形单元的集,壳体平面单元的应力状态是由平面应力和弯曲应力的叠加而成的，因此在构造壳体平面单元时，只要将第二章和第七章所讨论的相应单元进行简单的组合就可以了。同样，前述二章所导出的刚度矩阵可作为建立壳体平面单元刚度矩阵的基础。,现在把平面单元的计算步骤归纳如下,1,.划分单元，选定整体坐标系,oxyz,，定出节点在整体坐标系中的坐标值。,2,.对于各个单元利用节点坐标值，建立一个局部坐标系,例如三角形单元123，可以选取节点1为局部坐标系的原点，并且以1-2边为,轴的正方向，如图8-3所示。于是，方向的单位,e,1,求得是,(8-1),壳体平面单元的应力状态是由平面应力和弯曲应力,式中,是矢量,12,的长度。取单元的外法线方向作为 轴的正方向，于是它的单位矢量,图8-3 三角形单元局部坐标系,式中,容易看出，矢量,12,和,13,的矢性积的模等于三角形面积,的一倍，即|,12,13,|=2,。,最后，按右手定则可以决定,y,轴的正方向,它的单位矢量,e,2,是,e,2,=,e,3,e,1,(8-3),利用上述方法确定的局部坐标系，三角形单元123是在,平面内，它的三个角点的局部坐标值是很容易确定的。,对于柱面上的矩形单元，局部坐标的原点,选在矩形的形心，通常选 轴和,x,轴均沿柱面母线方向。如图8-4中所示，由矢量,12,确定单位矢量,e,1,再由矢量,14,确定单位矢量,e,2,，于是,e,3,=,e,1,e,2,。,容易看出，矢量12和13的矢性积的模等于三角形面积的一倍，,3,.对于各个单元，确立在局部坐标系,中的结点载荷列阵,。壳体载荷可以分解成二组：一组作用在平面内，另一组垂直于平面。为此，在计算各个单元的结点载荷列阵 （包括等效结点力）可以直接引用第二章和第七章中所叙述的载荷计算的相应公式。,各个单元的结点载荷列阵,图8-4 矩形单元局部坐标,求得后，建立变换矩阵公式，从而把 转换到整体坐标系中去求出在整体坐标下的单元节点载荷列阵，然后经各单元的简单叠加可以求出结构在整体坐标下的节点载荷列阵。,3.对于各个单元，确立在局部坐标系 中的,显然，平面单元在局部坐标系中，结点,i,有五个广义位移：即,，其中前两个对应于平面应力问题，后三个对应于平板弯曲问题。类似地，所对应的结点力列阵,显然，在上式中,实际上总是等于零的。,。为了经坐标变换后不影响在整体坐标系中对各特征量的计算，我们引进,(8-4),显然，平面单元在局部坐标系中，结点i有五个广,容易看出，把以上结点位移和结点力变换到整体坐标中后，他的结点位移和结点力列阵具有如下形式,上式右端的前三项分别表示位移和力，后三项分别表示转角和力矩，它们都是有明显物理意义的矢量。因此，（8-4）式和（8-5）式之间的坐标变换公式是,式中,而,(b),(8-5),(8-6),(a),容易看出，把以上结点位移和结点力变换到整体坐,于是，壳体单元,e,在局部坐标下的结点位移列阵是,或,而所对应的单元节点力（包括等效节点力）列阵是,或,式中,n,=3是对应于三角形单元；,n,=4对应于四边形单元。本节以下的,n,所指的意义均是如此，不再重复说明。,(c),(d),(e),(f),于是，壳体单元e在局部坐标下的结点位移列阵是或,4,.建立局部坐标系中的单元刚度矩阵,，从而求出整体坐标系中的单元刚度矩阵。如果将单元刚度矩阵 和 对应于单元节点划分为,n,n,个子矩阵，每个子矩阵都是66的，于是 的子矩阵有如下形式,式中,和 分别是平面应力问题和平面弯曲问题的相应子矩阵，它们是22和33矩阵。图8-5 示出了在局部坐标系中三角形壳体单元刚度矩阵用平面应力和平面弯曲刚度矩阵的构成方法。,(8-7),4.建立局部坐标系中的单元刚度矩阵,图8-5 三角形壳体单元刚度矩阵用平面应力和平板弯曲刚度矩阵的构成方法,图8-5 三角形壳体单元刚度矩阵用平面应力和平板弯曲刚度矩阵,单元,e,中任意结点,i,的平衡方程，在两个坐标系中分别为,，,式中,是刚度矩阵 的子矩阵。而对于局部坐标和整体坐标之间的变换公式是,把(h)式代入(g)式得,(g),(h),单元e中任意结点i的平衡方程，在两个坐标系中分别为，,将公式(g)中的第一式左乘矩阵,，并且同上式进行比较，可以得到,由于,是正交阵，容易证明 也是正交阵，即 。这样就得到关于矩阵 的转换公式,(8-8),将公式(g)中的第一式左乘矩阵,5,.集和单元刚度矩阵及等效结点力。线作简单求和,然后将它们放入整体刚度矩阵,K,和等效结点荷载列阵,的相应位置上去。,5.集和单元刚度矩阵及等效结点力。线作简单求,6,.修改整体刚度矩阵，然后求解平衡方程,式中,是总的结点位移列阵。特别值得注意，在局部坐标系中单元刚度矩阵,K,对于三角形单元它的第6、12、及18行和列全是零元素，对于四边形单元它的第6、12、18及24行和列全是零元素，其原因是转角 ，并不包含在平面应力的单元结点位移列阵中。当所有在一个结点相连接的单元共面时，壳体结构的刚度矩阵将是奇异的。避免这个奇异性的一个办法是引入关于壳体法线的转动为零的附加条件。,6.修改整体刚度矩阵，然后求解平衡方程式中,7,计算应力。首先是按照公式,求出局部坐标系中的结点位移，再按第二章中所给出的公式计算应力 、和,；通过第七章所给出的公式计,算,、和 进而求得应力,和,。,、和 。于是，壳体应力可以由简单的叠加求得；即,，,7计算应力。首先是按照公式,对于一个壳体结构，如果采用上节所述的平面单元，将会引起几何上的离散误差。人们希望采用曲面单元来描述壳体的真正几何形状，使之用不太多的单元来替代复杂形状的壳体，并得到具有一定精度的解答。另外，在薄壳理论中都是用中面位移来表示中面转动。正如在第七章中所述，这将要求在单元交界面上有横向位移及其一阶导数的连续性，于是增加了选择位移模式的困难。如果考虑横向剪切变形的影响就可以认为中面转动是独立变量而不依赖于位移的一阶导数。因此，只要利用单元交界面上位移函数的连续性就可以了，并不要求其一阶导数的连续性。,现在我们来论述一个考虑横向剪切影响的曲面单元，称为八结点四十个自由度的一般壳体单元，如图8-6所示。,8-2 考虑横向剪切变形影响的壳体单元,对于一个壳体结构，如果采用上节所述的平面单元,在图8-6中所示的壳单元，象空间等参数单元一样引进一个自然坐标系,。命 为壳体中面上的曲线坐标；对应于 的表面称为顶面（或上表面），对应于 的表面称为底面（或下表面）。在单元的中面上选取八个点称为结点，过各结点,i,(,i,=1,2,8)作中面的法线，交顶面和底面的点称为结点,i,的对点。结点,i,相对应的对点，它的整体坐标值分别记作,图8-6 八结点四十个自由度,的一般壳体单元,一.单元几何形状的确定,图8-6 八结点四十个自由度一.单元几何形状的确定,于是，中面上的结点,i,的整体坐标值是,显然，结点,i,处的中面法线方向可以由下列单位矢量所确定,(8-10),(8-9),于是，中面上的结点i的整体坐标值是显然，结点i处的中面法线方,式中,l,3,i,、,m,3,i,和,n,3,i,是结点,i,处中面法线方向对于整体坐标轴,oxyz,的方向余弦，而,h,i,是结点,i,处的壳体厚度，即,(a),结点,i,处法线上任意点的整体坐标值，可以通过矢量相加得到（图8-7），即,(b),式中l3i、m3i和n3i是结点i处中面法线方向对于整体坐标,于是，单元内任意点的坐标值可以通过形函数,的插值得到，即,式中形函数,由公式（5-1）表示。,这样，我们就可以通过八对点的整体坐标值，按照（8-11）式近似地确定了单元的形状。当,时，分别确定上下表面各点；当 时，确定中面各点；而单元的侧表面是由中面法线（或近似的中面法线）所构成。,(8-11),于是，单元内任意点的坐标值可以通过形函数,假设中面法线变形后仍为直线，但是不再是变形后的中面法线。也就是说，中面法线有绕两不同轴的转动，这两不同轴分别垂直于法线。设,V,1,i,和,V,2,i,分别是这两轴上的单位矢量，显然它们具有一定的任意性。为了唯一地确定起见，不妨设,式中,(c),(8-12),二.位移模式,式中(c)(8-12)二.位移模式,于是有,而,另一轴的单位矢量可以由,V,3,i,和,V,1,i,的矢量积求得,注意，在（8-12）式中如果,i,V,3,i,=0，则可用,j,V,3,i,来代替，而,(j),(d),(e),(8-13),于是有而另一轴的单位矢量可以由V3i和V1i的矢量积求得注意,设结点,i,处的中面法线,V,3,i,绕,V,1,i,和,V,2,i,两轴的转角分别为,和 ，,如图8-8所示。于是，转动矢量可以写成,图8-8 结点的自由度,(8-14),设结点i处的中面法线V3i绕V1i和V2i两轴,现在我们来计算结点,i,处法线上任意点的位移值。这可以通过结点,i,在整体坐标中的位移,u,i,、,v,i,、,w,i,以及法线绕,i,点转动而产生的位移相加而得到；关于转动引起的位移可以根据运动学中的公式计算。于是，在结点,i,处法线上各点在整体坐标系中的位移是,显然,。于是，上式就可以写成,(h),(g),现在我们来计算结点i处法线上任意点的位移值。,这样，利用插值方法可以得到单元内任意结点的位移列阵是,式中,(8-16),(8-15),这样，利用插值方法可以得到单元内任意结点的位移列阵是式中(,可以将（8-15）式写成标准形式,式中,I,是三阶单位阵。,(8-17),(i),可以将（8-15）式写成标准形式式中I是三阶单位阵。(8,在整体坐标系中，利用几何方程和（8-17）式可以将应变列阵写成标准形式是,（8-18）,三.应变计算,（8-18）三.应变计算,式中,而,(j),(8-19),式中而(j)(8-19),因此，为了计算应变需要算出六个偏导数,、和 、,前三个偏导数可以按照公式（5-54）进行计算，后三个则按定义分别是,J,-1,的第三列的三个元素，即,由（5-58）式得到,(8-20),因此，为了计算应变需要算出六个偏导数 、,出现在所提到的各式中的,S,、,T,和,V,分别按公式（5-55）所定义。将（8-11）式代入（5-55）式得,(8-21),出现在所提到的各式中的S、T和V分别按公式（5-55）所定义,我们还假设中面法线上的线段不伸长也不缩短。如果在单元中面各点处建立局部的直角坐标系,