动力学问题的有限元法

第七章动力学问题的有限元法结构动力学是研究动载荷作用下结构动力反应规律的学科，讨论结构在动力荷载作用下反应的分析方法，寻找结构固有动力特性、动力荷载和结构反应三者间的相互关系。研究结构在动力荷载作用下的反应规律，能够为结构的动力可靠性（安全、舒适）设计提供依据。前面介绍的静力学问题的研究对象是受不随时间变化的载荷作用。而动力学问题的对象受随时间而变的载荷的作用，从而使在结构中产生的位移、速度、应力和应变都随时间而变。当结构受随时间变化的载荷作用，且这种载荷的作用对结构的变形和应力的产生起主要作用，以致影响设备的安全性，或舒适性。这时就要进行动力学分析，充分认识其规律性，从设计阶段就抑制这种不利状况的发生。例如，有时虽然动载荷不大，但结构在交变力的作用下，其某些固有频率与激励力的作用频率相接近时，就会引起很大的振动、变形或应力，这时，就必须对结构作动力学分析。又如，要利用结构在周期性作用力驱动下的定向振动，例如利用这种运动输送产品，这时，就必须巧妙地设计结构，使其具有某些与激励频率一致的固有频率，并且使结构对激励具有适当的响应能力。总之，不管是利用振动，还是抑制振动，都需要进行结构动力学分析。当前结构动力学的研究内容有三类。第一类问题：反应分析（结构动力计算），第二类问题：参数（或称系统）识别，第三类问题：荷载识别。第一类问题是已知系统动态特性和动载荷作用部位及大小，求出系统的响应——随时间变化的位移，速度，加速度和应力等。第二类问题是已知系统的输入输出特性，分析系统固有的动态特性，结构模态分析就属于这一类问题。第三类问题是在已知系统动态特性的条件下,通过测量系统的响应,或由响应准则预先给出响应要求,以此识别对响应的外载荷。三类结构动力学研究内容的载荷、结构和响应之间的关系如图7-1所示。动载荷种类大致分类如图7-2所示。图7-1结构动力学研究的内容图7-2动载荷种类1本章主要介绍结构动力学分析的基础知识，并主要介绍系统固有特性的有限元分析方法——有限元模态分析。主要知识点：1.导出有限元的结构动力学方程和模态分析有限元方程。2.介绍模态分析方程中质量矩阵与阻尼矩阵的计算。3.把计算结构的固有频率与振型的问题归纳为一个求特征值的问题。7.1．结构动力学方程及有限元方程系统的受迫振动的微分方程：m()xcxkxRt++=其中：——质量，——阻尼系数，——弹性系数。mck,,xxx分别为加速度，速度和位移。1.单元的动力方程：（用“虚功原理”求解）结构动力学有限元分析与静力学的有限元分析一样，首先要建立单元的有限元方程。为求解起见，要引入达朗贝尔原理。质点的达朗贝尔原理可表述为：当非自由质点运动时，主动力，约束反力和惯性力构成一个动平衡力系。惯性力是一种虚构的力，由质点的加速度引起。mFNS即：0fns++=现考察单元的动力学方程的建立。1).准备工作：将位移函数{}δ表达成近似函数。单元内任一点位移{}[]{}[]{}(,,)()eNNxyzδδδ==et，是单元内位移的近似函数,这儿[]N只是位置的函数，与时间无关。而{}eδ与时间相关。[][]{}{}(,,)()eeNNxyztδδ⎧=⎪⎨=⎪⎩{}[]{}{}[]{}eeNNδδδδ==2).分析单元上的作用力。图7-3单元分析主动力：集中力，面力，体积力：{}psRRRRV=++阻尼力：与速度相关，并与速度方向相反，两者的比例关系为阻尼系数：{}{}[]{}ceRCCNδδ=−=−2惯性力：与加速度，质量相关：单位体积惯性力{}{}[]{}meRNρδρδ=−=−,阻尼力、惯性力与主动力方向相反。ρ——密度，单位体积质量。3).建立动力学平衡方程：利用虚功原理和达朗贝尔原理，系统处于动平衡，首先将载荷移置到节点上去。利用虚功原理。外力做的虚功：{}{}{}[]{}{}{}{}[]{}{}{}{}[]{}{}{}{}[]{}{}{}{}[]{}11****11******TeTTppTeTTSsSsTeTTVVVVTeTTVCVCTeTTVmVmRNRRdSNRdSWRdVNRdVNRdVRdVNRdVδδδδδδδδδδδ⎧=⎪⎪=⎪⎪⎪==⎨⎪⎪=⎪⎪=⎪⎩∫∫∑∫∫∫∫∫∫集中力：面力：体积力：阻尼力：惯性力：RdV内力做的虚功：{}{}{}[][][]{}**TeTTeVVUdVBDBδεσδδ==∫∫dV0虚功原理：()WUWUδδδ=⇒−=，最小势能原理根据虚功原理，应有WUδδ={}[]{}{}[]{}{}[]{}{}[]{}{}[]{}{}[][][]{}1***1***eTeTeTTTpSsVVeTeTeTTTTeVCVmVNRNRdSNRdVNRdVNRdVBDBdVδδδδδδ+++++=∫∫∫∫∫δ由{}*eTδ的任意性，可约去，且由阻尼力，惯性力的定义可知：{}{}[]{}{}{}[]{}cmeeRCCNRNδδρδρδ⎧=−=−⎪⎨⎪=−=−⎩，代入上式，有[]{}[]{}[]{}[][][]{}[][]{}[][]{}11TTTpSsVVeeTTTeVVVNRNRdSNRdVBDBdVCNNdVNNdVδδρ++=++∫∫∫∫∫δ=eRi令：[]{}[]{}[]{}{}{}[][][]{}[][]{}[][]11TTTepSsVVTeVTeVTeVNRNRdSNRdVRKBDBdVcCNNdVmNNdVρ⎧++⎪⎪=⎪⎨=⎪⎪⎪=⎩∫∫∫∫∫上式可表达为：{}，这就是单元的动力学有限元方程。{}{}{}{}{}{}eeemcKδδδ++=2.结构整体动力学有限元方程与静力学的集成原理相同：根据①局部单元节点编号与整体的节点编号的对应关系，利用{}，即一个节点上的节点力是该节点所在的所有单元的相应节点（同这个节点）上的节点力的总和来集成。{}eeiiRR=∑3②位移唯一性。得{}{}{}{}{}{}{}MCKδδδ++=R{}{}eMm=∑，{}，{}{}eCc=∑{}eKk=∑其中，当{}0C=，称为无阻尼受迫振动，此时{}{}{}{}{}MKRδδ+=当{}{}0,0CR==，称为无阻尼自由振动，此时{}{}{}{}0MKδδ+=7.2单元质量矩阵和阻尼矩阵的表达式在动力学有限元方程中除刚度矩阵之外，还有两个重要的系数矩阵：即质量矩阵和阻尼矩阵。这一节中将介绍单元质量矩阵和单元阻尼矩阵的表达式。一.单元质量矩阵单元质量矩阵有二种表达形式：一致质量矩阵和集中（堆积）质量矩阵。1.一致质量矩阵：用公式{}计算的单元质量矩阵，称为一致质量矩阵，是因其中的[[][]TeVmNNρ=∫dV]N就是位移模式形函数[]N而得名。或者说，推导质量矩阵采用的形函数与采用推导单元刚度矩阵时采用的形函数是一致的。以平面应力问题为例：三角单元，[]ijmNNININI⎡⎤=⎣⎦m（位移模式中的[]000000ijmijNNNNNNN⎡⎤=⎢⎥⎢⎥⎣⎦[]N）图7.3三节点三角形单元则三角形单元的一致质量矩阵[]em：（以单位厚度计）ii22jj2j0N0N0N00N00N00jmijmijmeTijNNNNNNNNmNNdxdyNNρρ==∫∫j2mm2mmNN0N000N0N0mijmijmdxdyNNNNNNN⎤⎢⎥⎥⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦∫∫由于m，利用数值积分：2iN⎡ii00N0N⎢⎢[],,iijjmNLNLNL===4!!!2(2)abcijmabcLLLdxdyabc∴=!+++∫∫Δ，有22!24!62(112)!12iijNdxdyNNdxdyΔΔ⎧==⎪⎪⎨ΔΔ⎪==⎪++⎩∫∫∫∫从而得到[]20210212010210102010102emρ⎡⎤⎢⎥⎢⎥⎢⎥Δ=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦对称2.集中质量矩阵：将单元质量假想地集中到节点上去，也就是说只在节点上有质量。这样，一个节点的加速度就不会影响其它节点的惯性力（相互之间不受影响）具体处理方法是将质量平均分到所有节点上去。例如三角形三个节点单元，每个节点分到1/3的质量，即1m3itρ=⋅Δ⋅,t为厚度，Δ为三角形面积。[]111111111331111ewmtρ⎡⎤⎡⎢⎥⎢⎢⎥⎢⎥==Δ⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎤⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎦资料表明，集中质量矩阵计算出来的频率与一致质量矩阵的计算结果相差无几。在单元相同的情况下，集中质量矩阵计算的结果低于一致质量矩阵的计算频率，而振型则是一致质量矩阵的计算结果精确。由于集中质量矩阵式对角线矩阵，所以计算中常采用。3.附加质量当有附加的质量在节点上时，要在计算整体质量矩阵时，加上这些附加质量[]cm即[][][]ecMmm=+∑二.单元阻尼矩阵根据阻尼的成因可以分为粘滞阻尼和结构阻尼。1．粘滞阻尼（比例阻尼）前面得到的阻尼矩阵[]是将阻尼看作是正比于节点的运动速度，称为粘滞阻尼。它消耗振动的动能（能量）。从[]可知，这个阻尼矩阵与质量矩阵成正比。[][]TeCCNNd=∫VdVeTmNNρ=∫[][]eeCCmα=，/Ccαρ=如果,cρ均为常数的话。所以称为比例阻尼。（相对于质量矩阵而言）[][][]NTeCCNd=∫V2．结构阻尼比例于结构的变形速度{}ε，由于材料内部的摩擦引起的，称为结构阻尼。这时的阻尼力可以简5化为{}{}cRDμε=，这样可得到单元的结构阻尼矩阵使其具有[]的形式，所以结构阻尼正比于刚度矩阵：[][][]BeTkCCBDdε=∫V[][]eekCkα=3．Rayleigh阻尼在实际分析中，精确决定阻尼矩阵是相当困难的，通常允许将整体阻尼矩阵简化为质量矩阵[]M和刚度矩阵[]K的线性组合，这样可用来解耦有阻尼的振动方程。[][][]M+KCαβ=，称为Rayleigh阻尼。其中j22i222(-)2(-)ijiijjijjijiξωξωαωωωωξωξωβωω⎧=⎪−⎪⎨⎪=⎪−⎩其中22iiiαβωξω+=称为结构的第i阶振型的振型阻尼比。2iαβω+称为结构的第i个校准振型的振型阻尼常数。工程中一般取iξ=0.02~0.24，它与结构的类型，材料的能量消耗特性及振型有关。可通过实测，从系统的固有频率和辐频特性曲线获得。质疑：相关内容的矛盾描述：在王勖成的书P476有如下字样：“,αβ是不依赖于频率的常数”。但是，在同一本书的P487，又出现下列字样：“如果根据实验或相近似结构的资料，已知二个振型的阻尼比iξ，jξ，可以得到j22i222(-)2(-)iiiijjijjijiξωξωαωωωωξωξωβωω⎧=⎪−⎪⎨⎪=⎪−⎩”7.3．结构动力学特性的计算求解结构的固有频率和对应于每一阶固有频率的振型，这是动力分析的基本内容，称为模态分析。振型是指在某一个自振（固有）频率下，结构各节点振动值的集合，即节点振幅的相对大小。由于阻尼对结构的自振频率和振型的影响不大，所以，一般只考虑无阻尼的自由振动。一.模态分析的有限元方程在[]{}[]{}[]{}{}MPCKδδδ++=中，令{}[]P0,0C==，得模态分析的有限元方程：6[]{}[]{}M0Kδδ+=（注：联系到温度场的瞬态热传导方程：[]{}[]{}{}3PKtKt+=，{}P0=——无内热源。有人作比对分析，对温度场的分析采用所谓的“热模态分析”，但要注意，两个方程只是“形似”，它们的解有区别，热是一阶方程，振动是二阶方程（相对于时间）热是指数型解，振动是正弦函数分布。）二。数学预备知识求解固有频率和振型，数学上涉及到矩阵的特征值问题。这儿做一个有关矩阵的特征值，特征向量和特征方程方面的知识的扼要回顾，以便将我们的主要精力放到模态分析的物理意义上去。1．特征值和特征向量对n阶方程：和n维非零列向量，如果有一个数1111nnnnnnaaAaa×⎛⎞⎜=⎜⎜⎟⎝⎠…⎟⎟12(,,)Tnαααα=λ使得矩阵nnA×和维非零列向量nα有Aααλ=，则称λ为矩阵A的特征值（特征根），向量α为矩阵A的特征值λ所对应的特征向量。()0AAIαλλα=⇒−=，∵α为非零列向量，∴若等式成立，应有0AIλ−=α2.特征矩阵11121212211