“有限元法基础及应用”补充讲义

有限元基础

有限元法基础及应用

补充讲义

2010年2月

有限元基础

“有限元法基础及应用”补充讲义

一、弹簧单元与弹簧系统

1、 弹簧单元分析

1）单元描述

弹簧系统受力平衡时，从中隔离出一个典型弹簧单元进行分析。

图 1-1

弹簧单元的端点i,j设置为单元节点。 基本未知量为节点位移：ui,uj

单元节点力（单元在节点处受到的作用力）：fi,fj 已知弹簧的物理特性：F k

其中：k为弹簧刚度, uj ui为弹簧伸长量，F为弹簧力（拉伸为正） 2）建立弹簧单元的单元刚度方程

考虑弹簧单元在系统中变形平衡时的条件：力平衡条件和弹簧物理特性，得到下列方程：

fi F k(uj ui) kui kujfj F k(uj ui) kui kuj

写成矩阵形式：

（1-1）

fi k fj k k ui

k uj

（1-2）

上式的矩阵符号形式为：

f kd

（1-3）

方程（1-2）或（1-3）称为弹簧单元的刚度方程，反映了单元的力学特性，

即节点力~节点位移之间的关系。

式（1-3）中：

有限元基础

k k k ,称为单元刚度矩阵

kk ui

d ,称为单元节点位移列阵

uj fi

f ,称为单元节点力列阵

fj

3）弹簧单元刚度方程的讨论

a. k有何特点？

对称、奇异、主对角元素恒正。 b. k中元素的物理意义是什么？ 刚度矩阵元素的大小等于弹簧刚度。

从对方程（1-2）分析的分析可以看出，矩阵中某列的各元素代表列序号对应节点有单位位移，其它节点位移为零时，单元各节点上的节点力；某行的各元素分别是单元各节点的位移对行序号对应节点的节点力贡献系数。

因此，矩阵中任意一个元素kij的物理意义是：j节点的位移对i节点的节点力贡献系数，或者j节点有单位位移，其他节点位移为零时，i节点上的节点力。

c. 单元刚度方程可以求解吗？为什么？

不可以。单元刚度方程仅仅表征一个单元的力学特性，单元水平上无法确定单元节点位移。只有把系统中所有单元特性集成后，在系统水平上才可能求出所有未知位移和反力。单元水平上，若已知单元的节点位移，可由刚度方程求出所有单元节点力分量。若节点力已知，单元节点位移不能确定，单元可作刚体运动。这也是单元刚度矩阵奇异性的物理解释。

2、弹簧系统整体分析求解

1）建立系统节点平衡方程

以右图的一个弹簧系统为例，研究如何由单元特性集成系统特性并建立对系统进行求解的控制方程。

由前面得到的弹簧单元的刚度方程公式（1-2），分别写出2个弹簧单元的特性方程如下： 单元1

单元2

图 1-2

（1-4）

（1-5）

（注：右端节点力分量的下标为单元节点的局部编号，上标是单元编号）

下面用两种方法装配单元特性、建立系统控制方程，并在特定条件下求解。

有限元基础

（1）由系统中节点平衡条件导出：

系统处于平衡时，考虑各节点（1，2，3节点）的平衡条件：

由于节点受到的外载荷F1,F2,F3与节点受到与其连接的所有单元对其作用： 力（单元节点力的反作用力）之和等于零。因此有下列（节点）平衡方程（组）

F1 f11

F2 f21 f12 F3 f22

把单元特性（1-4），（1-5）代入（1-6）得到：

（1-6）

F1 k1u1 k1u2

F2 k1u1 (k1 k2)u2 k2u3 F3 k2u2 k2u3

写成矩阵形式：

（1-7）

（1-8）

或矩阵符号形式：

KD F

与节点位移之间的关系，是求解节点位移的控制方程。

方程（1-9）中：

（1-9）

方程（1-8），（1-9）是系统节点平衡方程，该方程建立了离散系统的外载荷

K——弹簧系统的结构总刚度矩阵 D——系统节点位移列阵 F——系统节点载荷列阵 讨论：

a．K 有那些特点和性质？ b．上述方程能求解吗？

有限元基础

（2）由单元刚度方程叠加导出

将单元1，2的刚度方程（1-4），

（1-5）扩大到系统规模：

（1-10）

（1-11）

注意：1）对单元刚度方程扩大规模并不改变其表达的力学关系。

2）扩大后的单元刚度方程采用整体节点位移列阵。 3）扩大后的方程中矩阵元素按对应的整体节点序号排列！

将上述两个方程叠加，得到：

（1-12）

将系统中节点平衡条件（1-6）代入上式，就得到与（1-8）相同的系统节点平衡方程。上述两种建立系统平衡方程的方法都考虑了1）单元特性集成；2）系统中节点外载荷与系统的节点力（系统节点内力）的平衡。因此方程（1-8）的本质是系统中所有节点的力平衡关系，其左边是由节点位移表示的系统节点力，右边是节点所受外载荷。不难发现，系统总刚度矩阵可以直接由单元刚度矩阵扩大后叠加而得到。总刚度矩阵元素的含义可以由方程（1-12）分析出。 2）系统平衡方程求解

假如边界条件为：

u1 0 F F P3 2

则节点平衡方程（1-8）化为：

（1-13）

（1-14）

有限元基础

将该方程展开为两部分。第2，3

个方程变化为：

（1-15）

第1个方程变化为：

（1-16）

先后解方程（1-15）、

（1-16）得：

（1-17）

（1-18）

从而解出了系统的未知位移和未知反力，并可以进一步求弹簧力。

3、例题

例1 图示一个3弹簧系统。

k1 100N/mm,k2 200N/mm,k3 100N/mm,P 500N,u1 u4 0

求： （a） 系统总刚度矩阵 （b） 节点2，3的位移 （c） 节点1、4的反力 （d） 弹簧2中的力

图 1-3

解：

（a）分别写出各单元刚度矩阵：

有限元基础

参照方程（1-10）、（1-11）中单元刚度矩阵的扩大，用叠加法直接得到系统

总刚度矩阵：

该总刚度矩阵的特点：对称性、奇异性、稀疏、非零元素沿主对角线呈带状分布。

（b）参考方程（1-8）

，利用求出的总刚度矩阵，写出系统节点平衡方程：

（1-19）

考虑到位移边界条件：u1 u4 0 则平衡方程组（1-19）第2，3

方程化为：

求解上式得：

（c）由（1-19）的方程1，4得1，4节点的反力：

（d）弹簧2内力为：

F2 k2 2 k2(u3 u2) 200 3 2 200(N)

（拉力）

有限元基础

4、练习题

对图1-4所示弹簧系统，用扩并根据大-叠加法求其总刚度矩阵。节点平衡方程的物理意义，直接写出总刚度矩阵。

图 1-4

二、杆单元

目标：通过杆单元特性方程的建立，初步掌握有限元法单元分析的过程和原理。

了解杆系结构有限元分析的原理。

1、等截面杆单元及其刚度矩阵

1）单元描述

L— 杆长 A— 截面积

E— 弹性模量

单元上的力学量和基本关系如下：

u u(x) (x) (x)

——杆单元沿轴向位移分布 ——杆单元应变分布 ——杆单元应力分布

图 2-1

应变—

（2-1）

应力—应变关系： E （2-2）

单元节点位移： d

ui

u j

fi f 单元节点力：

fj

有限元基础

2)单元特性方程（刚度方程）

（1）直接法导出杆单元特性

采用材料力学基本知识对单元进行力学分析。

杆单元伸长量： uj ui （2-3）

（2-4） （2-5） （2-6）

（2-7） 由于轴向变形模式下，杆单元的行为与弹簧单元相同，因此可比照弹簧单元

的刚度方程（1-2），考虑到式（2-7），直接写出杆单元的刚度方程：

（2-8） 写成符号形式：

f kd

因此杆单元刚度矩阵为：

（2-9）

（2-10）

（2）公式法导出杆单元特性

a.在单元上假设近似位移场

对图2-1所示杆单元，首先用函数插值法构造以单元节点位移未知量为待定参量的简单多项式函数，作为单元上的假设位移分布函数。插值过程如下。

考虑到杆单元只有2个沿轴向的未知节点位移分量，因此假设单元上位移函数为一次多项式：

u(x) a0 a1x （2-11） 将单元两个节点的坐标值0，L分别带入上式得到：

ui u(0) a0

uj u(L) a0 a1L

对上面2个方程联立求解，得到用节点位移表示的多项式系数：

有限元基础

a0 ui

a1 (uj ui)/L

上两式代入式（2-11），并整理得：

u(x) Ni(x)ui Nj(x)uj （2-12） xx

,Nj(x) 分别为节点i,j的插值基函数，有限元法中LL

称为形状函数，简称“形函数”。

上式中：Ni(x) 1

式（2-12）就是插值形式的单元假设近似位移场，也称为单元位移模式，这里是线性位移模式。

单元位移模式（2-12）写成矩阵形式为：

u Ni

ui

Nj Nd uj

（2-13）

式（2-12）或（2-13）是有限元法中最重要的关系式之一，通过该式把单元上假设的近似位移分布函数用节点位移来表示，为进行单元层次上的分析打下了基础。

上式中N称为单元的形函数矩阵。 b.单元应变和单元应力公式

可由直杆的应变——位移方程

（2-1）和单元位移模式（2-13）求出单元的应变分布和节点位移的关系：

（2-14）

式中：

（2-15） B称为单元的位移——应变转换矩阵，简称应变矩阵。

由杆的应力——应变关系（2-2） ，得单元应力分布和单元节点位移的关系：

E EBd （2-16）

这里可以看出，有限元法中单元上位移、应变、应力等力学量的连续场函数都用节点上的未知量（位移）来表达！

c.用虚功原理导出杆单元刚度方程

变形体的虚位移：

假想在变形体上发生的，满足位移许可条件（内部连续，边界协调）的微小、

任意位移场。

有限元基础

可以理解为某个位移场的微小扰动（变分）。 虚位移的特征：

1）假想的，与真实位移无关； 2）几何上是许可的：连续、协调； 3）微小、任意大小。

虚功原理：

变形体受力平衡时，若发生虚位移，则外力虚功等于变形体内的虚应变能（应

力在虚应变上做的虚功）。下面把虚位移原理应用在所分析的杆单元上。 定义杆单元的虚位移：节点虚位移 单元虚位移 单元虚应变

ui

节点虚位移： d

uj

单元虚位移：

u N d

那么，单元节点力（外力）虚功为： dTf

T 单元虚应变能： dV dBEBddV d BEBdV d

VV V

T

T

T

T

根据虚功原理，上述节点力（外力）虚功等于虚应变能，因此有下列关系：

T

dTf dT B EBdV d （2-17）

V

考虑到 d的任意性，从上式可以得到：

T

d kd （2-18） f BBdVE

V

上式就是杆单元的刚度方程，杆单元的刚度矩阵为：

k BTEBdV （2-19）

V

其导出原理和计算方法可推广到其他类型的实体单元。 具体计算式如下：

（2-20）

显然，与前面直接法得到的单元刚度矩阵（2-10）式相同。

有限元基础

3）杆单元讨论

a. 单元自由度

只有拉伸、压缩变形的杆单元在局部坐标系下是一维问题，其2节点单元只有2个节点位移分量——单元有2个自由度，单元刚度方程、刚度矩阵为2阶。

b. 单元刚度矩阵元素的物理意义推导： 单元刚度方程为：

fi k11 fj k21

令：

ui 1 uj 0

k12 ui

k22 uj

（2-21）

（2-22）

上式代入（2-21）得到：

fi k11

fj k21

（2-23）

上式表明，单元刚度矩阵第一列元素就是当单元节点位移满足式（2-22）时的单元节点力分量。如果能设法求出此时的节点力，就得到第一列的刚度元素。

一般地，单元刚度矩阵的第i列元素表示当维持单元的第i个自由度位移为１，其它自由度位移为0时，单元上的节点力分量。

可以用此概念和推导过程直接求出杆单元的刚度矩阵元素。 c. 单元刚度矩阵性质

单元刚度矩阵对称、奇异、主元恒正。 4）例题

例1求图2-2杆中的应力。 解：结构划分为２个杆单元。根据杆单元刚度矩阵公式（2-10）分别写出两个

单元的刚度矩阵为：

图 2-2

有限元基础

参照前面弹簧系统的分析方法，装配系统的有限元方程（节点平衡方程）如

下：

（2-24）

考虑结构的约束条件：u1 u3 0和载荷F2 P，方程（2-24）化为：

上式的第2个方程为：

（2-25）

EA

3 u2 P （2-26）

L

求解该方程后得到系统的位移解：

计算单元应力：

单元1

（2-27）

单元

2

几点提示：

（1） 例中单元应力的计算采用了材料力学中的方法,与采用有限元单元应

力公式 E EBd的结果相同，请验证。

（2） 锥形杆，单元截面积可以用平均值，从而转化为类似本题的问题求解。 （3） 求应力之前需要先求出节点位移，因此本方法称为有限元位移法。 （4） 如果杆上受连续分布的轴向载荷或节点之间受轴向集中载荷，分析时

可以按照虚功相等的原则先把单元上的载荷等效移置到节点上。

有限元基础

例2 如图2-3所示直杆结构有限元模型。 已知：

求：杆两端的支反力。

图 2-3

解：单元、节点的定义如图2-3。须检查杆右端（节点3）与墙壁是否接触。先计算右端的自由伸长：

根据上式可以判断：结构受力平衡时，右端间隙将闭合，节点3与刚性墙壁接触。

参考前面的讨论，直接写出系统有限元平衡方程：

引入载荷和位移边界条件：

（2-28）

F2 P 6.0 104Nu1 0u3

1.2mm

则有限元平衡方程（2-28）成为：

分离出第二个方程：

即：

解得：

（2-31） （2-30） （2-29）

有限元基础

全部位移解为：

根据上式位移解，从系统平衡方程（2-28）的第1，3个方程分别求出支反力如下：

解毕。

（2-32）

2、 二维空间杆单元（平面桁架单元）

二维空间中建立杆单元刚度方程的基本思路是根据前面在杆的一维局部坐标系下建立的单元刚度方程通过坐标变换，转换为二维总体坐标系下的方程，同时得到坐标变换后的单元刚度矩阵。而系统整体分析的原理和方法与一维情况相同。

1）单元向量坐标变换

图 2-4

上图为一个杆单元及其二维局部坐标系与二维总体坐标的关系。

现在将节点位移分量和单元节点力分量在二维局部坐标系x-y下描述，则节点i形式上具有2个自由度：

有限元基础

位移分量为ui'，vi'；

''

节点力分量为fxi，fyi

其中只有x方向的位移分量和节点力分量用来描述单元特性。

i节点的位移和节点力向量在二维局部坐标系与总体坐标系下的变换如下：

ui

ui' uicos visin lm

vi

ui

vi' uisin vicos ml

vi

上式中 l cos ,m sin

上述变换的矩阵形式：

m ui ui' l

（2-33） ' vi ml vi 上式的矩阵符号形式：

~d Tdi （2-34） i

其中：

m ~ lT ml

（2-35）

~~

上式中T称为向量的坐标变换矩阵。显然T是正交矩阵，即：

~~T 1 TT

（2-36）

因此，由（2-33）可得单元节点位移向量的坐标变换式如下：

（2-37）

上式的矩阵符号形式为：

d Td

其中：

（2-38）

~

T0 T ~

0T

同理也可得到单元节点力向量的坐标变换式：

（2-39）

有限元基础

f Tf （2-40）

2）二维空间刚度矩阵

可导出二维空间杆单元的单元刚度方程根据上述讨论的单元向量坐标变换，和单元刚度矩阵。

已知杆单元在一维局部坐标系下的刚度方程为：

（2-41）

把该方程扩充到二维局部坐标系x-y下的4阶形式：

（2-42）

符号形式：

k d f

引入单元向量变换式（2-38），（2-40）得：

（2-43） （2-44） （2-45） （2-46） （2-47）

k Td Tf

考虑到变换矩阵T的正交性，得到： TTk Td f 或写成：

kd f

其中：

k TTk T

式（2-46）就是二维总体坐标系下杆单元刚度方程，k是二维空间杆单元刚度矩阵，由（2-47）可计算出：

（2-48）

3）单元应力计算

由单元应力计算公式（2-16）和位移向量变换式（2-37）得：

有限元基础

即：

（2-49）

4）例题

平面桁架由2根相同的杆组成（E，A，L）。求： （1）节点2位移 （2）每根杆应力

解：先求出每个单元在总体坐标下的刚度矩阵： 单元1，

（i,j）=(1,2)：

图 2-5

k1

T1k 1T1

T

（2-50）

（i,j）=(2,3): 单元2，

有限元基础

k2 T2k

2T2

T

（2-51）

刚度矩阵k1,k2也可按（2-48）直接计算。

将上述总体坐标系下的单元刚度矩阵扩充到整体规模，装配得到系统有限元

方程：

（2-52）

引入边界约束和载荷：

则系统有限元平衡方程凝聚为：

解得：

（2-53）

按公式（2-44）分别求得单元1，2的应力：

（2-54）

（2-55）

（2-56）

本文来源：https://www.bwwdw.com/article/qxcj.html