第9章(动力学fem)

4.1公式号、图号第九章机械结构的动力有限元本章介绍介绍机械结构动力学的有限元建模方法，然后利用有限元求解机械结构系统的固有动特性，对典型工况下的动态响应进行分析，并介绍了几个典型应用实例。8.1动力学方程的建立在动力学问题中，位移、速度、应变、应力和载荷都是时变的。弹性结构的动力学问题的基本方程如下平衡方程,,,0ijjiittitfuu(4.70)几何方程,,12ijijjiuu(4.71)物理方程ijijklklD(4.72)边界条件iiuu,(在Su边界上)、ijjinT,在S边界上(4.73)以及初始条件,,,,,0,,,,,,0,,iiitituxyzuxyzuxyzuxyz。式中，是密度，是阻尼系数，,,,ittituu分别是iu对t的二阶导数和一阶导数，即加速度和速度。,,,ittituu分别表示惯性力和阻尼力，作为体积力的一部分出现在平衡方程中。只对空间域进行离散，单元内位移,,uvw的插值函数分别为111,,,,,,,,,,,,,,,niiiniiiniiiuxyztNxyzutvxyztNxyzvtwxyztNxyzwt(4.74)机械产品动态优化设计理论与方法4.2其中，u,,,,,,,,,uxyztvxyztwxyzt。平衡方程式(4.71)以及力边界条件的等效积分形式的Galerkin提法如下0,,,d()d0iijjiiuitiijjiVSufuuVunTS(4.75)对上式第一项,diijjVuV进行分部积分，并引入物理方程，则由上式可以得到0,,dddijijklkliiuiitViiiiVsDuuuuVufVuTS(4.76)将位移空间离散后的表达式(4.74)代入上式，注意这里的321,,uuu相当于u,v,w，并根据变分原理，最终得到结构系统的动力学方程()()()()ttttMa+Ca+ka=Q(4.77)其中，()ta和()ta分别是结点加速度向量和速度向量，,,MCK和()tQ分别是结构系统的质量矩阵、阻尼矩阵、刚度矩阵和结点载荷向量，分别由各自的单元矩阵和向量集成eeMM，eeKK,eeCC,eeQQ(4.78)其中，deeTVVMNN,deeTVVKBDB,deeTVVCNN,ddeeeTTVSVSQNfNT分别是单元的质量矩阵、刚度矩阵、阻尼矩阵和单元载荷向量。4.3.2质量矩阵和阻尼矩阵式(4.79)所表达的单元质量矩阵deeTVVMNN（4.79）称为协调质量矩阵或一致质量矩阵，这是因为公式导出时采用了和导出刚度矩阵一致的位移插值函数，其质量分布是按照实际分布情况考虑的。式中的形函数矩阵N随单元类型而异。例如，对于平面三角形常应变单元，其位移插值函数为123[]NNNNI,其中，I是22单位矩阵，()/2,1,2,3iiiiNabxcyAi，系数,,iiiabc见前文，A是三角形单元面积。利用上式可以算得平面三角形常应变单元的协调质量矩阵具体表达式为4.320101002010110201001020112101020010102etAM如果是等参元，设其形函数矩阵为(,,)N，单元协调质量矩阵为111111||dddeTMNNJ(4.80)为计算方便，还经常采用集中质量矩阵，即假设单元的质量集中在结点上，即把每个单元的分布质量按静力学中的平行分解原理，平均分配到每个结点上，形成一个阶数等于单元自由度数的对角线质量矩阵，而非主角线元素均为0。例如，对于平面三角形常应变单元，其集中质量矩阵为1000000100000010000001003000010000001etAM(4.81)对于式(4.78)所表达的单元阻尼矩阵deeTVVCNN（4.82）称为协调阻尼矩阵，它是假设阻尼力正比于质点运动速度的结果。通常将介质阻尼简化为这种情况，这时单元阻尼矩阵比例于单元质量矩阵。除此之外，还有比例于应变速度的阻尼，例如由于材料内摩擦引起的结构阻尼通常可简化为这种情况，这时的阻尼力可表示成Dε，可以得到单位阻尼矩阵deeTVVCBDB（4.83）此单元阻尼矩阵比例于单元刚度矩阵。由于系统的固有振型对于M和K具有正交性，因此固有振型对于比例于M和K的阻尼矩阵C也具有正交性。所以这种阻尼矩阵称为比例阻尼矩阵或振型阻尼矩阵。通常允许将实际结构的阻尼矩阵简化为M和K的线性组合。即CMK（4.84）其中,是不依赖于频率的常数。这种振型阻尼称为Rayleigh阻尼。机械产品动态优化设计理论与方法4.48.2利用有限元分析机械结构的固有动特性4.4.1机械结构的固有频率和固有振型不考虑阻尼影响的机械结构的自由振动方程是0ttΜaKa(4.85)它的解可以假设为以下形式0sin()ttaφ(4.86)其中,φ是n阶向量，是向量φ振动的频率，t是时间变量，0t是由初始条件确定的时间常数。将式(4.86)代入式(4.85)，就得到一广义特征值问题20KφMφ(4.87)求解以上方程可以确定φ和，结果得到n个特征解221122,,,φφ…，2,nnφ，其中特征值12,，…，n代表系统的n个固有频率，并有120n.对于结构的每个固有频率，由式(4.87)可以确定出一组各节点的振幅值，它们互相之间保持固定的比值，但绝对值可任意变化，所构成的向量称为特征向量，在工程上通常称为结构的固有振型。设特征向量12,φφ，…，nφ代表结构的n个固有振型，它们的幅度可按以下要求规定11,2,,TiiinφMφ(4.88)这样规定的振型又称为正则振型，即所谓的固有振型。固有振型具有如下重要性质。将特征解22,,,iijjφφ代回方程式(4.86)，得到22,iiijjjKφMφKφMφ(4.89)上式前一式两端前乘以Tjφ，后一式两端前乘以Tiφ，并由K和M的对称性推知TTjiijφKφφKφ(4.90)所以可以得到220TijjiφMφ(4.91)由上式可见，当ij时，必有4.50TjiφMφ(4.92)上式表明固有振型对于矩阵M是正交的。和(4.88)式在一起，可将固有振型对于M的正则正交性质表示为1,0,TijijijφMφ(4.93)进而可得2,0,iTijijijφKφ(4.94)定义固有振型矩阵12nΦφφφ，则222212diag()nΩ(4.95)特征解的性质还可表示成TΦMΦI2TΦKΦΩ(4.96)式中，Φ和2Ω分别称为固有振型矩阵和固有频率矩阵。由此，原特征值问题还可以表示成2KΦMΦΩ(4.97)4.4.2固有频率和固有振型的求解方法结构振动的固有频率和固有振型求解是振动模态分析的关键问题。求解固有频率和振型的方法有振型截断法、矩阵逆迭代法、李兹法、广义雅可比法等。对于一个连续体的结构，其固有频率有无限多阶，在有限元中，结构被离散成细小的单元，对于大型复杂的结构，单元的数目可能数以万计，由这些单元形成的振动方程组规模是很庞大的，由此而得到的特征方程，其矩阵的阶次通常很高。有限元中经常用来求解结构的低阶模态。另外，同样规模的特征值问题,其计算量比静力问题的计算量要高出几倍，因此，如何降低特征值问题的计算规模、减少计算量也是一个重要的课题。如下介绍振型截断法即所谓的Guyan缩聚的基本原理。对于一个机械结构，设其静力问题的总体刚度方程为()(1)(1)rrrrKqP（4.98）其中,r为总自由度。将其写成分块矩阵的形式，有111211212211mmmmmsmmsssmssssqPKKKKqP（4.99）其中,rsm。上面的分块是根据节点位置的重要程度来划分的，一般情况下，将对机械产品动态优化设计理论与方法4.6应于结构关键位置的节点位移划分为mq，叫做“主自由度”（masterDOF），而剩下的节点位移划分为Sq，叫做“从自由度”（slaveDOF）。假定进行分块时考虑到0SP的特征（比如结构内部无载荷，或对应于内部自由度的载荷为零），同时还可以定义矩阵对角线的相对刚度系数来确定“主自由度”，即计算iiiiimk（4.100）其中iik，iim分别为系统的刚度矩阵和质量矩阵主对角线上的元素。当相对刚度系数超过某一临界值时，所对应的节点自由度选为“主自由度”，即cri（4.101）其中，cr为临界相对刚度系数值，可以根据需要缩聚的自由度数来确定。方程（4.99）可以写为111221220msmmsKqKqPKqKq(4.102)由上式中的第二式，有12221smqKKq(4.103)将q写成1111122211mmmmmmrmrmmssmsqIqqTqKKq(4.104)其中转换矩阵T为12221mmrmsmITKK(4.105)将该转换关系（4.105）代入到求动力学问题的虚功方程中，可以得到动力学系统的方程1T1T1T1TrmrmmmmrrrrmmmmrrrrmmmmrrrrmPTqTKTqTCTqTMT(4.106)即有mmmmmmmPqKqCqM~~~~(4.107)其中,~TmrrrrmmmmTKTK,~TmrrrrmmmmTMTM,~TmrrrrmmmmTCTC1T~rmrmmmmPTP(4.108)考虑无阻尼自由振动，则4.70~~11mmmmmmmmmmqKqM(4.109)该方程的解为ˆitmmteqq(4.110)其中为系统的固有频率，将式（4.110）代入到式（4.109）中，有0qMKmmmwˆ~~2(4.111)方程（4.111）就是缩聚后的自由振动方程。求出“主自由度”mqˆ后，由（4.104）求出所有节点的位移。以上的缩聚是基于静力问题中对应于“从自由度”上无外载的情形下推导的，即得到缩聚关系式（4.103）。对于动力学方程，同样可以写出完全的分块矩阵方程，即()()()()(1)2()()()()(1)11121112021222122ˆˆmmmsmmmsmsmsssmsssmsqkkMMkkMMq(4.112)对于（4.112）中的第二组方程，有inessmPqKqKˆˆˆ2221（4.113）其中ine22122ˆ()mssPMqMq为惯性力。与式(4.102)中的第二个方程进行比较，式（4.112）中由于存在惯性力inesPˆ，因为它不为零，由式（4.113）得到的sqˆ为mqKMMKqsˆ)()(ˆ21212122222(4.114)在获得“主自由度”的振型mqˆ后，由（4.113）式来求取关于“从自由度”的振型sqˆ和总的振型为TTˆ{}Tmsqqq。由于此时的转换关系（4.114）式与原转换关系（4.103）式有差别，计算出的