第五章-线性三角形单元.ppt

三角形单元*2引 言杆梁结构：由于有自然的连接关系，可以凭一种直觉将其进行自然的离散。连续体：它的内部没有自然的连接节点，必须完全通过人工的方法进行离散。三维问题平面问题平面应力平面应变平面问题平面应力平面应变离散*3三节点平面三角形单元节点1的位移节点2的位移节点3的位移三节点三角形单元的位移函数可假设为：“位移函数”也称“位移模式”，是单元内部位移变化的数学表达式，是坐标的函数。有限元分析必须事先给出（设定）位移函数。一般而论，位移函数选取会影响甚至严重影响计算结果的精度。弹性力学中，恰当选取位移函数不是一件容易的事情。有限单元法中当单元划分得足够小时，把位移函数设定为简单的多项式也可得到相当精确的结果。这正是有限单元法具有的重要优势之一。引入位移函数的概念：*4平面三角形单元显然，三角形三个节点的的位移可由下列方程给出，在各节点上的水平位移方程为：u1=1+2 x1+3 y1 u2=1+2 x2+3 y2 u3=1+2 x3+3 y3解出*5平面三角形单元假设求得其中A是三角形的面积*6平面三角形单元式中N1，N2和N3是坐标的函数，反映了单元内近似解的形态，称为单元的形函数，数学上它反应了由节点的场量对单元内任意一点场量的插值，也叫做插值函数。三个函数其实描述的就是单元上近似解的插值关系，它决定了近似解在单元上分布的形状，所以称它为形函数（shape function）。这里值得注意一下的是近似解，前面我们说过，假设位移模式是线性变化的，实际情况并不一定是线性变化的，所以我们通过所做假设得到的结果只能说是近似解，而不能说是精确解。为什么叫形函数?同理*7平面三角形单元其中ijki=1,2,3j=2,3,1k=3,1,2ijk三角形的形函数可统一表示为：*8形函数的性质在单元任一点上三个形函数之和等于1(单位分解性)1.三个形函数只有两个是独立的2.当三角形单元的三个结点的位移相等第一列与它的代数余子式乘积之和第一列与第二列的代数余子式乘积之和第一列与第三列的代数余子式乘积之和2A00*9形函数Ni 在节点i 上的值等于1，在其它节点上的值等于0。Ni=1ijmNj=1ijmNm=1ijm形函数的性质*10在三角形单元的边界ij上任一点（x，y），有:形函数的性质xxixjxyNi(xi，yi)j(xj，yj)m(xm，ym)Ni(x、y)1证ij方程*11形函数的性质相邻单元的位移在公共边上是连续的ijpm形函数在单元上的面积分和边界上的线积分公式为 式中 为 边的长度。xxixjxyNi(xi，yi)j(xj，yj)m(xm，ym)Ni(x、y)1Ni=1ijm*12形函数的性质完备性包含常应变项和刚体位移项如果在势能泛函中所出现的位移函数的最高阶导数是m阶，则选取的位移函数至少是m阶完全多项式。协调性相邻单元公共边界保持位移连续如果在势能泛函中所出现的位移函数的最高阶导数是m阶，则位移函数在单元交界面上必须具有直至(m-1)阶的连续导数，即Cm-1连续性。如果在单元交界面上位移不连续，表现为当结构变形时将在相邻单元间产生缝隙或重叠，这意味着将引起无限大的应变，这时必然会发生交界面上的附加应变能补充到系统的应变能中去，有限元解就不可能收敛于真正解。收敛单元尺寸趋于零时，有限元解趋于真解*13形函数的性质当单元的位移函数满足完备性要求时，称单元是完备的（通常较容易满足）。当单元的位移函数满足协调性要求时，称单元是协调的。当势能泛函中位移函数的导数是2阶时，要求位移函数在单元的交界面上具有C1或更高的连续性，这时构造单元的插值函数往往比较困难。在某些情况下，可以放松对协调性的要求，只要单元能够通过分片试验(Patch test)，有限元分析的解答仍然可以收敛于正确的解。这种单元称为非协调单元。分片试验由B.M.Irons首先提出，已经证明它给出了收敛性的充分条件。*14单元应变和应力矩阵应变矩阵*15单元应变和应力矩阵 由于 与x、y无关，都是常量，因此B矩阵也是常量。单元中任一点的应变分量是B矩阵与单元节点位移的乘积，因而也都是常量。因此，这种单元被称为常应变单元。*16单元应变和应力矩阵平面应力：应力矩阵平面应变：用平面应变弹性矩阵代入得到类似结果。*17单元应变和应力矩阵由于同一单元中的D、B矩阵都是常数矩阵，所以S矩阵也是常数矩阵。也就是说，三角形三节点单元内的应力分量也是常量。当然，相邻单元的E，A和bi、ci(i，j，m)一般不完全相同，因而具有不同的应力，这就造成在相邻单元的公共边上存在着应力突变现象。但是随着网格的细分，这种突变将会迅速减小。*18单元分析几何关系位移函数本构关系平衡关系单元刚度矩阵*19单元应变能单元应变能 U为：ijmxyh注意到弹性矩阵D的对称性*20刚度矩阵引入刚度矩阵K：则：注意：hdxdy的实质是任意的微体积dv，于是得 Ke的一般式：*21单元外力功单元受到的外力一般包括体积力、表面力和集中力。自重属于体积力范畴。表面力指作用在单元表面的分布载荷，如风力、压力，以及相邻单元互相作用的内力等。ijmxyqVijmxyqsijmxyfc*22单元外力功（1）体积力所做的外力功ijmxyqV*23单元外力功（2）面力所做的外力功ijmxyqsqs*24单元外力功（3）集中力所做的外力功ijmxyfc当结构受到集中力时，通常在划分单元网格时就把集中力的作用点设置为节点。于是单元集中力fc的势能Vc为综合以上诸式，单元外力的总外力功V为*25系统势能扩充叠加扩充叠加系统势能*26单元刚度矩阵的扩充叠加mijmij单元编号 ijm*27单元等效节点载荷列阵的扩充叠加mij单元编号 ijm*28能量原理和系统平衡方程系统势能根据弹性力学能量原理：结构处于稳定平衡的必要和充分条件是总势能有极小值。上式是从能量原理导出的系统平衡方程。这个方程表达了节点力与节点位移之间的关系。*29刚度矩阵单元刚度矩阵：1D:2D:系统刚度矩阵：弹性矩阵D的对称性Ke对称K 对称*30刚度矩阵刚度矩阵K的详细内容为：i、j是行列号，Ns为系统自由度数。*31刚度矩阵(1)刚度矩阵中每个元素有明确的物理意义。例如，Kij表示当节点位移中第j个元素为1(dj=1)其余元素为零时，引起的单元力中的第i个节点力fi。把平衡方程写开主对角线上元素Kii(i=1,Ns)恒为正值：位移和作用力同向*32刚度矩阵(2)K的每一行或每一列元素之和为零以上式中第i行为例:f2i-1=0 f2i=0f2j-1=0 f2j=0f2m-1=0 f2m=0ijmxyijm11当所有节点沿x向或y向都产生单位位移时，单元作平动运动，无应变，也无应力，因而单元结点力为零（不含初应力）。所以有即，K的每一行元素之和为零。由于对称性，每一列元素之和也为零。*33刚度矩阵(3)系统刚度矩阵是奇异矩阵（即K的行列式为零）(4)系统刚度矩阵是常量矩阵系统的节点力和节点位移成线性关系是基于弹性理论的结果。刚度矩阵是在系统处于平衡状态的前提下得出的。作用在它上面的外力必定是平衡力系。然而，研究系统平衡时没有引入约束。承受平衡力系作用的无约束系统，其变形是确定的，但位移不是确定的。所以出现性质（2）中的“平动问题”，即可以发生任意的刚体运动。从数学上讲，系统平衡方程的解不是唯一的或不能确定的。由此，系统刚度矩阵一定是奇异的。单元刚度矩阵也一定是奇异的。*34位移边界条件的处理系统刚度矩阵是奇异矩阵，其物理原因是结构缺少刚性位移的约束，实际的工程结构都受有足够的支承约束，排除了发生任何刚体位移的可能性，因此，必须引入位移约束。有限元中，位移约束都设置在节点处。这里，只讨论刚性约束情况，即被约束的位移分量为零。设讨论的结构有Nn个节点，每个节点有ndf个自由度。则系统的总自由度为Ns，且节点总位移列向量d中共包含Ns个分量*35为了引入位移约束，把节点总位移列向量d分成两部分。一部分是不受约束的位移分量，记为df。另一部分是受刚性约束的位移分量，记为dr。不失一般性，设1N号位移分量是不受约束的；N+1N+Nr共Nr个分量是受刚性约束的。即：位移边界条件的处理*36位移边界条件的处理显然不受约束的节点位移的总数N为 N=Ns-Nr 对方程中的刚度矩阵K和节点荷载向量列阵f也作相应分割，则得到式中，ff是已知力边界，fr是约束反力。*37位移边界条件的处理按矩阵乘法规则得每个受刚性支承约束的位移分量都等于零，即从而得到*38位移边界条件的处理 Kff 引入约束后的约化的系统刚度矩阵。这是一个非奇异 矩阵，它的逆矩阵Kff-1是存在的。引入约束后的约化的系统平衡方程在分析计算时，从无约束的系统刚度矩阵K中删去与受约束位移号对应的行和列，再将矩阵压缩排列成NN阶方阵，即为约化后的结构刚度矩阵。*39位移边界条件的处理00010置一法置一法置一法置一法显然*40位移边界条件的处理乘大数法乘大数法乘大数法乘大数法显然*41节点位移和约束反力通过求解平衡方程即可解出全部未知的节点位移：约束反力约束反力约束反力约束反力把解出的d代入未经修改的平衡方程，即可得到约束反力：关于上述方程的解算方法，一般不采用求逆的方法求解，而是直接采用高斯消元法等求解线性方程组的方法求解求解。施加边界条件后，得到修改后的平衡方程(未约化的)(约化的)或节点位移节点位移节点位移节点位移或*42单元应变和应力根据三角形节点的位移，求出单元应力应变为ijk如何求系统应变能和节点应力？*43有限元解的收敛性由于在有限元计算中引入了结构离散和位移模式，导致有限元计算结果和真实解的偏差。单元划分越小、位移模式取得越接近真实变形，解答越收敛于真实解。当单元的位移模式采用解析的位移解时，有限元的计算结果和解析结果是相同的。然而，许多情况无真实的位移模式可以借用，只能寻求其近似函数，不可避免带来计算精度问题。实践证明：只要位移模式满足单元的完备性准则和协调性条件，就保证了有限元的解答收敛于真实解。系统表现过刚系统表现过刚系统表现过刚系统表现过刚*44有限元计算过程框图剖分结构为有限个单元，对节点、单元编号构建单元刚度矩阵和单元等价节点力向量组装系统刚度矩阵并引入约束组装整体等价节点荷载向量从节点平衡方程解未知节点位移计算结构应变、应力*45解综合方程Kd=f得结构节点位移d从d中找单元位移de用公式用公式=Bde和和 =D，计算应力，计算应力应变应变把单元刚度矩阵组装成系统刚度矩阵K离散结构为若干单元建立单元刚度矩阵Ke形成等价节点荷载f形成单元等价节点力有限元计算流程图*46关于三角形单元形函数的一点补充 2-3-P:同样,3-1-PA2 1-2-PA3 面积坐标面积坐标面积坐标面积坐标 面积坐标的定义 在三角形内任意一点P定义 *47关于三角形单元形函数的一点补充面积坐标与形函数的关系 面积坐标与直角坐标的关系*q考考虑虑一一个个平平面面应应力力问问题题如如图图所所示示，假假设设厚厚度度h=1，材材料料为为各各项项同同性性，杨杨氏氏模模量量为为E=1，泊泊松松比比为为=0，相相关关力力和和位位移移边边界界条条件件如如图图中中所所示示，问问题题左左端端为为固固定定约约束束。试试用用两两个个三三角角形形单单元元分分析析此此问问题题，三三角角形形单单元元的的网网格格划划分分如如图图所所示示。试试求求问问题题各各节点位移节点位移u、v和应力和应力x，y和和xy。例题*对于三角形单元，其对于三角形单元，其B B矩阵的表达式为：矩阵的表达式为：123For e=1:1,(2)3,(1)2,(3)1,(2)2,(3)3,(1)例题*123For e=2:1,(3)3,(1)2,(4)1,(3)2,(4)3,(1)对于三角形单元，其对于三角形单元，其B B矩阵的表达式为：矩阵的表达式为：例题*组装刚度矩阵 1,(2)3+3,(1)2+1,(3)2,(4)3+3,(1)1,(2)2+1,(3)2,(4)例题*载荷向量123123F=11,(2)3,(1)2,(3)1,(3)3,(1)2,(4)2+1,(3)2,(4)3+3,(1)1,(2)例题*系统方程可表示为：系统方程可表示为：例题*施加 BCs:例题*55例题三角形单元为常应变单元，应力结果在一个单元内也为常数。对于单元1，*56例题节点2仅在第1个单元内，其应力结果直接为节点4仅在第2个单元内，其应力结果直接为节点1和节点3为两个单元共有，其应力结果为两个单元内应力的平均，为57三角形单元的不足l3节点三角形单元是常应变（常应力）单元。在应变梯度较大的部位（亦即应力梯度较大的部位），单元划分应适当密集，否则不能反映真实的应变变化而导致较大的误差。l提高计算精度的其它措施采用高精度三角形单元（2次单元、3次单元）采用四边形单元（1次单元、2次单元）