第8章空间实体单元

第8章空间实体单元

8.1概述许多工程实际问题属于空间问题。用有限元法分析空间问题和分析平面问题在原理、思路和解题方法完全相同，基本未知量仍然是节点位移。不同的是单元具有三维特点。节点位移在x、y、z三个坐标轴方向都有分量：u、v、w。它的基本方程比平面问题要多，有3个平衡方程，6个几何方程，6个物理方程。分析方法仍然是先进行单元分析，再进行系统分析，最后求解系统的节点平衡方程，解算内力或应力。空间离散化后的单元模型主要有：四面体单元、长方体单元、直边六面体单元、曲边六面体单元，如图8-1所示。(a)(b)(c)(d)图8-1空间实体单元模型(a)图为4节点四面体单元，是空间问题最简单的单元，也是常应变、常应力单元。类似平面问题三节点三角形单元进行分析。？(b)图w为长方体单元，可以类似平面四节点矩形单元进行分析。

(c)图为任意八节点六面体单元，可以类似平面四节点任意四边形等参元进行分析。(d)图为20节点曲边六面体单元，可以类似平面八节点曲边四边形等参元进行分析。8.24节点四面体常应变单元1、位移模式

如图8-2所示，取四面体的4个顶点i,j,m,n为节点。每一个节点有3个位移分量，即（8-1）单元节点位移向量为（8-2）与平面问题式（2-12）类似，假定单元内一点的位移分量为坐标的线性函数（8-3）将式（8-3）的第1式应用于4个结点，则ijmnxyz图8-2（8-4）由此可解出a1～a4，再代回到式（8-3）的第1式，与式（2-19）的第1式类似，有（8-5）式中形函数具有与式（2-18）类似的形式：（8-6）其中（8-7）（8-8）在式（8-8）中，V为四面体的体积。为使其计算值不为负，单元的节点（i,j,m,n）编号次序应遵循右手法则。(p4)采用同样的方法，可得（8-9）（8-10）将式（8-5）、（8-9）（8-10）统一用矩阵式表示，可得与平面问题式（2-20）类似的公式（8-11）式中[N]为单元形函数矩阵，其维数为3×12。进一步可写为与平面问题式（2-21）、（2-22）类似的子块形式（8-12）其中，子矩阵（8-13）式中，I为3阶单位矩阵。2、应变矩阵

在空间问题中，每点有6个应变分量。几何方程为：（8-14）将式（8-11）～（8-13）和（8-6）代入上式，得（8-15）式中（8-16）上式（8-15）、（8-16）与平面问题式（2-24）～（2-26）类似。与平面问题三节点三角形单元相同，在四节点四面体单元中，[B]的元素都是常量，因此是常应变单元。2、应力矩阵

三维问题的应力应变关系也可写为式（2-8）的矩阵形式（8-17）与平面问题不同，这里和分别由6个分量组成，弹性矩阵[D]是一个6×6的矩阵：（8-18）（8-19）（8-20）将式（8-20）所表示的[D]和式（8-15）、（8-16）所表示的[B]代入式（8-22），并将[S]定成分块矩阵的形式，有将式（8-15）代入式（8-17），得（8-21）应力矩阵[S]为（8-22）由于[D]、[B]都是常数矩阵，因此应力矩阵[S]也是常数矩阵。也就是说，单元中的应力分量也是常数。（8-23）式中（8-24）其中（8-25）3、单元刚度矩阵

仿照平面问题中的推导，可得单元平衡方程（8-26）单元刚度矩阵具有与式（2-33）类似的形式（8-27）式中，[k]是一个12×12的矩阵。由于[B]、[D]都是常数矩阵，所以[k]也是一个常量矩阵。并且（8-28）写成分块矩阵的形式，有（8-29）式中子矩阵[krs]为3×3的矩阵（8-30）4、等价节点力向量

式（8-26）中的单元等价节点力也包括体积力、表面力、集中力几部分。体积力与表面力的计算公式与平面三角形单元公式（2-36）、（2-37）类似：（8-31）（8-32）对于简单情形，也可采用静力等效原则简化计算。进一步的整体平衡方程的建立（即结构刚度矩阵、结构等价节点力列阵的组集）、位移约束条件的引入、线性方程组的求解等，和平面问题有限元法一样，不再赘述。8.3二十结点六面体等参数单元由于精度高，容易适应不同边界，在平面问题中常选用了八节点四边形等参数单元。与此类似，在三维问题中,常选用二十节点六面体等参数单元。如图8-3所示，在整体坐标系的二十节点六面体的实际单元与中心在局部坐标的原点、边长为2的立方体基本单元相对应。1234567891011121314151617181920••••••••••••••••••••xyz6121812345789101113141516171920••••••••••••••••••••=-1=1=1=-1=1图8-3实际单元基本单元1、位移模式、形函数和坐标变换式=-1（8-33）（8-34）形函数的表达式如下：

（8-35）式中位移模式和坐标变换式可写为如下形式：根据几何方程，单元中的应变为2、应变矩阵（8-36）其中（8-37）应变矩阵的子矩阵[Bi]为：（8-38）根据复合函数求导规则，有（8-39）式中[J]-1为雅可比矩阵[J]的逆矩阵。[J]的表达式为（8-40）3、应力矩阵单元中的应力为单元应力矩阵[S]为（8-41）4、单元刚度矩阵单元刚度矩阵可写成（8-42）[k]是一个60×60的矩阵，式（8-42）通常采用高斯法进行积分。5、等价节点力矩阵等价节点力计算公式如下：（1）体积力设单位体积力是，则等价节点力为（8-43）（2）表面力设某边界面上作用表面力，则等价节点力为（8-44）设该边界面对应基本单元=±1的面。由数学公式，结构坐标下曲面微元dS对应于单元坐标下的微元面积式为（8-45）则式（8-44）可写为单元坐标系下的积分公式（8-46）以上为=±1的表面力计算公式。对于其他表面力，可类似处理。式（8-46）通常也采用高斯法进行积分计算。8.4空间轴对称单元许多工程构件，其几何形状、约束条件及所受的荷载都对称于某一轴，因而所有的位移、应变和应力分量也都对称于该轴。这类问题称空间轴对称问题。对空间轴对称问题，采用圆柱坐标系，r表示径向坐标，z表示轴向坐标，任一对称面为rz面。在有限元分析时，采用绕对称轴旋转一周的轴对称环形单元。把轴对称的工程构件与rz坐标平面正交的截面划分为一些三角形，例如其中一个是ijm（图8-4）。由这些三角形旋转一周就构成上述环形单元。图8-4轴对称构件及三角形环形单元r(u)z(w)ij当然，也可以划分为矩形，构成矩形截面的环形单元。总之，单元截面形状，可采用三角形、四边形等平面有限元法所用单元形状。本节只讨论三角形截面的环形单元，其它截面形状环形单元可参照本节方法进行分析。各个环形单元（后简称单元）间用位于三角形顶点的环形铰联系起来，就成了离散结构模型。我们把环形铰在rz平面的交点（即三角形顶点）称为单元的节点（如i,j,m）。因此对轴对称问题进行有限元分析只需在rz平面划分网格，就像平面问题在xy平面上划分网格一样。照此思路，轴对称空间问题被大大简化。1、位移模式对于如图8-4所示轴对称三角形环形单元，是由如图8-5所示平面上的三角形绕对称轴轴回旋一周得到的。考虑到问题的轴对称特点，可以借助于平面问题三节点三角形单元位移模式进行单元分析。图8-5rzijmrirjrm由于轴对称，环向位移恒等于零。只有径（r）向位移和轴（z）向位移，它们恰在rz坐标平面上。设径向位移为u，轴向位移为w。对于图8-5所示情形，借助平面问题的三角形单元，取位移模式为？（8-47）代入节点位移后，可解出a1～a6，再代入上式，得（8-48）其中，形函数（8-49）式中A,ai,bi,ci，与平面问题三角形单元的对应公式（2-15）、（2-17）一致，区别仅仅是将那里的x,y换成这里的r,z。式（8-48）也可写为式（2-20）同样的形式（8-50）式中，[Ni]=NiI(i,j,m)，其中I为2阶单位矩阵。2、应变矩阵根据弹性力学理论，空间轴对称问题的几何方程为（8-51）与平面问题中只含有3个应变分量不同，这里含有4个应变分量。将u,w的表达代入式（8-51），得（8-52）式中（8-53）其中（8-54）由式（8-52）～（8-54）可见，矩阵中含有变量r,z，因此它不是常数矩阵。即轴对称问题的三角形环形单元不是常应变单元。3、应力矩阵根据弹性力学理论，空间轴对称问题的应力-应变关系为（8-55）式中[D]是轴对称问题的弹性矩阵（8-56）将式（8-52）代入式（8-55），得（8-57）应力矩阵[S]=[D][B]显然，除rz外，单元中其他应力分量不是常数。4、单元刚度矩阵轴对称问题的单元刚度矩阵可由式（8-27）计算由于被积函数与无关，故在三角形截面的环单元的积分可简化为在三角形截面上的积分：（8-58）式（8-58）中，[B]不是常数，而由式（8-53）、（8-54）确定。式（8-53）、（8-54）中，ai,bi,ci是（P32）（1）一般公式（2）近似计算可见，单元刚度矩阵计算比平面三角形单元麻烦。须对其进行数值计算或近似计算。现在，讨论面积分∬Ag(r,z)drdz的计算问题。虽为常数，由式（2-17）计算，但fi是r,z的函数。式（8-58）被积函数[B]T[D][B]r是r,z的函数，暂简写为g(r,z)。式（8-58）可写为图8-6rzijm对于一个典型三角形截面ijm（图8-6），过j点作垂直于r轴的直线，将三角形ijm，分成两个三角形：A1，A2。于是，有A2A1（8-59）完成上述积分可采用数值方法，有专门的二重数值积分程序可供引用。也可采用近似方法完成：当单元较小时，常把各个单元中的r,z近似看作常数，并且分别等于各单元形心的坐标，即ijmA2A1rzdefZ1(r)Z2(r)Z3(r)则式（8-54）成为常数：（8-60）这样，就可把各个单元近似地当做常应变单元，式（8-58）变成（8-61）式中[B]是从式（8-52）～（8-54）将[B]中的r,z用r,z代替后得出的。若将单元刚度矩阵[k]写为子块形式其中，子矩阵可近似为（8-62）写成显式（8-63）式中A1,A2,A3由式（8-25）确定。5、等价节点力轴对称问题的等价节点力可用类似平面问题的式（2-36）、（2-37）写出。对于作用于三角形环单元上的体积力、表面力的等效结点力为：（8-64）（1）体积力的等价节点力①一般公式体积力等价节点力的一般计算式为②重力设物体的容重为(=g)，则由式（8-64），有（8-65）注意到结构坐标与面积坐标的关系式：（8-66）以及Ni与Li等价的关系，利用积分公式（7-27），有（8-67）式中将式（8-67）代入式（8-65），得到（8-68）（2）表面力的等价节点力①一般公式由式（8-69）计算等价节点力时，其中r由（8-66）式确定。并注意到在ij边上有Lm=0。利用式（7-28）进行积分，最后得（8-69）表面力等