二、-弹性力学有限元法基本原理(一)

第二单元弹性力学有限元法基本原理（一）第一节里兹法的有限元形式由于需要在整个求解区域上假设试探函数，经典里兹法在解决实际问题时，尤其是几何形状复杂的二、三维问题，具有局限性。解决上述问题的办法是在求解区域上分片假设试探函数。下面以一维直杆的分析为例子，研究基于里兹法的有限元位移法基本原理和求解过程。离散化考虑图2-1（a）所示的一维直杆结构，杆的总长度为3L。载荷与图1-4相同。杆分为三个子区域，称为单元。单元之间的连接点称为节点。这一步骤称为离散化。节点位移是问题的基本未知量。（a）一维直杆的分域（截面积A，弹性模量E，轴向受力）（b）杆的有限单元图2-1分片假设单元上位移试探函数xbbu21Lx0xbbu43LxL2xbbu65LxL32是待定常数。ib首先按如下形式假定分片位移场为了使得上述假设位移场是“许可位移”，上述多项式待定系数必须满足一定约束关系（3个），显然，该问题的独立参量（广义坐标）只有3个。把上述假设的分段线性位移场代入势能泛函：分段计算上式积分，再应用驻值条件，可立即求出待定系数，位移场就完全确定，进而可求出各单元应力。LxLpqudxAdxuE302,302上述是原理上的里兹法有限元形式求解过程，其关键是采用分片多项式拟合全域上的可能位移场——分片试探函数。这个做法正体现了有限元法的实质。上面形式的分片位移试探函数有下列缺点：1)必须对它进行调整，使其满足连续条件和边界约束条件；2)多项式系数作为广义坐标缺乏明显的物理意义。因此，上述不是通常意义上标准的有限元形式，仍然具有局限性，如对于二维以上的问题使各单元之间分片多项式保持连续性很难处理。下面用节点位移未知量作为待定参数（广义坐标），得到其标准有限元形式。重新构造单元内位移试探函数dNu矩阵形式为：LsLsLN是插值基函数矩阵，称为“形函数”矩阵。21uud称为单元1的节点位移列阵。21uLsuLsLu对单元1有：其中：离散结构中，节点位移分量是问题的基本未知量。在每个单元内通过对节点位移插值，分片建立位移试探函数：其它两个单元也有同样的插值位移试探函数：单元2：dNu32,uud单元3：dNu43,uudN同上每个单元中，位移试探函数是位置坐标的简单函数（线性），任意一点函数值取决于单元节点位移，且在单元节点上试探函数值等于待定节点位移。单元上的位移试探函数又称为单元位移模式。显然，在使用插值试探函数情况下，整个杆上，由各单元位移函数拼接而成的试探位移场是连续的。只要我们记住，得到的就是可能位移场。这样的位移试探场以未知节点位移作为待定参量（广义坐标）。下面进一步实施里兹法求解。01u计算离散系统的势能泛函单元应变表达式：sxxuu,,得到单元上应变表达式：dBdNdsdxLxLpqudxAdxE302302由应变~位移关系：上式就是用节点位移表达单元应变的公式。LLB11称为单元应变矩阵其中ddsBAEBddsAEAdsEULTTLxTxxL0020212121或dkdUT21计算一个单元内应变能(单元1）：其中：——称为杆单元刚度矩阵。dBdNdsdx11110LAEdsBAEBkTL由于单元的几何、物理参数相同，上述单元应变能表达式对于3个单元相同，只是节点位移列阵的分量不同。载荷与上节例题相同：cxq载荷在三个单元内局部坐标下的表达式分别为：，，csq)(sLcq)2(sLcq按外力功的积分式，分别计算三个单元的外力功第1单元外力功为：2160200TLTTLLTdcLcsdsNdqdsuquds同理，第2、3单元外力功分别为：和5462TdcL8762TdcL计算单元上外力的功：为了适应矩阵形式运算，将单元势能矩阵表达式中的用整体节点位移向量代替，同时单元刚度矩阵扩展成总体规模（4×4）。则各单元相加后系统总势能为：dTuuuuD432187006054060021611001100000000000000011001100000000000000011001121222cLcLcLDDLAELAELAEDTTpdkdUT21系统总势能等于三个单元应变能之和减去三个单元外力功。应用势能驻值条件：0DpRDK812616110012100121001124321cLuuuuLAE01uRDDKDTTp21即划去第一个方程，解出其余三个方程得到：即：引入约束条件：简写为：得到有限元求解方程——系统平衡方程：单元应力由公式得到。27231333432AEcLuuudBEExx结合单元位移模式就得到整体上近似位移场。dNu图2-2受轴向力杆的精确结果和有限元结果本问题有限元解与精确解的比较如图2-2所示。有限元解与经典里兹解对比第二节常应变三角形单元解平面问题上一节以受轴向力的弹性杆为例讨论了有限元位移法的基本原理和步骤，揭示了有限元法的本质特征。本节讨论将该方法推广到解决弹性力学平面问题。弹性区域离散化采用3节点三角形单元（T3单元）。1、结构离散化连续求解区域划分为有限个三角形区域（单元）；或者连续体离散为有限个小单元的结合体，单元之间在节点处连接。问题的未知量转化为节点位移。图2-3二维区域离散化2、单元位移模式及插值函数从离散结构中取出一个典型的分片区域（单元），研究二维区域上分片多项式假设位移场（单元位移模式）。单元如图2-4所示。单元节点采用局部编号i，j，m（逆时针旋转），每个节点两个自由度（位移分量）：Tmmjjiimjievuvuvuaaaa图2-4三节点三角形单元——称为单元节点位移列阵。iiivua（i,j,m）则每个单元六个节点位移分量：在单元的三角形区域构造简单位移试探函数。采用x，y的一次多项式：yxvyxu654321~16是待定系数，称为广义坐标。显然，用这种形式的单元位移模式构造整个求解区域上的许可位移场将非常困难。解决办法：将上述多项式系数广义坐标代换为单元节点位移广义坐标（插值法）。将三个节点坐标分别代入上述位移多项式：iiiyxu321jjjyxu321mmmyxu321解上述方程，用节点位移表达多项式系数：AyxyxyxDmmjjii2111（A为三角形单元的面积。）yxu321定义上面方程的系数行列式为D：用公式法（克莱姆法则）解出广义坐标得：~132)(211111mmjjiimmjjiiubububAyuyuyuD3)(211111mmjjiimmjjiiucucucAuxuxuxD1)(211mmjjiimmmjjjiiiuauauaAyxuyxuyxuD上式中，系数分别是行列式D第i,j,m行第1，2，3列元素的代数余子式，其值取决于节点坐标。),,(,,mjicbaiiiiiiyxu321jjjyxu321mmmyxu321AyxyxyxDmmjjii2111mmjjiiuNuNuNummjjiivNvNvNv)(21ycxbaANiiii将求得的广义坐标代回前面位移多项式，整理得到以单元节点位移分量表示的位移模式：~16其中（i，j，m）iN称为单元的插值基函数或形函数。4)(21mmjjiivavavaA5)(21mmjjiivbvbvbA6)(21mmjjiivcvcvcA同理，对y方向位移函数作代换，可求得广义坐标~：46每个节点对应一个形函数。yxvyxu654321mmjjiimjimjivuvuvuNNNNNNvuu000000eeNaaNNNmji位移模式的矩阵形式为：式中：mjiNNNN称为单元的形函数矩阵mmjjiiuNuNuNummjjiivNvNvNv1mjiNNN3、形函数的性质和几何意义1）相应于某节点i的形函数在i节点上值为1，在j，m节点上值为0。2）单元中任一点各形函数之和等于1，即：提示：三个形函数只有两个独立。该性质保证各节点函数值相同时，插值得到单元上任意一点函数值亦相同。3）推论：对三节点单元，形函数在单元内部和边界上均为线性函数，形函数在单元一条边上的值，只跟边端点位置有关，与第三点位置无关；某节点的形函数在对边上恒为0。ijijyxNijjji当当01),(),,(mji4）形函数的几何意义图2-5形函数的几何意义根据三角形单元形函数的上述性质也可以推断，单元的假定位移在内部和边界上都线性分布，边界上的位移只跟两端节点位移有关，保证了整个求解区域上位移的连续性。4、用节点位移表达单元应变和应力矩阵称为单元应变矩阵，其子块为：单元位移模式确定后，很容易得到用节点位移表示单元应变和应力的表达式。应用弹性力学平面问题几何方程的矩阵形式，得到：eemjiemjiexyyxvuxyyxBaaBBBaNNNLLNaLu00B将形函数分别代入上式，最后求得应变矩阵如下：（i,j,m）xNyNyNxNNNxyyxiiiiiiii000000LNBmmjjiimjimjibcbcbccccbbbA00000021B该应变矩阵的元素取决于节点坐标，是常数矩阵。因此，单元内应变、应力是常数。该单元称为常应变单元。单元应力根据平面问题的物理方程得到：其中：称为应力矩阵。将平面应力或平面应变问题的弹性矩阵代入，就可以具体计算出应力矩阵。至此，已完成在二维弹性体区域上构造位移试探函数，并做好计算系统总势能的准备。eexyyxSaDBaDmjimjiSSSBBBDDBSS5、利用最小势能原理建立有限元求解方程弹性力学平面问题的总势能泛函表达式如下：STTTptdstdxdytdxdyTufuD21其中：t——弹性体厚度；f——平面问题体积力向量；T——边界面力向量。应用前面在整个求解区域上分片假设的位移场（位移模式），在有限元离散模型上对上述总势能泛函进行分片计算（对各单元区域积分）并求和：esTeTeTptdStdxdytdxdyeeeTNafNaaDBBaTeTeeTe21——单元刚度矩阵；——单元体力等效节点力列阵；——单元面力等效节点力列阵；eTtdxdyeKDBBebTtdxdy