计算力学课程设计模板改副本概要Word下载.docx

上传人:b****3 文档编号:18202446 上传时间:2022-12-14 格式:DOCX 页数:13 大小:144.61KB
下载 相关 举报
计算力学课程设计模板改副本概要Word下载.docx_第1页
第1页 / 共13页
计算力学课程设计模板改副本概要Word下载.docx_第2页
第2页 / 共13页
计算力学课程设计模板改副本概要Word下载.docx_第3页
第3页 / 共13页
计算力学课程设计模板改副本概要Word下载.docx_第4页
第4页 / 共13页
计算力学课程设计模板改副本概要Word下载.docx_第5页
第5页 / 共13页
点击查看更多>>
下载资源
资源描述

计算力学课程设计模板改副本概要Word下载.docx

《计算力学课程设计模板改副本概要Word下载.docx》由会员分享,可在线阅读,更多相关《计算力学课程设计模板改副本概要Word下载.docx(13页珍藏版)》请在冰豆网上搜索。

计算力学课程设计模板改副本概要Word下载.docx

2.基本力学模型

如图所示的一个块体,在右端面上端点受集中力F作用。

基于MATLAB平台,计算各个节点位移,支座反力以及单元的应力。

取相关的参数E=1.5e10Pa,u=0.25,F=1.5e5N。

对该问题进行有限元分析的过程如下:

图一

(1)结构的离散化与编号

将结构离散为5个4节点四面体单元,单元编号和坐标如上图所示,连接关系见表一,节点的坐标见表二。

表一单元连接关系

单元号

节点号

1

1426

2

1437

3

6751

4

6784

5

1467

表二节点的坐标

节点

节点坐标/m

x

y

z

6

7

8

0.4

2.0

1.0

节点位移列阵

节点外载列阵

其中

约束的支反力列阵

总的节点载荷列阵

(2)计算各单元的刚度矩阵(以国际标准单位)

首先在MATLAB环境下,输入弹性模量E、泊松比NU,然后针对单元1和单元2,分别5次调用函数Tetrahedron3D4Node_Stiffness,就可以得到单元的刚度矩阵k1(6×

6)~k5(6×

6)。

>

E=2.5e10;

NU=0.25;

k1=Tetrahedron3D4Node_Stiffness(E,NU,0,0,0,0.4,2.0,0,0.4,0,0,0.4,0,1.0);

k2=Tetrahedron3D4Node_Stiffness(E,NU,0,0,0,0.4,2.0,0,0,2.0,0,0,2.0,1.0);

k3=Tetrahedron3D4Node_Stiffness(E,NU,0.4,0,1.0,0,2.0,1.0,0,0,1.0,0,0,0);

k4=Tetrahedron3D4Node_Stiffness(E,NU,0.4,0,1.0,0,2.0,1.0,0.4,2.0,1.0,0.4,2.0,0);

k5=Tetrahedron3D4Node_Stiffness(E,NU,0,0,0,0.4,2.0,0,0.4,0,1.0,0,2.0,1.0);

(3)建立整体刚度方程

由于该结构共有8个节点,则总共的自由度数为24,因此,结构总的刚度矩阵为KK(24×

24),先对KK清零,然后5次调用函数Tetrahedron3D4Node_Assembly进行刚度矩阵的组装。

KK=zeros(24);

KK=Tetrahedron3D4Node_Assembly(KK,k1,1,4,2,6);

KK=Tetrahedron3D4Node_Assembly(KK,k2,1,4,3,7);

KK=Tetrahedron3D4Node_Assembly(KK,k3,6,7,5,1);

KK=Tetrahedron3D4Node_Assembly(KK,k4,6,7,8,4);

KK=Tetrahedron3D4Node_Assembly(KK,k5,1,4,6,7);

(4)边界条件的处理及刚度方程求解

由图4-22可以看出,节点1,2,5和6上3个方向的位移将为零,即

因此,将针对节点3,4,7和8的位移进行求解,节点1,2,5和6的位移将对应KK矩阵中的第1~6行,第13~18行和第1~6列,第13~18列,需从KK(24×

24)中提出,置给k,然后生成对应的载荷列阵p,再采用高斯消去法进行求解。

注意:

MATLAB中的反斜线符号“\”就是采用高斯消去法。

k=KK([7:

12,19:

24],[7:

24]);

p=[0,0,0,0,0,0,0,0,-2e5,0,0,-2e5]'

u=k\p

u=1.0e-003*

0.0958-0.0295-0.27690.0994-0.0430-0.2761[将列排成行排量[]___________________________________________________________________________________________________________________________]

0.09870.0500-0.28890.10020.0380-0.3000[将列排成行排量[]___________________________________________________________________________________________________________________________]

所求得的位移结果见表三。

表三空间块体的节点位移计算结果

=0.0958×

10-3

=0.0987×

=-0.0295×

=0.0500×

=-0.2769×

=-0.2889×

=0.0994×

=0.1002×

=-0.0430×

=0.0380×

=-0.2761×

=-0.3000×

(5)支反力的计算

在得到整个结构的节点位移后,由原整体刚度方程就可以计算出对应的支反力;

先将上面得到的位移结果与位移边界条件的节点位移进行组合(注意位置关系),可以得到整体的位移列阵U(24×

1),再代回原整体刚度方程,计算出所有的节点力P(24×

1),按式(4-192)的对应关系就可以找到对应的支反力。

U=[zeros(6,1);

u([1:

6]);

zeros(6,1);

u(7:

12)];

P=KK*U

P=1.0e+005*

0.92304.07170.3908-1.04813.92831.2071[将列排成行排量[]___________________________________________________________________________________________________________________________]

-0.0000-0.00000.00000.00000-0.0000[将列排成行排量[]___________________________________________________________________________________________________________________________]

-1.1625-4.07171.29651.2876-3.92831.1056[将列排成行排量[]___________________________________________________________________________________________________________________________]

-0.00000.0000-2.00000.0000-0.0000-2.0000[将列排成行排量[]___________________________________________________________________________________________________________________________]

由式(4-193)的对应关系,可以得到对应的支反力见表四。

表四空间块体的支反力计算结果

(6)各单元的应力计算

先从整体位移列阵U(24×

1)中提取出单元的位移列阵,然后,调用计算单元应力的函数Tetrahedron3D4Node_Stress,就可以得到各个单元的应力分量。

u1=[U(1:

3);

U(10:

12);

U(4:

6);

U(16:

18)];

stress1=Tetrahedron3D4Node_Stress(E,NU,0,0,0,0.4,2.0,0,0.4,0,0,0.4,0,1.0,u1)

stress1=1.0e+006*

-0.2150-0.6449-0.21500.4972-1.38070[将列排成行排量[]___________________________________________________________________________________________________________________________]

u2=[U(1:

U(7:

9);

U(19:

21)];

stress2=Tetrahedron3D4Node_Stress(E,NU,0,0,0,0.4,2.0,0,0,2.0,0,0,2.0,1.0,u2

stress2=1.0e+005*

0.0885-4.7076-4.16401.4161-5.89410.4868[将列排成行排量[]___________________________________________________________________________________________________________________________]

u3=[U(16:

21);

U(13:

15);

U(1:

3)];

stress3=Tetrahedron3D4Node_Stress(E,NU,0.4,0,1.0,0,2.0,1.0,0,0,1.0,0,0,0,u3)

stress3=1.0e+006*

0.25000.75010.25000.4936-1.44470[将列排成行排量[]___________________________________________________________________________________________________________________________]

u4=[U(16:

U(22:

24);

stress4=Tetrahedron3D4Node_Stress(E,NU,0.4,0,1.0,0,2.0,1.0,0.4,2.0,1.0,0.4,2.0,0,u4)

stress4=1.0e+005*

0.66673.7041-4.86282.0178-6.8964-2.6756[将列排成行排量[]___________________________________________________________________________________________________________________________]

u5=[U(1:

stress5=Tetrahedron3D4Node_Stress(E,NU,0,0,0,0.4,2.0,0,0.4,0,1.0,0,2.0,1.0,u5)

stress5=1.0e+005*

-0.1936-0.0239-1.6540-6.6710-9.47741.5635[将列排成行排量[]___________________________________________________________________________________________________________________________]

各个单元应力分量的计算结果列在表五中。

表五空间块体的各个单元应力分量的计算结果

1号单元

2号单元

3号单元

4号单元

5号单元

3计算结果分析

在如图一的悬臂式梁结构模型中,当受到如图的集中荷载作用时,1,2,5,6节点的y方向上的支反力较大,这符合竖向荷载作用下的支承力实际情况;

同时,通过表五对空间块体的各个单元应力分量的计算结果可以看出5号单元的剪切力最大,容易发生剪切破坏,因此需要对其进行加固处理。

4.参考文献

[1]曾攀.有限元基础教程[M].北京:

高等教育出版社,2009,7

[2]王勖成.有限单元法[M].北京:

清华大学出版社,2003,

[3]张国瑞有限元法北京机械工业出版社1991

[4]杨帆宋小军悬臂梁的有限元分析江西理工大学2008

[5]江见鲸陆新征混凝土结构有限元分析北京:

清华大学出版社2005

附录(程序)

1.Tetrahedron3D4Node_Stiffness函数

functionk=Tetrahedron3D4Node_Stiffness(E,NU,x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4)

%该函数计算单元的刚度矩阵

%输入弹性模量E,泊松比NU

%输入4个节点i、j、m、n的坐标xi,yi,zi,xj,yj,zj,xm,ym,zm,xn,yn,zn

%输出单元刚度矩阵k(12X12)

%---------------------------------------------------------------

%以下数组与书中的对应关系

xyz=[1x1y1z1;

1x2y2z2;

1x3y3z3;

1x4y4z4];

V=det(xyz)/6;

mbeta1=[1y2z2;

1y3z3;

1y4z4];

mbeta2=[1y1z1;

mbeta3=[1y1z1;

1y2z2;

mbeta4=[1y1z1;

1y3z3];

mgamma1=[1x2z2;

1x3z3;

1x4z4];

mgamma2=[1x1z1;

mgamma3=[1x1z1;

1x2z2;

mgamma4=[1x1z1;

1x3z3];

mdelta1=[1x2y2;

1x3y3;

1x4y4];

mdelta2=[1x1y1;

mdelta3=[1x1y1;

1x2y2;

mdelta4=[1x1y1;

1x3y3];

beta1=-1*det(mbeta1);

beta2=det(mbeta2);

beta3=-1*det(mbeta3);

beta4=det(mbeta4);

gamma1=det(mgamma1);

gamma2=-1*det(mgamma2);

gamma3=det(mgamma3);

gamma4=-1*det(mgamma4);

delta1=-1*det(mdelta1);

delta2=det(mdelta2);

delta3=-1*det(mdelta3);

delta4=det(mdelta4);

B1=[beta100;

0gamma10;

00delta1;

gamma1beta10;

0delta1gamma1;

delta10beta1];

B2=[beta200;

0gamma20;

00delta2;

gamma2beta20;

0delta2gamma2;

delta20beta2];

B3=[beta300;

0gamma30;

00delta3;

gamma3beta30;

0delta3gamma3;

delta30beta3];

B4=[beta400;

0gamma40;

00delta4;

gamma4beta40;

0delta4gamma4;

delta40beta4];

B=[B1B2B3B4]/(6*V);

D=(E/((1+NU)*(1-2*NU)))*[1-NUNUNU000;

NU1-NUNU000;

NUNU1-NU000;

000(1-2*NU)/200;

0000(1-2*NU)/20;

00000(1-2*NU)/2];

k=abs(V)*B'

*D*B;

2.Tetrahedron3D4Node_Assembly函数

functiony=Tetrahedron3D4Node_Assembly(KK,k,i,j,m,n)

%该函数进行单元刚度矩阵的组装

%输入单元刚度矩阵k

%输入单元的节点编号i、j、m、n

%输出整体刚度矩阵KK

DOF=[3*i-2:

3*i,3*j-2:

3*j,3*m-2:

3*m,3*n-2:

3*n];

forn1=1:

12

forn2=1:

KK(DOF(n1),DOF(n2))=KK(DOF(n1),DOF(n2))+k(n1,n2);

end

end

y=KK;

3.Tetrahedron3D4Node_Stress函数

functiony=Tetrahedron3D4Node_Stress(E,NU,x1,y1,z1,x2,y2,z2,x3,y3,z3,x4,y4,z4,u)

%该函数计算单元的应力

%输入单元的位移列阵u(12X1)

%输出单元的应力stress(6X1)

%由于它为常应力单元,则单元的应力分量为Sx,Sy,Sz,Sxy,Syz,Szx

delta4=det(mdelta4

展开阅读全文
相关资源
猜你喜欢
相关搜索

当前位置:首页 > 工程科技 > 材料科学

copyright@ 2008-2022 冰豆网网站版权所有

经营许可证编号:鄂ICP备2022015515号-1