2000.4 哈尔滨建筑大学 王焕定教授制作 1
第五章 大变形问题的有限单元法
1,弹性大变形问题的有限元法
2,弹性分支点稳定问题有限元分析
3,物质描述大变形增量问题的 T.L,
U.L法
2000.4 哈尔滨建筑大学 王焕定教授制作 2
1,弹性大变形问题的有限元法
弹性大变形问题,需要考虑变形的非线性项
和变形对平衡的影响。
若以初始自然平衡状态作初始位形,则物质
描述的格林应变为
L
ijij
j
k
i
k
j
i
i
j
ij X
u
X
u
X
u
X
u
E ?? ??
?
?
?
?
?
?
?
?
?
?
? )(
2
1
式中
)(
2
1
j
i
i
j
ij X
u
X
u
?
?
?
?
?
??
j
k
i
kL
ij X
u
X
u
?
?
?
?
?
2
1?
线性部分 非线性部分
2000.4 哈尔滨建筑大学 王焕定教授制作 3
为便于计算机编程,将张量转换为矩阵:
格林应变矩阵和张量的分量间有如下关系
? ? ? ? T312312332211 222 EEEEEEE ?
对应的克希荷夫应力矩阵和张量分量间关系为
? ? ? ? T312312332211 222 ??????? ?
? ? ? ? T312312332211 222 LLLLLLL ??????? ?
? ? ? ? T312312332211 SSSSSSS ?引入两个算子矩阵
? ?
T
?
?
?
?
?
?
?
?
?
?
?
123
312
321
000
000
000
XXX
XXX
XXX
A
d,d,d,
d,d,d,
d,d,d,
2000.4 哈尔滨建筑大学 王焕定教授制作 4
式中
iX
d,
?
??
iX
? ? uIII
j
i
T
??
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
3
3
2
3
1
3 XXXX
u
?
? ?
T
123
312
321
000
000
000
?
?
?
?
?
?
?
?
?
?
?
XXX
XXX
XXX
L
uuu
uuu
uuu
A
,,,
,,,
,,,
iX
,
?
?? uu
iX
? ? T321 uuu?u
再引入位移梯度向量的记号
? ?H
3阶单位矩阵
2000.4 哈尔滨建筑大学 王焕定教授制作 5
在上述符号基础上,格林应变由位移表为
? ? ? ? ? ? ? ? ? ? ? ???? LL AuAE
2
1????
则单元格林应变为 eNu ??
? ? eeLe BBBE ??? ~???
其中线性应变矩阵 B和有限元 (I)一样
设单元位移场和有限元 (I)一样为
ANB ?
GAHNAB LLL
2
1
2
1 ??
非线性部分“应变矩阵”为
LBBB ??~
2000.4 哈尔滨建筑大学 王焕定教授制作 6
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
3
3
3
3
2
3
3
1
3
2
3
2
2
3
2
1
3
1
3
1
2
3
1
1
III
III
III
HNG
m
m
m
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
??
??
??式中 G为如下 9× 3m的矩阵
eG ?? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 7
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
T
1
T
3
T
2
T
3
T
1
T
2
T
3
T
2
T
1
0
0
0
00
00
00
??
??
??
?
?
?
L
A
式中 AL为如下 6× 9的矩阵
i
i
u
X?
???
2000.4 哈尔滨建筑大学 王焕定教授制作 8
由 (AL)可见,格林应变 -位移关系是非线性的。
非线性部分因为 ? ? eδBε ??δ
? ? ? ? θAθA LL ?? ?
所以
为用虚位移原理建立单元特性方程,还得建
立应变增量和位移增量间的关系。对线性部分
综上所述,格林应变增量为
? ? eLLe δGABθAδBE ???? )( ????
验证
? ? ? ? ? ? ? ? eLLLL δGAθAθAε ???? ??? )
2
1(
如果记 )2( LBBB ??? ? ?
eδBE ?? ??
,则 。
2000.4 哈尔滨建筑大学 王焕定教授制作 9
对弹性问题,在物质描述下本构关系为
在上述基础上,由虚位移原理可得
e
V
e
V
FFBSES
e e
δ)(dδd TeEejTT ???? ? ???
0 0
00 V
?V
式中 为单元结点力矩阵,为单元等效结
点荷载矩阵 eFj
eFE
pqk l p qkl EDS o?
j
q
i
p
mn i j
n
l
m
k
k l p q x
X
x
X
D
x
X
x
X
JD
?
?
?
?
?
?
?
?
?o
? ? ? ?? ?EDS o?
由于应变、应力以矩阵表示,因此弹性矩阵
应按下式并考虑应变矩阵定义来建立 oD
2000.4 哈尔滨建筑大学 王焕定教授制作 10
Pd和 PE分别为直接和等效结点荷载矩阵,R为
综合等效结点荷载矩阵。上式也可写为
RPPSB
e V e
???? ? EdT d
0
0V
?
将克希荷夫应力表达式代入,可得
?? ??
ee SV
e φNFNF
?
000
0
SV dd TTE
按集成规则集装后可得
? ? ?????
e
T )(d)(
eV
e RUUKRBDBU
0
000 V~? ??
? ? ???
e
T d)(
eV
RSBU
0
00V??
2000.4 哈尔滨建筑大学 王焕定教授制作 11
? ?
e V e
BDB
0
0V
~? d0T
根据本构方程,则有
式中 K(U)是非对称的,为
? ? ???
e
TT )d()dd(d)(d
eV
T UKSBSBU
0
0V
???
eddd ?BDEDS TT ???
对非线性弹性问题, 和 都是位移的函数。B? 0D B~
根据非线性方程切线刚度矩阵的定义,可得
又因
GAB LL
2
1? GABBBB LL ???? )2(?
B,G为已知矩阵
2000.4 哈尔滨建筑大学 王焕定教授制作 12
式中
? ?? ? ?
e Ve V ee
SAGSB
00
00 VV
? dddd TLTT
?
?
?
?
?
?
?
?
?
?
?
123
312
321
000
000
000
???
???
???
ddd
ddd
ddd
d
TL
A
所以
又因,所以? ? T
312312332211 SSSSSS?S
?
???
???
???
d
ddd
ddd
ddd
d
T
MSA ?
?
?
?
?
?
?
?
?
?
?
??
??
??
?
311232333
233121222
313122111
SSS
SSS
SSS
L
2000.4 哈尔滨建筑大学 王焕定教授制作 13
由此可得
式中
?
?
?
?
?
?
?
?
?
?
?
333323331
323322312
331312311
III
III
III
M
SSS
SSS
SSS
ed
d
d
d
d ?
?
?
?
? G?
?
?
?
?
?
?
?
?
?
?
?
3
2
1
e
e Ve V ee
MGGSB ?dddd TT ? ?? ? ?
00
00 VV
?
基于上述说明,可得
? ? ???
e
TT )dd(d)(d
eV
SBSBU
0
0V
???
UUKBDBMGG Te
e V
T
e
)d(d)d( TT ??? ? ? ?
0
0V
??
2000.4 哈尔滨建筑大学 王焕定教授制作 14
如果引入如下记号
则单元切线刚度矩阵为
? ???
e
eee
e )( LT kkkK ?
? ???
eV
LTLTLLT
e
L BDBBDBBDBk
0
0V)d(
TTT
??
eV
e MGGk
0
0Vd
T
?
GABB LLL ?? 2 LBBB ???
??
eV
T
e BDBk
0
0Vd
T
e
初应力或几何刚度矩阵 线弹性刚度矩阵
大位移刚度矩阵
eLeeeT kkkk ??? ?e
“结构”切线刚度矩阵为
建立了切线
刚度矩阵,
用牛顿法等
即可求解。
2000.4 哈尔滨建筑大学 王焕定教授制作 15
本节的讨论没涉及具体单元,因此具有普
遍性。具体单元分析时,因形函数一般是自
然坐标的函数,故需作坐标变换后代入相关
公式,从而建立具体单元的切线刚度矩阵。
需要指出的是
建议自行对各种单元自行推导切线矩阵。
本节只讨论了全量形式的弹性大变形分析,
具体求解步骤如讲义所示。为了保证收敛,
拟用增量迭代法。对第二类稳定问题(极值
点失稳问题)、弹塑性问题等,必须用 3,所
介绍的增量形式来解决。
To
26
在我们、王勋成、谢贻权、徐次
达等的有限元教材中,都有一些
具体单元的切线刚度矩阵,需要
时可供参考。
2000.4 哈尔滨建筑大学 王焕定教授制作 16
2,弹性分支点稳定问题有限元分析
对分支点稳定问题,关键是建立几何刚度矩
阵。此时,以失稳前的平衡位置作初始位形,
以失稳形态作现时位形。
和前述弹性大变形不同的是:
大位移矩阵可忽略。
应变仅包含失稳位移的非线性项。
分支点处相应失稳位移的综合荷载为零。
注意到上述差异,即可用上节结果解决分支
点稳定的有限元分析。
失稳前变形是微小的。
2000.4 哈尔滨建筑大学 王焕定教授制作 17
设变形前单元长度为 l,截面积为 A,弹性模
量为 E。单元杆端位移矩阵为 δe
式中
2.1 桁架单元
? ? T22122111 uuuu?e? eNu ??
。单元位移为
? ? 为二阶单位阵,)1( 222 IIIN ????
式中 为失稳位移。基于此,单元格林应变
为 2iu
2
1
2
12111 2
1)(
2
1 ?? ????,u,uE
eδ??
??
?
? ?? 0101
1 LL? ??
??
?
? ?? 0101
LL
B
2000.4 哈尔滨建筑大学 王焕定教授制作 18
克希荷夫应力为
因分支点稳定关心的临界荷载,临界荷载时
平衡有两重性,故临界状态应力等于失稳前状
态的应力,也即
基于上述结果,单元几何刚度矩阵为
)
2
1( 2
111 ?? ??? EEES
101 SS ?
eδ??
??
?
? ??
LL
1010? ?
?
??
?
? ??
LL
1010G
,失稳前应力为
1110 ?? AENES ?? 内力为
?? ?? L
V
e GGMGGk
e
0 1
T
0
T dd
0
xASV?
2000.4 哈尔滨建筑大学 王焕定教授制作 19
设变形前单元长度为 l,截面积和惯性矩为 A、
I,弹性模量为 E。单元杆端位移矩阵为 δe
单元位移为
式中
2.2 梁单元
? ? T12221111,vvu,vvu?e? i)
d
d(
1 X
v,v
i ?
e
N
N
N
u ?
??
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
2
2
1
d
d
d
d
X
Y
v
X
v
Yu
? ?0000 411 NN?N
? ?65322 00 NNNN?N
2000.4 哈尔滨建筑大学 王焕定教授制作 20
式中 Ni和有限元 (I)一样。梁单元应变为
其中第一项为有限元 (I)里的线性项,非线性的
第二项为
266 ?? ???a
? ? ee δδG caba ??? 00?
2
1
2
2
2
1 2
1)
d
d(
2
1
d
d
d
d ?? ?????
X
v
X
vY
X
uE
)34(1 2?? ??? lb )32( 2?? ??? lc
基于上述结果,单元几何刚度矩阵为
??? 101 SS
?? ?? L
V
e GGMGGk
e
0 1
T
0
T dd
0
xASV?
象桁架单元说明一样,单元应力为
为初始状态的轴向应力 ?
同结构
力学
2000.4 哈尔滨建筑大学 王焕定教授制作 21
有限元 (I)里的二维问题单元位移为
2.3 二维单元 (薄板稳定 )
? ? Twvu?u
? ?Tvu?u
但失稳时的位形,将有出平面的位移 w,和杆
单元一样,应变需考虑 w及其非线性项
2
2
2
11 )(2
1
X
w
X
wZ
X
uE
?
??
?
??
?
??
2
2
2
22 )(2
1
Y
w
Y
wZ
Y
uE
?
??
?
??
?
??
Y
w
X
w
YX
wZ
X
v
Y
uEE
?
?
?
??
??
??
?
??
?
??? 2
1221 2
非
线
性
项
2000.4 哈尔滨建筑大学 王焕定教授制作 22
设出平面位移为
面内弹性应力为 ewN ??w
? ? T122211 ????ζ
??
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
22
2
2
1
11
2
1
1
X
N
X
N
X
N
X
N
X
N
X
N
w
m
ww
w
m
ww
HNG
?
?
对出平面位移,其 G矩阵为
对应的 M矩阵为
??
?
?
??
?
?
?
2212
1211
??
??
M
2000.4 哈尔滨建筑大学 王焕定教授制作 23
对各种具体单元,将 Nw的具体形函数代入 G矩
阵的表达式,即可积分得到具体单元的几何刚
度矩阵。
基于上述结果,单元几何刚度矩阵为
??
eV
e MGGk
0
0
T d At
?
对于分支点稳定问题,结构力学已经指出,
在比例加载下,最终归结为一个特征值问题
0?? UKK )( ??
解得特征值后,即可得到临界荷载,一般只
关心最小临界荷载。
2000.4 哈尔滨建筑大学 王焕定教授制作 24
3.1 物质描述大变形增量问题的 T.L法
对弹塑性、粘性 -蠕变和施工力学等问题,介
质的反应和变形的历史有关。对随时间变化的
荷载,需要将时间变量离散成序列,以求解各
时刻的响应。为此,都需要用增量法来解决。
从 t到 t+Δt的增量期间进行物质描述求解时,
一般可选两种参考位形:初始和 t时刻的位形。
前者称为全拉格朗日 (T.L)表述,后者称为修正
拉格朗日 (U.L)表述。
设从 0到 t时刻的全部反应、位形均已求得,
现在的问题是,如何求 t+Δt时刻的响应。
2000.4 哈尔滨建筑大学 王焕定教授制作 25
未
知
设 t0,t和 t+Δt的物理量分别用如下符号标记
对 T.L法,介质位移是初始位形坐标的函数
iii x,x,X
坐标,密度 ???,,
0
面积和体积 V,V,VA,A,A
00 和
iiiiii XxuXxu - - ?? 和
)(-)(- iiiiiii XxXxuuu ???
设有限元分析时单元形状描述为
eNxx ?eNXX ? exNx ?又设有限元分析时单元位移场为
eNu ??eNu ?? eNu ??? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 26
i
k
i
k N
X
N
??
??
?
? ? 1J
式中雅可比矩阵 J为
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
???
???
???
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
3
3
2
3
1
3
3
2
2
2
1
2
3
1
2
1
1
1
???
???
???
J
象有限元 (I)等参元分析一样,由于形函数一
般是对自然坐标 定义的,因此有限元分析中
的对坐标求导等,应象有限元 (I)一样进行转换。i?
To
15
To
36
2000.4 哈尔滨建筑大学 王焕定教授制作 27
在上述记号下,格林应变为
时刻 t +Δt和 t的应变增量为
)(
2
1
j
k
i
k
i
j
j
i
ij X
u
X
u
X
u
X
u
E
?
?
?
?
?
?
?
?
?
?
?
nijijijijijij EEEEEE ???? ????? 10式中
)(
2
10
i
j
j
i
ij X
u
X
u
E
?
?
?
?
?
?
??
?
)(
2
1
j
k
i
k
i
j
j
i
ij X
u
X
u
X
u
X
u
E
?
?
?
?
?
?
?
?
?
?
?
j
k
i
kn
ij X
u
X
uE
?
?
?
?? ???
2
1
)(
2
11
j
k
i
k
j
k
i
k
ij X
u
X
u
X
u
X
uE
?
?
?
??
?
?
?
?? ???
2000.4 哈尔滨建筑大学 王焕定教授制作 28
要强调指出的是,u在增量步内已知,因此
同大变形有限元,将张量转换为矩阵,则
成线性关系与和 10 iijij uEE ???
引入大变形所用算子记号,则有
uAE ?? ?0
成非线性关系与 inij uE ??
? ? ? ? T312312332211 222 EEEEEE?E
? ? ? ? T312312332211 222 EEEEEE?E
nEEEE ???? ??? 10
EEE ???
??? LAE ?1 ???? L
n AE 2
1?
eGuH ?? ??? ??
2000.4 哈尔滨建筑大学 王焕定教授制作 29
由于 u在增量步内已知,因此
B和 BL都是已知的。又若记 e
BE ??? ?0
eenL BBBBE ?? ??? ???? )(
en
L
n BAE ?? ???? ?? 2
1
eLL BAE ?? ??? ??1
GABB Lnn ??? 2则有 GAB Ln ??2
因为,因此? ? ? ? θAθA LL ?????? ?
)()()( enLn δBθAE ??????? ??综上可得
eenL δBδBBBE ??????
????? )()(
2000.4 哈尔滨建筑大学 王焕定教授制作 30
为进行有限元列式,还需讨论克希荷夫应力。
设 t和 t+Δt时刻的应力分量分别为
kl
l
j
k
i
ij x
X
x
XJS ?
?
?
?
??
SSS ???
基于上述分析,利用 t+Δt时刻的虚位移原理
虚功方程
kl
l
j
k
i
ij x
X
x
XJS ?
?
?
?
??
则象应变分析一样,可将 分成 。同
样换为矩阵表示,则有 ijS ijij SS ??
? ? 等 T312312332211 SSSSSS?S
?? ? ??
ee e SV V
uφuFES
00 0
0
T
00
T
00
T ddd
?
??? SVV
2000.4 哈尔滨建筑大学 王焕定教授制作 31
再次强调,t时刻及其前的量都是已知的,因
此变分为零。基于此
将此结果代入虚功方程,可得单元刚度方程
)()( eδNuu ????? ??
? ??
eV
e
ee
je δFFδBS
0
T
E0
T )(d ????V?
式中 和 是对初始位形定义的,t+Δt时刻的
体积力和表面力,它们是已知的。 0F 0?
)()( eδBEE ????? ???
式中
?? ??
ee SV
e φNFNF
00
00
T
00
T
E dd
?
SV
2000.4 哈尔滨建筑大学 王焕定教授制作 32
t时刻应力引起的
等效结点力矩阵
t+Δt时
刻荷载
引起的
等效结
点力矩
阵或
? ? ???
e V e
RPPSB
0
0Ed0
T )(d V?
将其按集成规则集装后可得
再引入如下记号
? ? ??
e V
L
e
VSBBR
0
0
T d)(
?
? ? ???
e V
e
e
RSBδ
0
000T d)( V???
? ?? ???
e SV ee
NFNPR )dd(
00
00
T
00
T
d0
?
? SV
将 和 的表达式代入,可得B? S
? ? ?????
e V
nL
e
RSSBBB
0
000T )d()( V?
2000.4 哈尔滨建筑大学 王焕定教授制作 33
几何或非线性应变
增量刚度矩阵
或
利用这些关系,非线性平衡方程可写为
? ? ? ??
e V V
L
n
e e
SAGSB
0 0
00 VV dd
TTT ?
? ??
e V e
MGGK
0
0
T d V
?
? ? ?? ??
e e V
e
V ee
MGGMG
00
00 ?? ?? VV dd
TT
? ? ?????
e V e
RRUKSB
0
000T d ?? ?? V?
? ? ???
e V
e
e
RSB
0
00
T d)( V????
为求解上述方程,尚需解决如下两方面问题
2000.4 哈尔滨建筑大学 王焕定教授制作 34
首先假设
?
?
?
EE
E
T EDS
?
? d
然后将 ΔS和 ΔE的关系线性化。根据本构关系
EDS tT ?? ?
? ?? ? ??
e V
L
e V ee
SBBSB
00
0
T
0
T d)(d VV ???
EDS T dd ?
则有
?
?
?
?
ep
D
D
D TT
弹塑性
非线性弹性
为使其线性化,设 ( t时刻的材料性质矩阵 )t
TD
2000.4 哈尔滨建筑大学 王焕定教授制作 35
在做了上述两方面处理后,可得
由此出发,用非线性方程的相关解法,即可解
决大变形非线性(材料非线性)问题。
? ?? ? ??
e V
L
e V ee
SBBSB
00
0
T
0
T d)(d VV ???
将其代回非线性平衡方程,可得
? ? ??
e V
L
e
EDBB
0
0T
T d)( V?
UKBBDBB Ne
e V
LL
e
?? ???? ? ? ?
0
0T
T )d()( V
0????? 0)()( RRUKKUΨ N ?? ??
讲义上给出了 T.L法的求解步骤,可供大家编
程序参考。
2000.4 哈尔滨建筑大学 王焕定教授制作 36
3.2 物质描述大变形增量问题的 U.L法
因 t+Δt的位移是用 t时刻位形为基准度量的,
因此
在 [t,t+Δt]间隔内,以 t时刻位形为参考位形,
其增量位移为
iii xxu ???
ii uu ?? )( ii uu ??? ?象 T.L法一样,设单元和位移的描述为
eNxx ? exNx ? eNu ??? ?但需指出的是,式中形函数 N是 t时刻单元自然
坐标的函数。在计算 等导数时,要先作坐
标变换 (Xi应换为 xi)。 jxi,N
To
26
2000.4 哈尔滨建筑大学 王焕定教授制作 37
类似地,用矩阵来表示则有
在时刻 t和 t+Δt的格林应变是以 t时刻位形定义
的,因而它们分别为
ij
j
k
i
k
j
i
i
j
ij Ex
u
x
u
x
u
x
u
E ?
????
?
?
?
?
?
?
?
?
?
?
?
? )(
2
1
NL EEE ??? ?? eL BE ??? ?
式中算子矩阵象 T.L法一样,但应将 Xi换为 xi。
0?ijE
ene
LL
N BGAAE ??? ?????? ??? 2
1
2
1
enBBE ??? )( ??
2000.4 哈尔滨建筑大学 王焕定教授制作 38
再次强调,式中算子符号象 T.L法一样,但应将
Xi换为 xi。
基于上述说明,象 T.L法一样可导得
)()(2)( enenN δBδBE ?????? ??
)()( eL δBE ???? ?
ijijS ??
)()()()( een δBδBBE ?????? ????
关于应力的处理也和 T.L法一样,对 t 时位形
定义的 t和 t+Δt时刻的克希荷夫应力分别为
kl
l
j
k
i
ij x
x
x
xJS ?
?
?
?
??
SSSS ?? ???? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 39
象 T.L法一样由虚位移原理虚功方程可导得
象 T.L法一样推导,引入如下符号定义
? ??
eV
e
ee
je δFFδBS ????
T
E
T )(d V?
式中
?? ??
ee SV
e φNFNF
?
SV dd TTE
其中 分别为 t时刻位形定义的单
元体积、应力表面、体力和表面力。 φFSV
ee,、,?
几何或非线性应变
增量刚度矩阵? ??
e V e
MGGK 0T d V?
荷载 引起的等效结点力矩阵
? ?? ???
e SV ee
φNFNPR )dd( TTd
?
SV
2000.4 哈尔滨建筑大学 王焕定教授制作 40
则可得
t时刻应力引起的等效结点力矩阵
? ?? ? ??
e Ve V ee
VBVSBR dd TT ??
? ? ?????
e V e
RRUKSBUΨ 0?? ??? Vd)( T?
t+ Δt时刻的非线性平衡方程
象 T.L法一样,为求解上述方程也需解决线性
化问题。首先讨论 ΔS的计算。因为
ij
k
j
l
i
lkij x
X
x
XJS ?? ?
?
?
?
??
2000.4 哈尔滨建筑大学 王焕定教授制作 41
m
m
klkl
l
j
k
i
ij x
v
x
X
x
XJS
?
??
?
?
?
?? ?? ??(
)
m
k
lm
m
l
km x
v
x
v
?
??
?
?? ??
因此
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
可改写为
由第四章已知 Chap4 64
klijk lij VD
ep???
式中
kijkkjikijij ?????? ???
??
Chap4 44
)(,tVE ijtij ?
Chap4 26
ij
k
j
l
i
lkij x
X
x
XJS ?? ?
?
?
?
??
2000.4 哈尔滨建筑大学 王焕定教授制作 42
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
可得 kijkkjikijij ?????? ??? ?
?
由如下两式消去 并利用,且
注意到 Vij对称,Ωij反对称
?ij?
ijijxi Vv j ??? ?,
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ??
GoTo 49
klijk lij VD
ep???又由于, 和,因
此 ijij VE ?
?
klkl
k
k V
x
v ??
?
?
?? ? klijk lij EDS
GoTo 50
klijlijkljiki j kli j kl SSSDD ??? ???? ep
2000.4 哈尔滨建筑大学 王焕定教授制作 43
上式最后一项将使本构张量不对称,对金属
类不可压缩介质,这一项可略去,也即
在有限的克希荷夫应力和格林应变增量之间仍
认为
这就是 U.L法的本构关系线性化。
)( ep kijkkjikkli j klij ESESEDS ???? ???
lijkljiki j kli j kl SSDD ?? ??? ep
与 T.L法一样,除本构关系线性化外,还需解
决几何方面的线性化。因为
? ?? ? ?
e Ve V
e
ee
SESBδ VV d)(d)( TTT ?????? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 44
为对上述三项积分作几何方面的线性化,可设
这一关系可有两种矩阵表达方式:其一是
eL BEE ???? ??其另一方式是
? ?? ? ???
e V
kli j k lij
e V
ijij
ee
VEDEVSE d)((d)( ep ??????
)d)(d)( ?? ??
ee V
kijkij
V
kjikij VESEVESE ??????
)(
2
1
j
i
i
j
ij x
u
x
u
E
?
?
?
?
?
?
??
?
? ? T333231131211,,,,,,EEEEEEL ???????? ??? EE
2000.4 哈尔滨建筑大学 王焕定教授制作 45
eeL BNAEE ?? ???? ???
引入讲义上式 (5,102)算子矩阵,则有A
? ?
e V e
SB VdT ??
基于上述讨论,积分 时,作如
下处理
?? ? ?? ?
e Ve V ee
BDBSB VV d(d epTT ??
UBMBBMB
ee VV
?)dd TT ?? ?? VV
式中 为由 组成的 6× 6矩阵,矩阵为
epD epijklD
M
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
S
S
S
M
00
00
00
?
?
?
?
?
?
?
?
?
?
?
?
333231
232221
131211
SSS
SSS
SSS
S
2000.4 哈尔滨建筑大学 王焕定教授制作 46
结
束
引入如下两个矩阵
? ??
e V
L
e
BDBK VdepT
? ? ??
e V
N
e
BMBK Vd2 T
? ? ?????
e V e
RRUKSBUΨ 0?? ??? Vd)( T?
则 t+ Δt时刻的非线性平衡方程
0?????? RRUKKKUΨ NL ?? ?? )()(
可改写成
需要指出的是,要提高精度可减小 Δt或采用
适当的修正技术。
2000.4 哈尔滨建筑大学 王焕定教授制作 47
验证 δALθ=ALδθ
? ? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
3
2
1
T
1
T
3
T
2
T
3
T
1
T
2
T
3
T
2
T
1
0
0
0
00
00
00
θ
θ
θ
θθ
θθ
θθ
θ
θ
θ
θA
L
??
??
??
?
?
?
?
2000.4 哈尔滨建筑大学 王焕定教授制作 48
? ? θA
θθθθ
θθθθ
θθθθ
θθ
θθ
θθ
θθθθ
θθθθ
θθθθ
θθ
θθ
θθ
L
?
??
??
??
?
?
?
??
??
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
3
T
11
T
3
3
T
22
T
3
2
T
11
T
2
3
T
3
2
T
2
1
T
1
3
T
11
T
3
3
T
22
T
3
2
T
11
T
2
3
T
3
2
T
2
1
T
1
纯量的转置不变
返
回
2000.4 哈尔滨建筑大学 王焕定教授制作 49
推证式( A)
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
kijkkjikijij ????? ? ???
??
k
k
ijkijkkjikijij x
vS
?
????? ?? ?????
?
)()( ikikjkjkjkik VV ?????? ????
k
k
ijkijkkjikij x
v
?
????? ? ?????
?
)()( kikijkkjkjik VV ???? ????
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ?? 返
回
2000.4 哈尔滨建筑大学 王焕定教授制作 50
式( 5,96)的推证
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ??
kijkkjikklklijkli j kl VVVVD ???? ???? ep
kijkkjikklklijkli j kl VSVSVSVD ???? ?ep
kllijkklljikklklijkli j kl VSVSVSVD ??? ???? ep
kllijkljikkliji j kl VSSSD )( ep ??? ????
????? ijlijkljikkliji j kl ESSSD )( ep ???
?? ?
kli jk lij EDS lijkljikkliji j kli j kl SSSDD ??? ????
ep 返
回
第五章 大变形问题的有限单元法
1,弹性大变形问题的有限元法
2,弹性分支点稳定问题有限元分析
3,物质描述大变形增量问题的 T.L,
U.L法
2000.4 哈尔滨建筑大学 王焕定教授制作 2
1,弹性大变形问题的有限元法
弹性大变形问题,需要考虑变形的非线性项
和变形对平衡的影响。
若以初始自然平衡状态作初始位形,则物质
描述的格林应变为
L
ijij
j
k
i
k
j
i
i
j
ij X
u
X
u
X
u
X
u
E ?? ??
?
?
?
?
?
?
?
?
?
?
? )(
2
1
式中
)(
2
1
j
i
i
j
ij X
u
X
u
?
?
?
?
?
??
j
k
i
kL
ij X
u
X
u
?
?
?
?
?
2
1?
线性部分 非线性部分
2000.4 哈尔滨建筑大学 王焕定教授制作 3
为便于计算机编程,将张量转换为矩阵:
格林应变矩阵和张量的分量间有如下关系
? ? ? ? T312312332211 222 EEEEEEE ?
对应的克希荷夫应力矩阵和张量分量间关系为
? ? ? ? T312312332211 222 ??????? ?
? ? ? ? T312312332211 222 LLLLLLL ??????? ?
? ? ? ? T312312332211 SSSSSSS ?引入两个算子矩阵
? ?
T
?
?
?
?
?
?
?
?
?
?
?
123
312
321
000
000
000
XXX
XXX
XXX
A
d,d,d,
d,d,d,
d,d,d,
2000.4 哈尔滨建筑大学 王焕定教授制作 4
式中
iX
d,
?
??
iX
? ? uIII
j
i
T
??
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
3
3
2
3
1
3 XXXX
u
?
? ?
T
123
312
321
000
000
000
?
?
?
?
?
?
?
?
?
?
?
XXX
XXX
XXX
L
uuu
uuu
uuu
A
,,,
,,,
,,,
iX
,
?
?? uu
iX
? ? T321 uuu?u
再引入位移梯度向量的记号
? ?H
3阶单位矩阵
2000.4 哈尔滨建筑大学 王焕定教授制作 5
在上述符号基础上,格林应变由位移表为
? ? ? ? ? ? ? ? ? ? ? ???? LL AuAE
2
1????
则单元格林应变为 eNu ??
? ? eeLe BBBE ??? ~???
其中线性应变矩阵 B和有限元 (I)一样
设单元位移场和有限元 (I)一样为
ANB ?
GAHNAB LLL
2
1
2
1 ??
非线性部分“应变矩阵”为
LBBB ??~
2000.4 哈尔滨建筑大学 王焕定教授制作 6
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
3
3
3
3
2
3
3
1
3
2
3
2
2
3
2
1
3
1
3
1
2
3
1
1
III
III
III
HNG
m
m
m
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
??
??
??式中 G为如下 9× 3m的矩阵
eG ?? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 7
? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
T
1
T
3
T
2
T
3
T
1
T
2
T
3
T
2
T
1
0
0
0
00
00
00
??
??
??
?
?
?
L
A
式中 AL为如下 6× 9的矩阵
i
i
u
X?
???
2000.4 哈尔滨建筑大学 王焕定教授制作 8
由 (AL)可见,格林应变 -位移关系是非线性的。
非线性部分因为 ? ? eδBε ??δ
? ? ? ? θAθA LL ?? ?
所以
为用虚位移原理建立单元特性方程,还得建
立应变增量和位移增量间的关系。对线性部分
综上所述,格林应变增量为
? ? eLLe δGABθAδBE ???? )( ????
验证
? ? ? ? ? ? ? ? eLLLL δGAθAθAε ???? ??? )
2
1(
如果记 )2( LBBB ??? ? ?
eδBE ?? ??
,则 。
2000.4 哈尔滨建筑大学 王焕定教授制作 9
对弹性问题,在物质描述下本构关系为
在上述基础上,由虚位移原理可得
e
V
e
V
FFBSES
e e
δ)(dδd TeEejTT ???? ? ???
0 0
00 V
?V
式中 为单元结点力矩阵,为单元等效结
点荷载矩阵 eFj
eFE
pqk l p qkl EDS o?
j
q
i
p
mn i j
n
l
m
k
k l p q x
X
x
X
D
x
X
x
X
JD
?
?
?
?
?
?
?
?
?o
? ? ? ?? ?EDS o?
由于应变、应力以矩阵表示,因此弹性矩阵
应按下式并考虑应变矩阵定义来建立 oD
2000.4 哈尔滨建筑大学 王焕定教授制作 10
Pd和 PE分别为直接和等效结点荷载矩阵,R为
综合等效结点荷载矩阵。上式也可写为
RPPSB
e V e
???? ? EdT d
0
0V
?
将克希荷夫应力表达式代入,可得
?? ??
ee SV
e φNFNF
?
000
0
SV dd TTE
按集成规则集装后可得
? ? ?????
e
T )(d)(
eV
e RUUKRBDBU
0
000 V~? ??
? ? ???
e
T d)(
eV
RSBU
0
00V??
2000.4 哈尔滨建筑大学 王焕定教授制作 11
? ?
e V e
BDB
0
0V
~? d0T
根据本构方程,则有
式中 K(U)是非对称的,为
? ? ???
e
TT )d()dd(d)(d
eV
T UKSBSBU
0
0V
???
eddd ?BDEDS TT ???
对非线性弹性问题, 和 都是位移的函数。B? 0D B~
根据非线性方程切线刚度矩阵的定义,可得
又因
GAB LL
2
1? GABBBB LL ???? )2(?
B,G为已知矩阵
2000.4 哈尔滨建筑大学 王焕定教授制作 12
式中
? ?? ? ?
e Ve V ee
SAGSB
00
00 VV
? dddd TLTT
?
?
?
?
?
?
?
?
?
?
?
123
312
321
000
000
000
???
???
???
ddd
ddd
ddd
d
TL
A
所以
又因,所以? ? T
312312332211 SSSSSS?S
?
???
???
???
d
ddd
ddd
ddd
d
T
MSA ?
?
?
?
?
?
?
?
?
?
?
??
??
??
?
311232333
233121222
313122111
SSS
SSS
SSS
L
2000.4 哈尔滨建筑大学 王焕定教授制作 13
由此可得
式中
?
?
?
?
?
?
?
?
?
?
?
333323331
323322312
331312311
III
III
III
M
SSS
SSS
SSS
ed
d
d
d
d ?
?
?
?
? G?
?
?
?
?
?
?
?
?
?
?
?
3
2
1
e
e Ve V ee
MGGSB ?dddd TT ? ?? ? ?
00
00 VV
?
基于上述说明,可得
? ? ???
e
TT )dd(d)(d
eV
SBSBU
0
0V
???
UUKBDBMGG Te
e V
T
e
)d(d)d( TT ??? ? ? ?
0
0V
??
2000.4 哈尔滨建筑大学 王焕定教授制作 14
如果引入如下记号
则单元切线刚度矩阵为
? ???
e
eee
e )( LT kkkK ?
? ???
eV
LTLTLLT
e
L BDBBDBBDBk
0
0V)d(
TTT
??
eV
e MGGk
0
0Vd
T
?
GABB LLL ?? 2 LBBB ???
??
eV
T
e BDBk
0
0Vd
T
e
初应力或几何刚度矩阵 线弹性刚度矩阵
大位移刚度矩阵
eLeeeT kkkk ??? ?e
“结构”切线刚度矩阵为
建立了切线
刚度矩阵,
用牛顿法等
即可求解。
2000.4 哈尔滨建筑大学 王焕定教授制作 15
本节的讨论没涉及具体单元,因此具有普
遍性。具体单元分析时,因形函数一般是自
然坐标的函数,故需作坐标变换后代入相关
公式,从而建立具体单元的切线刚度矩阵。
需要指出的是
建议自行对各种单元自行推导切线矩阵。
本节只讨论了全量形式的弹性大变形分析,
具体求解步骤如讲义所示。为了保证收敛,
拟用增量迭代法。对第二类稳定问题(极值
点失稳问题)、弹塑性问题等,必须用 3,所
介绍的增量形式来解决。
To
26
在我们、王勋成、谢贻权、徐次
达等的有限元教材中,都有一些
具体单元的切线刚度矩阵,需要
时可供参考。
2000.4 哈尔滨建筑大学 王焕定教授制作 16
2,弹性分支点稳定问题有限元分析
对分支点稳定问题,关键是建立几何刚度矩
阵。此时,以失稳前的平衡位置作初始位形,
以失稳形态作现时位形。
和前述弹性大变形不同的是:
大位移矩阵可忽略。
应变仅包含失稳位移的非线性项。
分支点处相应失稳位移的综合荷载为零。
注意到上述差异,即可用上节结果解决分支
点稳定的有限元分析。
失稳前变形是微小的。
2000.4 哈尔滨建筑大学 王焕定教授制作 17
设变形前单元长度为 l,截面积为 A,弹性模
量为 E。单元杆端位移矩阵为 δe
式中
2.1 桁架单元
? ? T22122111 uuuu?e? eNu ??
。单元位移为
? ? 为二阶单位阵,)1( 222 IIIN ????
式中 为失稳位移。基于此,单元格林应变
为 2iu
2
1
2
12111 2
1)(
2
1 ?? ????,u,uE
eδ??
??
?
? ?? 0101
1 LL? ??
??
?
? ?? 0101
LL
B
2000.4 哈尔滨建筑大学 王焕定教授制作 18
克希荷夫应力为
因分支点稳定关心的临界荷载,临界荷载时
平衡有两重性,故临界状态应力等于失稳前状
态的应力,也即
基于上述结果,单元几何刚度矩阵为
)
2
1( 2
111 ?? ??? EEES
101 SS ?
eδ??
??
?
? ??
LL
1010? ?
?
??
?
? ??
LL
1010G
,失稳前应力为
1110 ?? AENES ?? 内力为
?? ?? L
V
e GGMGGk
e
0 1
T
0
T dd
0
xASV?
2000.4 哈尔滨建筑大学 王焕定教授制作 19
设变形前单元长度为 l,截面积和惯性矩为 A、
I,弹性模量为 E。单元杆端位移矩阵为 δe
单元位移为
式中
2.2 梁单元
? ? T12221111,vvu,vvu?e? i)
d
d(
1 X
v,v
i ?
e
N
N
N
u ?
??
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
2
2
1
d
d
d
d
X
Y
v
X
v
Yu
? ?0000 411 NN?N
? ?65322 00 NNNN?N
2000.4 哈尔滨建筑大学 王焕定教授制作 20
式中 Ni和有限元 (I)一样。梁单元应变为
其中第一项为有限元 (I)里的线性项,非线性的
第二项为
266 ?? ???a
? ? ee δδG caba ??? 00?
2
1
2
2
2
1 2
1)
d
d(
2
1
d
d
d
d ?? ?????
X
v
X
vY
X
uE
)34(1 2?? ??? lb )32( 2?? ??? lc
基于上述结果,单元几何刚度矩阵为
??? 101 SS
?? ?? L
V
e GGMGGk
e
0 1
T
0
T dd
0
xASV?
象桁架单元说明一样,单元应力为
为初始状态的轴向应力 ?
同结构
力学
2000.4 哈尔滨建筑大学 王焕定教授制作 21
有限元 (I)里的二维问题单元位移为
2.3 二维单元 (薄板稳定 )
? ? Twvu?u
? ?Tvu?u
但失稳时的位形,将有出平面的位移 w,和杆
单元一样,应变需考虑 w及其非线性项
2
2
2
11 )(2
1
X
w
X
wZ
X
uE
?
??
?
??
?
??
2
2
2
22 )(2
1
Y
w
Y
wZ
Y
uE
?
??
?
??
?
??
Y
w
X
w
YX
wZ
X
v
Y
uEE
?
?
?
??
??
??
?
??
?
??? 2
1221 2
非
线
性
项
2000.4 哈尔滨建筑大学 王焕定教授制作 22
设出平面位移为
面内弹性应力为 ewN ??w
? ? T122211 ????ζ
??
?
?
?
?
?
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
??
22
2
2
1
11
2
1
1
X
N
X
N
X
N
X
N
X
N
X
N
w
m
ww
w
m
ww
HNG
?
?
对出平面位移,其 G矩阵为
对应的 M矩阵为
??
?
?
??
?
?
?
2212
1211
??
??
M
2000.4 哈尔滨建筑大学 王焕定教授制作 23
对各种具体单元,将 Nw的具体形函数代入 G矩
阵的表达式,即可积分得到具体单元的几何刚
度矩阵。
基于上述结果,单元几何刚度矩阵为
??
eV
e MGGk
0
0
T d At
?
对于分支点稳定问题,结构力学已经指出,
在比例加载下,最终归结为一个特征值问题
0?? UKK )( ??
解得特征值后,即可得到临界荷载,一般只
关心最小临界荷载。
2000.4 哈尔滨建筑大学 王焕定教授制作 24
3.1 物质描述大变形增量问题的 T.L法
对弹塑性、粘性 -蠕变和施工力学等问题,介
质的反应和变形的历史有关。对随时间变化的
荷载,需要将时间变量离散成序列,以求解各
时刻的响应。为此,都需要用增量法来解决。
从 t到 t+Δt的增量期间进行物质描述求解时,
一般可选两种参考位形:初始和 t时刻的位形。
前者称为全拉格朗日 (T.L)表述,后者称为修正
拉格朗日 (U.L)表述。
设从 0到 t时刻的全部反应、位形均已求得,
现在的问题是,如何求 t+Δt时刻的响应。
2000.4 哈尔滨建筑大学 王焕定教授制作 25
未
知
设 t0,t和 t+Δt的物理量分别用如下符号标记
对 T.L法,介质位移是初始位形坐标的函数
iii x,x,X
坐标,密度 ???,,
0
面积和体积 V,V,VA,A,A
00 和
iiiiii XxuXxu - - ?? 和
)(-)(- iiiiiii XxXxuuu ???
设有限元分析时单元形状描述为
eNxx ?eNXX ? exNx ?又设有限元分析时单元位移场为
eNu ??eNu ?? eNu ??? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 26
i
k
i
k N
X
N
??
??
?
? ? 1J
式中雅可比矩阵 J为
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
???
???
???
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
m
ii
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
X
N
3
3
2
3
1
3
3
2
2
2
1
2
3
1
2
1
1
1
???
???
???
J
象有限元 (I)等参元分析一样,由于形函数一
般是对自然坐标 定义的,因此有限元分析中
的对坐标求导等,应象有限元 (I)一样进行转换。i?
To
15
To
36
2000.4 哈尔滨建筑大学 王焕定教授制作 27
在上述记号下,格林应变为
时刻 t +Δt和 t的应变增量为
)(
2
1
j
k
i
k
i
j
j
i
ij X
u
X
u
X
u
X
u
E
?
?
?
?
?
?
?
?
?
?
?
nijijijijijij EEEEEE ???? ????? 10式中
)(
2
10
i
j
j
i
ij X
u
X
u
E
?
?
?
?
?
?
??
?
)(
2
1
j
k
i
k
i
j
j
i
ij X
u
X
u
X
u
X
u
E
?
?
?
?
?
?
?
?
?
?
?
j
k
i
kn
ij X
u
X
uE
?
?
?
?? ???
2
1
)(
2
11
j
k
i
k
j
k
i
k
ij X
u
X
u
X
u
X
uE
?
?
?
??
?
?
?
?? ???
2000.4 哈尔滨建筑大学 王焕定教授制作 28
要强调指出的是,u在增量步内已知,因此
同大变形有限元,将张量转换为矩阵,则
成线性关系与和 10 iijij uEE ???
引入大变形所用算子记号,则有
uAE ?? ?0
成非线性关系与 inij uE ??
? ? ? ? T312312332211 222 EEEEEE?E
? ? ? ? T312312332211 222 EEEEEE?E
nEEEE ???? ??? 10
EEE ???
??? LAE ?1 ???? L
n AE 2
1?
eGuH ?? ??? ??
2000.4 哈尔滨建筑大学 王焕定教授制作 29
由于 u在增量步内已知,因此
B和 BL都是已知的。又若记 e
BE ??? ?0
eenL BBBBE ?? ??? ???? )(
en
L
n BAE ?? ???? ?? 2
1
eLL BAE ?? ??? ??1
GABB Lnn ??? 2则有 GAB Ln ??2
因为,因此? ? ? ? θAθA LL ?????? ?
)()()( enLn δBθAE ??????? ??综上可得
eenL δBδBBBE ??????
????? )()(
2000.4 哈尔滨建筑大学 王焕定教授制作 30
为进行有限元列式,还需讨论克希荷夫应力。
设 t和 t+Δt时刻的应力分量分别为
kl
l
j
k
i
ij x
X
x
XJS ?
?
?
?
??
SSS ???
基于上述分析,利用 t+Δt时刻的虚位移原理
虚功方程
kl
l
j
k
i
ij x
X
x
XJS ?
?
?
?
??
则象应变分析一样,可将 分成 。同
样换为矩阵表示,则有 ijS ijij SS ??
? ? 等 T312312332211 SSSSSS?S
?? ? ??
ee e SV V
uφuFES
00 0
0
T
00
T
00
T ddd
?
??? SVV
2000.4 哈尔滨建筑大学 王焕定教授制作 31
再次强调,t时刻及其前的量都是已知的,因
此变分为零。基于此
将此结果代入虚功方程,可得单元刚度方程
)()( eδNuu ????? ??
? ??
eV
e
ee
je δFFδBS
0
T
E0
T )(d ????V?
式中 和 是对初始位形定义的,t+Δt时刻的
体积力和表面力,它们是已知的。 0F 0?
)()( eδBEE ????? ???
式中
?? ??
ee SV
e φNFNF
00
00
T
00
T
E dd
?
SV
2000.4 哈尔滨建筑大学 王焕定教授制作 32
t时刻应力引起的
等效结点力矩阵
t+Δt时
刻荷载
引起的
等效结
点力矩
阵或
? ? ???
e V e
RPPSB
0
0Ed0
T )(d V?
将其按集成规则集装后可得
再引入如下记号
? ? ??
e V
L
e
VSBBR
0
0
T d)(
?
? ? ???
e V
e
e
RSBδ
0
000T d)( V???
? ?? ???
e SV ee
NFNPR )dd(
00
00
T
00
T
d0
?
? SV
将 和 的表达式代入,可得B? S
? ? ?????
e V
nL
e
RSSBBB
0
000T )d()( V?
2000.4 哈尔滨建筑大学 王焕定教授制作 33
几何或非线性应变
增量刚度矩阵
或
利用这些关系,非线性平衡方程可写为
? ? ? ??
e V V
L
n
e e
SAGSB
0 0
00 VV dd
TTT ?
? ??
e V e
MGGK
0
0
T d V
?
? ? ?? ??
e e V
e
V ee
MGGMG
00
00 ?? ?? VV dd
TT
? ? ?????
e V e
RRUKSB
0
000T d ?? ?? V?
? ? ???
e V
e
e
RSB
0
00
T d)( V????
为求解上述方程,尚需解决如下两方面问题
2000.4 哈尔滨建筑大学 王焕定教授制作 34
首先假设
?
?
?
EE
E
T EDS
?
? d
然后将 ΔS和 ΔE的关系线性化。根据本构关系
EDS tT ?? ?
? ?? ? ??
e V
L
e V ee
SBBSB
00
0
T
0
T d)(d VV ???
EDS T dd ?
则有
?
?
?
?
ep
D
D
D TT
弹塑性
非线性弹性
为使其线性化,设 ( t时刻的材料性质矩阵 )t
TD
2000.4 哈尔滨建筑大学 王焕定教授制作 35
在做了上述两方面处理后,可得
由此出发,用非线性方程的相关解法,即可解
决大变形非线性(材料非线性)问题。
? ?? ? ??
e V
L
e V ee
SBBSB
00
0
T
0
T d)(d VV ???
将其代回非线性平衡方程,可得
? ? ??
e V
L
e
EDBB
0
0T
T d)( V?
UKBBDBB Ne
e V
LL
e
?? ???? ? ? ?
0
0T
T )d()( V
0????? 0)()( RRUKKUΨ N ?? ??
讲义上给出了 T.L法的求解步骤,可供大家编
程序参考。
2000.4 哈尔滨建筑大学 王焕定教授制作 36
3.2 物质描述大变形增量问题的 U.L法
因 t+Δt的位移是用 t时刻位形为基准度量的,
因此
在 [t,t+Δt]间隔内,以 t时刻位形为参考位形,
其增量位移为
iii xxu ???
ii uu ?? )( ii uu ??? ?象 T.L法一样,设单元和位移的描述为
eNxx ? exNx ? eNu ??? ?但需指出的是,式中形函数 N是 t时刻单元自然
坐标的函数。在计算 等导数时,要先作坐
标变换 (Xi应换为 xi)。 jxi,N
To
26
2000.4 哈尔滨建筑大学 王焕定教授制作 37
类似地,用矩阵来表示则有
在时刻 t和 t+Δt的格林应变是以 t时刻位形定义
的,因而它们分别为
ij
j
k
i
k
j
i
i
j
ij Ex
u
x
u
x
u
x
u
E ?
????
?
?
?
?
?
?
?
?
?
?
?
? )(
2
1
NL EEE ??? ?? eL BE ??? ?
式中算子矩阵象 T.L法一样,但应将 Xi换为 xi。
0?ijE
ene
LL
N BGAAE ??? ?????? ??? 2
1
2
1
enBBE ??? )( ??
2000.4 哈尔滨建筑大学 王焕定教授制作 38
再次强调,式中算子符号象 T.L法一样,但应将
Xi换为 xi。
基于上述说明,象 T.L法一样可导得
)()(2)( enenN δBδBE ?????? ??
)()( eL δBE ???? ?
ijijS ??
)()()()( een δBδBBE ?????? ????
关于应力的处理也和 T.L法一样,对 t 时位形
定义的 t和 t+Δt时刻的克希荷夫应力分别为
kl
l
j
k
i
ij x
x
x
xJS ?
?
?
?
??
SSSS ?? ???? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 39
象 T.L法一样由虚位移原理虚功方程可导得
象 T.L法一样推导,引入如下符号定义
? ??
eV
e
ee
je δFFδBS ????
T
E
T )(d V?
式中
?? ??
ee SV
e φNFNF
?
SV dd TTE
其中 分别为 t时刻位形定义的单
元体积、应力表面、体力和表面力。 φFSV
ee,、,?
几何或非线性应变
增量刚度矩阵? ??
e V e
MGGK 0T d V?
荷载 引起的等效结点力矩阵
? ?? ???
e SV ee
φNFNPR )dd( TTd
?
SV
2000.4 哈尔滨建筑大学 王焕定教授制作 40
则可得
t时刻应力引起的等效结点力矩阵
? ?? ? ??
e Ve V ee
VBVSBR dd TT ??
? ? ?????
e V e
RRUKSBUΨ 0?? ??? Vd)( T?
t+ Δt时刻的非线性平衡方程
象 T.L法一样,为求解上述方程也需解决线性
化问题。首先讨论 ΔS的计算。因为
ij
k
j
l
i
lkij x
X
x
XJS ?? ?
?
?
?
??
2000.4 哈尔滨建筑大学 王焕定教授制作 41
m
m
klkl
l
j
k
i
ij x
v
x
X
x
XJS
?
??
?
?
?
?? ?? ??(
)
m
k
lm
m
l
km x
v
x
v
?
??
?
?? ??
因此
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
可改写为
由第四章已知 Chap4 64
klijk lij VD
ep???
式中
kijkkjikijij ?????? ???
??
Chap4 44
)(,tVE ijtij ?
Chap4 26
ij
k
j
l
i
lkij x
X
x
XJS ?? ?
?
?
?
??
2000.4 哈尔滨建筑大学 王焕定教授制作 42
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
可得 kijkkjikijij ?????? ??? ?
?
由如下两式消去 并利用,且
注意到 Vij对称,Ωij反对称
?ij?
ijijxi Vv j ??? ?,
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ??
GoTo 49
klijk lij VD
ep???又由于, 和,因
此 ijij VE ?
?
klkl
k
k V
x
v ??
?
?
?? ? klijk lij EDS
GoTo 50
klijlijkljiki j kli j kl SSSDD ??? ???? ep
2000.4 哈尔滨建筑大学 王焕定教授制作 43
上式最后一项将使本构张量不对称,对金属
类不可压缩介质,这一项可略去,也即
在有限的克希荷夫应力和格林应变增量之间仍
认为
这就是 U.L法的本构关系线性化。
)( ep kijkkjikkli j klij ESESEDS ???? ???
lijkljiki j kli j kl SSDD ?? ??? ep
与 T.L法一样,除本构关系线性化外,还需解
决几何方面的线性化。因为
? ?? ? ?
e Ve V
e
ee
SESBδ VV d)(d)( TTT ?????? ?
2000.4 哈尔滨建筑大学 王焕定教授制作 44
为对上述三项积分作几何方面的线性化,可设
这一关系可有两种矩阵表达方式:其一是
eL BEE ???? ??其另一方式是
? ?? ? ???
e V
kli j k lij
e V
ijij
ee
VEDEVSE d)((d)( ep ??????
)d)(d)( ?? ??
ee V
kijkij
V
kjikij VESEVESE ??????
)(
2
1
j
i
i
j
ij x
u
x
u
E
?
?
?
?
?
?
??
?
? ? T333231131211,,,,,,EEEEEEL ???????? ??? EE
2000.4 哈尔滨建筑大学 王焕定教授制作 45
eeL BNAEE ?? ???? ???
引入讲义上式 (5,102)算子矩阵,则有A
? ?
e V e
SB VdT ??
基于上述讨论,积分 时,作如
下处理
?? ? ?? ?
e Ve V ee
BDBSB VV d(d epTT ??
UBMBBMB
ee VV
?)dd TT ?? ?? VV
式中 为由 组成的 6× 6矩阵,矩阵为
epD epijklD
M
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
S
S
S
M
00
00
00
?
?
?
?
?
?
?
?
?
?
?
?
333231
232221
131211
SSS
SSS
SSS
S
2000.4 哈尔滨建筑大学 王焕定教授制作 46
结
束
引入如下两个矩阵
? ??
e V
L
e
BDBK VdepT
? ? ??
e V
N
e
BMBK Vd2 T
? ? ?????
e V e
RRUKSBUΨ 0?? ??? Vd)( T?
则 t+ Δt时刻的非线性平衡方程
0?????? RRUKKKUΨ NL ?? ?? )()(
可改写成
需要指出的是,要提高精度可减小 Δt或采用
适当的修正技术。
2000.4 哈尔滨建筑大学 王焕定教授制作 47
验证 δALθ=ALδθ
? ? ?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
3
2
1
T
1
T
3
T
2
T
3
T
1
T
2
T
3
T
2
T
1
0
0
0
00
00
00
θ
θ
θ
θθ
θθ
θθ
θ
θ
θ
θA
L
??
??
??
?
?
?
?
2000.4 哈尔滨建筑大学 王焕定教授制作 48
? ? θA
θθθθ
θθθθ
θθθθ
θθ
θθ
θθ
θθθθ
θθθθ
θθθθ
θθ
θθ
θθ
L
?
??
??
??
?
?
?
??
??
??
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
?
3
T
11
T
3
3
T
22
T
3
2
T
11
T
2
3
T
3
2
T
2
1
T
1
3
T
11
T
3
3
T
22
T
3
2
T
11
T
2
3
T
3
2
T
2
1
T
1
纯量的转置不变
返
回
2000.4 哈尔滨建筑大学 王焕定教授制作 49
推证式( A)
k
k
ijijij x
vS
?
??? ?? ??
k
i
jk
k
j
ik x
v
x
v
?
??
?
?
? ??
kijkkjikijij ????? ? ???
??
k
k
ijkijkkjikijij x
vS
?
????? ?? ?????
?
)()( ikikjkjkjkik VV ?????? ????
k
k
ijkijkkjikij x
v
?
????? ? ?????
?
)()( kikijkkjkjik VV ???? ????
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ?? 返
回
2000.4 哈尔滨建筑大学 王焕定教授制作 50
式( 5,96)的推证
kijkkjik
k
k
ijijij VVx
vS ???
? ???
??? ??
kijkkjikklklijkli j kl VVVVD ???? ???? ep
kijkkjikklklijkli j kl VSVSVSVD ???? ?ep
kllijkklljikklklijkli j kl VSVSVSVD ??? ???? ep
kllijkljikkliji j kl VSSSD )( ep ??? ????
????? ijlijkljikkliji j kl ESSSD )( ep ???
?? ?
kli jk lij EDS lijkljikkliji j kli j kl SSSDD ??? ????
ep 返
回