有限元-计算结构力学-大作业

SHANGHAIJIAOTONGUNIVERSITY平面应力问题解的Matlab实现姓名:heiya168学号:帆哥班级：指导老师：目录1绪论.........................................................................................................................12平面问题的四节点四边形单元.............................................................................22.1单元的构造...............................................................................................22.2等参变换...................................................................................................52.3边界条件的处理——置“1”法.............................................................83有限元分析流程...................................................................................................103.1程序原理和流程.....................................................................................103.2使用的函数.............................................................................................113.3文件管理.................................................................................................113.4数据文件格式.........................................................................................114算例——开方孔的矩形板拉伸分析...................................................................134.1问题的具体参数与载荷.........................................................................134.2Matlab程序计算....................................................................................134.3ANSYS建模计算......................................................................................154.4误差分析.................................................................................................175总结.......................................................................................................................18参考文献......................................................................................................................18附录..............................................................................................................................1911绪论有限元方法(finiteelementmethod)，是求取复杂微分方程近似解的一种非常有效的工具，是现代数字化科技的一种重要基础性原理。将它用于在科学研究中，可成为探究物质客观规律的先进手段。将它应用于工程技术中，可成为工程设计和分析的可靠工具。弹性体在载荷作用下，其基本方程可写成以下的三类方程和两种边界条件。平衡方程——应力与外载荷的关系；几何方程——应变位移关系；物理方程——应力应变关系；力的边界条件；几何边界条件。应用最小位能原理，并利用上述关系，最终建立由刚度方程，节点位移和等效节点载荷所构成的求解方程。带入边界条件求解方程，就可以得出弹性力学问题的一般性解答。本次大作业基于有限元方法的基本原理，使用Matlab这一平台，针对平面应力问题，采用四节点四边形单元编写了求解单元节点位移的程序。主要内容包括：1）介绍有限元的基本原理；2）编程基本思路及流程介绍；3）程序原理及说明；4）具体算例这四个部分。22平面问题的四节点四边形单元2.1单元的构造（1）单元的几何和节点描述平面4节点矩形单元如图4-6所示，单元的节点位移共有8个自由度。节点的编号为1、2、3、4，各自的位置坐标为(xi,yi),i=1,2,3,4，各个节点的位移(分别沿x方向和y方向)为(ui,vi),i=1,2,3,4。图2.1平面4节点矩形单元若采用无量纲坐标,=xyab（2.1）则单元4个节点的几何位置为112233441=11=11=11=-1（2.2）将所有节点上的位移组成一个列阵，记作eq；同样，将所有节点上的各个力也组成一个列阵，记作eP，那么311223344(81)11223344(81)[vvuvuv][]eTeTxyxyxyxyquuPPPPPPPPP（2.3）若该单元承受分布外载，可以将其等效到节点上，也可以表示为如式(2.3)所示的节点力。利用函数插值、几何方程、物理方程以及势能计算公式，可以将单元的所有力学参量用节点位移列阵及相关的插值函数eq来表示；下面进行具体的推导。（2）单元位移场的表达从图2-1可以看出，节点条件共有8个，即x方向4个1234(,,,)uuuu，y方向4个1234(,,,)vvvv，因此，x和y方向的位移场可以各有4个待定系数，即取以下多项式作为单元的位移场模式01230123(,)(,)uxyaaxayaxyvxybbxbybxy（2.4）它们是具有完全一次项的非完全二次项，以上两式中右端的第四项是考虑到x方向和y方向的对称性而取的，除此外xy项还有个重要特点，就是“双线性”，当x或y不变时，沿y或x方向位移函数呈线性变化，这与前面的线性项最为相容，而2x或2y项是二次曲线变化的。因此，未选2x或2y项。由节点条件，在x=xi，y=yi处，有(,)(,)iiiiiiuxyuvxyv1,2,3,4i（2.5）将式(2.5)代入式(2.4)中，可以求解出待定系数a0，…，a3和b0，…，b3，然后代回式(4-52)中，经整理后有1122334411223344(,)(,)(,)(,)(,)(,)(,)(,)(,)(,)uxyNxyuNxyuNxyuNxyuvxyNxyvNxyvNxyvNxyv（2.6）其中412341(,)(1)(1)41(,)(1)(1)41(,)(1)(1)41(,)(1)(1)4xyNxyabxyNxyabxyNxyabxyNxyab（2.7）如以无量纲坐标系(2.1)来表达，式(2.7)可以写成1(1)(1)4iiiN1,2,3,4i（2.8）将式(2.6)写成矩阵形式，有11212342(21)(28)(81)12343344(,)u(,)=(,)euvuNNNNvuxyxyNqNNNNuvxyvuv（2.9）其中，(,)Nxy为该单元的形状函数矩阵。（3）单元应变场的表达由弹性力学平面问题的几何方程(矩阵形式)，有单元应变的表达(38)(31)(21)(28)(81)(81)(32)(32)(,)uxxeeyyxyxyNqBq（2.10）其中几何矩阵B(x,y)为1234(38)(28)(32)1234(,)xNNNNBxyNNNNNyyx5=1234(32)(32)(32)(32)BBBB（2.11）式(4-59)中的子矩阵iB为(32)1,2,3,4iiiiiNxNBiyNNyx（2.12）（4）单元应力场的表达由弹性力学中平面问题的物理方程，可得到单元的应力表达式(33)(33)(38)(31)(31)(38)(81)(81)eeDDBqSq（2.13）（5）单元应力场的表达以上已将单元的三大基本变量(,,)u用基于节点位移列阵来进行表达，见式(2.9)、式(2.10)及式(2.13)；将其代入单元的势能表达式中，有12eeTeeeTeqKqPq（2.14）其中eK是4节点矩形单元的刚度矩阵。将单元的势能对节点位移eq取一阶极值，可得到单元的刚度方程(88)(81)(81)eeeKqP（2.15）2.2等参变换由前面的单元构造过程可以看出，一个单元的关键就是计算它的刚度矩阵，而由刚度矩阵的构成可知要实现两个坐标系中单元刚度矩阵的变换，必须计算两个坐标系之间的三种映射关系：6坐标映射(,)(,)xy（2.16）偏导数映射(,)(,)xy（2.17）面积映射eeAAdxdydd（2.18）图2.2矩形单元映射为任意四边形单元（1）两个坐标系之间的函数映射设如图4-17所示的两个坐标系的坐标映射关系为(,)(,)xxyy（2.19）x和y方向上可以分别写出各包含有4个待定系数的多项式，即01230123(,)(,)xaaaaybbbb（2.20）其中待定系数a0，…，a3和b0，…，b3可由节点映射条件(2.19)来唯一确定。对照前面4节点矩形单元的单元位移函数式(2.5)，映射函数式(2.20)具有完全相同的形式，同样，将求出的待定系数再代回式(2.20)中，重写该式为1122334411223344(,)(,)(,)(,)(,)(,)(,)(,)(,)(,)xNxNxNxNxyNyNyNyNy（2.21）7其中1(1)(1)i=1,2,3,44iiiN（2.22）11212342(21)(28)(81)12343344(,)x(,)=(,)(,)exyxNNNNyxNqNNNNxyyxy（2.23）这就可以实现两个坐标系间的映射。（2）两个坐标系之间的偏导数映射对物理坐标系(x,y)中的任意一个函数Φ(