《弹性力学与有限元》第5章平面问题有限元分析Word格式.docx

1、在平面应力问题中,最简单、最常用的单元类型是三角形三节点单元。如果我们选择这种类型的单元,就可以用一组网格线把求解区域图5-2实线部分划分为若干个三角形单元,如图5-3所示。网格的交点或单元的角点就是节点。每个单元有三个节点,把每个节点沿坐标轴的位移取为基本未知点数。由于是平面问题,所以每个节点有两个自由度,所有节点都图5-2结构的简化图5-3有限元剖分图可以取为铰接。对于非节点上的外载荷按静力等效的原则移置到节点上去,成为节点载荷。这样就把平面连续体（或结构离散成为由有限个单元和节点所组成的等效集合体。有了计算模型,就可以用处理杆系问题的基本思想来分析单元和结构特性。3.结构离散化之后,对单

2、元和节点分别进行编号将结构离散化之后,把所有的单元及节点按一定顺序分别进行编号,以便进行计算。编号时应注意:单元号及节点号均不能有错漏或重复。需要指出,编号的顺序不影响计算结果,原则上可以任意编排。但是为了节省计算机内容,减少计算时间,单元中每个节点的编号与周围节点的编号应尽可能接近,即同一个单元内的节点编号差值的最大值尽量的小,这样就可以使在求解线性方程组时,使系数矩阵中非零元素的宽度（即带宽为最小。4.定义单元单元和节点分别编号后,还要定义单元（即定义每个单元的i,j,m,且单元一经定义后,整个计算中不允许再改变了。最好把求解区域分成内部单元和边界单元。对于内部单元,三个节点i,j,m的位

3、置可以是任意的。对于边界单元,为了使计算公式简化,格式统一起见,只准有一条单元边界处于求解区域边界上,并且规定节点j,m处于边界上。5.1.3 离散化应注意的问题划分有限单元网格时,必须注意以下几点:1.单元数目的确定应兼顾精度、经济性和计算机容量。从有限元法本身来讲,单元划分得越细,节点布置得越多,越接近真实求解区域,因而计算精度越高。但是随之而来的就是计算时间、计算费用和计算机内存的增加。所以要考虑各因素来恰当地确定单元的数目。其原则是,在满足工程精度要求的前提下,单元数目应尽量划分得少一些。2.在初步分析的基础上,再离散化。首先根据求解区域的形状、载荷分布情况和边界条件大致分析,在此基础

4、上进行有限元分割（或剖分,划分。3.如果求解区域具有不同的厚度或者由两种以上的材料组成时,不要把厚度不同或材料不同的区域划分在一个单元中。4.单元的边长尽量接近,以便提高计算精度。对于三角形单元,不应出现过大的钝角或过小的锐角。对于四边形单元,长度和宽度也不宜相差太大。5.任意一个单元的角点（顶点或节点必须同时也是相邻单元的角点,而不能是相邻单元边上的内点。5.2 单元位移模式的插值函数常应变三角形单元是最早提出、最简单的一种单元,对于从杆系结构转入连续体问题,理解有限单元法的原理、方法和步骤,从而更好地掌握有限单元法是十分必要的。这种单元虽然简单,但它也可用来拟合复杂边界体形。虽然边界为曲线

5、时,存在以三角形直边来代替而带来离散误差,但比用矩形单元进行离散,离散误差要小得多。5.2.1 单元结点位移和结点力图5-4为任一典型单元。单元局部结点编号记作1、2、3,逆时针进行标记,其对应的整体结点编号记作i 、j 、m 。由于平面问题局部坐标和整体坐标是一致的,因此没有坐标转换等问题,故也可只标记整体编号以便形成定位向量。如图所示,单元每个结点有两个位移分量,称为结点位移,记为:i i i u v =将三个结点位移按结点编号排在一起称为单元结点位移,记为:（ Te T T Ti j m =与结点位移相对应,每个结点上受有2个其它单元对它作用的力,称为结点力,记为ix i iy F F

6、F =将三个结点力按结点编号排在一起称为单元结点力,记为（T eTT T i jmF F F F=单元上作用的体积力记为bx by ebF F F =若单元的边界是物体边界,并且该边界有表面力作用的话,该表面力记为jx图5-4 典型单元sx sy e s F F F =体积力和表面力表达式中的矩阵元素均是沿坐标方向的分布荷载集度。5.2.2 插值函数的建立如图5-5所示是一个典型三节点三角形平面单元,三个节点为m j i 、,按逆时针方向排列。节点位移分量用一个向量表示为:=n m j j ii m j i ev u v u v u 只知道节点位移,并不能直接求出单元内的应变和应力。因此,必须

7、假定单元内任一点的位移为坐标的某种函数。有限元中,单元的位移模式一般采用多项式,即yx y x v y x y x u 654321,（,（+=+=将三个节点的坐标及位移代入上式得到六个多项式:m m j j j i i i m m m j j j i i i y x v y x v y x v y x u y x u y x u 654654654321321321+=+=+=+=+=+=,按照线性代数中的克莱姆法则,从上面式中解出系数i ,得到单元位移函数如下:XY图5-5 三节点三角形平面单元（m m m m j j j j i i i i m m m m j j j j i i i i

8、 v y c x b a v y c x b a v y c x b a Av u y c x b a u y c x b a u y c x b a A u +=+=21其中ij m j i m i j j i m m i j i m j m i i m j j m i m j i j m m j i x x c y y b y x y x a x x c y y b y x y x a x x c y y b y x y x a =, mj j i i yxy x y x A 11121=按照解析几何的理论,A 等于三角形ijm 的面积,为了保证面积A 为正,所以节点m j i 、的次序必

9、须是逆时针方向。为了简化位移函数的表达式,记Ac x b a N A y c x b a N A y c x b a N m m m m j j j j i i i i 222+=+=+=,代入上面式子,得到位移函数的简易表达式如下:j j m m i i j j m m i i v N v N v N v u N u N u N u +=+=, （5-2写成矩阵形式如下:e e m ji m m j j i i m jm j iN IN IN IN v u v u v u N N N N N N v u u =0000（5-3N 是坐标的函数,它们反映了单元的位移形态,称为单元的形函数或插值

10、函数。上式按虚线所示分块,其各子块分别记作22222200010;00001ji i j m j i m N N N N I N I N I N N N N = 式中I 为二阶单位矩阵。由式（5-2可知,当1i i u v =,而其它节点位移为零时,单元内任何一点处的位移为i iu N v N =就是说节点I 发生位移时,函数i N 就表示单元内部的位移分布形状,故称i N ,j N ,m N 为位移的形函数（或形状函数。形函数是定义于单元内部的坐标的连续函数。对于三角形三节点单元,其形函数是坐标的线性函数,因此位移函数式（5-2也是坐标的线性函数。线性位移函数的一个特点是:单元上的任意一条直

11、线在变形后,仍然是一条直线。所以相邻单元的公共边线,在变形后仍是一条直线,且只要相邻两个单元在公共节点处保持位移相等,则公共边线在变形后仍然保持密合。因此,对于三角形单元来说,选择线性位移函数是能够满足协调性要求的。5.2.3 插值函数的性质（1在单元节点上,插值函数的值为=j i ji y x N j j i 当当01,（; （i ,j ,m 即有（,1i i i N x y =,（,0i j j i m m N x y N x y = 。（2在单元中任意一点各插值函数之和应等于1,即i j m N N N +=（3.2.12（3插值函数是坐标x 、y 的线性函数,因此单元内及边界上的位移也

12、是线性变化的。5.3 单元应变和应力平面问题单元内具有三个应变分量xy y x 、,可用矩阵表示如下:m m j j i i m jm j i e v u v u v u N N N N N N x y y x N x y y x u L =00000em m j j i i m m j j i i m j im j i B v u v u v u x N yN xN yN y N y N y N x N x N x N =000000（5-3 其中=m mji m jj ib c b c b c c c c b b b A B B B B 00000021 （m j i b c c b AN

13、 L B i ii i ,0021,= （5-4其余两式脚码用j,m 代换。由于单元面积A 以及系数i i c b 、都是常量,所以矩阵B的元素都是常量,因此应变的元素也是常量,也就是说,在每个单元中,三个应变分量都是常量,这种单元称为常应变单元。在求出单元应变之后,根据前面对平面问题的研究不难求出单元应力如下:e e S B D D = （5-5其中D 称为弹性矩阵,S 称为应力矩阵,由弹性矩阵和应变矩阵相乘而得。 对于平面应力问题=2100010112E D （5-6 210000100002（11002ij m im i j m ii j j mm b b b ES s s s c c

14、c A c b c b c b =（5-72,（,2（11122i i i i i i i b c ES b c i j m A c b = （5-8平面应变问题 把1,21EE 代入上式（5-6,可的得到平面应变公式中的D : D =+1（2210001101121（1（1（E同理,可得到平面应变公式中的i S :+=i i ii i ii b c c b c b E S 1（2211（2211121（1（21（ 因为以上式子中i i c b E 、等都是常量,所以矩阵S的元素都是常量,也就是说,单元内的应力都是常量。同样,节点力eF 可由虚功方程导出:Te eem F F F F k =其

15、中:iiij im eji jj jm mi mjmm k k k k k k k k k k =+=s r s r sr s r s r s r s r s r rs b b c c cb bc b c c b c c b b A Et k 212121211（42（r,s=i,j,m 。 对于平面应变问题,E 做与上式相应的代换:1。 5.4 单元等效节点载荷列阵5.4.1平面问题节点力的一些约定:1、分别用U 和V 加节点号下角字表示某节点的水平与垂直节点力分量;2、节点力的方向,假设节点对单元的节点力为正,方向取为坐标轴正方向;单元对节点的节点力为负,方向取为坐标轴负方向。5.4.2关

16、于载荷移置载荷移置就是把弹性体所受的各种形式的载荷按等效原则移置到节点上成为节点载荷。如果弹性体所受的载荷全都是集中力,那么把所有集中力作用点都取为节点,就不存在移置问题,集中力就是节点载荷,把它们分解到x 轴和y 轴即可。但实际工程问题往往受有分布的面力和体力作用,在这种情况下就必须把面力和体力向节点移置。就是对集中力,若其作用点未被取为节点,也要向节点移置。单元载荷的移置是载荷移置的第一步也是最主要的一步,只要有了单元载荷的移置结果,总体节点载荷就容易计算了。要把载荷移置到节点上,首先要遵循静力等效原则,即原载荷与节点载荷在任意轴上的投影之和以及对任一轴的力矩之和都相等。对连续弹性体来讲,

17、静力等效原则实质上就是虚功等效原则,即单元的原载荷与移置后的节点载荷在任何虚位移上的虚功都相等。必须注意的是,这里的单元虚位移与所设选的单元位移函数要采用同一位移模式。我们前面所选的单元位移函数是线性模式,所以单元虚位移也都采取线性模式。设单元所受载荷有:体力f ,面力T ;移置后的节点载荷为eP ;在外载荷作用下节点产生的任意虚位移为e*,则单元内各点的任意虚位移为eN u *=。根据静力等效原则有:+=SV T edS T u dV f u P *即+=SVdS T N dV f N P *由于虚位移e*是任意的,所以有:+=S TV T my mx jy jx iy ix m j i e

18、dS T N dV f N p p p p p p p p p P （5-9 5.4.3 载荷移置举例1、单元某点处的集中载荷设在三角形单元的（x,y 点受到集中载荷y xP P P =,则:=y m x m y j x j y i x i my mx jy jx iy ix eP N P N P N P N P N P N p p p p p p P 2、均质等厚单元的自重设单元厚度为t ,材料的比重为,则单元自重的合力为S t G =,并作用在单元形心（即三角形的重心c 上,如图5-6所示,我们先求节点载荷iy p 。为了简化计算,我们假设单元发生这样的虚位移:节点i 只沿y 方向移动一个

19、单位,即01*=i i u v ,而其余两个节点都不移动,0*=m m j j v u v u 。由于采取的是线性位移模式,因此单元形心将产生1/3的虚位移,所以由静力等效原则可知:31=G p iy故得mxP3St G p iy = 按照上述分析方法同样可得。0,3=mx jx ix my jy p p p St p p 综上,写成矩阵形式为=1010103S t p p p p p p P my mx jy jx iy ix e （5-10 由此可见,像重力这种在单元内均匀分布的体力,把它平均地移置到单元的三个节点上即可。3、单元一边上受有均布侧压设单元在jm 边上受有均匀分布的垂直压力,

20、jm 边与x 轴的夹角为,如图5-7所示,求其节点载荷。我们先来求节点载荷分量mx p 。为了简化计算同样假设01*=m m v u ,其余两个节点都不移动,即0*=j j i i v u v u 。 由图可知,原载荷的两个分量为sin cos j m x jm mj y jm y y q q ql x x q q q l =（5-11其中jm l 为单元该边的长度。设s 是自节点m 到jm 边上任一点的距离,由于单元虚位移也采取线性模式,所以当节点m 发生水平单位虚位移时,s 点的虚位移是图5-7根据静力等效原则有+=jmy x xm tds s v q s u q p （1i m j l

21、jm jmm j l jmjm xmx tb q y y t q ds s l l y y qttds l s l q p jmjm（2（0= 按照上述类似的分析方法,可得其它节点力如下:综上,写成矩阵形式为=i i i i my mx jy jx iy ix ec b c b qt p p p p p p P 002 （5-12由此可见,对于面力的移置,在采取线性位移模式的情况下,只要把原载荷静力等效地移置到受力边的两个节点上即可,均布压力如此,对于非均布压力也是这样。 3、单元一边上受有x 方向的均布载荷设单元的ij 边上受有沿x 方向的集中力P ,其作用点离i 及j 点的距离分别为i l

22、 及j l ,则:（0jm jm l s u s l v s =（22（22（2200my m j i jx j m ijy m j iix iy q qp x x tc q q p y y tb q q p x x tc p p =0000l ll l P p p p p p p P i j my mx jy jx iy ix e其中,l 为ij 边长度。由此推导:设单元在ij 边上受有均匀分布的x 向压力,如图5-8所示,求其节点载荷。边界上的作用力写成:Q=0qqlt tds s u q p ij x ix 21（1=0=jm im jy iy p p p p由以上同样的分析方法可以得出等效节点载荷为:=0001012qtl p p p p p p P my mx jy jx iy ix e（5-134、单元一边上受有x