有限元第5章-等参数单元

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,5,等参数单元,5-1,等参数单元的引入,三角形单元内的应力为常量，不同单元的应力互不相同，提高精度的方法：,（,1,）减小单元尺寸；,（,2,）提高单元插值函数的阶次。,为了适应不规则边界，要求用曲边单元。,基于以上原因，引入等参数单元。,5-2,四节点四边形等参数单元,四节点四边形单元的位移插值函数可以写成（以,x,方向的位移插值函数为例）,对于边界来讲，将,带入上式，经简化可得,上式中有三个待定系数，由所在单元的节点场变量值确定，但是不能由这个单元的这条边界的两个节点的场变量值唯一确定，因此相邻两单元在同一边界上的位移表达式并不一致，使相容性条件不能得到满足。,这种情况该怎样处理？,我们知道，矩形单元满足相容性条件。,以图示为例。它有四个节点，各条边与总体坐标轴平行。单元内任意一点的位移插值函数可以包含四个待定系数,在矩形单元的任意一条边上，把该边的方程,或,带入上式，总可以得到,或,上式只含有两个未知参数，由边界上的两个节点的位移值唯一确定。,可见矩形单元的特点：,（,1,） 矩形单元满足相容性条件。,（,2,）含有一次项和常数项，故也满足收敛性条件。,（,3,）单元插值函数含有交叉项,xy,，,比三节点三角形单元的阶次要高。,如果通过坐标变换，将任意四边形单元变换成矩形单元，只要在坐标变换中，任意四边形单元与矩形单元之间的点是一一对应的（称为坐标变换的几何相容性），而变换后的位移插值函数又满足解的收敛性条件，这两条合在一起，就能保证任意四边形在原坐标系中满足解的收敛性条件。,即：使四节点四边形单元满足解的收敛性的途径是,（,1,）将四边形通过坐标变换，转化为矩形单元；（几何相容）,（,2,）以四边形节点位移值作为矩形单元的节点位移值。（收敛性要求）,以上两条结合，即可保证四节点四边形单元的几何相容性和有限元解的连续性。,在建立四边形和矩形单元的坐标变换关系时应注意：四边形单元定义在总体坐标系中，而矩形单元定义在局部坐标系中。坐标系的变换是一个四边形单元到一个矩形单元的变换。矩形单元的局部坐标系，仅仅适用于每个要变换的单元。,为此，首先讨论局部坐标系下的位移插值函数、形状函数和收敛性条件，然后再讨论具体的坐标变换。,根据前述，矩形单元四个节点的位移值，就是原四节点四边形单元的节点处的位移值。因此，局部坐标系下的矩形单元内任意一点的位移可以表示为,或者,将,和节点坐标 带入位移插值函数表达式，可得,其中,解上面的方程,从而,其中,或,其中,写成统一的形式,形函数的性质：,（,1,） 保证位移在节点连续。又因为是双线性单元，故也保证在边界连续。,（,2,） 保证单元包含刚体位移。,这两条性质，保证了解的收敛性。,下面讨论坐标变换。,可以证明：视整体坐标系下的四节点四边形单元的节点坐标值为“位移值”，采用与矩形单元内任意一点的插值函数完全相同的插值方式，就可以满足坐标变换的相容性（几何相容性），即,现证明如下：,从四边形到矩形的坐标变换是点点对应，并能保证相邻单元的几何相容（前面的位移插值可以看成是位移相容）。所谓几何相容，即是指总体坐标系下的两四边形单元在转换到局部坐标系下的矩形单元后,：（,1,）相邻单元的公共节点位置重合；（,2,）相邻单元的公共边界不开裂，不重叠，反之亦然。,关于（,1,）因为 ，所以相邻单元的公共节点位置重合；,关于（,2,）：局部坐标系下的矩形单元边界上的 或 保持常数，转换到总体坐标系下后，,边界为线性函数，该线性函数可以由边界上的两个节点坐标完全确定。因此，保证了相邻单元的公共边界既不开裂，也不重叠。,对于矩形单元中的点也可同样证明（提示：用通过矩形单元中任意点的水平或垂直的直线在总体坐标系和局部坐标系中的对应关系来证明）。,我们看到：矩形单元的插值函数对于场变量和坐标变换完全一样，故称之为的等参数单元。如果两者不一样，就称为超参元或亚参元。在此不予介绍。,5-3,等参数单元平面问题的有限元格式,前述有限元求解的七个步骤中,第,13,步：形成插值函数；,第,46,步：求出单元刚度矩阵，并集成求解；,第,7,步：用已知节点位移计算应力。,对于等参元，已经得到四边形四节点的等参数单元的形状函数。下面主要讨论单元刚度矩阵的形成，即上述中的,46,步。,一 等参数单元刚度矩阵,第,4,步：单元应变,单元位移,节点位移的关系,由平面问题几何方程和位移插值函数，有,第,5,步：单元应力,应变,节点位移的关系,由平面问题的物理方程，有,第,6,步：节点力,节点位移间的关系,由虚功原理，可得节点力于节点位移间的关系式,对于平面问题有,其中,积分区域为四节点四边形单元四条边所围成的区域。,需要积分。遇到的问题：,（,1,） 不是常量，不能提到积分号外面来；,（,2,） 是基于局部坐标而建立的，而 涉及,。,（,3,）上述积分区域是总体坐标系下的，不是局部坐标系下的。,为了能够计算上式，可以,（,1,）从局部坐标系向总体坐标系转换；,也可以,（,2,）从总体坐标系向局部坐标系转换。,统一之后，再积分。,比较方便的做法：把 的积分从总体坐标系向局部坐标系转换。,下面讨论坐标变换。,二 等参数坐标变换,现在的关键是（,1,）求解 在局部坐标系中的表达式；（,2,）将 转换成 。,已知：,根据复合函数求导法则,写成矩阵形式,其中,称为雅可比矩阵。,由此可得,其中,雅可比逆矩阵；,雅可比行列式。,利用 的表达式，可以将形状函数 对整体坐标,x,y,的偏导数，转换成对局部坐标 的偏导数。,例如,其中,此外，整体坐标系与局部坐标系的面积微分之间的关系为,到此，完成 计算中，从总体坐标向局部坐标的变换。,积分区域：从总体坐标系下的四节点四边形转换为局部坐标系下的正方形（边长为,2,）。,即,上式积分域简单，但是函数形式复杂，常采用数值积分的方法。,三 能进行等参数变换的条件,在 的计算中，涉及：,（,1,） 计算中 ，否则雅可比矩阵的逆矩阵不存在；,（,2,） 中，也要求 。,因此， 是 中将总体坐标转换成局部坐标的必要条件。,问题：有使 的四边形吗？如果有，什么样的四边形会导致 ？怎样避免？,为此，展开雅可比矩阵,其中,由此可得雅可比行列式为,可见， 是 的线性函数。,要使,在单元上不等于,0,，则只需 在四个节点上具有相同的正负即可，否则， 在单元上某些地方一定等于零。,以节点,1,为例，将局部坐标 和整体坐标 带入 中，可得,其中， 为总体坐标系中节点,1,、,2,边和节点,1,、,4,边的长度， 是这两条边的夹角。,同理，在节点,2,、,3,、,4,， 的值分别为,由于四边形的内角之和为,所以要使 全部为负是不可能的。全部为正是可能的。部分为正、部分为负是可能的，但是，是不能接受的（将使 ）,大于,0,的条件是：,上式意味着四边形为凸四边形，而不能有一个内角大于或等于 。,可见，任意四边形的任意不包括：任意两条边通过延伸在单元内出现交点。,为了保证计算精度，在划分单元时，应尽量使四边形单元的形状与正方形相差不远。,四 单元负荷向量的形成,单元上的外力不一定是集中力。如果是体力和面力，就涉及到把外载荷等效到节点上去的问题。即使是集中力，也有可能不落在节点上。,等效方法：虚功原理。即体力和面力对单元的外力功，等于等效的节点外力对单元的功。,上式是等效的节点外力对单元的虚功等于单元内储存的虚应变能。,1,等效体力负荷节点向量 的形成,设单位体积所受到的体积力，即体积力密度（二维）,则该单元上的体积力所做的虚功为,其中 分别为单元节点虚位移向量和单元形状函数。,是单元内任意一点的位移,而等效的节点外力的虚功为,与体积力等效的节点外力。,由此,2,与表面力等效的节点向量 的形成,设单元表面所受到的表面力密度（二维为单位长度上）为,则该单元上的表面力所做的外力虚功为,同样，等效的节点外力的虚功为,与表面力等效的节点外力。,从而,其中,比如，在单元 和 边上， 则,比如，在单元 和 边上， 则,3,没有落在单元节点上的外力集中力向单元节点的等效，即 的形成,设外力没有作用在节点上，而是作用在单元边界上某点,C,处，力的大小为 ，该力在虚位移上所做的功为,式中 分别是外力集中力作用点处的虚位移和形函数在该点的值。,从而,这样，单元节点的外力由等效的体积力节点向量、等效的表面力节点向量和等效的集中力节点向量合并成为,得到了单元刚度矩阵和单元节点负荷向量，就可以按节点号叠加总体刚度矩阵和总体节点负荷向量，最后得到方程组,置入边界条件后，即可求解。,5-4,八节点曲边四边形等参数单元,四节点四边形等参数单元的不足,：（,1,）阶次仍嫌不高；（,2,）直边对曲边拟合精度不高。,故提出八节点曲边四边形单元。,一 平面八节点等参数单元位移插值函数,还是完成两件事情：（,1,）以总体坐标系的节点位移，作为局部坐标系下的节点位移，建立插值函数；,（,2,）以（,1,）的插值函数，建立总体坐标系下的节点坐标和局部坐标系下的节点坐标的关系。,先讨论第一个。,局部坐标系的单元和节点配置：边长为,2,的正方形，节点位于正方形的角点和各边中点上。节点坐标为,位移插值函数取为下列多项式,注意到在单元边界上，最高次数为二次。这样的插值函数称为双二次插值函数，相应的插值称为双二次插值。显然，它比双线性插值的阶次要高。,不出现边界上的三次项是为了相容性要求。边界上的二次函数，可以由边界上的三个节点的位移值唯一确定。因为相邻单元边界上的节点是公共节点，其位移唯一，所以相邻单元边界上任意一点的位移都相等。,下面来讨论第（,2,）个问题。,如果用上述插值函数联系总体坐标系下的单元坐标和局部坐标系下的单元坐标，那么，总体坐标系中相邻单元的边界和节点值也是相等的，不会开裂，也不会重叠，这样，几何协调性（相容性）也得到了保证。,二 八节点等参数单元的形状函数,记节点处的位移值为,则单元内任意一点的位移可以写为,为形状函数，它可以利用形状函数的性质而得到。,其中 是节点,i,的局部坐标。,是双二次函数。因为形状函数的值是唯一的，所以只要求出来，就是结果。,以 为例。,（,1,） 在节点,28,其值为,0,；,（,2,） 在节点,1,的值为,1,。,节点,28,在直线,35,、,57,和,28,上，这三条直线的方程为,由此，可以构造 的表达式为,它是双二次函数。再利用 ，即,可以求得,同理可得,再求中点。以,2,为例，它在节点,1,和节点,38,为,0,。注意到三条直线,17,、,35,和,57,通过这些点,于是可以构造,它也是双二次函数。再利用 ，即,可以求得,同理可得,将,上述式子写成统一形式,三 等参数变换及其实现条件,局部坐标到总体坐标的坐标变换式,上式可以满足坐标变换的几何相容性。,雅可比矩阵,相应地可以写出雅可比矩阵的逆矩阵 和雅可比行列式 。,从局部坐标系到总体坐标系的变换，局部坐标系下的直边，转换成了总体系下的二次曲线（特殊情况下成为直线）。由此，可以提高单元对于曲边边界的拟合精度。,为了计算 ，同样要保证 。类似于,4,节点直边四边形单元，对于,8,节点曲边四边形的节点，也要作一些限制。这些限制是,（,1,）对边延长不能相交；,（,2,）中间点在边界的中间部分。如果位于,1/3,分点，计算出现较大偏差；如果位于,1/4,分点，则计算结果会完全不正确，甚至方程出现奇异性，计算无法进行下去。,单元划分时应做到,（,1,）尽量接近正方形；,（,2,）中间节点尽量位于每边的,1/2,分点处。,四 等参数单元的特点,计算精度高。有资料表明，分别采用三节点三角形线性单元和,8,节点等参数单元求解，节点总数相同（,200,），三角形单元数,300,个左右，等参元只用了,50,个。三角形在应力集中处计算结果不理想，等参元与实际情况基本相同。因此，现在有限元分析大多采用等参元计算。,缺点,：（,1,）等参数变换，使程序复杂；,（,2,）数值积分，使计算时间增加。,如果认为精度还不够高，还可以在每边再加一个节点，成为,12,节点等参元。实际情况是,8,节点等参数单元求解平面问题，精度通常足够了。,除了上述四节点或,8,节点等参数单元外，还可以在单元的四条边上配置不同数目的节点，形成其它类型的等参数单元。这样的单元可以实现单元分划时的疏密过渡。,5-5,数值积分及其应用,计算单元刚度矩阵,被积函数 复杂，需要使用数值积分。,采用常规数值积分，要提高精度，就要增加计算时间，使用中必须考虑精度和计算开销的统一。,为简单计，我们从一维函数积分 开始讨论数值积分方法。之后，再将其推广到二维和三维情况。,一般来讲，有两种数值积分方法。,一 牛顿,柯特斯求积法（,Newton-Cotes,）,基本想法：在,-1,1,区间取,n,各分点：,计算各点处的函数值 。用这些函数值去构造一个,m,次,多项式，并用这个多项式的积分代替对 的积分。,实用中，将,m,次多项式积分后，带入相应的节点坐标和节点处的函数即可。,比如，,m,=0,时，为矩形积分公式,m,=1,时，为梯形积分公式,m,=2,时，为辛普森积分公式,按照这样的办法，可以将类似的公式一直写下去。统一的形式为,求积系数或积分系数。,因为,n,个函数值最高可以确定一个,(,n,-1),次,多项式，所以牛顿,柯特斯积分公式具有,(,n,-1),次的,代数精度，即公式的误差是 。,二 高斯求积法,如果不指定积分点的位置，只给定积分点的数目,n,，,这样，需要求解,2,n,个未知量，即,n,个积分点的位置和这,n,个积分点位置处的函数值,2,n,个未知量求出后，可以构造一个最高次数达,(2,n,-1),次的多项式，其误差是 。,这样做可以提高所构造的求积公式的精度。,例如，仍采用公式,按照上面这样的想法得到的数值积分公式，称为高斯求积公式。,高斯积分中， 的求解比较复杂。,例如，积分点应是勒让德多项式 的,n,个不同的实根，而积分系数 也与勒让德多项式的导数有关。,工程中是利用已被求出的 所列成的表格来查取使用的。,积分点坐标,k,n,积分系数,H,k,2,3,4,5,0.577 350 269 2,0.774 596 669 2,0,0.861 136 311 6,0.339 981 043 6,0.906 179 845 9,0.538 469 310 1,0,1,0.555 555 556,0.888 888 889,0.347 854 845 1,0.652 145 154 9,0.236 926 885 1,0.478 628 670 5,0.568 888 889,比较牛顿,柯特斯法可见，高斯法求积精度高（积分点数目相同），或者积分点数少（精度相同）。,对于等参元的有限元分析，被积函数的计算十分复杂，因此，计算次数少的高斯积分法更为适用。所以，等参元积分中的数值积分都采用高斯求积法。,三 高维情况的高斯求积法,对于二维情况，将二重积分化为二次积分即可。,其中， 为高斯积分点，,为相应的积分系数。此时，平面上的积分点数有 个。,对于三维情况也是类似的,此时，空间区域有积分点 个。,四 合适的高斯积分阶数,数值积分很费时间，高精度数值积分尤其如此。,例如，对于平面问题，当,n,=2,时，每个单元有,4,个积分点；当,n,=3,时，每单元有,9,个积分点；当,n,=4,时，每单元有,16,个积分点。,二维是平方，三维是立方。所以，为了提高计算速度，在保证收敛和计算精度的允许范围内，应尽量降低数值积分的阶次。,经验表明：二维情况下，四节点单元,n,=2,，,八节点单元,n,=3,较好。有时，八节点单元,n,=2,也会有较好的计算精度。,小结：,（,1,）等参元是指单元内的插值函数的公式和坐标变换的公式一致。等参元的种类很多，可以根据实际需要选用。,（,2,）以局部坐标系中的规则的正方形构造插值函数，继而完成从局部坐标系到总体坐标系的坐标变换和单元刚度矩阵从总体坐标系到局部坐标系的转换。整个讨论和计算都是在局部坐标系中的规则的正方形单元内进行，最后在总体坐标下叠加各单元刚度矩阵，然后求解。,（,3,）等参元的优点：,1,）单元选择的自由度大；,2,）能很好地模拟曲线边界；,3,）计算精度高；,4,）输入数据少。这一点对于复杂区域的求解特别突出。,（,4,）使用等参元：,1,）要坐标变换；,2,）数值积分；,3,）程序编制复杂；,4,）形成刚度矩阵时间长。但不对计算带来本质的影响。,（,5,）为了保证计算精度和使计算切实可行，单元分割时要满足单元不过分歪斜的要求，中间节点尽量等分边长。,（,6,）等参数单元法只能满足在单元边界上函数值的相容性条件，其一阶导数在边界上不连续。因此，原则上只能适用于二阶微分方程所描述的问题，如应力分析，稳态温度场，电磁场分析等。,对于板弯曲问题所必须满足的四阶微分方程，上述等参元方法一般不适用。,