桥梁结构有限元

第三章 桥梁结构分析有限元方法

第一节 平面杆系有限元

有限元分析结构的工作为两大部分： 元素分析（单元刚度矩阵）：探讨单个元素的力学特性。 结构分析（总刚度矩阵）：用众多元素集合成整个结构。

一、单元刚度矩阵（局部坐标）

如图所示为一个等截面的直杆e 。设杆长为l ，截面面积为A ，截面惯性矩为I ，弹性模量为E 。

图1

如图所示一个杆件在平面内的变形，可用杆件两端的结点位移来描述。则一个杆件单元有六个结点位移和六个杆端力分量。用矩阵公式表示为：

T

j j j i i i e

T j j j i i i e M Y X M Y X F v u v u ]

,,,,,[}{],,,,,[}{==??δ (1)

在分析计算中，还必须对这些分量的方向做规定：如图2所示的位移、力分量方向为正方向。

图2

下面寻找结点位移和杆端力的关系—刚度矩阵：

??

??

???

?

????????????????

???????????

????????????????

?--------=??????????????????????j j j i i i

j j j i i i V U V U l EI l EI l

EI l EI

l EI l EI l EI l EI l EA l

EA l EI l EI l EI l

EI

l EI l EI l EI

l EI l EA l

EA

M Y X M Y X ??460

260

6120612000026046061206120

0002

22323222323

(2) e e e K F }]{[}{δ= (3)

1 用结构力学直接法推导

运用叠加原理，将六个结点位移各自发生单位位移时所产生的杆端力分量与对应位移分量相乘后叠加。

2 用弹性力学虚功原理推导 （1）选择位移函数

342321x a x a x a a v +++=

利用边界条件，可获得用位移向量e }{δ表示的矩阵形式为： e N v }]{[δ= (4) 所以[N]称之为形函数。

（2）位移与应变的关系（几何条件）

22dx

v

d =κ y κε=

对形函数求导，可获得：e

e B y dx

N d }]{[}{][}{2

2δδε=-= (5)

（3）位移与应力的关系（物理条件）

εσE =

将上式改写为矩阵形式，即：

e B D D }]{][[}]{[}{δεσ== (6) （4）虚功原理求单元刚度矩阵

虚应变能：?=V

e T e dV U }{)}({**σεδ (7)

?=V

e T T e dV B D B U }]{][[][)}({**δδδ

e V

T T e dV B D B U }{]][[][)}({**δδδ?=

外力虚功：e e T e e T e K F W }]{[)}({}{)}({**δδδδ== (8) 虚功原理有：**W U δδ=

则单元刚度矩阵为：?=V

T e dV B D B K ]][[][][ (9)

二、整体坐标刚度矩阵 1 坐标转换矩阵

图3

将整体坐标的杆端力投影到局部坐标系上，有：

?????

?

???????????????????????????

???????

??--=??????????????????????j j j i i i j j j i i i M Y X M Y X M Y X M Y X 10

0cos sin 0sin cos 1000cos sin 0sin cos ααααααα

α (10) e e F T F }]{[}{= (11)

同理，可得到整体坐标位移向量投影到局部坐标系上，得

e e T }]{[}{δδ= (12)

2 整体坐标单元刚度矩阵

由局部坐标的结点位移和杆端力的关系式：

e e e K F }]{[}{δ= 有：e e e T K F T }]{][[}]{[δ=，则：

e e e e T e K T K T F }]{[}]{][[][}{δδ== 即：]][[][][T K T K e T e = (13)

为整体坐标下的单元刚度矩阵。

三、总刚度矩阵

若将单元杆端力向量和位移向量分解为单元杆端i 的力向量和杆端j 的力向量以及杆端i 的位移向量和杆端j 的位移向量，记为：

T e j e i e F F F )}{,}({}{=， T e j e i e )}{,}({}{δδδ= (14) 那么刚度矩阵可写成如下分块矩阵：

[]

?

??????

?=e jj e ji e

ij

e ii e

k k k k k (15) 则有：

??

?

???????????????=??????????e j e

i e jj e ji e ij e

ii

e j e i k k k k F F }{}{}{}{δδ 展开为：

e j e

ij e i e ii e i k k F }{}{}{δδ+= 和 e j e jj e i e ji e j k k F }{

}{}{δδ+=

(16) 图4 按如图所示结构，对于单元①，有：

121121111111}{}{}{δδk k F += 和 121

221112112}{}{}{δδk k F +=

对于单元②，也有：

232232222222}{}{}{δδk k F += 和 232332223223}{}{}{δδk k F +=

对于单元③：

343343333333}{}{}{δδk k F += 和 343443334334}{}{}{δδk k F +=

对于结点2，其平衡方程与变形协调为： }{}{}{22212P F F =+，}{}{}{22212δδδ== 同理，对于结点3，有：

}{}{}{33323P F F =+，}{}{}{33323δδδ== 整理可得到：

??????????????=????????????????????????????++}{}{}{}{}{}{}{}{0

00

0043214321344343

334333

2332

32223

222

1221211

121

11P P P P k k k k k k k k k k k k δδδδ (17)

上式又可写为：

}{}]{[P K =δ （18） 因此，总刚度矩阵一般形式为：

?????

???

?

??

?

?????????????

?=∑∑∑∑====m e e

nn e

ni

e n e n e

in m

e e

ii

e i e i e

n e i

m

e e e e n e i

e m e e k k k k k k

k k k k k

k k k k k K 12

11

2

1221

22

211111

111...

...

......

......

...............

............

............

][ （19）

分析总刚度矩阵的特点。

四、支承条件的引入

上面的总刚度矩阵是在不考虑结构约束条件的情况下，结构可以自由移动，所以总刚度矩阵为奇异矩阵。在求解总刚度矩阵方程之前，必须加入约束条件。

对上例进行计算，引入的支承条件是结点1和结点4的固端约束，其结点位移为零，即：

{δ1}={0} {δ4}={0} 代入方程（17），有：

??????????????=???????????????????????????

?++}{}{}{}{}0{}{}{}0{0

000

00432132344343

334333

2332

32223

222

1221211

12111P P P P k k k k k k k k k k k k δδ 求解结点位移只需中间两个方程，即

?

??

???=????????????++}{}{}{}{3232333233232

2

23

222122P P k k k k k k δδ （19） 而剩下的两个方程则可以用来求约束反力，即

}{}{21

121δk P =

}{}{33434δk P = （20）

从上可以得知：支承条件的引入的问题归结为划去总刚度矩阵中位移为零所对应的行和列的问题。

五、有限元计算简单程序 1 程序流程

2 程序设计

目前计算程序设计均采用模块（子程序）式建造方法。上述框图中的每个方框均编制为子程序模块，而主程序仅负责模块数据间的传递和数据的输入输出（I/O ）。