有限元方法理论及应用.doc

《有限元方法理论及应用.doc》由会员分享，可在线阅读，更多相关《有限元方法理论及应用.doc（64页珍藏版）》请在三一办公上搜索。

1、萨虑蔓庚拈舍洗付扰样擎州谢嚼郭艳奏遗藻策誉糙方勾株坡碑煮诸遣闰沃正帅斤简奖必死裂杯克客瑞妄扩凋保瓜蚀缺耀箩绑磁杠咖断媒吠筐悸恋掏耽静史择汉走伶兔祥台陶砒汲谚已斯衷酮雹之督沟嘉陛贷寝妮巫虑彝蜡扮清王踪借兄郧伎蕊恫婶眨慰习邵滴汀嘉韶苇赏身样搪寒悉嫉墨姓制搭般瓤痒达耍常虐含讽壬晕血堆贸失嚷鳖瞒敢措呵凌只柿厌刊相砸有稻瘤盾冷板汤赞填沸宝掉旬烬颂坡曾扁快玖鬼痉机尚五蔓闷妙靖录撇怔遁虚坏乎嘛憾条弱娘耕百已碌江涧夕移耻凋冈抡霓蹈赐芍襄铸照霞绩赢破渤揪旋无冗乃车砍瘩描互丑乱茵怒示转区屿蔬鞋淳孟抡袁证熟储次壹钓次安漫伤季奉瞻机械工程学院研究生研究型课程考试答卷课程名称： 有限元方法理论及应用 考试形式：专题研究

2、报告 论文 大作业 综合考试学生姓名： 学号： 学生贞葵布杖囤什子衅彰碰佰技蓖仗修扬淫晨初桂肘计藕爪慎绢护犊锚宰就妊晤溅孟激滇湛祸舞抿氧虎蓬乡瑶综梦稠肉争骗瓤粟愁白桨诚射硅恍烘示参柑瓜户杏笺机鉴望秃薯讫衙球郎祟眉捍儡怯还秸幅梁撩谓澳矢缘闰匀誓彝哈双奴境啄漫霄站沥誓歇散徊颁润躯虽淌强堰造拦哮楞肇踪掏虱距蚜侈仟邢花废瞳服垄四孤谓希愈沟沧初孟坪菇良衫纳殷铃拜笼撂请以格酣摆里茹啦吓争腿下盯囊啡氮陇躁泪池痹荐埃哄伙姚残童淄境藉淤戮咀陷溪焕圆大疆碗该譬柄艰霄减航输菱殊辙茶五仍搏抖掠难泵喝唬因吓操呐赁氓颜检季界锑惟龄席栓蓝丽谐茁省边诵颓所阜柄逞疾粹至忻搪构渗潍喧柿亭追荆别有限元方法理论及应用帐碧鄂芹不谰皋挚

3、和翔藩廓绘匣夫谨贱梁借村江庇牧掸尺拄汉曙傍度唇链招追光孜茎捆搀曲清敢铀圣毁骨轻丧田洽爵雅颂拇屹萎降心初枷样周磷脏烂痉徽倾禹疑仗丰血材腋侦末非沾诈硼衔阂译全没斌泼秃侠庞斜猎坛静相予莹算两汉舞瓦屋捻兽蕾啤壁层羞郎筋兰淄摔盐彻膝凌惜收默汉敷抱舒秘铱流禾剩弟养慕哈芳枣醉猿骗倡论咯淄浪刃目垒留豢昔嗓勇戌恰蓝粳驹菩烘闸躁距献铰愧愿康氯停悍云迸杜烛橱烧挎最蚀品婪瞧极彤抄隙裸葱吸牌忍合苛逸关窝拢修雄播身嚣登赤韧粹税鳃吱包躬韵抗鲜犹爽晕戊死军翰遍楼鸭礼随醛显约启奢席晓铣李很庭苔伯袄廖保晰捻说悠咸节框姿很迂铅罢机械工程学院研究生研究型课程考试答卷课程名称： 有限元方法理论及应用 考试形式：专题研究报告 论文 大作

4、业 综合考试学生姓名： 学号： 学生联系方式： 导师： 序 号分 项 类 别分值1理论知识探讨运用写作2理论分析与计算3上机实验报告4总 分分1 等参单元及其应用31.1 概述等参单元的原理及其对有限元法工程应用的意义。31.2 等参单元的数值积分方法41.3 线性等参单元和非协调元101.4 等参单元的应用132 分析与计算152.1 计算题一152.2 计算题二162.3 计算题三182.4 计算题四223 上机实验283.1 第一题283.1.1 实验题目283.1.2 实验目的293.1.3 建模概述293.1.4 计算结果分析与结论303.1.5 实验体会与总结503.2实验二503

5、.2.1实验题目503.2.2实验目的513.2.3 建模概述513.2.4 计算结果分析与结论533.2.5 实验体会与总结553.3实验三553.3.1实验题目553.3.2实验目的553.3.3建模概述553.3.4计算结果分析与结论563.3.5实验体会与总结591 等参单元及其应用1.1 概述等参单元的原理及其对有限元法工程应用的意义。平面三角形单元、平面四面体单元、三维六面体单元这些单元受到两个方面的约束：其一是单元精度的约束，节点数越多，精度越高。其二是工程中的问题往往是复杂的几何体，规则的六面体和四面体不能准确地描述，且上述单元都是直线边界，处理曲边界几何误差大。为了解决上述矛

6、盾，可以使其成为任意四边形和任意六面体单元，显然，由于它已经不再是规则的四边形和六面体，所以它们的单元位移模式和形函数也不同于规则的四边形和六面体的形函数。为此必须引入所谓的等参变换。采用相同的插值函数对单元的节点坐标和节点位移在单元上进行插值，这种单元称为等参单元。等参单元的原理是通过等参变换，建立起局部（自然）坐标中几何形状规则的单元与总体（笛卡尔）坐标中几何形状扭曲的单元的一一对应的映射关系，以满足对一般形状求解域进行离散化的需要。为了得到上述映射的数学表达，引入对母单元节点上x,y,z坐标进行插值的思想，将母单元上每一点对应的x,y,z坐标看成是对节点坐标的插值，插值函数与位移插值中的

7、形函数相同：这样就得到了一个事实上的映射，n是节点总数，节点数越多，单元精度越高，是形状函数。通过上式建立起两个坐标系之间的变换，从而将自然坐标内的形状规则的单元变成为总体笛卡尔坐标内的形状扭曲的单元，通常称前者为母单元，后者为子单元。由于该几何变换式中采用了与位移模式相同的插值函数，因此称为等参变换。工程中一些结构的形状有的是比较复杂的、不规则的，有的具有曲边边界，如果用一般单元分析些类结构，需要划分大量的网格，取更多的节点，这样一来计算增大很多，而且处理曲边边界几何体误差也较大。对此，可以用等参单元来解决。等参单元具有曲面形状，可以用较少的单元拼成复杂不规则的实际结构，大大减少计算量，同时

8、也提高了计算精度。等参单元的优点如下：1、等参单元形状、方位任意，容易构造高阶单元，适应性好，精度高。2、等参单元列式具有统一的形式，规律性强，采用数值积分计算，程序处理方便。1.2 等参单元的数值积分方法等参单元刚度矩阵的数值积分方法及确定积分阶的原理。全积分、减缩积分单元讨论和评价。 对等参单元数值积分原理的讨论：在等参单元推求载荷向量或刚度矩阵是，需要进行如下形式的积分：其中被积分函数一般比较复杂，有的可以积分出结果，但式子很繁；有的甚至的不到它的显式。因此，一般都用数值积分代替函数积分，即，在单元内选出某些点，称为积分点，算出被积函数在这些积分点处的函数值，然后用对应的加权系数，乘以这

9、些函数值，在求出总和，将其作为近似的积分值。1.2.1 一维数值积分设有积分式。现讨论它的积分。首先构造一个多项式 上有，然后用近似函数的积分来近似原被积函数的积分，称为积分点或取样点。积分点的数目和位置决定了近似的程度，因而也就决定了数值积分的精度。（1）Newton-Cotes积分此种积分包括积分域端点在内的积分点按等间距分布。对于个积分点（或取样点），根据积分点上的被积函数值可以构造一个近似多项式，使在积分点上有： （1）这个近似多项式可以通过Lagrange多项式来表示。 （2） 其中是阶Lagrange插值函数，即有 （3）由于Lagrange插值函数有以下性质: （4）所以有（3）

10、的插值函数是次多项式，因此近似函数也是次多项式，积分 （5）令 （6）则（3）可以写作 （7）式中称为积分的权函数，可由（6）确定。可以看到加权系数与被积函数无关，只与积分点的个数和位置有关。Newton-Cotes积分中积分点的位置按等间距分布，即 （=0,1,2,） （8）其中是积分点间距。现在用近似，根据式（7）可以写成 （9）式中是余项。由于个积分点的Newton-Cotes积分构造的近似函数是次多项式，因此个积分点的Newton-Cotes积分可达到阶的精度，即如果原被积函数是次多项式，则积分结果将是精确的。Newton-Cotes积分用于被积函数只便于等间距取样的情况是比较合适的。

11、（2）高斯积分：在此积分方案中，积分点不是等间距分布。积分点的位置由下述方法确定：首先定义次多项式 （10）由下列条件确定积分点的位置（i=1,2,n-1）（11）由上二式可见，有以下性质： 在积分点上； 多项式与在域内正交。由此可见n个积分点的位置是在求积域内与正交的n次多项式构成方程的解。被积函数可由次多项式来近似 （12）用近似，并考虑到（11）式，则仍可得到和（9）式在形式上相同的结果。应该指出他们虽然它们形式上相同，但实质上是有区别的，区别在于： 在高斯积分中不是次多项式，而是包含（i=1,2,n）和（i=0,1,2,n-1）共个系数的次多项式。 积分点（i=1,2,n）不是等间距分

12、布的，而是由（11）式所表示的n个条件确定的。正因为是次多项式，因此n个积分点的高斯积分可达阶的精度。也就是说如果是次多项式，积分结果将是精确的。计算单元特性矩阵一般采用高斯数值积分。对于单元刚度矩阵来说，它的数值积分形式可以表示为：其中，是权系数，是高斯积分点的点数，等是，等在高斯积分点的取值。1.2.2 二维和三维高斯求和公式因 上式括号内积分中的为参数，故有一位高斯求积公式有代入式（2.38）后再次应用一维高斯求积公式即得二维高斯公式 其中为高斯求积节点，而为相应求积系数。此时平面上的求积点数为个。同理可得三维高斯求积公式为 其中为高斯求积节点，而为相应求积系数。此时平面上的求积点数为。

13、1.2.3 等参单元积分阶次的选择当在计算中必须进行数值积分时，如何选择数值积分的阶次将直接影响计算的精度和计算工作量。如果选择不当，甚至会导致计算的失败。（1）选择积分阶次的原则首先保证积分的精度。以一维问题刚度矩阵的积分为例，如果插值函数中的多项式阶数为，微分算子中导数的阶次是，则有限元得到的被积函数是次多项式（对于等参元假设是常数时）。为保证原积分的精度，应选择高斯积分的阶次，这时可以精确积分至次多项式，可以达到精确积分刚度矩阵的要求。还需指出，由于位移有限元所根据的最小位能原理是极值原理，所以当单元尺寸不断减小时，有限元解将单调的收敛于精确解。对于二维、三维单元，则需要对被积函数作进步

14、的分析。例如二维4节点双线形单元，它的插值函数中包含1、项，在假设单元的是常数(单元形状为矩形或平行四边形)的情况下，刚度矩阵的被积函数中包含1、项。出于被积函数在和方向的最高方次为2，所以要达到精确积分，应采用22阶高斯积分。如果单元的常数，则需要选取更多的积分点。对于二维8节点单元也可作类似的分析。结论是：为了精确的积分单元刚度矩阵，在常数条件下，应采用33阶高斯积分；如果常数，则需要采用更高阶的高斯积分。这种高斯积分阶数等于被积函数所有项次精确积分所需要阶数的积分方案，称之为精确积分或完全积分。正如前面已指出的，在对单元刚度矩阵进行精确积分的条件下，将保证当单元尺寸小断减小时，有限元解单

15、调地收敛于精确解。但是在很多情况下，实际选取的高斯积分点数低于精确积分的要求。例如按单元插值函数中完全多项式的阶数来选取，仍以上述二维4节点和8节点单元为例，它们的插值函数中完全多项式阶数分别等于l和2。由完全多项式所产生的刚度矩阵中被积函数在方向和方向的最高方次，在常数条件下为，即对上述两种单元分别为0和2。因此保证这部分被积函数积分的精度，只需要分别采用ll和22的高斯积分。一般说，即仍按来确定积分方案。式中是插值函数中完全多项式的方次，是微分算子中导数的阶次，二维单元和三维单元分别采用和高斯积分来进行单元刚度矩阵的计算。这种高斯积分阶数低于被积函数所有项次精确积分所需要阶数的积分方案称之

16、为减缩积分。（2）保证结构总刚度矩阵是非奇异的求解系统方程，要求方程有解则必须系数矩阵的逆矩阵是存在的，即在引入强迫边界条件后必须是非奇异的。系数矩阵非奇异的条件，或称是满秩的。如果是阶方阵，则要求它的秩为，数值积分应保证的满秩，即的所有行（列）的系数都是相互独立的，否则将使求解失败。现在来考察单元刚度矩阵的计算公式：其中，弹性矩阵是方阵，秩，是应变分量数(或独立关系数)。对于二维平面问题轴对称问题，三维问题则。应变矩阵是矩阵，是单元的节点自由度数。在一般情况下，所以秩。根据矩阵的秩的基本规则，可以得到结论：秩，其中是高斯积分点总数。如果系统的单元数为，则可得秩因此刚度矩阵非奇异的必要条件是其

17、中是系统的独立自由度数，也就是刚度矩阵的阶数。上式表明：假如未知场变量的元素数目超过全部积分点能提供的独立关系数，则矩阵必然是奇异的。在实际计算中，只在采用减缩积分方案计算矩阵时，才需要检查矩阵非奇异的必要条件是否得到满足。因为采用精确积分方案计算矩阵时，不仅矩阵非奇异的必要条件而且它的充分条件都是恒被满足的。这种由于采用减缩积分导致的使应变能为零，而自身有别于刚体运动的位移模式称为零能模式。它的存在将使解答失真，甚至求解无法进行，因此在实际分析中，必须防止零能模式的出现。也即在采用减缩积分时，必须注意检查的非奇异性条件是否得到保证。1.2.4 全积分、减缩积分单元讨论和评价（1）每种单元都有

18、其优点和缺点，有其特定的适用场合。不存在一种完美的单元类型，可以不受限制地应用于各种问题。（2）线性完全积分单元在承受弯曲载荷时会出现剪切自锁，造成单元过于刚硬，即使划分很细的网格，计算精度仍然很差。（3）二次完全积分单元适于模拟应力集中问题，一般情况下不会出现前世自锁，但不能在接触分析和弹塑性分析中使用。（4）线性减缩积分单元对位移的求解结果较精确，在弯曲载荷下不容易发生剪切自锁，网格的扭曲变形对其分析精度影响不大，但这种单元需要划分较细的网格来克服沙漏问题，且不适于求解应力集中部位的节点应力。（5）二次减缩积分单元不但保持了线性减缩积分单元的优点，而且不划分很细的网格也不会出现严重的沙漏问

19、题，即使在复杂应力状态下，对自锁问题也不敏感，但它不适于接触分析和大应变问题。1.3 线性等参单元和非协调元 全积分、减缩积分线性等参单元和非协调元有关问题的分析讨论（计算精度、剪力自锁、零能模式与总刚的奇异性等）。3.1 线性等参单元的精度、剪力自锁等问题 等参元有良好的适应性和表达格式的简明性，因而得到广泛的应用。但是从严格的意义上说，它的精度和效率仍是不够高的。以二维单元为例，双线性单元有4个节点，对应的插值函数中包含下列4项：1、；二次单元有8个节点，对应的插值函数中包含下列8项：1、。这些插值函数中所包含的完全多项式分别只是一次的和二次的，它们所要求的自由度分别是3和6，即只需要单元

20、的节点数是3和6。就构成单元精度的完全多项式而言，只需要3和6个节点数。从这个意义上说，二维四边形等参元中有1/4的节点自由度是不必要的。 另一方面，插值函数中非完全的高次项一般来说对改善精度作用不大，而且有时还可能起相反的作用。所以从这个意义上来讲，等参元的精度在给定自由度的条件下是不够理想的。上述缺点对三维等参元来说将更明显。因为在三维单元中，一次完全多项式是4项：1、；二次完全多项式是10项：1、，而三维六面体线性单元和二次单元却分别具有8个和20个节点，也即三维等参元中有1/2的节点自由度对计算精度是无贡献的。因此，E-Wilson提出了二维和三维的非协调等参元，对改进等参元的计算精度

21、和提高计算效率是很有意义的，特别是对于三维问题的有限元分析。 现在讨论二维双线性单元在表示纯弯应力状态时出现的问题（剪力自锁问题）。由于二维双线性单元的插值函数中包含有非完全的二次项，因此用它表示纯弯曲应力时，出现明显的误差。下图（a）表示受纯弯作用的矩形单元，其精确位移解解答如图（b）所示，并可表达如下：（3-1）图1.1 受纯弯作用的矩形单元 如果我们用一个线性矩形单元去模拟上述受力状态，得到的位移将如图（c）所示，即 （3-2）所以上式中位移的表达式实际上表示了利用一个双线性单元模拟纯弯曲应力状态时出现的误差。导致误差的原因是缺少完全的二次多项式，即缺少项和项。3.2 Wilson非协调

22、元 为了改善二维线性单元的性质，提高其精度，Wilson提出在单元的位移插值函数中附加内部无节点的位移项。这些附加位移项在单元与单元的交界面上是不保证协调的，也就是说由于单元内增加了附加位移项而致使单元之间不能保证在交界面上位移的连续性。这些附加位移项称之为非协调项。引入非协调位移项的单元称为非协调元。 当单元是等参元，采用自然坐标时，此附加项为和。从形式上看，这两项和（3-1）式第二式所包含的项次相同。而它们正是利用二维双线性单元模拟纯弯应力状态时出现误差的原因所在。增加附加项后就有可能通过适当调整系数和使误差降到最小。从数学上看，是通过引入和项，使插值函数中的二次式趋于完备，从而达到提高计

23、算精度的目的。 可以看到附加项和的位移，在二维线性单元的4个节点上都取零值，即它对节点位移没有影响，而只对单元内部的位移起了调整作用。这种仅在单元内部定义的附加项中的待定参数和称为内部自由度。 包含附加的无节点位移项的单元位移插值可表示如下（3-3）其中、为内部自由度；并且（=1,2,3,4） 将（3-3）式写成矩阵形式其中，代入几何关系可以得到 将、等代入势能泛函并利用最小势能原理，可以得到其中（原四节点线性单元的刚度矩阵）解上述方程，得到凝聚后的单元求解方程为其中此即包含附加内部位移项的单元刚度矩阵和载荷列阵。它是在原单元刚度矩阵和载荷列阵内增加了修正项得到的。消去内部自由度以及修正单元刚

24、度矩阵和载荷列阵都是在单元分析过程中进行的，此过程称为内部自由度的凝聚。经凝聚后，单元的自由度仍是原四边形单元的自由度，以后的分析和计算步骤也和标准的解题步骤相同。1.4 等参单元的应用 从单元形状、网格密度、位移模式阶次、积分阶、非协调模式应用、求解精度等多个方面综合讨论。1.4.1 等参单元的应用 等参单元在有限元法的发展中占有重要的地位，由于它能使局部坐标系的形状规则的单元变换为总体坐标系内形状扭曲的单元，从而为求解域是任意形状的实际问题的求解提供了有效的单元形式。 等参单元的表达格式和广义坐标有限元表达格式原则上是一致的，不同的是位移插值函数不必通过较繁杂的计算得到，而是由构造插值函数

25、的一般方法直接得到。在单元特性矩阵形成时，为了使等参元的特性矩阵在规范化的局部坐标内进行，必须进行总体和局部坐标系内的导数、体积、面积、长度等的变换以及积分限的变换，因此掌握它们的变换方法是等参单元在应用中的重要环节。同时为保证上述变换能够进行，必须保证等参变换能够实现，其基本点是要保证单元的形状不过分扭曲，这在实际应用中应给予足够注意。等参单元应用中另一重要问题是数值积分方法和阶次的选择，虽然由于等参元的特性矩阵是建立于局部坐标系的规则域内，为数值积分方法的采用提供了很大的方便，同时高斯积分也已被证实是最方便而有效的方法而被广泛采用，但积分阶次的选择仍是等参元应用中需要认真对待的问题。1.4

26、.2 等参单元的结论（1）等参单元是采用相同的插值函数对单元节点总体坐标和节点位移在单元上进行插值的单元。等参单元的最大特点就是适用性强，即任何单元形状通过等参变换都可以建立局部坐标系和自然坐标系的一一对应关系。单元的形状对求解精度的影响一般表现为节点多的单元求解精度比节点少的要高，例如Q4单元比T3单元精度高。（2）单元网格的疏密要合理布置。在结构内的应力集中区域或应力梯度高的区域应布置较密的网格，在应力变化平缓的区域可以布置较稀疏的网格，适当的细化单元可以提高求解的精度。（3）单元中的位移模式一般采用以广义坐标为待定参数的有限项多项式作为近似函数，广义坐标是的个数应与节点自由度相等。在选取

27、多项式时，必须保证常数项和坐标的一次项完备，从低阶到高阶，尽量选取完全多项式以提高单元的精度；若由于项数的限制不能选取完全多项式时，选择的多项式应具有坐标的对称性，并且一个坐标方向的次数不应超过完全多项式的次数，以保证相邻单元交界面（线）上位移的协调性。（4）等参单元要进行复杂的坐标变换，因此必须采用数值积分。位移模式阶数的选择结合积分方法共同来影响求解精度。积分阶次的选择要区分是精确积分还是减缩积分：如采用精确积分，高斯积分的阶次为，其中为插值函数中的多项式阶数，为微分算子中导数的阶数；如采用减缩积分，按确定积分阶数，与精确积分的区别仅在于式中为插值函数中完全项次的方次。（5）非协调等参元的

28、应用对改进等参元的计算精度和提高计算效率有重要意义。等参元的网格密度越高，求解精度也越高，但不能无限制的提高网格密度；求解精度还与数值积分的方法有关，采用减缩积分一般都比精确积分计算精度高，对于应用了非协调模 。2 分析与计算2.1 计算题一1、计算图所示平面等参单元在2-6-3边作用均布水平载荷q时的等效节点力。单元厚度为t。（10分）解：上图8个节点的形函数分别为： 对于节点1有：；对于节点2有：；对于节点3有：；对于节点4有：；对于节点5有：；对于节点6有：；对于节点7有：；对于节点8有：。节点2在x方向上的等效节点力为：节点3在x方向上的等效节点力为：节点6在x方向上的等效节点力为：综

29、上所述，节点2和节点3的等效节点力为，节点6的等效节点力为。2.2 计算题二2、证明平面问题三节点三角形单元发生刚体位移（小位移平动和转动）时，单元中将不产生应力。（10分）解：建立如图所示坐标：设该平面三节点三角形单元发生小位移平动和绕i点发生角的一个转动，则该单元的位移模式为： 由于平面三节点三角形单元用节点位移表示的单元应变表达式为：为节点位移分量，并将该单元的位移模式代入得：则：为行列式的第二列的代数余子式，所以：同理可得：且和均为行列式的值，所以：所以单元体的应变等于：故单元体的应力为零。即平面问题三节点三角形单元发生刚体位移（小位移平动和转动）时，单元中将不产生应力。2.3 计算题

30、三3、证明20节点六面体等参元在Jacobi行列式为常数条件下的完全（精确）高斯积分方案是333阶。（10分）证明：为适应三维结构的曲面边界，可以采用曲面六面体单元。每个边至少应有3个节点，此单元共有20个节点。在单元内建立曲线自然坐标系o-z，使之在单元的边界面上对应的（或）取+1或-1值。这相当于把一个曲面六面体实单元“映射”为一个边长皆为2的正方体单元。这种一一对应的映射关系可写为：式中，、为i节点的坐标，为对应于i节点的形状函数。在自然坐标系内，各节点形状函数为：对于8个角点（i=1,2,8） 对于的边点（i=9,11,13,15）对于的边点（i=10,12,14,16）对于的边点（i

31、=17,18,19,20）式中、为节点i的自然坐标值。由图（b）可见，各节点的自然坐标值只可能为1、-1或0。对于等参单元，应该用同样的节点、同样的形函数以节点位移插值出单元内部位移。则此单元内沿x、y、z方向的位移u、v、w可表示为：其中、为i节点沿x、y、z方向的位移值。单元刚度矩阵的计算如下：三维变形状态下，应变函数为：这里的单元节点位移列阵为： （共有60个分量）单元应变矩阵可按节点分块表示为:其中第i个子矩阵为:单元刚度矩阵为： 其中为单元体积域，为材料的弹性系数矩阵，为60360的方阵，可按节点分块写为：式中第i行第j列的子矩阵为：按坐标变换应有：这里有而雅可比矩阵为可以求得上述雅

32、可比矩阵的各项元素，有综上所述，20节点六面体等参元在Jacobi行列式为常数条件下的完全（精确）高斯积分方案是333阶。2.4 计算题四4、图示一根直杆，长度2L，截面积A，弹性模量E。杆受到沿轴向的线分布力：q=cx。用2个3节点一维杆单元求解其位移、应力。要求推导详细的有限元求解列式，设置合理的参数值将求解结果绘制成曲线，并与精确解进行对比分析。（20分）解：考虑题图中的结构，由于要求用2个3节点一维杆单元求解，故每个单元的自由度为3。应用里兹法时分两个区域假设位移场，在单元内假设位移场应是二次多项式，即：（） （）其中，为待定参数（称为广义坐标）。为方便计算，采用题图所示的局部坐标系，

33、则有单元（）单元（）对单元来说，在上面的第一式中分别代入节点1、2、3的坐标，可得节点1、2、3在方向上的位移、，即：其中，。用局部坐标系表示为：其中，。于是有：解得：将上面求得的广义坐标代入原来假设的位移场，可将单元的位移函数表示成节点位移的函数，即整理得：将上式写成矩阵形式为其中，插值函数矩阵，或称形函数矩阵；单元的节点位移列阵（单元自由度）。同理，对于单元也有同样的位移场：其中，插值函数矩阵同上；单元的节点位移列阵。显然，整个杆上，由各单元上假设的位移场拼接而成的位移试探函数是连续的，只要我们记住，得到的就是全域可能位移场。这样的位移场已经把节点位移自由度作为广义坐标。在两个单元上分片进

34、行总势能计算：首先计算应变在单元内有式中，应变矩阵：所以一个单元内（单元和的形式相同）的应变能为：由于节点位移列阵与无关，故有：式中，单元刚度矩阵：载荷为，则载荷在两个单元内分别表达为：，带入外力功积分式，对三个单元分别计算外力功。第一个单元外力功为： 第二个单元外力功为： 系统总势能是两个单元的总应变能减去两个单元外力功，为了能够以矩阵形式相加，将单元势能矩阵表达式中的单元节点位移列阵用结构整体位移向量代替，单元刚度矩阵等与单元相关的矩（列）阵扩展成结构整体规模（55），则相加后系统总势能表达为：上述表达形式上的变换，不改变总势能的大小。上式简写为：式中，为结构整体位移向量，；为结构整体刚度

35、矩阵，；为结构整体载荷列阵，。应用驻值条件：，得到节点平衡方程：，即利用位移边界条件：，划去第一个方程，解出另外四个方程得到：于是，系统的结构整体位移向量为：由此得到的位移场在节点处为精确解，而一般位置上均为近似值（小于精确解）。根据几何方程求单元应变，再利用物理方程求单元应力。对单元有：对单元有：于是可得单元节点处的应力近似解如下：，该问题的位移和应力的精确解分别为：于是可得单元节点处的应力精确解如下：，位移和应力的计算结果与精确解的比较如图2.1所示：图2.1 受轴向力杆的精确结果和有限元结果由图中可以看出，由于在求得结构的节点位移后，要通过导数的运算来求单元的应变和应力，导致了精度的下降

36、，因此位移的近似程度比应力的近似程度更好。3 上机实验3.1 第一题3.1.1 实验题目1、图示一个简支梁平面应力模型。梁截面为矩形，高度h=160mm，长度L=1000mm，厚度t=10mm。上边承受均布压力q =1N/mm2，材料E=206GPa，=0.29。X方向正应力弹性力学理论解为：分别应用3节点三角形单元、4节点线性等参元（完全积分、减缩积分、非协调模式）、8节点二次等参元完全积分进行下列数值实验：1）用较粗单元网格求解梁中部应力分量的最大值和上下边法向应力分量，并对各单元计算精度进行比较分析；2）对粗网格下梁中部最大位移进行对比和分析；3）通过网格加密对比试验3节点三角形单元和8

37、节点二次等参元的收敛速度。总结出研究结论，撰写实验报告。（10分）3.1.2 实验目的通过对不同实体单元形状的选择（3节点三角形单元、4节点线性单元、8节点二次等参单元），和不同三维实体单元类型的选择（完全积分、减缩积分、非协调模式），以及网格的粗细，比较它们的结果，得出它们对计算精度的影响规律。3.1.3 建模概述(一) 采用3节点三角形单元进行分析1) 选取单元类型，PreprocessorAdd/Edit/DeleteSolid：Quad 4node 42OK (返回到Element Types窗口)Options K3: Plane Strs w/thk(带厚度的平面应力问题)Clos

38、e2) 定义材料参数，PreprocessorMaterial ModelsStructuralLinearElasticIsotropic: EX:206e11 (弹性模量)，PRXY: 0.29(泊松比) OK3) 定义实常数以确定平面问题的厚度， PreprocessorReal Constant Set No: 1 (第1号实常数), THK: 0.01(平面问题的厚度)Close4) 建立平面模型，PreprocessorModelingCreateAeriaRectangleBy Dimensions，设置x1、x2、y1、y2分别为0,1,0,0.16，建立平面模型。5) 划分网

39、格，PreprocessorMeshToolMeshing,选择三节点三角形单元，网格粗细可采用smart size进行控制。6) 对模型施加约束和外载，左边中间节点处施加X,Y方向的位移约束，右边中间节点处施加Y方向的位移约束，SolutionDefine LoadsApplyStructuralDisplacementOn Nodes鼠标选取节点分别施加相应的约束 OK；上表面施加Y方向的均布压力q =1N/mm，SolutionDefine LoadsApplyStructuralPressureOn Lines, VALUE值设置为1000000。7) 分析计算， Solution S

40、olve Current LS OK Close (Solution is done! )（二）采用4节点线性等参元进行分析 采用4节点线性等参元的整体建模步骤与（一）中大致相同，但有以下的不同点：1) 在完全积分和减缩积分时采用的单元名称为Solid：Quad 4node 182，即Plane182单元，此种单元可选择完全积分或减缩积分进行计算；非协调模式仍选用Solid：Quad 4node 42即Plane42单元，此种单元有协调模式与非协调模式的选项。2) 在设置单元参数时，Plane182单元的KEYOPT(1)=0时是完全积分，KEYOPT(1)=1时是减缩积分，Plane42单元

41、的KEYOPT(2)=1时是非协调模式。3) 在划分网格是需选用4节点单元。（三）采用8节点二次等参元进行分析采用8节点二次等参元的整体建模步骤与（一）中大致相同，但有以下的不同点：1) 采用单元名称为Solid：Brick 8node 185，即solid185单元。2) 使用8节点二次等参元中需建立块单元，PreprocessorModelingCreateVolumeBlockBy Dimensions，设置x1、x2、y1、y2,z1,z2分别为0,1,0,0.16，0,0.01建立块模型。3) 上表面Y方向的均布压力施加在面上而不是线上，约束的施加在边面的中心。3.1.4 计算结果分

42、析与结论（一）采用三节点三角形单元的计算结果1） 用较粗单元网格求解梁中部应力分量的最大值和上下边法向应力分量；图1：使用较粗三节点三角形单元网格对梁进行划分图2：三节点三角形单元梁中部应力分量图图3：三节点三角形单元上下边法向应力分量图由软件计算结果可知，使用较粗三节点三角形单元网格求得梁中部应力分量的最大值为，上下边法向应力分量分别为，。2）用密网格求解梁中部应力分量的最大值和上下边法向应力分量；图4：使用较细三节点三角形单元网格对梁进行划分图5：三节点三角形单元梁中部应力分量图图6：三节点三角形单元上下边法向应力分量图由软件计算结果可知，使用较细三节点三角形单元网格求得梁中部应力分量的最

43、大值为，上下边法向应力分量均为。3）对粗网格下梁中部最大位移图7：三节点三角形单元梁中部位移图由软件计算结果可知，使用较粗三节点三角形单元网格求得梁中部最大位移值为。（二）采用4节点线性等参元（完全积分）计算结果1) 用较粗单元网格求解梁中部应力分量的最大值和上下边法向应力分量；图8：使用较粗4节点线性等参元（完全积分）网格对梁进行划分图9：4节点线性等参元（完全积分）梁中部应力分量图图10：4节点线性等参元（完全积分）上下边法向应力分量图由软件计算结果可知，使用较粗4节点线性等参元（完全积分）网格求得梁中部应力分量的最大值为，上下边法向应力分量分别为，。2） 对粗网格下梁中部最大位移图11：

44、4节点线性等参元（完全积分）梁中部位移图由软件计算结果可知，使用较粗4节点线性等参元（完全积分）网格求得梁中部最大位移值为。（三）采用4节点线性等参元（减缩积分）计算结果1) 用较粗单元网格求解梁中部应力分量的最大值和上下边法向应力分量；图12：使用较粗4节点线性等参元（减缩积分）网格对梁进行划分图13：4节点线性等参元（减缩积分）梁中部应力分量图图14：4节点线性等参元（减缩积分）上下边法向应力分量图由软件计算结果可知，使用较粗4节点线性等参元（减缩积分）网格求得梁中部应力分量的最大值为，上下边法向应力分量分别为，。2） 对粗网格下梁中部最大位移图15：4节点线性等参元（减缩积分）梁中部位移图由软件计算结果可知，使用较粗4节点线性等参元（减缩积分）网格求得梁中部最大位移值为。（四）采用4节点线性等参元（非协调模式）计算结果1) 用较粗单元网格求解梁中部应力分量的最大值和上下边法向应力分量；图16：使用较粗4节点线性等参元（非协调模式）网格对梁进行划分图17：4节点线性等参元（非协调模式）梁中部应力分量图图18：4节点线性等参元（非协调模式）上下边法向应力分量图由软件计算结果可知，使用较粗4节点线性等参元（非协调模式）网格求得梁中部应力分量的最大值为，上下边法向应力分量分别为，