计算土力学
主讲教师:张爱军
§ 4.3 形函数的选择
? 位移模式 是表征单元内部位移的形状的
函数,位移模式一般用多项式表示,因
为多项式可以无限逼近任意连续函数。
其形式为:
22
1 2 3 4 5 61 2 3
22
4 5 6 7 8 9 1 0 1 1 1 2
22
7 8 9 1 3 1 4 1 5 1 6 1 7 1 8
u a a x a y a x y a x a yu a a x a y
v a a x a y v a a x a y a x y a x a y
w a a x a y w a a x a y a x y a x a y
? ? ? ? ? ? ?? ? ??
??
? ? ? ? ? ? ? ? ???
??? ? ? ? ? ? ? ? ?
? ?
或 等 等
? 形函数是表征单元内部任意一点的位移
与结点位移之间关系的函数。形函数由
位移模式推导产生。其一般形式为:
? ? ? ? ? ?
? ?
? ?
[ ],
e
N
N
??
?
?
?
e
其 中,, 为 单 元 内 部 任 一 点 的 位 移 向 量
,为 单 元 结 点 上 的 位 移 向 量
为 形 函 数, 也 是 插 值 函 数
位移模式的选择
?包含常数项,反映刚体位移
?包括常应变项。单元应变包括两部分,一部
分与结点的位置有关,称为变量应变;一部
分与结点的位移无关,称为常应变。当单元
尺寸较小时,单元应变趋于均匀,其常应变
量成为应变的主要部分,位移模式中必须包
括这个部分。包括常应变的意思在于位移模
式中必须包括 x,y,z的一次项
?位移在单元内部连续,在边界上协调,相邻
单元在边界上不开裂也不重叠
?满足几何对称性,即:在各个方向上均形式
一样,只是参数不同
? 位移模式的选择一般采用帕斯卡三角形
22
3 2 2 3
4 3 2 2 3 4
1
xy
x x y y
x x y x y y
x x y x y x y y
1ua?
1 2 3u a a x a y? ? ?
221 2 3 4 5 6u a a x a y a x a x y a y? ? ? ? ? ?
22
1 2 3 4 5 6
3 2 2 3
7 8 9 1 0
u a a x a y a x a x y a y
a x a x y a x y a y
? ? ? ? ? ?
? ? ? ?
22
1 2 3 4 5 6
3 2 2 3
7 8 9 1 0
4 3 2 2 3 4
1 1 1 2 1 3 1 4 1 5
u a a x a y a x a x y a y
a x a x y a x y a y
a x a x y a x y a x y a y
? ? ? ? ? ?
? ? ? ?
? ? ? ? ?
形函数的建立
? 将各个结点的坐标和各个结点的位移值
(在这里认为是确定值)代入位移模式,
得到一组联立的代数方程组,求解这个
方程,得到由结点坐标和位移值表示的
待定参数值,即:
(,,)
,,:
i i i
i i i
a A x y u
x y u
?
A
其 中, a, 为 待 定 参 数
,表 示 一 种 表 达 式
表 示 各 个 结 点 的 坐 标 和 位 移
? 等参元的形函数有所不同
§ 4.4 常应变三角形单元
? 以平面应变问题的三角形单元有限元为
例,阐明有限元解题的过程,从位移模
式的建立,到形函数的形成以及单元刚
度矩阵、总体刚度矩阵,求解等全过程。
? 三角形单元是最简单的有限单元,但是
其代表性较强,便于理解。
? 单元位移模式、形函数的建立
设:一个三角形单元, 其结点编号为 1,2,3。
相应结点坐标为 {x1,y1,x2,y2,x3,y3},而结点
的位移为 {u1,v1,u2,v2,u3,v3},这样单元有 6
个自由度 。
取:单元的位移模式为 ( 常应变 ),
1
2
3
v
u(x,y)
x
y v2
u2
1 2 3
4 5 6
u a a x a y
v a a x a y
? ? ??
? ? ? ?
?
其中,x,y为单元内
部任一点坐标
u,v 为该点沿 x,y
方向的位移
1 1 1 2 2 3 3
1 1 1 2 3 3 2
2 1 1 2 2 3 3 2 2 1 2 3
3 3 1 2 3
3 1 1 2 2 3 3
1
2
1
11
1 1,2,3
22
1
1
2
a u u u
x y x y x y
a u u u x y y y
x y x x
a u u u
? ? ?
?
? ? ? ?
?
? ? ?
?
??
?
?
???
?
? ? ? ? ? ??
? ?
??
?
??
?
? ?
= ( )
= ( )
= ( )
?将三个结点的坐标与位移值代入位移模式中得:
1 1 2 1 3 1
2 1 2 2 3 2
3 1 2 3 3 3
u a a x a y
u a a x a y
u a a x a y
? ? ??
?
? ? ??
? ? ? ?
?
?得到:三个待定系数的值为:
只与结点坐标
有关
1 4 5 1 6 1
2 4 5 2 6 2
3 4 5 3 6 3
v a a x a y
v a a x a y
v a a x a y
? ? ??
?
? ? ??
? ? ? ?
?
?将待定参数代入位移模式:
1 2 3
4 5 6
u a a x a y
v a a x a y
? ? ??
? ? ? ?
?
4 1 1 2 2 3 3
1 1 1 2 3 3 2
5 1 1 2 2 3 3 2 2 1 2 3
3 3 1 2 3
6 1 1 2 2 3 3
1
2
1
11
1 1,2,3
22
1
1
2
a v v v
x y x y x y
a v v v x y y y
x y x x
a v v v
? ? ?
?
? ? ? ?
?
? ? ?
?
??
?
?
???
?
? ? ? ? ? ??
?
?
??
?
??
?
? ?
= ( )
= ( )
= ( )
? ? ? ? ? ?
? ? ? ?
? ?
1 1 2 2 3 3
1 1 2 2 3 3
1
1
2
2
3
3
1 2 3
1 2 3
0 0 0
0 0 0
e
e
e
u N u N u N u
v N v N v N v
N
u
v
uu
vv
u
v
N N N
N
N N N
??
??
? ? ??
?
? ? ?
?
?
??
??
??
????
??
?
? ? ? ?
?? ??
??
??
????
??
?
??
??
e
其中,=
可以得到:
1 1 1 1
1 2 3 3 2
2 2 2 2 1 2 3
1 2 3
3 3 3 3
1
2
1
1,2,3
2
1
2
N x y
x y x y
N x y y y
xx
N x y
? ? ?
?
? ? ? ?
?
? ? ?
?
??
?
?
???
?
? ? ? ??
?
?
??
?
??
?
? ?
其中,= ( )
= ( )
= ( )
1
1,2,3
2i i i i
N x y? ? ???
?
用简化的方式表示为:
= ( ) i=
?分析:
?位移模式确定后就可以得到形函数
?位移模式是表征单元内部任意点的
位移值的函数,存在许多待定参数
?形函数是根据位移模式,建立单元
内部任意点的位移值与为单元结点
位移值之间关系的矩阵
?形函数是单元内任一点坐标的函数,
并且与单元结点坐标有关。
? 单元几何方程的建立
?对于平面问题,其几何方程为:
? ?
1
1
1 2 3 2
1 2 3 2
3
3
0
0
0
0 0 0
0
0 0 0
e
ex
e
y
xy
u
xx
uv
vyy
uv
y x y x
u
v
x
N N N u
N N N vy
u
yx v
?
??
?
? ? ? ?
??
? ? ? ?
??
? ? ? ???
? ? ? ? ??????
? ? ?
? ? ? ? ? ???
??
??? ? ? ?
??
??
?? ??
? ? ? ?
???
??
? ? ? ?
? ? ? ?
????
?
????
? ??
??
???? ??
? ??
?
???? ??
?
?? ??
??
????
??
????
??
???? ??
B N δe
代入形函数
1
1
1 2 3 2
1 2 3 2
3
3
1
312
1
2
312
2
3331 1 2 2
3
0
0 0 0
0
0 0 0
0 0 0
0 0 0
u
v
x
N N N u
N N N vy
u
yx v
u
NNN
v
x x x
uNNN
vy y y
uNNN N N N
y x y x y y v
????
?
????
?
????
???? ??
? ??
?
???? ??
?
?? ??
??
????
??
????
??
????
??
????
???
???
? ? ?
???
???
??? ?
?
???
? ? ?
???
???
??? ? ? ?
???
? ? ? ? ? ?
???
?
?
?
?
?
?
?
?
?
?
?
B
?得到最后由单元结点位移表示的几何
方程如下:
? ? ? ?? ?e eB???
1
1
1 2 3
2
1 2 3
2
1 1 2 2 3 3
3
3
000
1
0 0 0
2
u
v
u
v
u
v
? ? ?
? ? ?
? ? ? ? ? ?
??
??
??
??
????
??
? ??
?
??
??
??
??
??
????
1,2,3
1
2i i i iN x y? ? ???? (i = )= ( ) 求导
B
? 单元物理方程的建立
?对于平面应变弹性问题,有:
? ? ? ? ? ?
10
1
( 1 )
10
( 1 ) ( 1 2 ) 1
12
00
2( 1 )
ee
e
x
y
xy
D
E
??
?
?
?
??
?
? ? ?
?
?
?
?
??
??
?
?? ??
??? ??
? ??
??
? ? ?
????
??
?? ?
??
???
? 单元刚度矩阵形成
?由单元刚度矩阵的计算公式得到:
? ? ? ? ? ? ? ?
? ? ? ? ? ?
? ? ? ? ? ?
1
e
Te e
v
T
T
k B D B d V
B D B
B D B
?
? ? ?
??
? 单位宽度
2
11
11
1 2 3
22
1 2 3
22
1 1 2 2 3 3
36
33
3333
63
11 12 13 14 15 16
22 23 24 25 2
12
2 ( 1 )
1
1
( 1 ) 1
()
( 1 ) ( 1 2 ) 2
0
10
0
000
0
1 0 0 0 0
0
0
00
0
E
k k k k k k
k k k k k
?
?
?
?
?
?
?
??
??
??
? ? ?
??
? ? ?
??
? ? ? ? ? ?
??
??
?
?
?
?
?
?
?
?
??
? ? ?
??
??
??
??
??
?? ??
??
?? ??
?
??
?? ??
??
?? ??
??
??
??
??
??
????
??
?
? ? ? ? ? ?
? ? ? ? ? ?
? ? ? ? ? ?
11 12 13
21 22 23
31 31 33
6
22
33 34 35 36
44 45 46
55 56
66
66
1k k k
k k k
k k k
k k k k
k k k
kk
k
?
?
??
??
??
??
??
??
?
??
??
??
??
??
??
??
??
??
结 点
结 点 2
结 点 3


两个自由度
?单元刚度矩阵的性质
?对称性:由弹性力学中弹性体的功
互等性质决定。并且行数与列数相

?奇异性:其每一行的和为 0,由于在
无约束情况下,单元位移不定决定
的。 —— 刚体运动。
?单元刚度矩阵的行、列数 =单元节点
数 × 单元结点自由度数 = 3× 2= 6
? 自由度就是未知量的个数。
? 整体分析
?用一个有 4个单元的例子分析
1
2
4
3
5 6




y
x
体力 p
单元与结点编
号如右图,在
单元上作用有
体力 p,在单元
2上作用有面力
q。
q
?每个单元可以形成一个单元刚度矩阵,
将其叠加就得到了总体刚度矩阵。
?总体刚度矩阵的行数和列数= 总结点
数 × 单元结点自由度数 的方阵,其排
列的顺序按照结点标号依次从小到大
排列,每个结点的按照自由度的顺序
排列。上例中总刚度矩阵排列见后,
这样排列是为了与单元结点位移向量
对应。
? ?
1
1
2
2
3
3
4
4
5
5
6
6
1 1 2 2 3 3 4 4 5 5 6 6
...
u
v
u
v
u
v
K
u
v
u
v
u
v
u v u v u v u v u v u v
? ? ? ? ? ? ? ? ? ? ? ???
??
? ? ? ? ? ? ? ? ? ? ?
???
?
??
?
??
?
?
??
?
???
?
??
???
??
??
?? ?
??


结点 1
结点 1
?单元刚度矩阵组合时,就是“对号入
座”,即按照单元在整体结构上的总
体编号(有别于单元自身的结点编号)
找相应的位置叠加,也就是将相应的
单元刚度矩阵上的值 kij,加到上面即
可。
?也可以理解为将所有单元刚度矩阵也
按照整体刚度矩阵的行列数及位置,
扩大 到整体刚度矩阵的大小,用矩阵
加法得到整体刚度矩阵。
?以上题为例说明:
? ?
1 1 1
1 2 3 1 3 2 2 3
1 3 4 3 4 4
22
2 3 4 4
4
1000
20
30
40
5
6
1 2 3 4 5 6
K
? ? ? ?
? ? ?
??
??? ? ?
??
? ? ? ?
?? ? ? ?
?
????
??
??
???


1
2
4
3
5 6




?1+2+3的意思是将单元 ①②③ 单元刚度矩阵上的相应
元素相加 。
?为了方便, 每个结点只列一项, 不再列出各结点的各
个自由度相应项
?相关结点才有叠加, 不相关结点相应值为 0。
?整体刚度矩阵的特点
?对称性
?奇异性
?稀疏性( 0多)
?带状性(集中在主元周围,两边为 0)
?最大半带宽=(结点编号最大差值
+ 1) × 结点自由度数=( 6- 1+ 1)
× 2= 12
因此需要按结点差值最小编整体
结点号,方能够保证整刚最小。
?整刚变带宽一维存储
? ?
1 1 1
1 2 3 1 3 2 2 3
1 3 4 3 4 4
22
2 3 4 4
4
1000
20
30
40
5
6
1 2 3 4 5 6
K
? ? ? ?
? ? ?
??
??? ? ?
??
? ? ? ?
?? ? ? ?
?
????
??
??
???


? ? 1 1 1 1 2 3 1 1 6K ????? ? ? ? ??? ×
在我书上有详细说明, 大家看一下 。
? 等效结点矩阵的形成
?由变分可以推出:
? ? ? ?{ } [ ] [ ]
e
e T T
SR N p d V N q d S? ??eV +
?本例为:
? ?
1
1
1
2
22
3
3
3
{ } [ ]
0
0
0
00
0
0
0
0
0
e
eT
p
e
v
e
e
R N p dV
N
N pdV
N
N
dV
N pdVN p
N
N
N pdV
?
??
?? ??
?? ??
?? ??
?? ?? ??
???? ? ? ? ?
???? ??
?? ??
?? ??
?? ????
??
?
??
??
??
??
e
V
?将各个单元等效结点向量,用与整刚集合的
方法进行叠加,就得到相应的整体等效结点
向量。
? 约束处理
?将整刚中,存在支承结点的相应支承自由度
的行和列均划去,使该列行所有元素均为 0;
?将相应等效结点向量相应行元素也充 0
?这样就将支撑条件引入了方程中。
?对于已知结点位移的结点,可在相应整刚主
元中乘以大数,将相应等效结点向量元素等
于给定位移乘以该大数即可。
? 方程求解
?高斯消元法
?波前法
?三角分解法
?迭代法等等
?在线性代数中讲过,这里不再讲授,
在我书里有三角分解法的详细说明。
? 应力求解
?整体方程求解出来以后,可以得到各
个结点的各个自由度上的位移值
?结点应力可以根据结点位移值得到
{ } [ ] { } [ ] [ ] { } [ ] { }eeD D B S? ? ? ?? ? ?
三角形单元单元内部任意一点的应变相
等,也就是应力相等,由于各个单元可
能应力不等,在单元结点上的应力不连
续,因此一般不计算结点应力,而认为
单元应力就是单元形心处的应力。
Any questions?