有限单元法基础第2版 课件 第7、8章 杆件结构力学问题、板壳结构力学问题

第7章杆件结构力学问题工程中的杆件结构重庆东水门大桥广州塔输电杆塔卢浮宫金字塔7.1引言1结构单元2杆件单元杆单元（桁架单元），梁单元，板单元，壳单元，索单元，膜单元等。轴力单元扭转单元梁弯曲单元拉伸压缩变形仅受拉伸索单元扭转变形弯曲变形3杆件结构单元桁架结构单元（平面和三维空间）轴力单元平面框架结构单元轴力单元+梁弯曲单元空间框架结构单元轴力单元+扭转单元+梁弯曲单元7.2轴力单元1基本方程受轴向载荷作用的等截面直杆几何关系本构关系平衡方程（A：杆的截面积）端部位移边界条件沿杆轴向的分布力b(x)载荷杆件的应变能势能（泛函）为应变能和外力势之和uk是作用有集中力2势能截面上的轴力的点的位移位移插值：xξξ=+1ξ=-1ξ=0x1xnxC坐标变换：Lagrange插值函数轴向应变3有限元平衡方程均布力等效结点载荷：集中力直接加在结点上（一般在集中力作用处设置结点）单元势能（泛函）：2结点Lagrange单元：f(x)=const坐标变换式等参变换等参单元受扭矩作用的等截面直杆几何关系（α

：扭率）本构关系（Iρ：截面极惯性矩）平衡方程或端部位移边界条件1基本方程7.3扭转单元载荷mt(x)扭转轴的应变能或作用有集中扭矩结点的扭转角度2势能势能单元的势能：单元的扭角插值单元的扭率3有限元平衡方程注意：前述扭转方程为自由扭转方程，对于非园截面杆，扭转后截面不再保持为平面，即发生翘曲，需应用约束扭转理论。计算截面任意点的剪应力2结点一次扭转单元的刚度矩阵7.4Euler-Bernoulli梁单元1Euler-Bernoulli平面梁基本理论Euler-Bernoulli梁变形特征变形前垂直于中面的截面变形后仍然垂直于中面，且保持为平面。忽略层间剪切变形影响。Euler-Bernoulli梁微元体及应力分布端部位移边界条件载荷q(x)m(x)几何关系（κ

：梁中面变形后的曲率）应力应变关系（I：截面惯性矩）平衡方程或端部位移边界条件基本方程应变能和势能载荷q(x)m(x)22结点Euler-Bernoulli梁单元2结点Euler-Bernoulli梁单元Hermite插值函数有限元离散

（Hermite插值函数）广义应变矩阵

由单元的势能泛函单元刚度矩阵计算单元中任意一点的值7.5Timoshenko梁单元1Timoshenko平面梁基本理论Timoshinko梁截面的转动关系梁的高度相对于跨度不太小，层间剪切变形不能忽略，横向剪切力产生的剪切变形引起的附加挠度不能忽略，因此原垂直于中面的截面变形后不再垂直于中面，但假设该截面仍为平面。正应力分布剪应力实际分布剪应力等效均匀分布截面上均匀分布的剪应变，等于中面处的剪应变剪力

等效均布剪应力校正因子实际上，剪应变和剪应力在截面上按抛物线分布，变形后截面不再是平面，引入校正因子k

。按实际剪应力分布的剪应变能截面的等效剪切应变能矩形截面：k=5/6圆截面：k=9/10令梁的应变能：弯曲应变能+剪切应变能梁的势能泛函等截面梁忽略位能泛函中载荷贡项，并除以，得到如下泛函

该式反映了弯曲和剪切变形对刚度的相对贡献当梁的高度很小时会很大如：矩形截面梁（罚参数）讨论：梁的高度很小时剪切变形的约束条件得以满足。Euler-Bernoulli梁2在Euler-Bernoulli梁单元中引入剪切应变的梁单元1）位移插值法向位移：弯曲引起的法向位移剪切变形引起的附加法向位移2结点单元的单元结点位移与不考虑剪切变形梁相同的Hermite插值2结点Lagrange插值Hermite插值函数Lagrange插值函数2）有限元平衡方程

利用位移插值关系，由可得（与不考虑剪切变形梁弯曲问题有限元方程一致）单元的势能单元平衡方程可以进一步简化，由弹性关系可得由平衡条件有几何关系由上两个方程可以得到代入Euler-Bernoulli梁单元的平衡方程剪切变形的影响可以忽略单元刚度矩阵简化为Euler梁单元刚度矩阵剪切变形的影响通过系数β反映。当梁的高度很小时，忽略横向剪切变形的影响，可以理解为剪切模量G趋于无穷大。讨论：0挠度和界面转动独立插值Ni：Lagrange插值函数。C0连续单元，便于推广到板壳。32结点Timoshenko梁单元单元结点处的集中力和集中力矩列向量由得若q和m为均匀分布当梁的高度变小时，剪切应变为零常数项和一次项应分别为零单元内：梁不能发生弯曲变形——单元剪切锁死（shearlocking）原因：w，θ

独立插值，同阶，导致不同阶。解决途径：计算剪切应变时，使预先保持同阶。剪切自锁现象具体方案：i）减缩积分（reducedintegration）采用比精确积分要求少的积分点数。如：2结点单元，精确积分：2点，减缩积分：1点。采用减缩积分，以积分的一个θ值代替单元内的线性变化，同阶，消除自锁现象。ii）假设剪切应变（assumedshearstrains）收敛性问题2结点单元不能描述纯弯曲状态。故推荐采用3结点或4结点Timoshenko梁单元。描述常弯曲状态的位移模式：对于剪切变形可以忽略的情况，应尽量采用Euler-Bernoulli梁单元。7.6平面桁架结构二维桁架结构每一根杆都是二力杆每一根杆件的取向不同局部坐标下的轴力单元局部坐标和整体坐标下的典型平面桁架单元局部坐标下单元的结点位移列向量整体坐标系下单元结点位移列向量结点位移变换1结点位移坐标变换2刚度矩阵坐标变换局部坐标系整体坐标系应变能是标量，与坐标系无关轴向分布载荷分解到整体坐标系的两个方向结点集中载荷分解到整体坐标系的两个方向单元方程的组集单元等效结点载荷单元结点载荷施加在整体平衡方程右端列向量中对应的自由度方向，无需进行坐标转换。3刚度矩阵及结点载荷的组集4单元应变和应力计算算例E=29.5X106Pa，A=1.0m2。（1）计算每一个单元的刚度矩阵；（2）组装得到结构的刚度矩阵；（3）计算结点位移；（4）计算每个单元的应力；（5）计算支座反力。结点坐标结点结点坐标（m）xy1002400340304030单元定义单元局部结点1局部结点2112232313443单元长度和方向余弦单元le/m方向余弦1401.00.02300.0-1.03500.80.64401.00.0整体自由度组集结构得到整体刚度矩阵注意各单元结点自由度在整体坐标系中的对应关系结构的右端载荷列向量引入位移边界条件，高斯消去法求解得到结点位移：计算单元应力利用整体平衡方程可以求支座的反力：注意计算过程中单位的统一。局部坐标和整体坐标下的典型三维桁架单元局部坐标下单元的结点位移列向量整体坐标系下单元结点位移列向量结点位移变换单元刚度矩阵变换7.7三维桁架结构7.8平面框架结构平面框架杆件受力特点：杆件受轴力和弯矩共同作用局部坐标系内梁单元结点位移：单元刚度矩阵：轴力单元刚度子矩阵弯曲单元刚度子矩阵2结点Euler梁单元刚度矩阵1单元矩阵及其变换局部坐标系下的二维梁单元轴力结点载荷子矩阵弯曲结点载荷子矩阵结点等效载荷矩阵转换局部坐标系下结点位移列向量整体坐标系下结点位移列向量局部坐标和整体坐标下的梁单元线位移和转角的变换关系矩阵变换关系2铰结点的处理铰结点处（a）各杆线位移相同（b）截面转动不相同（c）铰结点不承受弯矩单元②的铰接端，转动属于内部自由度，可以在单元层次上凝聚掉。

——自由度释放铰结点由第二个方程代入第一个方程为编程方便，K*保留原来的阶次，将对应于释放自由度的相应刚度系数置为零。实现方法：单元内需凝聚掉的自由度3单元截面应力计算轴向位移轴力单元几何关系计算轴向应变轴向应力轴向变形弯曲变形截面上的正应力需要将轴向力引起的正应力和弯曲引起的正应力进行叠加空间杆件系统受力特点：杆件受轴力、弯矩和扭矩共同作用。7.9三维框架结构三维梁变形及受力特征三维框架结构三维空间2结点梁单元力和线位移力矩和截面转动1单元矩阵及其变换局部坐标系内空间杆单元的特性矩阵矩阵转换空间梁局部坐标和整体坐标的转换2单元截面应力计算轴向位移轴力单元几何关系计算轴向应变轴向应力轴向变形（与平面框架类似）弯曲变形（需计算两个方向弯矩对应的应力）结点位移扭转变形截面上的应力需要将轴向变形、弯曲和扭转引起的应力进行叠加。7.10杆件结构弹性稳定性分析1梁单元的几何刚度矩阵梁单元梁的微元变形受轴向力和弯矩作用的梁单元梁的微元，其在变形前的长度为dx，发生弯曲变形后的长度改变为利用二项式定理展开忽略高阶项，得到挠度在中性层上引起的附加应变由于轴向力作用，对应于轴向位移u的应变梁厚度方向上任意一点（坐标为z）的弯曲应变梁单元内任意一点的应变（三项之和）梁单元的应变能轴力与轴向位移之间存在如下关系忽略高阶项梁的应变能可分为轴向应变能和弯曲应变能杆件同时承受轴向力和弯曲载荷作用时，单元的结点位移包括轴向位移、挠度和转角。轴向位移和弯曲位移可以分别插值。以2结点单元为例梁单元的弯曲刚度矩阵几何刚度矩阵几何刚度矩阵是由于单元中轴力Fx的作用引起的弯曲刚度的增加。当杆件承受轴向拉力作用时，单元刚度增加，受压缩时刚度减小，且其大小与材料的物理常数无关，与单元的几何尺寸和内力有关。也称为初应力刚度矩阵。2杆件结构失稳临界载荷杆件系统的弯曲应变能为所有Me个单元应变能之和在横向力p的作用下系统的势能为单元几何刚度矩阵与单元的轴向力成正比，如果轴向载荷按比例增大λ倍如果轴力很小，可以忽略，结点位移与p成正比。当p=0时，结点位移也为零。当轴向载荷达到临界载荷时，系统发生弹性失稳。此时即使横向载荷P=0，结点位移也不为零有非零解的条件求解该特征方程得到其最小特征值λcr结构失稳的临界载荷将特征值λcr代入计算得到对应的特征向量，即系统的失稳模态。失稳模态反映了结构发生失稳时的变形形态稳定性分析分两步：（1）计算结构的内力分布（2）计算结构的失稳临界载荷7.11算例：受压框架强度及稳定性计算工字型梁：截面的宽度b=0.2m，高度h=0.3m，厚t=0.02m。材料为Q420，其杨氏弹性模量为200GPa，泊松比为0.3，屈服极限为420MPa。采用空间梁单元模拟计算该结构的变形和应力，并分析其强度。受压均布压力的框架空间梁有限元模型变形及位移分布梁截面最大Mises应力分布受压框架结构静力分析受集中力作用的框架结构第1阶屈曲模态，特征值：3.041×106工字型梁：截面的宽度b=20mm，高度h=30mm，厚t=2mm。材料为Q420，其杨氏弹性模量为200GPa，泊松比为0.3，屈服极限威420MPa。计算其失稳载荷稳定性分析时假设载荷的大小为F0=1.0N7.12小结轴力单元Lagrange插值（2结点Lagrange单元）扭转单元Lagrange插值（2结点一次扭转单元）Euler-Bernoulli梁单元2结点Euler-Bernoulli梁单元Hermite插值函数Timoshenko梁单元在Euler-Bernoulli梁单元中引入剪切应变的2结点梁单元Hermite插值Lagrange插值2结点Timoshenko梁单元平面桁架结构坐标变换结点位移单元刚度单元应变和应力计算三维桁架结构分解到整体坐标系的两个方向平面框架结构单元矩阵及其变换（2结点Euler梁单元）铰结点的处理：自由度凝聚单元截面应力计算三维框架结构单元矩阵及其变换单元截面应力计算杆件结构弹性稳定性分析几何刚度矩阵杆件结构失稳临界载荷初应力刚度矩阵第8章板壳结构力学问题超参数单元薄板单元厚板单元平面壳体单元8.1引言轴对称壳单元薄板结构弹性稳定性8.2薄板单元1薄板弯曲基本理论受横向力和弯曲载荷作用的平板基本假设中面（z=0）上的点面内位移为零垂直于中面的法线上的点具有相同的竖向位移垂直于中面方向的正应力可以忽略不计变形前正交于板中面的直线变形后仍正交于中面且长度不变板中面外的其它点的面内位移直法线假设薄板中任意一点的位移薄板的变形特征板中任意一点的应变广义应变（曲率和扭率表达）横向剪切应变为零面内应力分量薄板微元体及应力分布广义应力应变关系（本构关系）板的抗弯刚度板的抗弯刚度矩阵平衡方程基本方程几何方程本构方程平衡方程板弯曲问题的三种边界条件外力边界n:法向；s:切向位移边界混合边界条件或弹性薄板的弯曲应变能板弯曲时为平面应力状态，其应力和应变关系为系统的势能泛函离散化：2基于薄板理论的平板单元结点变量插值函数绕y轴的转角：绕x的转角：w结点位移：对板的中面进行离散，中面以外其它点的位移由中面挠度计算得到。矩形4结点单元单元位移插值模式：代入结点位移和坐标：利用坐标变换xo和yo为单元中心点的整体坐标仅作用有竖向均布载荷q收敛性检查a）完备性检查刚体位移：z向刚体位移绕y

轴刚体转动绕x

轴刚体转动常应变：满足完备性要求。b）连续性检查w连续；一般不连续。非协调单元讨论：a）为非协调单元。b）能通过分片试验，能收敛于精确解，但收敛不是单调的。c）不能推广到一般的四边形单元。三角形3结点单元1233结点，9个参数每个结点自由度：若用完备三次项，应包含10项用面积坐标构造位移函数可以证明：当c=1/2时，满足常应变要求；w协调；不协调。非协调单元8.3厚板单元1厚板弯曲基本理论基本假设变形前正交于板中面的直线变形后仍为直线，其它三个假设与薄板相同。厚板界面的转动关系中面位移：中面外其它点位移中面法线的转角板中任意一点的应变层间剪切应变曲率和扭率（广义应变）板的面内应力状态为平面应力状态厚板微元体及应力分布平衡方程k为截面剪切校正因子抗弯刚度通过修正的剪应变能等于按实际剪应力和剪应变计算得到的应变能相等，可得到矩形截面的k=5/6。抛物线分布广义应力应变关系（本构方程）平衡方程板的应变能弯曲应变能剪切应变能板的势能泛函第一项：弯曲应变能第二项：剪切应变能后几项：外力势其中挠度和转角均为独立变量，由此构造的有限元格式适用于较厚的平板弯曲问题。该类单元的位移和转动各自独立插值，称为

Mindlin板单元。板很薄时，忽略层间剪切变形可以理解为剪切刚度罚函数法第二项可以理解为利用罚函数法引入薄板约束条件插值：C0连续函数，二维C0型插值函数均可。2Mindlin板单元由驻值条件剪切锁死和零能模式当板很薄时，为避免剪切锁死，Ks需奇异。采用减缩积分可能致使K奇异，使解答中包含除刚体运动之外的且对变形能无贡献的变形模式——零能模式a）Midlin板单元保证K非奇异的必要条件Me：单元数；：高斯积分点数；：弯曲应变和剪切应变分量数；N：系统独立自由度数。b）保证Ks奇异性的充分条件积分方案的选择分别采用不同阶次的积分方案。应用中的处理实际应用和通用程序中，普遍采用减缩积分方案。零能位移模式的出现：板较厚且单元较少，边界约束较少时容易出现。较厚板：精确积分，不会出现零能模式和锁死。应用中薄板分析时，推荐采用4结点和9结点Lagrange单元。8.4平面壳体单元1弹性薄壳基本假设（Kirchhoff-Love假设）直法线假设在变形前正交于壳体中面的直线，变形后仍然为正交于中面的直法线段，且保持长度不变。切平面应力假设正应力比应力分量小得多，相对而言可以忽略。模拟一般三维空间壳结构的平面壳体单元用三角形单元组成曲壳用矩形单元组成曲壳用三角形和矩形组合代替曲壳壳体与平板的异同a）相同点a）厚度方向的尺度比其它两个方向小得多b）均使用Kirchhoff假设b）不同点平板：u，v位移[平面应力（薄膜应力状态）]与弯曲状态不耦合。壳体：u，v，w同时发生，弯曲状态与薄膜状态相互耦合。1局部坐标系下的平板单元平面壳元＝平面应力单元＋平板弯曲单元结点位移：结点力：1）平面应力状态（薄膜应力状态）三角形平面薄壳单元的结点位移和结点力2）平板弯曲状态在局部坐标系下在整体坐标系下2向整体坐标的转换局部坐标系下整体坐标系下1）结点位移变换L：正交矩阵平面壳单元局部坐标系和整体坐标系2）单元结点位移变换3）单元刚度矩阵变换4）单元载荷向量的变换组集刚度矩阵后特殊情况的处理：若汇交于结点i的各单元在同一平面内，处理方法：（1）在局部坐标系建立该点的平衡方程，删去对应的方程（2）在此点上给出任意的刚度系数，不会影响计算结果。8.5超参数壳体单元1几何形状典型的厚壳单元三维实体单元用于壳体结构时，考虑到壳体在厚度方向的尺寸相对于其它方向的尺寸很小，可以在厚度方向只设置两个结点。在引入壳体理论的基本假设之后，可以蜕化为超参数壳单元。壳体单元中任意一点的坐标上式可改写为：如果V3i与法线同向，i点处的壳厚度V3i的单位向量的方向余弦：单元结点数：2n底结点到顶结点的向量2位移函数单元内任意点的位移：是v3i绕v3i和v3i的旋转。与v3i垂直并正交的单位向量超参单元的局部和整体坐标及位移位移插值表达式单元结点位移数：nX5单元几何形状参数：nX2X3超参单元i为x坐标方向的单位矢量i=[100]T3应变和应力根据壳体理论，局部坐标系下计算壳体应变能涉及的应变为局部坐标系和整体坐标系下的应变转换关系k＝5/6,考虑剪切边沿厚度方向不均匀分布的影响而引入的修正。4单元平衡方程超参壳元在引入一定的几何假设后和位移及转动各自由度独立插值的壳元等价。在计算单元刚度阵时，同样需保正K非奇异，Ks奇异。可采用Mindlin板单元中讨论的各种方法解决“锁死”问题。实用中普遍采用减缩积分方案。单元刚度矩阵可以分解为薄膜-弯曲和横向剪切两部分1轴对称壳的基本公式1）壳体中任一点的应变a）中面广义应变8.6轴对称壳单元b）中面外任意一点的应变（z：离中面的距离）2）广义内力分量和应力分量（假设应力沿壳厚度方向线性变化）3）广义应力和广义应变之间的弹性关系将变形分解为薄膜状态和弯曲状态广义应力－应变关系4）载荷和边界条件对称时的简化周向位移u，w仅为s的函数，与θ

无关。壳体的应变能势能考虑剪切变形的轴对称壳体理论的基本公式应变能2基于薄壳理论的轴对称单元结点位移单元结点位移单元中任意一点在局部坐标系下的位移模式将结点1，2处的代入，可得le：截锥单元经线的长度。轴对称截锥壳单元位移表达式局部坐标系下位移与整体坐标系下位移的变换关系位移插值表达式2结点直线单元，代入利用最小势能原理或加权残值法可得平衡方程单元刚度矩阵：单元等效结点载荷：用一维高斯数值积分，通常3～4点积分可达足够精度。优点：表达式简单，精度较好。应用较广。缺点：用直线代替曲线，有时会产生较大误差。任有必要研究曲边单元。3位移和转动各自独立插值的轴对称壳单元截面转动β是对立函数，在整体坐标系下构造插值函数代入考虑剪切变形影响的轴对称壳单元代入泛函并利用极值条件可得可将Ke分解轴向分布载荷单元的结点载荷径向分布载荷分布力矩讨论：类似于Timoshenko梁和Mindlin板要求：K非奇异，奇异。可采用：减缩积分、选择积分、假设剪切应变法等。轴对称3结点二次曲边壳单元坐标和总体坐标系下的位移及转动采用相同的插值表达式轴对称二次曲边单元最后可得注意：（1）结点3的位