第六章 三维问题有限元分析,讲 授：陈得良 Tel：13327311008 ：416501065 Email：deliang_,1,四 教学基本内容,第五章 三维问题有限元分析 第一节 三维应力状态 第二节 4节点四面体单元 第三节 8节点六面体等参单元 第四节 20节点等参单元 第五节 ansys空间问题实例 第六节 空间轴对称问题有限元法 第七节 Ansys轴对称旋转问题实例,2,工程实际中的很多问题难于简化为平面问题，如受任意空间载荷作用的任意形状几何体，受对称于轴线载荷作用的回转体,这类问题经典弹性力学往往无能为力在FEM中，空间问题只要求0阶连续，因此构造单元方便,空间问题简介,3,空间问题的主要困难： （1）离散化不直观；（网格自动生成） （2）分割的单元数量多，未知量的数目剧增 （对某些问题简化） （轴对称问题）,空间分析的优点 精确,4,6.1 三维应力状态,工程结构一般都是空间的弹性体受力作用后，其内部各点将沿x、y、z坐标轴方向产生位移，是三维空间问题，其应力状态如图6-1所示图6-1 空间结构应力状态,各点沿x、y、z方向的位移以u、v、w表示，这些位移为各点坐标的函数，即：,u=u( x、y、z) v=v( x、y、z) w=w( x、y、z),5,由弹性力学知，应变与位移间的几何关系是,,(6-1),三维弹性体的应变分量，用矩阵表示为,(6-2),6,弹性体受力作用，内部任意一点的应力状态也是三维的，用列向量表示为,弹性范围内，应力与应变间的物理关系矩阵表达式为,对于各向同性弹性体，在三维应力状态下，弹性矩阵 的形式为,(6-3),(6-4),7,(1)空间问题常用单元：四面体单元、长方体单元、直边六面体单元、曲边六面体单元 、轴对称单元。
4结点四面体单元：是空间问题最简单的单元，也是常应变、常应力单元，可以类似平面问题三结点三角形单元进行分析 8结点长方体单元：可以类似平面四结点矩形单元进行分析 8结点直边六面体单元：可以类似平面四结点任意四边形等参元分析 20结点曲边六面体单元：等参单元,可以类似平面八结点曲边四边形等参元进行分析 轴对称单元：一平面单元绕一对称轴旋转形成的空间问题只需在rz平面划分网格，就像平面问题xy平面中的网格一样，这样这类空间问题可以得到简化 （环向位移等于零）,（2）结点位移3个分量 （3）基本方程比平面问题多3个平衡方程，6个几何方程，6个物理方程8,6.2 四节点四面体单元,图6-2表示任一简单四面体单元，其中四个结点编号设为 i、j、m、n (或1、2、3、4)单元变形时，各结点沿x、y、z方向上的位移，以列向量表示为,图6-2四面体单元,8,单元变形时，单元内各点也有沿x、y、z方向的位移u、v、w，一般应为坐标x、y、z的函数对于这种简单的四面体单元，其内部位移可假设为坐标的线性函数，为满足变形协调条件，取为,(6-5),式(6-5)含有12个待定系数a，可由单元的12项结点位移决定.将4个结点的坐标值代入式(6-5)的u式中。
i、j、m、n共4个结点，分别有,(6-6),1 单元形函数,9,其中,式中，V为四面体的体积，且有,(6-7),由式(6-6)求出 ，再代回式(6-5) 中，整理后得,10,为使四面体的体积V不为负值，在右手坐标系中，使右手旋转按着由i- j- m的转向转动时，且法向n方向前进用求位移u的同样方法，可求得,将位移的3个线性方程形成的线性方程组用矩阵表示为,(6-8),式中,(6-9),11,2 单元刚度矩阵,将式(6-8)代入几何方程式(6-2)，经过微分运算，可得单元内应变为,(6-10),式中,(6-11),简单四面体单元内，各点的应变都是一样的，这是一种常应变单元（是三维单元中精度最低的单元）这一点与平面问题的简单三角形单元相似，由于单元内位移都假定为线性变化的，因而由位移一阶导数组成的应变也为常量12,同样，用虚功原理建立结点力和结点位移间的关系式，从而得出简单四面体单元的刚度矩阵6-12),(6-13),按结点分块表示，此单元刚度矩阵可表示为,(6-14),13,(r=i、j、m、n， S=i、j、m、n ) (6-15),式中,,弹性体三维（空间）问题的原始平衡方程组，即,其中,其中任一子矩阵为,14,3 整体结构载荷列向量,整体结构的结点载荷列向量,(6-16),式中,单元上集中力等效结点载荷列向量； 单元上表面力等效结点载荷列向量； 单元上体积力等效结点载荷列向量； 单元结点载荷列向量。
等效结点力公式为,式中,15,6.3 8节点六面体等参单元,8 (x5,y5,z5),,7,如同二维等参单元一样，三维等参单元的有关公式的建立也是采用局部自然坐标（曲面坐标）可参考的母体单元则为一正六面体，图6-3则表示了任意六面体单元与母体单元、局部的三维自然坐标与整体的直角坐标系的几何关系,母体单元,任意六面体单元 （对面不全平行）,图6-3 8结点三维等参单元,16,用形函数表示的位移插值形式的位移模式，可以直接利用拉格朗日插值公式，得到单元位移函数为,根据等参单元的定义自然坐标与整体直角坐标之间的关系可以写为,其中形函数为,例,其中 为节点坐标值（角点）,17,由节点位移求单元应变时，他要求形函数在整体坐标下的导数，但形函数是建立在局部坐标下的，这就需要将局部坐标中的表达式转换到整体坐标系中，如同平面等参单元一样，需要通过雅克比矩阵来实现，由偏导法则,同理可得,写成矩阵,18,求单元刚度矩阵，尚需对积分的单元体积进行积分变换,反之,有了上式很容易得到单元的应变应力矩阵,19,6.4 20结点等参元,为适应三维结构的曲面边界，可以采用曲面六面体单元正方体基本单元内任一点与实际曲面单元内的点一一对应，结点也一一对应。
这里，实际单元边界线中间的结点9、10、、20，都“映射”成为正方体的棱边中点8结点单元是线性单元，其位移模式是三维线性的，在8结点单元的基础上每边增加一个中点作为节点就构成了20节点单元此时六面体单元每条边上有3个节点，他们既可以是直线的，也可以是曲线的，因此每个面也可以是平面的，也可以是曲面的1 形状函数,20,(a)直角坐标系与实际单元 (b) 自然坐标系与基本单元,图6-3 20结点三维等参单元,位移函数和几何坐标的变换式应取为相同的参数，其坐标变换关系可表示为,(6-17),则单元的位移函数可写成,21,(6-18),在自然坐标系（局部坐标系）中，各结点的形状函数可写成如下形式， 对于8个顶角结点（ i1，2，，8）,式中 xi、yi、zi结点i的坐标； ui、vi、wi结点i沿x、y、z方向的位移； Ni对应于i结点的形状函数22,对于 的边上点（i17，18，19，20）,(6-19),对于 的边上点（i9，11，13，15）,对于 的边上点（i10，12，14，16）,23,2 单元刚度矩阵,三维变形状态下，一点的应变与位移的几何关系为,(6-20),24,(6-22),为便于以下计算，弹性矩阵D可分块写为,(6-23),令 ,,则 ，,为6060的方阵，可按结点写为子块形式,25,式中第i行j列的子矩阵为,(6-24),将(6-20)、(6-22)分块式代入(6-23)，其被积函数可写为,(6-25),式中,与式(5-9)相似，按坐标变换式(6-17)，应有,26,同样可有,(6-26),三维六面体的雅可比矩阵为,(6-27),同理可采用三维高斯求积公式计算单元刚度矩阵。
即,27,式中，L，M，N为沿 、 、 方向的积分点数目，而积分点坐标 及权重 可由高斯积分表查得对于20节点的三维单元，通常可取积分点数目m=3(即3 3 3），查表可得对应的积分点坐标和权重为,28,参照方法，可以很方便的得到另外两个方向的积分点坐标值和权重29,6.5 ANSYS空间问题示例,1 问题描述,如图6-4所示，一个圆柱实体柱高0.2m，圆柱横截面直径为0.1m约束方式：底面全约束承受载荷： A点承受Z方向集中载荷Fz=5000N和Y方向集中载荷 Fy=-5000N；B点承受X方向集中载荷Fx=5000N；C点承受Z方向集中载荷Fz=-5000N；D点承受X方向集中载荷 Fx=-5000N弹性模量为EX=210GP，泊松比=0.3,2 ANSYS求解操作过程,30,（1）选择单元类型 运行PreprocessorElement TypeAdd/Edit/Delete，弹出Element Types对话框，如图6-5所示然后单击 Add，弹出Library of Element Types窗口，如图6-6所 示，选择SOLID45单元，单击OK图6-5 单元类型对话框,图6-6 单元类型库对话框,31,（2）设置材料属性 运行PreprocessorMaterial PropsMaterial Models，弹出如图6-7所示对话框。
双击Isotropic，弹出Linear Isotropic Properties for Material Number1对话框，如图6-8所示，在EX选项栏中设置数值2.1e11，在PRXY选项 栏中设置数值0.3设置完毕单击OK按钮图6-7 选择材料属性对话框,图6-8 设置材料属性对话框,32,（3）建立模型 运行PreprocessorModelingCreateAreasRectangleBy 2 Corners，弹出如图6-9所示对话框，在WP X选项栏中填写0，在WP Y选项栏中填写0，在Width选项栏中填写0.05，在Height选项栏中填写0.2，点击OK生成如图6-10所示图形图6-9 两点建立矩形对话框,图6-10 生成的长方形面,33,将长方形旋转成柱体，运行 PreprocessorModeling Operate ExtrudeAreasAbout Axis，弹出如图6-11所示拾取框选择图7中长方形后弹出单击OK，再选择长方形左上角和左下角结点后，单击OK.弹出如图6-12所示对话框在ARC选项栏中填入旋转角度360度，设置完毕单击OK按钮，生成如图6-13所示 圆柱体。
图6-11 拾取对称轴对话框,图6-12 设置绕轴旋转参数对话框,图6-13 圆柱模型,34,运行MeshingSize CntrlsManualSizeGlobalSize弹出如图6-14所示对话框，设置SIZE选项栏中的数据为0.01运行MeshingMeshVolumesFree自由划分网格后得到如图6-15所示图形图6-14 设置网格尺寸对话框,图6-15 圆柱有限元模型,（5）施加约束 运行SolutionDefine LoadsApplyDisplacementOn Areas，拾取圆柱的底面，施加全约束35,（6）施加载荷 显示图形的关键点，运行PlotCtrlsNumbering弹出如图6-16所示对话框，激活KP Numbers后面的选框，使它变成on形式选择菜单SolutionDefine LoadsApply,Structure Force/Moment On Keypoints，载荷分别如下：8点承受Z方向集中载荷Fz=5000N和Y方向集中载荷Fy=-5000N；10点承受X方向集中载荷Fx=5000N；3点承受Z方向集中载荷Fz=-5000N；6点承受X方向集中载荷Fx=-5000N。
施加载荷，图形如图6-17 所示图6-16 编号显示设置对话框,36,图6-17圆柱实体示意图,（7）求解 选择Solution SolveCurrent LS，开始计算，计算结束会弹出计算完毕对话框，单击Close关闭对话框计算完毕8）后处理 运行 General PostprocPlot ResultsContour Plot Nodal Solu 弹出如图6-18所示对话框，运行DOF SolutionDisplacement vector sum和Stressvon Mises stress，分别显示圆柱体的位移和应力云图37,图6-18 云图显示对话框,结果显示如图6-19和图6-20所示38,图6-19 位移云图 图6-20 应力云图,39,6.6空间轴对称问题的有限元法,对空间轴对称问题，常采用圆柱坐标系r表示径向坐标，z表示轴向坐标，任一对称面为rz面在有限元分析时，可采用轴对称的环形单元进行环形单元 可以是任何平面单元某一平面图形绕平面上某一轴旋转形成的回转体称为轴对称物体，此平面称为子午面在动力机械，特别是叶轮机械中，有很多零件都具有轴对称特性，比如轮盘、旋转轴、承力环等。
对于直齿圆柱齿轮，由于齿的存在，严格地说它并非轴对称物体如果忽略齿的部分(将齿用外载荷表示)，则所得到的齿根以内的旋转体部分为轴对称物体轴对称物体的变形及应力分布不一定是轴对称的，只有当其约束和载荷都对称于旋转轴时，轴对称物体的变形和应力分布才是轴对称的 轴对称物体+轴对称约束+轴对称载荷=轴对称系统 对轴对称系统的应力分析=轴对称物体,40,1）几何形状关于轴线对称； 2）作用于其上的载荷关于轴线对称 3）约束条件关于轴线对称因过z轴的任一子午面都是对称面，其 上任一点p只在该平面上发生位移，即弹性体内任一点的位移、应力与应变只与坐标r、z有关，与 无关从而，轴对称问题可转化为二维问题，但因与平面问题有区别，常称为二维半问题柱坐标系,41,注意：应变 虽然与 无关，但是周向应变 ，周向应力 ，由径向位移 引起，因为径向位移会导致周长的改变1、基本方程,位移分量,应力分量,应变分量,42,虚功方程,,应变分量,轴对称问题的弹性矩阵：,43,2、轴对称问题的离散化,对于轴对称问题，利用其轴对称特性，在对其进行网格划分时可知取任意通过Z轴的截面进行，类似平面问题的网格形式。
本节以三角形单元为例1、位移模式,轴对称问题的环向位移恒等于零，径向r位移与轴向z位移不等于零对于图示情形，依照平面问题的三角形单元分析，取位移模式为,,代入结点位移后，可解出a1-a6，再代入上式，得,,xr, yz,44,其中形函数：,,,；,,单元中位移,根据弹性力学理论，空间轴对称问题的几何方程为,,2、单元中应变,45,将u,w表达式代入上式，整理后,,46,式中,,,其中,B矩阵中含有变量r，z，因此它不是常数矩阵，即轴对称问题的三角形环形单元不是常应变单元47,3、单元中应力,根据弹性力学理论，空间轴对称问题的应力-应变关系为,,,弹性矩阵:,48,单元中任意一点的应力：,,4、单元刚度矩阵,,由于被积函数与无关，故在三角形截面的环单元的积分可简化为在三角形截面上的积分故有：,,49,,单元刚度矩阵的积分参照图示分区，按下式采用数值积分的方法进行,50,当单元较小时，可把各个单元中的r，z 近似看作常数，并且分别等于各单元形心的坐标，即,,,这样，就可把各个单元近似地当做常应变单元,,,,51,单元刚度矩阵k的分块形式,,,,其中的近似子矩阵为,,52,5、等效结点荷载,类似平面问题。
对于作用于三角形环单元上的体积力、表面力的等效结点力为：,,,体力,,,,53,面力：,（1）均布表面力 设单元ij边上作用均布表面力，其集度为,,,,,l,当 ri=rj 时，静力等效原则,54,(2)三角形分布表面力 沿单元ij边作用了三角形分布的表面力，表面力在i点集度为,,,,当 ri=rj 时，静力等效原则2/3集中在i点，1/3集中在j点55,56,圆筒直径0.4m，高度0.6m，壁厚0.005m；材料Q235，弹性模量E=2.1e11Pa，泊松比=0.3；约束：圆筒的下部在轴线方向固定，其它方向自由；载荷：顶部环线上承受轴向线压力P-200000N/m图6-4 圆筒示意图,图6-5单元类型对话框,1 问题描述,6.7 ANSYS轴对称旋转单元计算示例,57,（1）选择单元类型 运行PreprocessorElement TypeAdd/Edit/Delete，弹出Element Types对话框单击Add，弹出Library of Element Types对话框，如图7-6所示，选择SHELL51单元2 ANSYS求解操作过程,,图7-6 单元类型库对话框,图7-7 选择材料属性对话框,58,（2）设置材料属性 运行PreprocessorMaterial PropsMaterial Models，弹出Define Material Model Behavior对话框，如图7-7所示。
双击Isotropic选项，弹出Linear Isotropic Properties for Material Number1对话框，如图7-8所示图7-8 设置材料属性对话框,59,（3）定义单元实常数 选择Main MenuPreprocessor Real Constants Add/Edit/Delete，弹出如图7-9所示对话框，单击Add按钮弹出Element Type for Real Constants对话框，如图7-10所示，选择Type 1 SHELL51，单击OK，弹出Real Constant Set Number 1,for SHELL51对话框，如图7-11所示，在TK(I)项输入0.005，单击OK图7-9 实常数对话框图 7-10选择要设置实常数的 单元类型,60,图7-11设置SHELL51实常数对话框,（4）建立模型 首先生成关键点，运行主菜单Preprocessor ModelingCreateKeypointsIn Active CS，弹出如图 7-12所示对话框 创建关键点1（0.2,0,0），2（0.2,0.6,0）。
生成圆筒母线：运行Main MenuPreprocessorModeling CreateLinesLinesStraight Line，弹出拾取关键点对 话框，拾取关键点1、2，单击OK61,,图7-12创建关键点对话框,（5）设置单元属性 运行Main MenuPreprocessorMeshingMesh Tool，弹出Mesh Tool对话框，在Element Attributes下拉列表中选择Lines，然后单击其后的Set按钮弹出拾取线对话框，单击Pick All，弹出分配线单元属性对话框，将 MAT,TEAL,TYPE依次设置为1，1，1，单击OK6）划分网格,62,单击Mesh Tool中Lines后的Set按钮，弹出拾取线对话框，单击Pick All弹出控制线单元尺寸对话框，将NDIV设置为10，单击OK在Mesh Tool对话框中的Mesh下拉列表中选择Lines单击Mesh，弹出拾取线对话框，单击Pick All，划分网格完毕 运行Plot CtrlsStyleSize and Shape，弹出如图 7-13所示对话框在Display of element选项后面选择 on，单击OK按钮。
显示如图7-14所示图形图7-13 尺寸和形状显示设置对话框 图7-14 有限元模型,63,（7）施加约束 运行Main MenuPreprocessorLoadsDefine LoadsApplyStructuralDisplacementOn Keypoints，弹出拾取关键点对话框，拾取关键点1，约束其Y方向上的自由度 （8）施加载荷 运行Main MenuPreprocessorLoadsDefine LoadsApplyStructuralForce/MomentOn Keypoints在关键点2上施加竖直向下的集中载荷F-200000N，如图7-15所示图7-15 施加载荷后模型,（9）求解 运行Main Menu SolutionSolveCurrent LS, 开始计算，计算结束会弹出计算完毕对话框，单击Close，关闭对话框，计算完毕64,（10）扩展成圆筒 运行Main MenuGeneral PostprocRead Results Last Set运行Utility MenuPlotCtrlsStyleSymmetry Expansion2D Axi-Symmetric弹出轴对称扩展设置对话框，如图7-16所示，选择Full expansion，单击OK。
显 示圆筒图形如图7-17所示图7-16 2D轴对称扩展设置对话框,图7-17圆筒图形,65,（11）结果显示 运行 General PostprocPlot ResultsContour Plot Nodal Solu弹出如图7-18所示对话框，运行DOF SolutionDisplacement vector sum和Stressvon Mises stress，分别显示圆筒位移和应力云图结果显示如图7-19和图7-20所示图7-18 云图显示对话框,,,66,The End,。