基于MATLAB的有限元法分析平面应力应变问题刘刚.docx
《基于MATLAB的有限元法分析平面应力应变问题刘刚.docx》由会员分享,可在线阅读,更多相关《基于MATLAB的有限元法分析平面应力应变问题刘刚.docx(10页珍藏版)》请在冰豆网上搜索。
基于MATLAB的有限元法分析平面应力应变问题刘刚
姓名:
刘刚
学号:
15
平面应力应变分析有限元法
Abstruct:
本文通过对平面应力/应变问题的简要理论阐述,使读者对要分析的问题有大致的印象,然后结合两个实例,通过MATLAB软件的计算,将有限元分析平面应力/应变问题的过程形象的展示给读者,让人一口了然,快速了解有限元解决这类问题的方法和步骤!
1.基本理论
有限元法的基本思路和基本原则以结构力学中的位移法为基础,把复朵的结构或连续体看成有限个单元的组合,各单元彼此在节点出连接而组成整体。
把连续体分成有限个单元和节点,称为离散化。
先对单元进行特性分析,然后根据节点处的平衡和协调条件建立方程,综合后做整体分析。
这样一分一合,先离散再综合的过程,就是把复杂结构或连续体的计算问题转化简单单元分析与综合问题。
因此,一般的有限揭发包括三个主要步骤:
离散化单元分析整体分析。
2.用到的函数
1.LinearTriangleElementStiffness(E,XU,t,xi,yi,xj,yj,xm,ym,p)
2.LinearBarAssemble(KkIf)
3.LinearBarElementForces(ku)
4.LinearBarElementStresses(kuA)
5.LinearTriangleElementArea(EXUt)
三•实例
例1•考虑如图所示的受均布载荷作用的薄平板结构。
将平板离散化成两个线
性三角元,假定E=200GPa,v二0.3,t=0.025m,w=3000kN/m・
H-2示例11」的嚴平板结构
谢11・3用双线件三角形离散化薄平板
1.离散化
2•写出单元刚度矩阵
通过matlab的LinearTriangleElementStiffness函数,得到两个单元刚度
矩阵
和
每个矩阵都是6
6的。
»E=210e6
210000000
>>kl=LinearTriangleElementStiffness(E,XU,t,0,0,0.5,0.25,0,0.25,1)
kl二
1.0e+006*
Columns1through5
2.0192
-1.0096-2.0192
0.8654
5.7692-0.8654
1.0096-5.7692
0.8654-0.5048-1・8750
Column6
1.0096
-5.7692
0.8654
-0.5048
-1.8750
6.2740
»NU二0.3
NU二
0.3000
»t二0.025
t二
0.0250
>>k2=LinearTriangleElementStiffness(E,NU,t,0,0,0.5,0,0.5,0.25,1)k2二
1.0e+006*
Columns1through5
1.44230-1.44230.86540
00.50481.0096-0.5048-1.0096
-1.44231.00963.4615-1.8750-2.0192
0.8654-0.5048-1.87506.27401.0096
-0.8654
0-1.0096-2.01921.00962.0192
0.8654-5.7692
Column6
-0.8654
0
0.8654-5.7692
0
5.7692
3.集成整体刚度矩阵
8*8零矩阵
K二
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
0
>>K=LinearTriangleAssemble(K,kl,1,3,4)
1.0e+006*
Columns1through5
00.8654-5.7692
000
000
0.5048
1.0096
-0.5048
1.0096
3.4615
-1.8750
-0.5048
-1.8750
6.2740
>>K=LinearTriangleAssemble(K,kl,1,2,3)
K二
l・0e+007*
0-0.1010
0.0505
-0.2019
-0.1442
0.1010
0.4904
0.1875
-0・1442
0.0865
-0.5769
0.0865-0.0505
0.1875
0.6779
0.1010-0・0505
-0.2019
0.0865
0.14420.10100.3462-0・1875
0.0505-0・18750.6274
4.引入边界条件•用上一步得到的整体刚度矩阵,可以得到该结构的方程组如下形式
本题的边界条件:
3.46150・1.44230.86540・1.8750-2.01921.0096
5
06.27401.0096・0.5048-1.875000.8654-5.7692
%
町
-1.44231.00963.4615・1.8750-2.01920.865400
5
0.8654・0.5048-1.87506.27401.0096・5.769200
%
0-1.8750-2.01921.00963.46150-1.44230.8654
%
Pax
-1.875000.8654-5.769206.27401.0096-0.5048
5
片
-2.01920.865400-1.44231.00963.4615-1.8750
u4X
1.0096-5.7692000.8654-0.5048・1.875062740
.5.
F*
F2x=9.375如=0,心=9375,^=0
将边界条件带入,得到:
3.46150-1.44230.86540-1.8750-2.01921.0096
0
Fix
06.27401.0096-0.5048-1.875000.8654-5.7692
0
町
-1.44231.00963.4615-1.8750・2.01920.865400
5
9375
0.8654・0.5048・1.87506.27401.0096-5.769200
%
0
0-1.8750-2.01921.00963.46150-1.44230.8654
5
9.375
-1.875000.8654-5.769206.27401.0096・0.5048
0
-2.01920.865400•1.44231.00963.4615-1.8750
0
1.0096-5.7692000.8654-0.5048-1.87506.2740
0
[%
5.解方程
分解上述方程组,提取总体刚度矩阵K的第3-6行的第3-6列作为子矩阵
3.4615
-1.8750
-2.0192
0.8654
9・375「
-1.8750
6.2740
1.0096
•5.7692
0
-2.0192
1.0096
3.4615
0
9.375
0.8654
■
-5.7692
0
6.2740
%-
0
■■
Matlab命令
»k=K(3:
6,3:
6)
k-
1.0e+006
*
3.4615
-1.8750-2.0192
-1.8750
6.27401.
0096
-2.0192
1.00963.4615
0.8654
-5.7692
0
»f二[9・375;0;9・375;0]
0.8654
-5.7692
0
6.2740
f二
9.3750
0
9.3750
0
»u二k\f
u=
1.0e-005*
0.7111
0.1115
0.6531
0.0045
现在可以清楚的看出,节点2的水平位移和垂直位移分别是o.7111m和0-111施。
节点3的水平位移和垂直位移分别是0.6531m和0.0045m。