清华大学-弹性力学有限元大作业
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
弹性力学有限元大作业
一、模型信息：
已知：材料为铝合金。E=71GPa ，v=0.3.
矩形平板的几何参数：板长为480mm ，宽为360mm ，厚度为2mm ；图形如下图；
加肋平板：
二、matlab 编程实现
1、程序相关说明：
计算使用的软件为：matlab2010a 主函数：main.m 主要计算部分
子函数：Grids.m 生成网格，节点数为：+1*+1I J （）（）
、单元数： 2**I J AssembleK.m 将单元刚度矩阵组装成总刚度矩阵（叠加方法）
GenerateB.m 生成单元格e B 矩阵 GenerateS.m 生成单元格e S 矩阵 GenerateK.m 生成单元刚度矩阵
2、网格划分：
利用Grid.m 子函数，取2020I J ==、，即可以得到网格如下： 节点数为：441个，单元格数：800个
3、计算过程及结果 （1）、网格划分：通过Grid.m ，生成节点数为：441个、单元格数：800个的网格 （2）、生成总刚度矩阵K ：通过GenerateK.m 、AssembleK.m 生成总刚度矩阵 采用常应变三角单元，e e u N a =,易得=e e B LN
由平面应力问题，可以确定2101011002E D νννν⎡⎤
⎢⎥
⎢
⎥=⎢⎥-⎢⎥
-⎢⎥⎣⎦
即e e S DB =
单元刚度矩阵为：e eT e K AtB DB = 总刚度矩阵为：eT
e e e
K G
K G =
∑
（3）、求解过程：
系统平衡方程为：Ka P = 将方程进一步划分为：E EF E E E T F F EF
F K K d f r d f K K +⎡⎤⎡⎤⎡⎤
=⎢
⎥⎢⎥⎢⎥
⎣⎦⎣⎦
⎣⎦ 通过已知边界条件（位移、载荷），确定E E F d f f 、、 ，从而将K 矩阵划分为四个模
块：E EF T
EF F K K K K ⎡⎤⎢⎥⎣⎦
1
()
E E E E
F F E T
F F F EF E r K d K d f d K f K d -=+-=-支反力：部分位移：
即整体位移向量为：E F d a d ⎡⎤
=⎢⎥⎣⎦
整体力边界条件为：E E F f r P f +⎡⎤
=⎢
⎥⎣⎦
（4）后处理：（应力、应变、抹平）
a、单元应力、应变：
e e e e e e
S a
B a σ
ε
=
=
b、抹平得到节点应力、应变：将每个节点参与组成的单元应力、应变叠加，然后除以叠加的单元数，得到抹平后的节点应力、应变。
（5）计算结果：
由于K矩阵行、列数过多，故以附件形式放在excel表格中
具体数据见附件：计算结果.xlsx
题目1：
（1）位移场
节点位置图：蓝色为初始位置，红色为最终位置
x、y的位移场：
（2）、应力场
（3）、应变场
题目2
（1）位移场：
节点位置：蓝色为初始位置，红色为最终位置
x、y的位移场：
（2）、应力场：（分别为x y xy σσσ、、 ）
（3）、应变场（分别为x y xy εεε、、 ）
三、Abaqus 分析
1、平板：（文件为：plate1.cae ；plate2.cae ）
题目1：
（1）位移场
（2）应力场（分别为x y xy σσσ、、 ）
（3）应变场（分别为x y xy εεε、、 ）
题目2：
（1）位移场
X 方向： y 方向：
（2）应力场（分别为x y xy σσσ、、 ）
（3）应变场（分别为x y xy εεε、、 ）
2、加肋后的平板（文件为：plate3.cae ；plate4.cae ）
题目1：
（1）位移场
（2）应力场（分别为x y xy σσσ、、 ）
（3）应变场（分别为x y xy εεε、、 ）
题目2：
（1）、位移场：
（2）、应力场：（分别为x y xy σσσ、、 ）
（3）、应变场（分别为x y xy εεε、、 ）
四、误差分析
题目1：
1、位移场
Matlab计算结果：
Abaqus计算结果：
误差分析：由于采用的网格划分都是三角单元，并且节点都比较多，所以x、y方向的位移两者计算结果、变化趋势大致相似。matlab计算：x方向最大值为0.0017mm，y方向最小值为-0.0100mm，abaqus计算：x方向最大值为0.001758mm，y方向最小值为-0.01mm。x 方向的变形大致呈鼓形，y方向的变形大致呈线性（层状），可能由于边缘效应，两边会出现一定的角度。
2、应力场
Matlab计算结果：
Abaqus计算结果：
误差分析：两者变化趋势大致相同，但是计算结果误差较大，尤其是角点位置。以x方向为例，matlab：x方向最大应力值5.4140e+003，最小应力值-8.7261e+005，abaqus：x方向最大应力值为5.98e+003，最小应力值-12.42e+005，这可能跟角点边缘效应、抹平方法在角点效果不大以及abaqus固支定义不同，所以误差比较大。
3、应变场
Matlab计算结果：