弹性力学有限元位移法原理.docx

资源描述

弹性力学有限元位移法原理.docx

《弹性力学有限元位移法原理.docx》由会员分享，可在线阅读，更多相关《弹性力学有限元位移法原理.docx（85页珍藏版）》请在冰豆网上搜索。

弹性力学有限元位移法原理.docx

弹性力学有限元位移法原理

一、课程论文：

弹性力学有限元位移法原理（30分）主要针对一维（直杆）问题，撰写一篇论文，对有限元位移法的原理作一般性虬括和论述。

要求论文论及但不限于下列内容：

1）弹性力学有限元位移法的数学、力学基础和基本思想；2）有限元法求解的原理和过程，推导所有计算列式；对基本概念和符号进行解释和讨论；3）收敛性、收敛准则及其数学、力学意义的讨论。

弹性力学有限元位移法原理

一、有限单元法的起源

有限单元法的形成可以追溯到20世纪50年代甚至更早些时间，基本思路来源于固体力学中矩阵位移法的发展和工程师对结构相似性的直觉判断。

对不同结构的杆系、不同的载荷，用矩阵位移法求解都可以得到统一的公式。

在1952-1953年期间，

R-W-Clough和M-J-Turner在分析飞机三角翼振动问题时，

提出了把平面应力三角形或矩形板组合起来表达机翼刚度的方

法，当时被称为直接刚度法。

1956年M-J-Turn

er，R-W-Clough，H•C•Martin，L-J-Topp在纽约举行的航空学会年会上发表论文《Stiffnessanddeflectionanalysisofcomplexstructures》（复杂结构的刚度和变形分析）介绍了这种新的计算方法，从而将矩阵位移法推广推广到求解弹性力学平面应力问题。

它们把平面板壳结构划分为一个个三角形和矩形的“单元”，利用单元中近似位移函数，求得单元节点力与结点位移关系的单元刚度矩阵。

1960年,R-W-Clough在论文《Thefiniteelementinplanestressanalysis》（平

Finite

面应力分析的有限元法）中首次提出了有限单元（

Element）这一术语，他也因此被称为“有限单元之父”

二、有限元法的基本思想

有限元法是一种结构分析的方法，正如O-C-Zienkiewicz所

说的：

“人类思维的限制在于不能通过一步运算就掌握复杂环境和事物的行为。

因此，先把所有系统分解为它们的元件或单元，这些元件的行为已经被充分的了解，再把元件重新组装成原来的系统来研究系统的行为”。

可以看出有限元法的基本思想是将连续的求解区域离散为一组由有限个单元组成并按一定方式相互连接在一起的单元组合体来加以分析。

三、有限单元法的数学基础当有限单元法成功的应用于求解弹性力学平面问题之后，下一步要解决的问题就是能否把这种方法应用于求解其他连续介质问题。

在寻找连续介质问题近似算法的时候，数学家们发展了微分方程的近似解法，包括有限差分方法，变分原理和加权余量法。

四、有限元分析的基本步骤

⑴建立研究对象的近似模型

⑵将研究对象分割成有限数量的单元

⑶用标准方法对每一个单元提出一个近似解

⑷将所有单元按标准方法组合成一个与原有系统近似的系统

⑸用数值方法求解这个近似系统

⑹计算结果处理与结构验证五、一维杆的有限位移法分析

本文以一维直杆的分析为例子，研究有限元位移法基本原理和求解过程。

⑴虚位移原理推到一维直杆单元的刚度方程

如下图所示一维直杆，已知直杆杆长为L，横截面积为A，材料弹性模量为E,所受轴向分布载荷集度为q（x）。

杆端位移分别记为Ui，w，杆端力分别记为Si，Sj。

q（x）

设局部坐标系下杆中A点的坐标为Xa，因为只有两个边界条件Ui,Uj，

因此杆轴任意一点（例如A点）的位移可假设为

式中a，b为待定常数它们可由杆端位移条件来确定:

将式⑵代入式⑴可得:

U=（VL）Ui

XUi

若引入无量纲变量:

则式（3）可改写成：

fuil

u=Nm+NjUj=[NiNj]<>=NUe⑷

式中

Nj=1-Nj=

称为形函数，矩阵N称作形函数矩阵；矩阵Ue称为杆端位移矩阵或节点位移矩阵。

由式（4）可以看出，形函数具有如下性质：

1、本端为它端为零

叫（0）=1叫（0）=0

（1）=0Nj

（1）=1

2、任意一点总和为1

）叫（）「

现采用虚位移原理给出该杆单元的特性公式，设杆端i，j分别产生虚位移Uj，由此引起的单元内任意一点的虚位移为：

__T

u二Nue=NUjuj

dNi

dNj

-1

—

-丨

|Ue

=|一—

—

.dx

」

L一

二【Bi

B2】ue=Bu

式中B为应变矩阵。

由此可得「；=BUe又

=EBu「E

根据虚位移原理：

对任意虚位移,

外力所做的总虚功恒等于变形体所

接受的总虚变形功,

所以有

qxudx二

ue0qxNuedx

（S'qqxNdx）ue

=（SqqxNTdx）Tue

LTLTTTLT

W变=oAdx二oueTBtEABuedx二ueT（°BTEABdx）ue

由W外=“变

可得：

（S:

qxNTdx）Tue

ueT（BtEABdx）u

（S:

qxNTdx）T=u

（BtEABdx）

若记

LTe

qxNTdx=FeoE

BtEABdx二Ke

FEe称为该杆单元等效节点载荷；Ke局部坐标单元刚度矩阵。

所以可得单元刚度方程:

SFE二KeUe

式中单元刚度矩阵的显式为:

心宜1T

可见单元刚度矩阵具有对称性。

即单元刚度矩阵的每一个元素可写成

L^EA独dx

0dxdx

⑵将一维直杆离散为三个单元进行分析

现考虑下图所示一维直杆：

长度为L,分为三个单元，每个单元长度

NJNgNjN4

对于单元①，节点位移分别为ui，U2,对应形函数为N,Nb;由

LdN

0dx

dNi

得:

dATjdNt血+血ax

dN2dNx

—AE—-dx+血dx

学AE学血+dxdx

dNadN\

--^AE—1dx+axdx

dN2dN^

—AE—-血+dxdx

dN^dN2

—^AE—1血+dxdx

dN4dN2

—AE—-血血dx

笄3

dxdx

dx+血ax

dN4dNA

—/IE—血+dxdx

dNrdN\一/IE一必+dx

^AE

dxdx则4dNx~^AE~dxdN*dN\—AE—-

dxdxdN2dN殳~dxAE~dxdR4dN,盂ae石

dNAdN2

dxdx

帕晰

（krAEdx

dN4dN4

盂ae盂

叫+

dx+

rdX1dN\

—-^AE~^dxdxdx

dN2dN\

—AE—-drdxdx学业学血

dxdr

dgdN\

—-AE—-

dxdx

dN?

dN2盂AE石

dN^dN?

盂盂

dN4dN2

—AE—-dx血血

dN.A冲,

——AE——dxdrdx

dgdN*

-^AE-^dxdxdx

dN4dN4

—血.drdx

又对于单元②、③，形函数N=0;对单元③,

形函数N2=0;对单元①,

形函数N3=0;对于单元①、②，形函数N=0。

因此可得：

[dN\dNx

/<11=，4月一必；S=

Jiii血血

fdN.dN2

=I——AE——dx+

J偽必必fdNA血

酹=-r±AE—^dx+丿血dxax

fdN4dN^

/<44=-r^AE—^dx.

Jq3血血

dN「dNr

——AE——dx;Km=0;K]4=0£*dbcdx

「誓E学血；心=/学卫e学血；隔=oq2dxdxJq2axax

fdN^血V3fdN4dg

I-~^AE-^dx;=-~^AE--^dx

阳必必/sb血必

且对于单元①，

Ni=i——

N2=

对于单元②，

N2「-上

N3=

对于单元③，

N3「-仝

N4a

所以

式中宀,2,3分别对应于单元①，②，③

Kii

dN—八dNdx

Ki2

hiEA

tdx

2dx=

0h/

Ki3

0,Ki4

K22

dN2EAdN

2dx

adx

hiiih

—EA—dxhh.

--）EA^-）dxhh.

EAEA

dN3EAdN

厂丄EA（-丄）=

0h2h2

K24

K33

dNoAdN3

3EA3

匕dxdx

h211h3

EAdx

0h2h20

EAEA

hzh3

」dNoAdN3,

dx3EA3dx

门3dxdx

（）EA（）dx

馆h3

K34

dNAEA^dx

13dxdx

h31

EA（-丄）二h3h37

K44

dN±EAd%x=

3dxdx

h3EA,dx二

所以单元整体刚度矩阵

0h32

心EA

——+——

h1h2

整体单元刚度矩阵元素

h2h3

Kj的物理意义为：

移时，在节点i上需施加的节点力

节点在X方向产生单位位

又由单元等效节点载荷

LTL__T

F,=打q（x）Ndx=打q（x）lMNJdx

可得：

每一节点等效载荷分别为：

=/N】q血+/TVjqdx+/TVjqdx

JillJSisJ口3

h=IN2q血+/N?

q血+/Mqdx

JfljJOjJO3

k=IMq血+/yV3q血+/Ngdx

JQiJO2vQj

f^N4q（x）dxN4q（x）dxN4q（x）dxR

又对于单元②、③，形函数N=0;对单元③，形函数N2=0;对单元①,形函数N3=0；对于单元①、②，形函数N=0。

所以有：

A=/Mq血

/2=/N皇q血+/Ar2q血

h=[N：

tdx+fMqdx

JQaJfig

f4二N4q（x）dxR

O-3

以上四式可写为:

h=Njq血+

血+/血=於+斤

Njq血+/用q血=於+胃

N2q（x）dxR

式中：

N「，表示单元e的形函数Nk;在本例子中，e=1、2、3,

k=1、2；

而且有

’X2—X…X—X1

N；（x）=—;逍（x）=—

he表示单元长度。

所以有

J1=

aIx2—町hi

a/云—房h2

af時一咯h：

由单元刚度方程可得:

1•丄

hih2

衍（囲一叱）、a

+隔

乃（啜一述）\.a

工4（诚一错）

丄丄

h2h3

■u/l

a（X27

z22、33

a/X2（X2-X1）X2-洛、R（2

X1（X22-N2）、a/X3（X32-X22）

）（

X3X2

h/32'g'23）

33（22、,22.33

aM3—X2X2（X3—X2）、丄a,X4（X4—X3）X4—X3、

—（）十—（）

032佗2，

33「22、

a/X4-X3X3（X4-X3R

（—-

h33

例：

假设A,E,L,a,和R都等于1。

且g=h2=h3=h，则有h=1,X1=

0,X2=1/3,X3=2/3,X4=1,

_2-10_

0,037

-12一1

血

01）74

0-11

■■

0.383

■■

由整体单元刚度方程得:

■■

■0.0■

0.494

0.951

_^4_

L333

对于单元e（e二①、②、③）可写出其位移函数:

应变为:

殂==

BN；

■■

.dx

由「丁二B护

得每个单元的应力:

b」=1*48

er2=1.37

cf3=L15

下图1和图2分别表示有限元解和精确解的比较:

图1位移对比

图2应力对比

二、分析与计算（40分）

1、图示两个结构和单元相似，方位相同的平面应力有限元模型，两模型的单元厚度和材料相同。

两个模型右端单元边上受均匀剪切面力。

对于下列2种情况,

试根据有限元法和力学有关知识来分析两个模型求解后对应节点的位移值和对

应单元的应力值之间的关系：

1）两个模型面力的合力相等；2）两个模型面力的集度相等。

（10分）

解：

建立坐标系如图所示，对（a）图，各节点坐标

点号

x坐标

y坐标

点号

x坐标

y坐标

单元节点信息数组可记为

单元①，面积A1=0.52020=200，由式

ai=Xjym-Xmyj

bi=yj-ymi,j,m（a）

Ci=Xm-Xj

可得

=0&=20

b2=-20c2=-20

b3=20c3=0

由式

i,j,m（b）.

单元几何矩阵为

-20

-10

0〕

-20

0-1

-20

20一

-1-1

［」

弹性矩阵大为简化,

可设

卩=0且为单位厚度,

B1=

400

为了计算简便,

1」

，可得D=E0

00.5

-20

-2

-20

10一

一40

-2

也0

-10

-1

=Se6e（c）,得单元的应力矩阵

=tABeTDBe（d），

s1亠

400

由式Ke=tBDBd'.1

单元①的单元刚度矩阵为

K⑴二tAB⑴TS⑴

根据单元刚度矩阵的性质可得

■1

-1

—1

-2

-1

-2

1一

二K⑴

二K⑶二K⑷

（2）

对（b）各节点坐标

点号

x坐标

y坐标

点号

x坐标

y坐标

单元节点信息数组可记为：

单元①，面积A⑴二0.51010=50，由式（a）可计算出

bj=0G=10

匕2=-10C2=-10b3=10C3=0

由式（c）

由式（b）,单元几何矩阵为

-10

—

■0

-1

-10

10一

-1

'10

-10

-1

得单元的应力矩阵

-10

-2

-10

-2

5一

-5

-1

S1丄

100

同理由单元刚度矩阵的性质可得

由式（d）得单元①的单元刚度矩阵

-1

-2

—1

-2

—1

-2

—1

-2

-1

二tAB⑴TS⑴

（2）

占北⑷

综上可知，两个模型中的单元刚度矩阵均相同，所以它们的总刚度矩阵也相

当两个模型面力的合力相等,

同，即

-1

-2

0〕

-2

-1

-2

-1

-2

-1

-2

-1

-2

-1

-2

-1

-4

-1

-2

-1

-2

-1

-4

-1

-2

-1

-2

-1

-2

它们的载荷列阵都为

0000U5V5U6V6

位移列阵为

-lu^

U3V3

位移约束条件为u5

形成整体平衡方程二V5二5二V6=0，将此约束条件引入整体刚度方程，对

■3

-1

-2

0「5]

■01

-2

-1

P/2

-1

-2

-1

-2

-1

P/2

-2

-1

-2

-1

-4

-2

-1

—4

01一

•°一

由上式可知在这种情况下,

两种模型求解后对应节点的位移相等。

有整体节

点位移获取单元节点位移，所以对应的单元节点位移也相等。

以单元①为例，两种模型应力矩阵的关系s（（b）=2s；a），又由｛可=〔six｝可

得，模型（b）中单元①的应力是模型（a）的2倍，其它单元可得到类似的结论。

当两个模型面力的集度相等，可设模型（a）右端受剪力的合力为2P，模型

（b）右端受剪力的合力为P,则两种模型的载荷列阵分别为

F（a）j0POPOOOOU5V5U6VJ

00000U5V5U6V6

由上述内容可得F（a）=2F（b），进而可知（a）模型中对应节点的位移是（b）模型的2倍

在单元①中，醐=2S；a），所以由2｝=[S】｛3｝可知模型（a）和模型（b）

中单元①的应力相等，其他单元可类似说明

2、证明平面问题三节点三角形单元发生刚体位移（小位移平动和转动）时，单元中将不产生应力。

（10分）

证明：

三节点三角形单元位移模式选取一次多项式:

（1）

u=12X

V=5X6讨

在

（1）的1式中带入节点i的坐标（x，yi）得到节点i在x向的位移ui，

同理可得，

u二爲•：

2人:

3yi

（2）

比八1T~yj

Um='1：

2Xm：

3ym

解

（2）式可得到广义坐标由节点位移表示的表达式:

1=2A（aiU二丄（bu

r（Cu

1/（aiv

（bv

（CV

ajUj

展开阅读全文