资源描述:
《平面应力应变有限元分析程序设计课程设计22》由会员上传分享,免费在线阅读,更多相关内容在行业资料-天天文库。
1、平面应力/应变有限元分析程序设计指导教师李欣宇200807152032吴义军200807152015刘刚21-引言有限元法是求解微分方程边值问题的一种通用数值方法,该方法是一种基于变分法(或变分里兹法)而发展起来的求解微分方程的数值计算方法,以计算机为手段,采用分片近似,进而逼近整体的研究思想求解物理问题。摘要有限元法的基本思想是将物体(即连续的求解域)离散成有限个且按一定方式相互联结在一起的单元的组合,来模拟或逼近原来的物体,从而将一个连续的无限自由度问题简化为离散的有限自由度问题求解的一种数值分析方法。物体被离散后,通过对其中各个单元进行单元分析,最终得到对整个物体的分析。
2、网格划分中每一个小的块体称为单元。确定单元形状、单元之间相互联结的点称为节点。单元上节点处的结构内力为节点力,外力(有集中力、分布力等)为节点载荷。三角形常应变单元求解结构主程序l功能:运用有限元法中三角形常应变单元解平面问题的计算主程序。l基本思想:单元结点按右手法则顺序编号。l荷载类型:可计算结点荷载。l说明:主程序的作用是通过赋值语句、读取和写入文件、函数调用等完成算法的全过程,即实现程序流程图的程序表达。21-1.1程序流程图计算单刚并将单刚加入总纲引入边界条件生成载荷向量计算节点位移及支反力计算单元应力及主应力输出结果开始导入数据根据题意将得到的程序运行所需的初始数据
3、的input.txt导入。由MATLAB中的LinearTriangleElementStiffness计算单刚矩阵和LinearTriangleAssemble依次加入总刚矩阵。根据问题将约束信息引入计算得到载荷向量,由F=KU计算支反力,再次调用MATLAB中的LinearTriangleElementStresses和LinearTriangleElementPStresses计算得到到单元应力和主应力。详细说明见附录。21-1.2程序应用举例【例1】如图所示,受均匀分布载荷作用的薄平板结构。将平板离散化成两个线性三角单元,如图所示。假定E=210Gpa,,t=0.025
4、m,w=3000kn/m。解答如图所示单元连通性如下表所示单元编号节点i节点j节点m11342123由题意可得输入数据文件input.txt为:4224210e60.30.0251341230.00.00.50.00.50.250.00.2529.375039.3750111200300411说明:第一行:读入程序控制信息NPION=fscanf(FP1,'%d',1)%结点个数(结点编码总数)NELEM=fscanf(FP1,'%d',1)%单元个数(单元编码总数)NFORCE=fscanf(FP1,'%d',1)%结点荷载个数21-NVFIX=fscanf(FP1,'%d'
5、,1)%受约束边界点数YOUNG=fscanf(FP1,'%e',1)%弹性模量POISS=fscanf(FP1,'%f',1)%泊松比THICK=fscanf(FP1,'%d',1)%厚度第二、三行:读入单元连接信息:LNODS=fscanf(FP1,'%d',[3,NELEM])';%单元定义数组,单元结点号,逆时针输入第四、五、六、七行:读入结点坐标COORD=fscanf(FP1,'%f',[2,NPOIN])';%结点坐标值,第1列为x坐标值,第2列为y坐标值第八、九行:读入结点荷载信息FORCE=fscanf(FP1,'%f',[3,NFORCE])';%结点号,X
6、向结点荷载数值,Y向结点荷载数值(与坐标轴方向一致为正)第十、十一、十二、十三行:读入零位移信息FIXED=fscanf(FP1,'%d',[3,inf])';%结点号,X向约束,Y向约束下述两个例子思路相同不再赘述%-----------------------------------------------------------------------------------------------------运行该程序得到输出文件output.txt如下所示节点号x位移y位移10.00000000e+0000.00000000e+00027.11111747e-0061
7、.11517786e-00636.53122498e-0064.46071143e-00840.00000000e+0000.00000000e+000节点号FxFy1-9.37500000e+000-5.62950360e+00029.37500000e+0006.66133815e-01639.37500000e+0002.22044605e-0164-9.37500000e+0005.62950360e+000单元号单元应力σx单元应力σy单元应力τxy13.01441153e+003