弹性力学平面问题

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,四、弹性力学平面问题,的有限元分析及程序,1,四、弹性力学平面问题,的有限元分析及程序,引 言,常应变三角形单元,矩形双线性单元,平面问题程序(一),平面等参数单元,平面问题程序(二),Wilson 非协调元,2,杆系问题以结点作为分割单元的“结点”是很自然的，但对于平面问题，待分析物体是连续的，并不存在实际结点。要将物体“拆”成单元，必须用一些假想的线或面作人为地分割。实际计算时，可将连续体分成多种形状单元，为讨论简单，现暂时规定只用一种单元来分割。,平面问题有限单元法可用的单元很多，作为初学，先介绍两种最简单的单元：三角形和矩形。然后再介绍高级些的单元“等参数单元”。,将物体进行分割时，必须保证相邻单元具有公共边界。假定相邻单元仅在一些点（顶点或顶点加边中点）相连接。这些点即为“结点”。,4.1 引 言,3,4.2 常应变三角形单元,4.2.1 面积坐标,三角形单元中任一点,P,可用直角坐标 (,x,y,) 表示。,P,2,1,3,y,x,如图所示连,P,1、,P,2、,P,3，则可得三个小三角形。它们和大三角形,123的面积比，记作,L,i,（= ,P,jk/ 123），称为面积坐标。,三个面积坐标显然,L,1,+,L,2,+,L,3,= 1，只有两个是独立的。三角形中任一点,P,的,位置可用面积坐标,L,1,、,L,2,确定。,当,P,点在1时,L,2,=,L,3,= 0，,L,1,= 1。余类推。可见面积坐标具有“形函数”的性质。,4,4.2 常应变三角形单元,4.2.2 位移模式,由于面积坐标有形函数性质，因此根据试凑法可得,P,2,1,3,y,x,形函数=,N,i,=,L,i,= 面积坐标,1) 面积坐标和直角坐标关系,如果结点,i,位移为,u,i,、,v,i,，则单元位移模式（位移场）为,u,=,N,i,u,i,；,v,=,N,i,v,i,5,4.2 常应变三角形单元,2) 矩阵表达,P,2,1,3,y,x,6,4.2 常应变三角形单元,4.2.3 单元列式,1) 微分算子矩阵,2) 应变、应力矩阵,平面应力问题,式中,平面应变时,7,4.2 常应变三角形单元,由此可见，单元应变、应力都是常量。,当所分析的问题具有初应变时，单元的弹性应变为,e,=,-,0,，应力为,=,D,e,。,3) 单元应变能,将上述应变、应力代入,4) 单元外力势能,第一项体积力、第二项结点力、第三项表面力的外力势。,代入位移后，经整理可得,8,4.2 常应变三角形单元,5) 令总势能一阶变分等于零，推导单元刚度方程,当有初应变时推导结果如何？,6) 单元刚度矩阵、等效荷载矩阵,当有初应变时结果如何？,具体显式表达式见教材P。47 式（3,2-39）,9,4.2 常应变三角形单元,7) 关于等效结点荷载,等效结点荷载可用公式积分计算，但由于形函数的图形是一平面（边界处为一直线），因此可证明也可按杠杆原理通过静力等效来求。,如 P.48 图3-4所示。,4.2.4 解答的收敛性准则,1) 位移模式（也称位移函数）必须包含刚体位移。,2) 位移模式必须包含常应变位移。,3) 位移模式必须保证单元间位移协调。,1)、2) 对平面问题也即要求具有常数项和坐标一次项，这称作“完备性准则”。,3) 称作“协调性准则”。既完备又协调的单元一定是收敛的。但不等于说非协调单元一定不收敛。,10,4.3 矩形双线性单元,三角形单元划分灵活，能较好拟合边界复杂（如曲线边界）物体计算。但是，单元应变、应力是常量，对一般问题精度较低，要提高精度就的增加结点、增加未知量，为此讨论其他单元。其一为本节单元。,1) 自然坐标,2,a,2,b,2,2,图示矩形单元，设,=x/a,，,=y/b,，则转换成正则单元。,2) 形函数,2,2,2,3,4,1,由形函数的性质“本点1，它点零”，利用试凑法可设：,N,1,=a,（1-,)(1-,）它满足“它点零条件”。再令本点为1，可得,a,=1/4，代回可的形函数,N,1,。,同理可得：,N,i,=,1/4(1+,0,)(1+,0,) (,i,=1,2,3,4)。,式中,0,= ,i,；,0,=,i,。,请大家验证,N,i,是否满足形函数性质。,11,4.3 矩形双线性单元,3) 位移模式,u=,N,i,u,i,； v=,N,i,v,i,。,或以矩阵表示为,2,2,2,3,4,1,可以用势能原理，也可以用虚位移原理。一经建立单元位移模式后，剩下的工作和杆系、三角形单元类似，因此这里从略。,d,=,u,，,v,T,单元结点,位移矩阵,4) 关于单元列式,12,5) 关于单元特性结果请看 P. 53 式(3,6-1315)。,6) 关于计算结果的整理,里兹法已经知道：位移结果比应力、内力结果精度高。位移达到满意结果，有几何方程求应变，再由物理方程求应力，结果精度较差。上述三角形单元常应力，矩形单元应力线性变化，许多工程问题的应力是复杂的。为更好标征性，需要对计算结果进行整理。常用处理方法有两种。,4.3 矩形双线性单元,6-1) 绕结点平均法,以交于同一结点各单元此结点处某应力分量的代数平均值，作为此结点该实际应力的近似值。,对于边界处的结点，由内结点结果的外插得到。,6-2) 两单元平均法,三角形单元时,，以两相邻单元应力平均值作为边中点的应力近似值。,矩形单元时,，以两相邻单元公共边两端结点四个应力的平均值作为边中点的应力近似值。对于边界处的结点，同样由内结点结果的外插得到。,13,1) 程序功能,本程序可用三角形或矩形单元计算平面应力问题。当计算平面应变问题时需要自行转换弹性常数。,本程序为了减少计算数据的准备，对规则问题具有做、单元结点编码等自动生成功能。,本程序荷载生成功能较弱，请自行修改。,4.4 平面问题程序(一),本程序可以用来计算如墙梁、剪力墙（可以带孔洞）等结构。,2）程序数据文件说明,2-1）基本数据,结点位移数，单元结点数，结点总数，最大半带宽，总约束位移码数，单元总数，点的坐标数，规则标志，问题标志，弹性常数及厚度。,14,2-2) 结点坐标,如果不规则：,按结点号顺序读入全部结点的坐标值。,如果规则无孔：孔标志,X,方向单元数，,Y,方向单元数，,X,方向单元长度，,Y,方向单元长度。,如果规则有孔：,控制结点数，生成结点类数。,结点号，,X，Y,4.4 平面问题程序(一),起点号，终点号，生成的点数，相邻点号差值，“相邻两点间距”。,2-3）读入结点荷载值,有荷载的结点数,结点号，,X,方向荷载值,，Y,方向荷载值。,2-4）单元的整体结点码,如果规则无孔：,不需要输入,如果规则有孔：,待修改,15,如果不规则：,按单元类型读入单元整体位移码。,2-5) 读入全部零位移约束的位移码。,如果问题类型不等于零,2-6）读入第二种材料的弹性常数，厚度，首单元号，终单元号，循环步长。,如果问题类型等于零,没有第六组,4.4 平面问题程序(一),2-7）一算例数据,2,4,441,46,22,400,2,0,0,6.93e4,0.3,2.0e-2,1,20,20,0.05,0.1,21,421,0.0,0.05,422,0.0,0.1,423,0.0,0.1,424,0.0,0.1,425,0.0,0.1,426,0.0,0.1,427,0.0,0.1,428,0.0,0.1,16,429,0.0,0.1,430,0.0,0.1,431,0.0,0.1,432,0.0,0.1,结,433,0.0,0.1,点,434,0.0,0.1,荷,435,0.0,0.1,载,436,0.0,0.1,数,437,0.0,0.1,据,438,0.0,0.1,439,0.0,0.1,440,0.0,0.1,441,0.0,0.05,运行程序,查看计算结果,1,2,4,6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38,40,42,17,4.5 平面等参数单元,三角形和矩形单元是最简单的单元形式，前已提及，由于位移模式是线性和双线性的，精度较低。对于曲线边界问题，还有以直边代替曲边的离散误差。为此，介绍本节等参数单元族。,首先以四结点等参元为例进行介绍，然后再介绍其他等参元。,2,2,2,3,4,1,x,y,1,2,3,4,母单元,4.5.1 单元描述,为克服矩形单元不能拟合曲线边界，用图示任意四边形单元，但在直角坐标下要描绘属于单元的点比较困难。,自然坐标下的图示单元（母单元）形状都是规则图形，考虑到母单元形函数性质，并设,i,点的坐标为：(,x,i,y,i,)，则由,x=,N,i,x,i,；,y=,N,i,y,i,可以描绘子单元形状。,为什麽？,子单元（等参元）,18,4.5 平面等参数单元,x,y,1,2,3,4,母单元,利用上述转换公式可看出：,1) 母单元正交坐标线映射后成为图示子单元斜角坐标线。仍为直线。,2) 子单元有两套坐标系：整体,x,y,坐标，和局部,坐标。,3),坐标又是母单元正交坐标。,4) 根据子单元结点坐标情况，规则母单元可映射出任意四边形单元。,5) 相邻单元映射后仍然连续。,子单元（等参元）,2,2,2,3,4,1,4.5.2 单元位移模式,设,i,点的位移为：(,u,i,v,i,)，则单元位移场为,u=,N,i,u,i,；,v=,N,i,v,i,。,和矩形单元一样，可用矩阵表示为：,19,4.5 平面等参数单元,x,y,1,2,3,4,母单元,4.5.3 坐标系间的转换关系,2,2,2,3,4,1,由复合函数求导数的规则可得,用矩阵表示则为：,引入记号,称作雅可比矩阵,其逆矩阵为：,4.5.4 其他等参元母单元形函数,1) 八结点等参元（矩形族）,x,y,1,2,3,4,5,6,7,8,子单元（等参元）,2,2,2,3,4,1,6,7,8,5,用试凑法,N,1,可设为：,N,1,=,a,(1-,)(1-,)(1+,+,)；它能自动满足它点为零。本点为1得：,a=,-1/4,同理可得角结点,N,i,通式为：,N,i,=,-1/4(1+,0,)(1+,0,)(1-,0,-,0,) (,i,=1,2,3,4).式中,0,= ,I,;,0,=,i,.,用试凑法,N,5,可设为：,N,5,=,a,(1-,2,)(1-,) ；它能自动满足它点为零。本点为1得：,a=,1/2,同理可得边中点,N,i,通式为：,N,5,7,=1/2(1-,2,)(1+,0,),N,6,8,=1/2(1-,2,)(1+,0,),映射后,子单元可以,是曲边单元,20,4.5 平面等参数单元,2) 六结点(三角形族) 单元如图所示,3,1,2,1,1,4,5,6,L,1,=,L,2,=,母单元,子单元,3,1,2,4,5,6,x,y,设,N,1,=,a,(,-1/2)；它满足它点为零条件，为使本点为1，可得,a,=2。同理可得,N,i,=,L,i,(2,L,i,-1) (,i=,1,2,3),设,N,4,=,aL,1,L,3,；它满足它点为零条件，为使本点为1，可得,a,=4。由此可得,N,4,=4,L,1,L,3,N,5,=4,L,2,L,1,N,6,=4,L,3,L,2,映射后,子单元可以,是曲边单元,21,4.5 平面等参数单元,3) 常用的两族单元,矩形：四、八、十二结点等参元。,三角形：六、十结点等参元。,问题：试建立十、十二结点母单元形函数。,各种单元分别包含坐标几次完整多项式？,能否给出建立两族形函数的一般方法？,4.5.5 单元描述和位移模式,1) 单元描述,2) 位移模式,结点数,22,4.5 平面等参数单元,4.5.6 等参元单元分析,记微分符号:,记,1) 应变矩阵,23,4.5 平面等参数单元,2) 应力矩阵,3) 单元应变能,4) 单元外力势能,5) dA、ds 的计算,为用虚位移原理推导单元刚度方程，必须解决dA、ds 的计算。,母单元规则微元体d,d,映射后变成图示（曲边）四边形。,x,y,1,2,3,4,5,6,7,8,24,4.5 平面等参数单元,如图示，此微面积为,x,y,1,2,3,4,5,6,7,8,坐标的积分上下限均为-1，1。,沿边线的积分（,=,1为例）,一般情况见 P.72 式(3,6-23)。,25,4.5 平面等参数单元,有了上述结果，经虚位移原理或势能原理即可推得式(3,6-24)(3,6-26)单元刚度和等效荷载结果。,4.5.7 数值积分,三角形和矩形单元可以写出刚度显式表达式，但对于等参元，由于两套坐标的转换，导致刚度、荷载的被积表达式十分复杂，一般不可能积出显式结果。只能用数值积分由程序来得到。,目前常用的是高斯积分（矩形族）和哈默尔积分（三角形族）。它们的积分点位置、加权系数等见表3-1、3-2（P.7476)。,其积分公式见式（3,6-40）、（3,6-41）。,26,4.5 平面等参数单元,4.5.8 作等参元分析时应注意的问题,等参元分析中要用,det,J,-1，,可见雅可比行列式等于零将导致刚度矩阵等无法积分，使分析失效。因此要避免以下可能使,det,J,=0的情况：,1) 子单元边界不能过于扭曲。,2) 矩形子单元不能退化成三角形。,3) 子单元角顶处单元边线切线角角不能等于180,0,。,上述情况如 P.81 图3-39 示意。,此外,子单元边界上结点应尽可能是或接近等分点，避免产生奇异单元。,可能情况下应采用直边子单元，这样可使雅可比矩阵简单，提高计算效率。,27,4.5 平面等参数单元,4.5.9 离散化时应注意的问题,除对等参元应注意上述问题外，任何有限元分析都还应注意以下几点：,1) 相互邻接的单元大小应尽可能均匀。,2) 单元最大尺寸与最小尺寸之比应尽可能接近一，最多不应大于二。,3) 应合理编码，使单元结点间的整体编号差值最小。,4) 应尽可能使各界点的单元数目相同,如 P.82 图3,6-42左图示意.,28,4.6 平面问题程序(二),应用本程序时,数据文件的准备见 P.217。,运行程序,查看计算结果,29,4.7 Wilson非协调单元,至今为止所介绍的单元都是能保证收敛的协调单元。单元计算结果的精度，取决于位移模式中坐标完全多项式的次数。为改进精度Wilson提出非协调的单元，简单介绍如下。,在矩形双线性单元基础上，增加一位移修正项：,显然在四个结点处修正项等于零，因此它只影响单元内部位移，可见,a,e,是单元内部位移参数。所求得的结点位移将仍直接是实际问题的近似值。,30,4.7 Wilson非协调单元,根据这一位移模式，可得应变矩阵为,经虚位移原理做单元列式，可得,双线性,修正项,相关项,相关项,结点力,双线性等效,修正项等效,各项的具体计算公式见下页。,31,4.7 Wilson非协调单元,双线性,修正项等效,相关项,修正项,双线性等效,32,4.7 Wilson非协调单元,由于,a,e,是单元内部位移参数，与单元集合体（整体）的其他单元无关，因此在集成之前可做如下消去,单元刚度,单元等效荷载,33,4.7 Wilson非协调单元,由此可见做了上述消去（也称静力凝聚）后，单元刚度方程形式上和双线性单元完全一样。,问题：试增加 PSTE 程序的单元库，使其能用非协调元进行分析。,能否在四边形等参元基础上来作非协调元？,协调元是肯定收敛的，非协调元收敛性不能保证。下一章将介绍如何检查非协调元的收敛性。这里只给出结论：,1) 当单元网格为矩形或平行四边形时，分析将试收敛的。,2) 任意四边形时，以,=,=0点的,J,作为单元的各点的,J,，这样计算的修正单元刚度也能使分析收敛。,34,