FE-Ch04.1-5等参元与数值积分

形函数Ni (ξ,η)在单元的i结点上的值为1， 在其它结点上的值均为0。
坐标变换式采用如下相似的公式，
x = ∑ N i (ξ ,η ) xi
i =1 8 8
（4-12）
y = ∑ N i (ξ ,η ) y i
i =1
将ξ=1代入公式（4-12），可以得到单元345边在 ξ=1 4-12 345 整体坐标下的参数方程：
N3 =
1 (1 − ξ )(1 + η ) 4
（4-5）
1 N 4 = (1 + ξ )(1 + η ) 4
四个结点的坐标为
(ξ i ,ηi )
，定义新的变量，
(i=1，2，3，4)
ξ 0 = ξ iξ ， η 0 = ηiη
(4-6)
形态函数表示为，N i
1 = (1 + ξ 0 )(1 + η 0 ) (i=1，2，3，4) 4
如果单元不是等参的, 即坐标插值:
x = ∑ N i (ξ ,η ) xi
i =1 m
y = ∑ N i (ξ ,η ) yi
i =1
m
中的节点数m和插值函数Ni, 各自不等于未知函数 n ’: φ插值中的节点数n和插值函数Nk ϕ = N ' (ξ ,η )ϕ
∑
k =1
k
k
这时,可以分两种情况: (1)超参单元, 即坐标插值节点数m>未知函数φ插值节 点数n, 单元一般不满足完备性的要求 (2)次参单元, 即m<n, 从构造变节点单元的一般方法, 假定一2D等参单元在各节点有线性变化的场函数, Φ=a+bx+cy 其在各节点有对应的场函数值, Φi=a+bxi+cyi (i=1,2,…,n)
u = a1 + a 2ξ + a3η + a 4ξ 2 + a5ξη + a 6η 2 + a 7ξ 2η + a8ξη 2
（5-9）
v = b1 + b2ξ + b3η + b4ξ 2 + b5ξη + b6η 2 + b7ξ 2η + b8ξη 2
该位移模式实际上是一个双二次函数，待定系数由结 点位移分量确定。在单元的每条边上，局部坐标 ξ=±1或η=±1，位移是局部坐标ξ或η的二次函数， 完全由边上的三个结点的位移值确定，所以这个位移 模式满足位移连续性条件。实际单元内的位移用形函 数表示为，
[K ]
e
= ∫∫ [ B ] [ D ][ B ]tdxdy
T
e
单元的应变为，
{ε } = [ B ]{δ }
{δ }
e
单元的结点位移，
= [u1 v1 u 2 v 2 ... u8 v8 ]
T
将形函数代入后，可以得到应变的矩阵表达式，
∂N1 ∂x {ε } = 0 ∂N 1 ∂y 0 ∂N1 ∂y ∂N1 ∂x ∂N 2 ∂x ∂N 2 ∂y ∂N 2 ∂y ∂N 2 ∂x ∂N 8 ... ∂x ... ∂N 8 ... ∂y u1 v1 ∂N 8 . ∂y . ∂N 8 u 8 ∂x v 8
Chap 4 等参元与数值积分
§4.1 引言 本章包括以下内容： 本章包括以下内容： •等参单元的基本概念 等参单元的基本概念 •四边形八节点等参单元 四边形八节点等参单元 •等参单元的单元分析 等参单元的单元分析 •六面体等参单元 六面体等参单元
§4.2 等参单元的基本概念与刚度矩阵变换 一、 等参单元的基本概念 在进行有限元分析时，单元离散化会带来 计算误差，主要采用两种方法来降低单元离散 化产生的误差： 1） 提高 单元划 分的密 度，被 称为h方 法（ hmethod）； 2）提高单元位移函数多项式的阶次，被Байду номын сангаас为p 方法（p-method）。
u = ∑ N i (ξ ,η )u i
i =1
8
v = ∑ N i (ξ ,η )vi
i =1
8
（4-10）
其中的形函数为：
N1 = 1 (1 − ξ )( 1 − η )( − ξ − η − 1 ) 4
N N N N
2
1 = (1 − ξ 2 = = = 1 (1 − ξ 2 1 (1 − η 2 1 (1 − η 2
在平面问题的有限单元中，我们可以选择四 结点的矩形单元，如图所示，该矩形单元在x 及y方向的边长分别为2a和2b。
同第三章的方法类似，将单元的位移模式选为，
u = a1 + a 2 x + a3 y + a 4 xy
v = a5 + a6 x + a7 y + a8 xy
可得到，
u = N i ui + N j u j + N m u m + N p u p
i =1 i =1 k =1 k =1 n m n ' 注意：x=∑ N i xi = ∑ N k xk = a ∑ N i + bx + cy i =1 k =1 i =1 n n m m
n
m
只要满足
∑N
j=1
n
j
=1
单元一般满足完备性的要求
§4.4 等参单元的一般格式
一、平面问题的四边形八结点等参单元 在本节，以平面问题的四边形八结点等参单元为 例，介绍构造等参单元的单元刚度矩阵的基本过程。 弹性力学平面问题的单元刚度矩阵为，
v = N i vi + N j v j + N m vm + N p v p
形态函数为，
1 x y N i = (1 − )(1 − ) 4 a b
1 x y N j = (1 + )(1 − ) 4 a b
1 x y N m = (1 + )(1 + ) 4 a b
1 x y N p = (1 − )(1 + ) 4 a b
x = aη + bη + c
2
y = dη + eη + f
2
（4-13）
可见在整体坐标系中，单元的边是一条抛物线或 退化为一条直线。
图4-6 ANSYS提供的Plane82单元
如图4-6所示，ANSYS提供的PLANE82单元是一个四 边形八结点等参单元，局部坐标定义为s和t，如图4-7 所示。PLANE82单元可以退化为三角形六结点单元。
参照矩形单元，四结点正方形单元的位移模式为，
u = N1u1 + N 2 u 2 + N 3u3 + N 4 u 4
v = N1v1 + N 2 v 2 + N 3v3 + N 4 v4
（4-4）
其中，
N1 =
N2 =
1 (1 − ξ )(1 − η ) 4
1 (1 + ξ )(1 − η ) 4
图4-2任意四结点四边形单元
图4-3四结点正方形单元
在图4-2所示的任意四边形单元上，用等分四条边 的两族直线分割四边形，以两族直线的中心为原点，建 立局部坐标系(ξ,η)，沿 ξ 及 η 增大的方向作为ξ轴 和η轴，并令四条边上的ξ及η值分别为±1。为了求 出位移模式，以及局部坐标与整体坐标之间的变换式， 在局部坐标系中定义一个四结点正方形单元，如图4-3 所示。
按照第3章构造插值函 N i' = ∑ Cij N j x j = ∑ Cij xi j=1 i =1 数的一般方法可得到: 其中Cij是常系数. 表示单元内任一点的坐标可用指 定点坐标重新插值得到, 因此可变换成:
φ = ∑ N iφi ( 注意：φi =a + bxi + cyi ) = a ∑ N i + b∑ N k' xk + c ∑ N k' yk
上述单元位移模式满足位移模式选择的基本要求： 反映了单元的刚体位移和常应变，单元在公共边界上 位移连续。 在矩形单元的边界上，坐标x和y的其中一个取常 量，因此在边界上位移是线性分布的，由两个结点上 的位移确定。 与三结点三角形单元相比，四结点矩形单元的位 移模式是坐标的二次函数，能够提高计算精度，但也 有显著的缺点，两种单元的比较如下:
2
)( 1 − η ) )( 1 + η ) )( 1 + ξ ) )( 1 − ξ )
N3 = N5 =
N
1 (1 + ξ )( 1 − η )( ξ − η − 1 ) 4 1 (1 + ξ )( 1 + η )( ξ + η − 1 ) 4
2
6
2
4
7
1 = (1 − ξ )( 1 + η )( − ξ + η − 1 ) 4
（4-7）
把 ξ 及 η 作为任意四边形单元的局部坐标，把（4-4） 的位移模式和（4-7）的形态函数用于任意形状的四边单 元，可得： 1.在四个结点处可以得到结点的位移； 2.在单元的四条边上，位移线性变化，保证了单元公共边 界上位移的连续性。
因此给出任意四边形单元的结点位移就能得到整个 单元上的位移，（4-4）的位移模式就是所要找的正确 的位移模式。 把局部坐标与整体坐标的变换式也取为，
2
8
将形函数归纳为，
1 (i = 1,3,5,7) 4 (1 + ξ iξ )(1 + η iη )(ξ iξ + η iη ) 1 （4-11） N i (ξ ,η ) = (1 − ξ 2 )(1 + η iη ) (i = 2,6) 2 1 (1 − η 2 )(1 + ξ iξ ) (i = 4,8) 2
图4-7 Plane82的基本单元 ANSYS理论手册中给出的PLANE82单元的位移模式如图 4-8所示，位移模式与公式（4-10）展开后是一样的。
§4.3 等参变换的条件和等参单元的收敛性
一、等参变换的条件 （教材122-123），为保证等参变换的一一对应性质， 应当避免单元一边上的两点退化为一个节点。还要防 止单元的任意两边的夹角接近180o，更不允许夹角等 于、大于180o。
把图4-3中的局部坐标系中的正方形单元 正方形单元称为基本 正方形单元 基本 单元。 单元 把图4-2中的在整体坐标系中的任意四边形单元 任意四边形单元看 任意四边形单元 作由基本单元通过坐标变换得来的，称为实际单元 实际单元。 实际单元 单元几何形状和单元内的未知量采用相同数目的 结点参数以及相同的插值函数进行变换，称为等参变 换。采用等参变换的单元，称为等参单元 等参单元。 等参单元 由于形态函数，正好反映了单元形状的变化，也称为 形函数（Shape function）。 形函数 采用等参单元，使我们可以在局部坐标系中的规 则单元上进行单元分析，然后在映射到实际单元上。 等参单元同时具有计算精度高和适用性好的特点，是 有限元程序中主要采用的单元形式。
表4-1 三结点三角形单元与四结点矩形单元比较
单元类型 三结点三角 形单元 四结点矩形 单元 优点 适应复杂形 状，单元大 小过渡方便 单元内的应 力、应变是 线性变化的， 计算精度较 高 缺点 计算精度低 多项式阶数 一阶线性
不能适应曲 线边界和非 正交的直线 边界
二阶双线性
如果任意形状的四边形四结点单元采用矩形单元的位 移模式，则在公共边界上不满足位移连续性条件。为了既 能得到较高的计算精度，又能适应复杂的边界形状，可以 采用坐标变换。
正常
不正常
正常
不正常
二、等参单元的收敛性 在单元分析已经提出有限单元解的收敛性要求, 即, 单元必须是完备的和协调的。对于等参单元: 1.完备性:对于C0型单元,由于等参单元的形函数中包 含有常数项和线性项,满足完备性的要求。 2. 协调性:由于单元之间的公共边上有完全相同的节 点, 同时每一单元沿这些边的坐标和未知函数均采用 相同的插值函数加以确定。因此, 只要在划分网格时, 遵守单元选择和节点配置的要求, 则等参单元满足协 调性的要求。
x = N1 x1 + N 2 x 2 + N 3 x3 + N 4 x 4
y = N1 y1 + N 2 y 2 + N 3 y3 + N 4 y 4
（4-8）
将坐标变换式用于任意四边形单元，可得： (1)在四个结点处给出结点的整体坐标， (2)在四条边上的整体坐标是线性变化的。 只要给出任意四边形单元四个结点的整体坐标，用 （4-8）式就可以建立局部坐标系中的正方形单元和整 体坐标系中的任意四边形单元之间的坐标变换关系。
二、 等参单元的矩阵变换
为了更好地反映物体内的应力变化，适应曲线边界， 在弹性力学平面问题的分析中经常使用四边形八节点等 参单元。如图5-4所示，由于每条边上增加了一个结点， 单元的边是一条二次曲线，可以更好地适应曲线边界，
图4-4四边形八结点单元
图4-5 八结点基本单元
对于等参单元，先在图5-5所示的八结点基本单元上进 行分析。八结点单元一共有16个已知的结点位移分量， 基本单元中取如下的位移模式：