代入各节点位移及坐标，可整理成
其中形状函数
10
10
10
∑ ∑ ∑ u = ui ， v = vi ， w = wi
i =1
i =1
i =1
角点： N1 = (2L1 − 1)L1 ， N2 = (2L2 − 1)L2 ， N3 = (2L3 − 1)L3 ， N4 = (2L4 − 1)L4 ；
边中点： N5 = 4L1L2 , N 6 = 4L1L3 , N 7 = 4L1L4 , N8 = 4L3 L4 , N9 = 4L2 L3 , N10 = 4L2 L4 .
(1
+
ξiξ
)(1
+
ηiη
)(1
+
ζ
iζ
)
其中 (ξi ,ηi ,ζ i ) 是 i 节点的局部坐标。并且规定在局部坐标系中，面 2376 上ξ = 1，面 1485
上ξ = −1；面 3487 上η = 1，面 1265 上η = −1；面 5678 上 ς = 1，面 1234 上ς = −1 。
Li
=
Vi V
, Lj
=
Vj V
, Lm
=
Vm V
, Lp
=
Vp V
为体积坐标。其中 V 是四面体 ijmp 的体积， Vi、Vj、Vm、Vp 分别是四面体 ojmp、oimp、 oijp 和 oijm 的体积。
热动3同学修改
王正伟 13601363209
11 1 1 V = 1 xi x j xm xp ，
A
则单元节点力矩阵
{ } { } ∑ {P}e = { } ( Pp e + PF e + Pq e ) .
高次四面体单元 常应变四面体单元中的应力分量都是常量，难以适应急剧变化的应力场。为了保证计算 精度，一方面可以分割更加密集的常应变四面体单元，另一方面可以采用高次四面体单元。 体积坐标 用于空间问题的高次四面体单元通常采用 体积坐标。 如图 2，定义
x j xm xp
111
ai = y j ym yp ， bi = − y j ym yp ，
z j zm zp
z j zm zp
x j xm xp
x j xm xp
ci = − 1 1 1 ， di = − y j ym yp 。
z j zm zp
111
热动3同学修改
王正伟 13601363209
u = a1 + a2x + a3 y + a4z + a5x2 +a6y2 +a7z2 + a8xy + a9 yz + a10zx; v = b1 + b2x + b3 y + b4z + b5x2 +b6y2 +b7 z2 + b8xy + b9 yz + b10zx; w = c1 + c2x + c3 y + c4z + c5x2 +c6y2 +c7z2 + c8xy + c9 yz + c10zx.
2. 六面体单元
在平面问题中可以采用比三角形单元更复杂的矩形单 元一样，在空间问题中也可以采用比四面体更复杂六面体 单元。 图 5 所示为 8 节点正六面体单元，它是研究高次六 面体单元和六面体等参单元的基础。
单元位移模式取
u = a1 + a2 x + a3 y + a4 z + a5 xy + a6 yz + a7 zx + a8 xyz; v = b1 + b2 x + b3 y + b4 z + b5 xy +b6 yz +b7 zx + b8 xyz; w = c1 + c2 x + c3 y + c4 z + c5 xy +c6 yz +c7 zx + c8 xyz.
(3L4
− 1)(3L4
−
2)L4
；
边三等分点：
N5
=
9 2
L1L2 (3L1
− 1)
，
N6
=
9 2
L1L2 (3L2
− 1)
，
N7
=
9 2
L1L3 (3L1
− 1)
，
N8
=
9 2
L1L3 (3L3
− 1)
，
N9
=
9 2
L1L4 (3L1
− 1)
，
N10
=
9 2
L1L4 (3L4
− 1)
，
N11
=
9 2
热动3同学修改
王正伟 13601363209
形状函数
而且满足
∑8
u = N1u1 + N2u2 +" + N8u8 = Niui
⎫ ⎪
i =1
⎪
∑8
v = N1v1 + N2v2 +" + N8v8 = Nivi
⎪ ⎬
i =1
⎪
8
⎪
∑ w =
N1w1
+
N 2 w2
+" +
N8 w8
=
i =1
Ni
wi
二、单元应变矩阵 空间问题中，每一点有 6 个应力分量
( { }ε e = ε x ε y ε z γ xy γ yz
) γ T zx
= [L][N ]{δ }e
= [B]{δ }e ，
其中
⎜⎛ ⎜
∂ ∂x
⎜ ⎜
0
⎜
⎜0
[L] =
⎜ ⎜
∂
⎜ ∂y ⎜
⎜0
⎜
0
∂ ∂y
0
∂ ∂x ∂ ∂z
0 ⎟⎞
⎟
0
⎟ ⎟
[ ] 可见矩阵 Bi 中的元素都是常量，所以采用线性位移函数的四面体单元是常应变单元。
三、单元应力矩阵
{σ }e = [D]{ε}e = [D][B]{δ }e ，
弹性矩阵
热动3同学修改
王正伟 13601363209
⎜⎛ 1 a a 0 0 0⎟⎞
⎜a 1 a 0 0 0⎟
[D] =
E(1 − µ)
6 yi y j ym y p zi z j zm z p
11 1 1
11 1 1
Vi
=
1 6
x y
xj yj
xm ym
xp yp
，Vj
=
1 6
xi yi
x y
xm ym
xp ， yp
z z j zm zp
zi z zm z p
1 111
11 11
Vm
=
1 6
xi yi
xj yj
x y
xp yp
，Vp
可得到形状函数矩阵
[N
]
=
⎜⎛ ⎜
N1 0
0 N1
0 0
N2 0 0 N2
0 0
" N8 0 " 0 N8
0 ⎟⎞ 0⎟
⎜⎝ 0 0 N1 0 0 N2 " 0 0 N8 ⎟⎠
得到位移插值函数之后，利用应力应变关系得出单元刚度矩阵并合成整体刚度矩阵，列 出等效节点载荷列阵，给出边界条件，即可解空间有限元问题。分析方法同四面体单元，不 再详诉。
代入各节点位移及坐标，可整理成
20
20
20
∑ ∑ ∑ u = ui ， v = vi ， w = wi .
i =1
i =1
i =1
其中形状函数
角点：
N1
=
1 2
(3L1
− 1)(3L1
−
2)L1
，
N2
=
1 2
(3L2
− 1)(3L2
−
2)L2
，
N3
=
1 2
(3L3
− 1)(3L3
−
2)L3
，
N4
=
1 2
第 7 章 空间问题的有限元分析
空间问题的有限元分析方法与平面问题是一致的。首先是将结构离散化，分割成有限 个单元；然后进行单元特性计算，确定单元刚度矩阵；最后进行整体计算，将各单元的刚度 矩阵集成整体刚度矩阵，引入边界条件，解方程得到各节点位移及其它所求物理量。下面将 重点介绍单元特性分析过程。
1. 四面体单元
以下分析方法与常应变四面体单元相同。
热动3同学修改
王正伟 13601363209
20 节点四面体单元 分别取角点（4 节点）、边三等分点（12 节点）和面形心（4 节点）共 20 节点。如图 4 所示。 20 节点四面体单元的单元位移函数取完全三次多项式
u = a1 + a2 x + a3 y + a4z + a5x2 +a6y2 +a7z2 + a8xy + a9 yz + a10 zx + a11x3 + a12 y3 + a13z3 + a14 x2 y + a15x2 z + a16 y2 x + a17 y2z + a18z2x + a19 z2 y + a20 xyz
单元位移矩阵
其中形状函数矩阵
⎧u ⎫
{u}e
=
⎪ ⎨
v
⎪ ⎬
=
[N
]{δ
}e
.
⎪⎩w⎪⎭
[N