清华大学：《有限元分析》课程教学资源（电子讲义）第五章 .等参单元

5.等参单元 本章包括以下内容 51等参单元的基本概念 52四边形八节点等参单元 53等参单元的单元分析 54六面体等参单 51等参单元的基本概念 在进行有限元分析时,单元离散化会带来计算误差,主要采用两种方法来降低单元离散 化产生的误差:1)提高单元划分的密度,被称为h方法(h- method):2)提高单元位移函 数多项式的阶次,被称为p方法( p-method) 在平面问题的有限单元中,我们可以选择四结点的矩形单元,如图5-1所示,该矩形单 元在x及y方向的边长分别为2a和2b 图5-1四结点矩形单元 同第三章的方法类似,将单元的位移模式选为, u=a1+a2x+a3y+a4xy v=a5+a6xta7y+agx (5-1) 可得到, Ni ui+n V=Ni+NVj+Nm,m+Npp (5-2) 形态函数为 N=7(1+1-2

5．等参单元 本章包括以下内容： 5.1 等参单元的基本概念 5.2 四边形八节点等参单元 5.3 等参单元的单元分析 5.4 六面体等参单元 5.1 等参单元的基本概念 在进行有限元分析时，单元离散化会带来计算误差，主要采用两种方法来降低单元离散 化产生的误差：1）提高单元划分的密度，被称为 h 方法（h-method）；2）提高单元位移函 数多项式的阶次，被称为 p 方法（p-method）。 在平面问题的有限单元中，我们可以选择四结点的矩形单元，如图 5-1 所示，该矩形单 元在 x 及 y 方向的边长分别为 2a 和 2b。 图 5-1 四结点矩形单元 同第三章的方法类似，将单元的位移模式选为， u a a x a y a xy = 1 + 2 + 3 + 4 v a a x a y a xy = 5 + 6 + 7 + 8 （5-1） 可得到， u = Niui + N ju j + Nmum + Npu p i i j j m m p p v = N v + N v + N v + N v （5-2） 形态函数为， (1 )(1 ) 4 1 b y a x Ni = − − (1 )(1 ) 4 1 b y a x N j = + −

(1+-)1+2) P=4(1-M1+ (5-3) 上述单元位移模式满足位移模式选择的基本要求: 1)反映了单元的刚体位移和常应变, 2)单元在公共边界上位移连续 在矩形单元的边界上,坐标x和y的其中一个取常量,因此在边界上位移是线性分布的 由两个结点上的位移确定 与三结点三角形单元相比,四结点矩形单元的位移模式是坐标的二次函数,能够提高计 算精度,但也有显著的缺点,两种单元的比较如下, 表5-1三结点三角形单元与四结点矩形单元比较 单元类型 点 缺点 三结点三角形单元适应复杂形状, 计算精度低 单元大小过渡方便 四结点矩形单元单元内的应力、应变是线性变不能适应曲线边界和非正交 化的,计算精度较高 的直线边界 如果任意形状的四边形四结点单元采用矩形单元的位移模式,则在公共边界上不满足位 移连续性条件。为了既能得到较高的计算精度,又能适应复杂的边界形状,可以采用坐标变 换 4 r=1 2=1 55=1 7=-1 图5-2任意四结点四边形单元 图5-3四结点正方形单元 在图5-2所示的任意四边形单元上,用等分四条边的两族直线分割四边形,以两族直线 的中心为原点,建立局部坐标系(5,m),沿5及n增大的方向作为5轴和7轴,并令四条边 上的ξ及n值分别为±1。为了求出位移模式,以及局部坐标与整体坐标之间的变换式,在 局部坐标系中定义一个四结点正方形单元,如图5-3所示 参照矩形单元,四结点正方形单元的位移模式为 l=N141+N2a2+N3u3+N4u4 N2v2+N3V3+Nav (5-4)

(1 )(1 ) 4 1 b y a x Nm = + + (1 )(1 ) 4 1 b y a x N p = − + （5-3） 上述单元位移模式满足位移模式选择的基本要求： 1）反映了单元的刚体位移和常应变， 2）单元在公共边界上位移连续。 在矩形单元的边界上，坐标 x 和 y 的其中一个取常量，因此在边界上位移是线性分布的， 由两个结点上的位移确定。 与三结点三角形单元相比，四结点矩形单元的位移模式是坐标的二次函数，能够提高计 算精度，但也有显著的缺点，两种单元的比较如下。 表 5-1 三结点三角形单元与四结点矩形单元比较 单元类型 优点 缺点 三结点三角形单元 适应复杂形状， 单元大小过渡方便 计算精度低 四结点矩形单元 单元内的应力、应变是线性变 化的，计算精度较高 不能适应曲线边界和非正交 的直线边界 如果任意形状的四边形四结点单元采用矩形单元的位移模式，则在公共边界上不满足位 移连续性条件。为了既能得到较高的计算精度，又能适应复杂的边界形状，可以采用坐标变 换。 图 5-2 任意四结点四边形单元 图 5-3 四结点正方形单元 在图 5-2 所示的任意四边形单元上，用等分四条边的两族直线分割四边形，以两族直线 的中心为原点，建立局部坐标系 ( ,) ，沿  及  增大的方向作为  轴和  轴，并令四条边 上的  及  值分别为  1 。为了求出位移模式，以及局部坐标与整体坐标之间的变换式，在 局部坐标系中定义一个四结点正方形单元，如图 5-3 所示。 参照矩形单元，四结点正方形单元的位移模式为， u = N1u1 + N2u2 + N3u3 + N4u4 1 1 2 2 3 3 4 4 v = N v + N v + N v + N v （5-4）

其中, N1=(1-5)(1-m) (1+5)(1-n) N3=(1-51+n) N4=(1+5)(1+m) (5-5) 四个结点的坐标为(512n),定义新的变量 50=55,m=n(=1,2,3,4) (5-6) 形态函数表示为 N=(1+50)(1+70)(i=1,2,3,4)(5-7) 把5及7作为任意四边形单元的局部坐标,把(5-4)的位移模式和(5-7)的形态函数 用于任意形状的四边单元,可得: 在四个结点处可以得到结点的位移 2)在单元的四条边上,位移线性变化,保证了单元公共边界上位移的连续性。 因此给出任意四边形单元的结点位移就能得到整个单元上的位移,(5-4)的位移模式就 是所要找的正确的位移模式。 把局部坐标与整体坐标的变换式也取为, x=Nix1+N2x2+N3x3+N4x4 y=N1y1+N2y2+N3]3+Naya 将坐标变换式用于任意四边形单元,可得 1)在四个结点处给出结点的整体坐标, 2)在四条边上的整体坐标是线性变化的。 只要给出任意四边形单元四个结点的整体坐标,用(5-8)式就可以建立局部坐标系中 的正方形单元和整体坐标系中的任意四边形单元之间的坐标变换关系。 把图5-3中的局部坐标系中的正方形单元称为基本单元。 把图5-2中的在整体坐标系中的任意四边形单元看作由基本单元通过坐标变换得来的 称为实际单元 单元几何形状和单元内的未知量采用相同数目的结点参数以及相同的插值函数进行变 换,称为等参变换。采用等参变换的单元,称为等参单元。 由于形态函数M;,正好反映了单元形状的变化,也称为形函数( Shape function)o 采用等参单元,使我们可以在局部坐标系中的规则单元上进行单元分析,然后在映射到 实际单元上。等参单元同时具有计算精度高和适用性好的特点,是有限元程序中主要采用的 单元形式

其中， (1 )(1 ) 4 1 N1 = − − (1 )(1 ) 4 1 N2 = + − (1 )(1 ) 4 1 N3 = − + (1 )(1 ) 4 1 N4 = +  + （5-5） 四个结点的坐标为 ( , )  i i ，定义新的变量，  0 =  i ，0 =i (i=1，2，3，4) (5-6) 形态函数表示为， (1 )(1 ) 4 1 Ni = + 0 +0 （i=1，2，3，4） （5-7） 把  及  作为任意四边形单元的局部坐标，把（5-4）的位移模式和（5-7）的形态函数 用于任意形状的四边单元，可得： 1）在四个结点处可以得到结点的位移； 2）在单元的四条边上，位移线性变化，保证了单元公共边界上位移的连续性。 因此给出任意四边形单元的结点位移就能得到整个单元上的位移，（5-4）的位移模式就 是所要找的正确的位移模式。 把局部坐标与整体坐标的变换式也取为， 1 1 2 2 3 3 4 4 x = N x + N x + N x + N x 1 1 2 2 3 3 4 4 y = N y + N y + N y + N y （5-8） 将坐标变换式用于任意四边形单元，可得： 1）在四个结点处给出结点的整体坐标， 2）在四条边上的整体坐标是线性变化的。 只要给出任意四边形单元四个结点的整体坐标，用（5-8）式就可以建立局部坐标系中 的正方形单元和整体坐标系中的任意四边形单元之间的坐标变换关系。 把图 5-3 中的局部坐标系中的正方形单元称为基本单元。 把图 5-2 中的在整体坐标系中的任意四边形单元看作由基本单元通过坐标变换得来的， 称为实际单元。 单元几何形状和单元内的未知量采用相同数目的结点参数以及相同的插值函数进行变 换，称为等参变换。采用等参变换的单元，称为等参单元。 由于形态函数 Ni ，正好反映了单元形状的变化，也称为形函数（Shape function）。 采用等参单元，使我们可以在局部坐标系中的规则单元上进行单元分析，然后在映射到 实际单元上。等参单元同时具有计算精度高和适用性好的特点，是有限元程序中主要采用的 单元形式

52四边形八节点等参单元 为了更好地反映物体内的应力变化,适应曲线边界,在弹性力学平面问题的分析中经常 使用四边形八节点等参单元。如图54所示,由于每条边上增加了一个结点,单元的边是 条二次曲线,可以更好地适应曲线边界, 6 n 8 2=-1 =1 图54四边形八结点单元 图5-5八结点基本单元 对于等参单元,先在图5-5所示的八结点基本单元上进行分析。八结点单元一共有16 个已知的结点位移分量,基本单元中取如下的位移模式 l=a1+a25+a3n+a152+a3m+a4n2+a252n+a25n v=b1+b25+b+b2+b5n+b6n2+b,52+b5 (5-9) 该位移模式实际上是一个双二次函数,待定系数由结点位移分量确定。在单元的每条边 上,局部坐标ξ=±1或门=±1,位移是局部坐标ξ或刀的二次函数,完全由边上的三个结 点的位移值确定,所以这个位移模式满足位移连续性条件。实际单元内的位移用形函数表示 ∑N(5,n N,(5,n) 其中的形函数为 N1=,(1-5)(1-m)(-5-n-1) N3=-(1+5)(1-m)(-n-1)

5.2 四边形八节点等参单元 为了更好地反映物体内的应力变化，适应曲线边界，在弹性力学平面问题的分析中经常 使用四边形八节点等参单元。如图 5-4 所示，由于每条边上增加了一个结点，单元的边是一 条二次曲线，可以更好地适应曲线边界， 图 5-4 四边形八结点单元 图 5-5 八结点基本单元 对于等参单元，先在图 5-5 所示的八结点基本单元上进行分析。八结点单元一共有 16 个已知的结点位移分量，基本单元中取如下的位移模式： 2 8 2 7 2 5 6 2 u = a1 + a2 + a3 + a4 + a  + a  + a   + a  2 8 2 7 2 5 6 2 v = b1 + b2 + b3 + b4 + b  + b  + b   + b   （5-9） 该位移模式实际上是一个双二次函数，待定系数由结点位移分量确定。在单元的每条边 上，局部坐标  = 1 或  = 1 ，位移是局部坐标  或  的二次函数，完全由边上的三个结 点的位移值确定，所以这个位移模式满足位移连续性条件。实际单元内的位移用形函数表示 为， i i i u N ( , )u 8 1    = = i i i v N ( , )v 8 1    = = （5-10） 其中的形函数为： (1 )(1 )( 1) 4 1 N1 = − − − − − (1 )(1 )( 1) 4 1 N3 = + −  − −

N5=(1+(1+m)(+n-1) (1-5)(1+m)(-5+n-1) (1-2)(1-m) (1-5)(1+m) (1-n2)1 N3=(1-mn2)1-5) 将形函数归纳为 (1+5;5)(1+n)X55+n)(=13.5,7) N,(5,m) (1-22)1+m)(=2,6) (5-11) (1-n2)1+55)(i 形函数N(5,m)在单元的i结点上的值为1,在其它结点上的值均为0。 坐标变换式采用如下相似的公式, x=∑N(5,n)x (5-12) y=∑N(5,m)y 将ξ=1代入公式(5-12),可以得到单元345边在整体坐标下的参数方程: x=an+bn+c (5-13) y=dn+en+f 可见在整体坐标系中,单元的边是一条抛物线或退化为一条直线

(1 )(1 )( 1) 4 1 N5 = + +  + − (1 )(1 )( 1) 4 1 N7 = − + − + − (1 )(1 ) 2 1 2 N2 = − − (1 )(1 ) 2 1 2 N6 = − + (1 )(1 ) 2 1 2 N4 = − +  (1 )(1 ) 2 1 2 N8 = − − 将形函数归纳为，          − + = − + = + + + = = (1 )(1 ) ( 4,8) 2 1 (1 )(1 ) ( 2,6) 2 1 (1 )(1 )( ) ( 1,3,5,7) 4 1 ( , ) 2 2 i i i N i i i i i i i                 （5-11） 形函数 (,) Ni 在单元的 i 结点上的值为 1，在其它结点上的值均为 0。 坐标变换式采用如下相似的公式， i i i i i i y N y x N x ( , ) ( , ) 8 1 8 1       = = = = （5-12） 将  = 1 代入公式（5-12），可以得到单元 345 边在整体坐标下的参数方程： y d e f x a b c = + + = + +     2 2 （5-13） 可见在整体坐标系中，单元的边是一条抛物线或退化为一条直线

PLANES2 2-D8-Node structural Solid (or Axial Triangular Opti X(or Ra dial 图5-6 ANSYS提供的 Plane82单元 如图5-6所示, ANSYS提供的 PLANE82单元是一个四边形八结点等参单元,局部坐 标定义为s和t,如图5-7所示。 PLANE82单元可以退化为三角形六结点单元 X.RU 图5-7Pane82的基本单元 ANSYS理论手册中给出的 PLANE82单元的位移模式如图5-8所示,位移模式与公式 (5-10)展开后是一样的。 Equation 12.113 u=(u(1-s)(1-t)(s-t-1+u(1+)(1-t)(5-t-1 +uk(1+)(1+)(s+t-1)+uL(1-)(1+t)(-s+t-1 +r(M(1-s2)(1-t)+u(1+s)(1-t2) +u0(1-521+t)u(1-s1-t2) Equation 12. 114 a((1-)…( analogous to u 图5-8 Plane82单元位移模式

图 5-6 ANSYS 提供的 Plane82 单元 如图 5-6 所示，ANSYS 提供的 PLANE82 单元是一个四边形八结点等参单元，局部坐 标定义为 s 和 t，如图 5-7 所示。PLANE82 单元可以退化为三角形六结点单元。 图 5-7 Plane82 的基本单元 ANSYS 理论手册中给出的 PLANE82 单元的位移模式如图 5-8 所示，位移模式与公式 （5-10）展开后是一样的。 图 5-8 Plane82 单元位移模式

53等参单元的单元分析 在本节,以平面问题的四边形八结点等参单元为例,介绍构造等参单元的单元刚度矩阵 的基本过程。弹性力学平面问题的单元刚度矩阵为, IK]=JIB][DJ[B]dxdy 单元的应变为, {s}=[B]{o} 单元的结点位移, }=[a1n1t2n2 将形函数代入后,可以得到应变的矩阵表达式 ax VI 0 (5-14) aNaN an aNs lg 可得应变矩阵的分块矩阵, [B,] oN (5-15) ON. aN 由于等参单元的形函数是局部坐标(,m)的函数,因此应变矩阵[B]也是局部坐标(5,n) 的函数。形成等参单元的单元刚度矩阵需要在整体坐标系中对局部坐标的函数进行积分,包 括以下三个基本步骤:1)计算用局部坐标表示的形函数M(5,n)对整体坐标x、y的偏导 数:2)将整体坐标系中的面积积分转换为在局部坐标系中的面积积分;3)用数值积分计算 出单元刚度矩阵中的元素 (一)计算形函数对整体坐标x,y的偏导数 由于局部坐标与整体坐标之间存在坐标转换关系,因此形函数N是局部坐标的函数 同时也可以看作是整体坐标的函数。由复合函数求导法则可得

5.3 等参单元的单元分析 在本节，以平面问题的四边形八结点等参单元为例，介绍构造等参单元的单元刚度矩阵 的基本过程。弹性力学平面问题的单元刚度矩阵为，   =  K B D B tdxdy e T [ ] [ ][ ] 单元的应变为， e {} = [B]{} 单元的结点位移，   e T [u v u v ... u v ]  = 1 1 2 2 8 8 将形函数代入后，可以得到应变的矩阵表达式，                                                                 = 8 8 1 1 1 1 2 2 8 8 1 2 8 1 2 8 . . ... 0 ... 0 ... v u v u x N y N x N y N x N y N y N y N y N x N x N x N  （5-14） 可得应变矩阵的分块矩阵，                           = x N y N y N x N B i i i i i 0 0 [ ] （5-15） 由于等参单元的形函数是局部坐标 ( ,) 的函数，因此应变矩阵[B]也是局部坐标 ( ,) 的函数。形成等参单元的单元刚度矩阵需要在整体坐标系中对局部坐标的函数进行积分，包 括以下三个基本步骤：1）计算用局部坐标表示的形函数 (,) Ni 对整体坐标 x、y 的偏导 数；2）将整体坐标系中的面积积分转换为在局部坐标系中的面积积分；3）用数值积分计算 出单元刚度矩阵中的元素。 （一）计算形函数对整体坐标 x，y 的偏导数 由于局部坐标与整体坐标之间存在坐标转换关系，因此形函数 Ni 是局部坐标的函数， 同时也可以看作是整体坐标的函数。由复合函数求导法则可得：

ax aN, a ax an, ay an ax an ay an aNl「axoy as[as a4 ax ay ax 定义 (5-17) ax ay aN (5-18) 矩阵[J称为雅可比矩阵( Jacobian matrix),单元的整体坐标可以形函数来表示,因此用 坐标变换公式可以计算雅可比矩阵 aN。1xy V2 aM 由(5-18)可得 ON/=[]as (5-20) ay ay (5-21) an as 将(5-21)代入(5-20)得到, Oy an Oy an an as as an (5-22) Jl_Ox aN,+ Ox aN an asas an

          +     =       +     =   y y x N x N N y y x N x N N i i i i i i （5-16）                                       =                   y N x N x y x y N N i i i i       定义，                     =     x y x y [J ] （5-17）                   =                   y N x N J N N i i i i [ ]   （5-18） 矩阵[J]称为雅可比矩阵（Jacobian Matrix），单元的整体坐标可以形函数来表示，因此用 坐标变换公式可以计算雅可比矩阵。                                     = 8 8 2 2 1 1 1 2 8 1 2 8 ... ... ... ... [ ] x y x y x y N N N N N N J       （5-19） 由（5-18）可得，                   =                   −   i i i i N N J y N x N 1 [ ] （5-20）                 −   −   = = −     x x y y J J J J [ ] 1 [ ] * 1 （5-21） 将（5-21）代入（5-20）得到，                   +     −     −     =                           i i i i i i x N x N y N y N J y N x N 1 （5-22）

(二)将整体坐标系中的面积积分转换为在局部坐标系中的面积积分 K=IB]IDIB]dxdy (5-23) 在整体坐标系中,面积微元为x方向和y方向微矢量的叉乘的模量, dA=dx (5-24) as an 币=d5 dA ds +)x( Gds Needn as an an as 77 代入(5-23),得到单元刚度矩阵在局部坐标系中的积分公式: 于=,[BB1“BDBB2BW(52 单元刚度矩阵中的任意一个分块矩阵的积分公式为, IK]n -LL [B, IDIB, ] diDn (5-27) (三)用数值积分计算出单元刚度矩阵中的元素 等参单元刚度矩阵的每个元素都是局部坐标的函数,在有限元程序中不用解析的办法来 计算局部坐标系中的积分,而采用数值积分方法。通常采用高斯积分方法计算单元刚度矩阵 中的元素。 高斯积分方法预先定义了积分点和相应的加权系数,求出被积分的函数在指定积分点上 的数值,加权后求和,就得到了该函数的积分。这种方法具有比较高的计算精度。已经证明, 采用n个积分点的高斯积分可以达到2n-1阶的精度,也就是说,如果被积分的函数是2n-1 次多项式,用高斯积分可以得到精确的积分结果。例如, ANSYS提供的 PLANE82单元采 用的积分点数为2×2。关于高斯积分方法的更多内容参见,王勖成、邵敏编著“有限元方 法基本原理和数值方法”的4.5节数值积分方法。 将作用在单元上的外载荷同样表示为局部坐标的函数,就可以在局部坐标下完成单元的 载荷移置。体力移置的公式为 (R=L INT(P) dEdr (5-28) 面力移置的公式也类似,例如在5=1的边上受到面力作用, {R}°=[N={P=4(529)

（二）将整体坐标系中的面积积分转换为在局部坐标系中的面积积分   =  K B D B tdxdy e T [ ] [ ][ ] （5-23） 在整体坐标系中，面积微元为 x 方向和 y 方向微矢量的叉乘的模量， dA dx dy   =  （5-24）        d x d x dx   +   =        d y d y dy   +   =               d d x y x y d x d y d x d x dA ( ) ( ) ( )     −     =   +      +   =     dA = J dd （5-25） 代入（5-23），得到单元刚度矩阵在局部坐标系中的积分公式： K B B B D B B B t J dd e T [ ... ] [ ][ ... ] 1 2 8 1 2 8 1 1 1 −1 − = （5-26） 单元刚度矩阵中的任意一个分块矩阵的积分公式为， K B D Bs t J dd T rs r [ ] [ ] [ ][ ] 1 1 1 −1 − = （5-27） （三）用数值积分计算出单元刚度矩阵中的元素 等参单元刚度矩阵的每个元素都是局部坐标的函数，在有限元程序中不用解析的办法来 计算局部坐标系中的积分，而采用数值积分方法。通常采用高斯积分方法计算单元刚度矩阵 中的元素。 高斯积分方法预先定义了积分点和相应的加权系数，求出被积分的函数在指定积分点上 的数值，加权后求和，就得到了该函数的积分。这种方法具有比较高的计算精度。已经证明， 采用 n 个积分点的高斯积分可以达到 2n-1 阶的精度，也就是说，如果被积分的函数是 2n-1 次多项式，用高斯积分可以得到精确的积分结果。例如，ANSYS 提供的 PLANE82 单元采 用的积分点数为 2  2 。关于高斯积分方法的更多内容参见，王勖成、邵敏编著“有限元方 法基本原理和数值方法”的 4.5 节数值积分方法。 将作用在单元上的外载荷同样表示为局部坐标的函数，就可以在局部坐标下完成单元的 载荷移置。体力移置的公式为， R N p t J dd e T { } [ ] { } 1 1 1 −1 − = （5-28） 面力移置的公式也类似，例如在  = 1 的边上受到面力作用，    R N P t J d e T 1 1 1 1 { } [ ] { } − = = = （5-29）

在点(0,70)集中力移置的公式为, {R;=[M5m){P (5-30) 54六面体等参单元 多数弹性力学问题需要按照三维空间问题来求解。三维弹性力学问题的有限元法的基本 步骤与平面问题的步骤一样,包括单元离散化、选择单元位移模式、单元分析、整体分析和 方程求解。在分析三维问题时,所选择的单元主要为四面体单元和六面体单元。每个单元节 点上定义有三个位移分量u、V、w 三维问题有限元法有以下两个主要难点 1)单元划分比较复杂 无法采用人工方法完成复杂三维实体的单元划分,需要有功能强大的单元划分程序,从 CAD模型直接生成离散的单元网格。现在的有限元软件可以读入IGES、STL等格式的 图形交换文件。六面体单元的计算精度比较高,但是对于复杂三维实体无法实现六面体 单元的自动划分。采用四面体单元能够实现单元自动划分,但是四面体单元的计算精度 比较低 2)计算规模大 三维问题的单元数目大,节点自由度多,导致计算规模大,对计算机硬件的要求很高。 为缩短计算时间,有许多问题需要采用巨型计算机,如CRAY,或并行计算机 常用的三维等参单元有六面体八节点等参单元和六面体二十结点等参单元。等参单元的 位移模式和坐标变化式采用相同的形函数 =∑N(5,n.5)u1 v=∑N(5,, (5-31) N,(5,n,5)w x=∑N5,n5)x y=∑N(5,n5) (5-32) =∑N(5,n,5) ANSYS提供的Sold45单元就是六面体八节点等参单元,每个节点有代表x、y、z三 个方向位移的三个自由度(DOF, Degree of freedom),可以退化为五面体棱柱和四面体单 元,如图5-9所示。单元局部坐标为r、s、t,如图5-10所示,六面体八结点等参单元的基 本单元如图5-11所示

在点 ( , )  0 0 集中力移置的公式为， { } [ ] { } ( , ) 0 0 R N P e T =   （5-30） 5.4 六面体等参单元 多数弹性力学问题需要按照三维空间问题来求解。三维弹性力学问题的有限元法的基本 步骤与平面问题的步骤一样，包括单元离散化、选择单元位移模式、单元分析、整体分析和 方程求解。在分析三维问题时，所选择的单元主要为四面体单元和六面体单元。每个单元节 点上定义有三个位移分量 u、v、w。 三维问题有限元法有以下两个主要难点： 1）单元划分比较复杂 无法采用人工方法完成复杂三维实体的单元划分，需要有功能强大的单元划分程序，从 CAD 模型直接生成离散的单元网格。现在的有限元软件可以读入 IGES、STL 等格式的 图形交换文件。六面体单元的计算精度比较高，但是对于复杂三维实体无法实现六面体 单元的自动划分。采用四面体单元能够实现单元自动划分，但是四面体单元的计算精度 比较低。 2）计算规模大 三维问题的单元数目大，节点自由度多，导致计算规模大，对计算机硬件的要求很高。 为缩短计算时间，有许多问题需要采用巨型计算机，如 CRAY，或并行计算机。 常用的三维等参单元有六面体八节点等参单元和六面体二十结点等参单元。等参单元的 位移模式和坐标变化式采用相同的形函数， i i n i i i n i i i n i w N w v N v u N u ( , , ) ( , , ) ( , , ) 1 1 1             = = = = = = （5-31） i i n i i i n i i i n i z N z y N y x N x ( , , ) ( , , ) ( , , ) 1 1 1             = = = = = = （5-32） ANSYS 提供的 Solid45 单元就是六面体八节点等参单元，每个节点有代表 x、y、z 三 个方向位移的三个自由度（DOF，Degree of Freedom），可以退化为五面体棱柱和四面体单 元，如图 5-9 所示。单元局部坐标为 r、s、t，如图 5-10 所示，六面体八结点等参单元的基 本单元如图 5-11 所示