机械振动数值分析

机械振动数值分析第1页，共148页，2023年，2月20日，星期一参考书目有限单元法，王勖成，清华大学出版社计算动力学，王天舒，清华大学出版社有限元分析的概念与应用，RobetD.Cook，西安交通大学出版社第2页，共148页，2023年，2月20日，星期一求解域的离散——将连续体求解域离散为若干个子域，通过子域边界协调条件构造连续体；利用近似函数求解——每个子域内用近似函数分片表示全域内待定的未知变量；原问题的等效——利用变分或加权余量法，建立基本变量代数方程组或常微分方程组。概述有限元法要点第3页，共148页，2023年，2月20日，星期一复杂几何构型的适应性；各种物理问题的可应用性；建立于严格理论基础上的可靠性；适合计算机实现的高效性。有限元法特性第4页，共148页，2023年，2月20日，星期一单元类型和形式；有限元法的理论基础和离散格式；有限元方程的求解方法；有限元分析程序开发。有限元法的发展和现状第5页，共148页，2023年，2月20日，星期一标准化——任意复杂问题→模块化分解，单元建模→有限种类模块化单元规范化——几何建模→力学建模→求解→后处理分析通用化——形成标准模块化程序应用规模化、普及性——求解问题规模庞大，易于为工程技术人员掌握有限元法的优势第6页，共148页，2023年，2月20日，星期一连续性假设——变形体内部处处连续均匀性假设——变形体内部物质分配均匀各向同性假设——物质在各方向上特性相同线弹性假设——变形与外力作用的关系为线性小变形假设——变形量远小于物体本身尺寸基本假设第一章线弹性动力学变分原理第7页，共148页，2023年，2月20日，星期一第8页，共148页，2023年，2月20日，星期一运动微分方程应变方程本构方程（物理方程）第9页，共148页，2023年，2月20日，星期一边界条件初始条件第10页，共148页，2023年，2月20日，星期一变形体域内任意一点在任意时刻均满足运动微分方程。变形体边界上任意一点在任意时刻均满足边界条件。问题的精确解的特点加权余量法的特点变形体域内和边界上任意一点在任意时刻均近似满足运动微分方程。第11页，共148页，2023年，2月20日，星期一残余力方程加权余量法允许运动平衡方程和边界条件在若干局部存在残差，但要求残余力在域内和边界上的加权积分为零，即第12页，共148页，2023年，2月20日，星期一当上式对任意权函数均满足，则称为式（7）为微分方程（1）和边界条件（4）的等效积分形式。一般情况可选择近似解将式（8）代入式（7），通过确定系数强迫残余力在域内和边界上在某种平均意义下为零。下面的讨论中假设近似函数完全满足边界条件，只考虑域内残差问题。第13页，共148页，2023年，2月20日，星期一权函数可以选N个函数的线性组合，即将式（9）代入式（7），得第14页，共148页，2023年，2月20日，星期一配点法——取Dirac函数为权函数子域法——权函数在N个子域内取1，在子域外取零，即几类常用权函数第15页，共148页，2023年，2月20日，星期一最小二乘法——调整近似函数中的参数，使余量均方和最小，即几类常用权函数第16页，共148页，2023年，2月20日，星期一迦辽金法——取试探函数为权函数迦辽金法的特点几类常用权函数余量方程相当于虚功；求解方程系数矩阵有对称性；当存在泛函时与变分法有等效结果。第17页，共148页，2023年，2月20日，星期一例1用各种加权余量法计算图示弹性基础梁的挠度第18页，共148页，2023年，2月20日，星期一弹性基础梁的基本微分方程和边界条件为第19页，共148页，2023年，2月20日，星期一取试探解为无弹性基础时的精确解则近似解可表示为容易发现，式（16）严格满足边界条件第20页，共148页，2023年，2月20日，星期一残差方程可写为精确解配点法子域法伽辽金法最小二乘法1.00000.17880.17240.18380.17900.183210.0000.078360.067570.089290.078910.08304100.000.011340.009540.014530.011970.058181000.00.0010250.0009950.0015510.0012620.006068第21页，共148页，2023年，2月20日，星期一达朗伯—拉格朗日原理在式（7）中取权函数为真实位移的变分，可得与运动微分方程和边界条件的等效积分形式。对方程中的第一项进行分部积分，得式（19）为运动微分方程和边界条件等效积分的“弱”形式第22页，共148页，2023年，2月20日，星期一阶连续性函数n-1导数连续，且第n阶导数仅有有限个可积间断点的函数第23页，共148页，2023年，2月20日，星期一哈密顿原理对式（18）在任意时间间隔内积分对给定时刻，方程中的第一项可转化为第24页，共148页，2023年，2月20日，星期一将式（21）代入式（20），得普遍意义下的哈密顿原理式（22）说明，对真实运动，系统动能变分和内、外力虚功在任意时间间隔内对时间的积分为零。第25页，共148页，2023年，2月20日，星期一考虑粘滞力后式（23）中考虑了粘滞力的虚功。第26页，共148页，2023年，2月20日，星期一考虑到应力应变关系，外力虚功可改写为将式（25）代入式（22），可得式（26）说明，完整有势系统在任意时间间隔内满足几何关系和给定位移边界条件的所有可能运动中，真实运动使哈密顿作用量取驻值第27页，共148页，2023年，2月20日，星期一例2用哈密顿原理推导受均布动载荷的等截面悬臂梁的振动微分方程梁内任意一点位移可表示为任意一点的应变可表示为第28页，共148页，2023年，2月20日，星期一系统总势能可表示为系统总动能可表示为第29页，共148页，2023年，2月20日，星期一对系统势能取变分对系统动能取变分第30页，共148页，2023年，2月20日，星期一代入哈密顿方程得第31页，共148页，2023年，2月20日，星期一拉格朗日乘子法约束泛函罚函数法引入附加泛函后，原泛函的有附加条件的驻值问题转化为修正泛函无附加约束条件的驻值问题。第32页，共148页，2023年，2月20日，星期一拉格朗日乘子法修正泛函的变分以离散结构为例，需要满足位移边界条件修正泛函为驻值条件第33页，共148页，2023年，2月20日，星期一由此得拉格朗日乘子法修正泛函的变分拉格朗日乘子法的特点方程组的阶数增加；拉格朗日乘子有明确物理意义，相当于力；导出的系数矩阵存在零对角元。第34页，共148页，2023年，2月20日，星期一罚函数法修正泛函的变分仍以离散结构为例，修正泛函为驻值条件由此得第35页，共148页，2023年，2月20日，星期一罚函数法的特点附加条件近似满足；不增加方程阶数；罚参数相当于刚度系数；解的精度与罚参数有关。罚函数法修正泛函的变分第36页，共148页，2023年，2月20日，星期一例3如图所示弹簧系统，在第一个自由度上施加外力，在第二个自由度上施加强制位移。用拉格朗日乘子法和罚函数法计算第一个自由度的位移未施加强制位移时可写出系统平衡方程第37页，共148页，2023年，2月20日，星期一将强制位移边界代入，可解得精确解为满足强制约束条件施加在第二个自由度上的载荷拉格朗日乘子法第38页，共148页，2023年，2月20日，星期一可解得为约束反力的负值罚函数法第39页，共148页，2023年，2月20日，星期一广义变分原理将场函数应满足的附加条件引入泛函，将的变分原理转化为无附加约束条件的变分原理。利用拉格朗日乘子法将哈密顿原理附加条件引入泛函，得H-W变分原理泛函为其中第40页，共148页，2023年，2月20日，星期一广义变分原理由于所有变分项均独立，可推导出泛函取驻值条件第41页，共148页，2023年，2月20日，星期一广义变分原理拉格朗日乘子的物理意义为由此得H-W变分原理泛函为第42页，共148页，2023年，2月20日，星期一习题2如图所示弹簧系统，在第一个和第三个自由度上施加外力，在第二个自由度上施加强制位移。用拉格朗日乘子法和罚函数法计算第一个和第三个自由度的位移，并计算第二个自由度上的作用力第43页，共148页，2023年，2月20日，星期一几何关系——应力应变关系——杆单元的刚度矩阵和质量矩阵第二章有限元离散泛函第44页，共148页，2023年，2月20日，星期一杆单元的刚度矩阵和质量矩阵取位移插值函数，用单元结点位移将杆内任意点位移表示为第45页，共148页，2023年，2月20日，星期一将式（2）代入式（1），得取变分，得杆单元的刚度矩阵和质量矩阵第46页，共148页，2023年，2月20日，星期一记将式（8）代入式（7），结合哈密顿原理，并考虑到变分的任意性可解得杆单元的刚度矩阵和质量矩阵第47页，共148页，2023年，2月20日，星期一2结点杆单元的刚度矩阵和质量矩阵如式（10）所示杆单元的刚度矩阵和质量矩阵若采用集中质量矩阵，则为第48页，共148页，2023年，2月20日，星期一考虑在平面中变形的杆单元，如图所示平面杆单元的坐标转换总体坐标系内的杆单元第49页，共148页，2023年，2月20日，星期一则单元刚度矩阵、质量矩阵和载荷向量转换到总体坐标系中为平面杆单元的坐标转换对空间杆系，转换原理与平面杆系原理相同，转换矩阵可表示为第50页，共148页，2023年，2月20日，星期一考虑如图所示杆系结构平面杆单元例子平面杆系结构第51页，共148页，2023年，2月20日，星期一总体坐标系下的单元刚度矩阵和单元质量矩阵为平面杆单元例子第52页，共148页，2023年，2月20日，星期一平面杆单元例子第53页，共148页，2023年，2月20日，星期一有限元离散的一般流程几何体的离散逼近；自由度选择与形函数的确定；单元应变场的表示；单元应力场的表示；哈密顿原理或达朗伯-拉格朗日原理；系统泛函变分；单元内部积分；单元坐标系到总体坐标系的转换；总体刚度矩阵和质量矩阵的组装；约束条件的处理或引入。第54页，共148页，2023年，2月20日，星期一形函数的确定形函数的阶数与单元形式和自由度数有关以平面三结点三角形单元为例2、结点位移表示待定系数1、假设域内位移的插值函数第55页，共148页，2023年，2月20日，星期一形函数的确定3、解出系数，并合并关于结点位移的项次，得形函数其中注意下标轮换第56页，共148页，2023年，2月20日，星期一形函数的确定4、确定由结点位移表示的域内位移场，并表示为矩阵形式第57页，共148页，2023年，2月20日，星期一形函数的确定形函数的几个性质结点上的插值函数满足单元中任意点处形函数和应有选取多项式时，常数项和一次项要完备为表达统一起见，可采用参数变换第58页，共148页，2023年，2月20日，星期一单元内部积分单元刚度矩阵和单元质量矩阵的一般形式单元等效结点载荷的一般形式第59页，共148页，2023年，2月20日，星期一单元内部积分定积分的数值积分（Guass积分）核心思想——通过若干个点的函数值的加权组合计算单点积分两点积分第60页，共148页，2023年，2月20日，星期一单元内部积分Guass积分点位置及权系数第61页，共148页，2023年，2月20日，星期一单元内部积分二维和三维问题第62页，共148页，2023年，2月20日，星期一约束条件的处理直接代入法特点方程阶数降低第63页，共148页，2023年，2月20日，星期一约束条件的处理对角元素置1法特点适用于零位移情况第64页，共148页，2023年，2月20日，星期一约束条件的处理对角元素置大数法特点相当于罚函数法第65页，共148页，2023年，2月20日，星期一等参变换等参变换的基本思想就是局部系向总体系的映射局部坐标系向物理坐标系的映射类似于构造形函数，可构造物理坐标和局部坐标之间的关系。仍然以三角形平面单元为例第66页，共148页，2023年，2月20日，星期一等参变换需要注意的地方等参变换第67页，共148页，2023年，2月20日，星期一等参变换对面积分来说三维问题？？？等参变换第68页，共148页，2023年，2月20日，星期一等参变换后的单元刚度矩阵、单元质量矩阵格式等参变换单元等效结点载荷的一般形式第69页，共148页，2023年，2月20日，星期一第三章大型系统特征值问题系统运动方程对特征值问题，可建立系统运动方程如下自由振动解将式（2）代入式（1），得第70页，共148页，2023年，2月20日，星期一系统运动方程由式（3）得特征方程解得n个特征根及相应的特征向量n次代数方程另一种形式第71页，共148页，2023年，2月20日，星期一特征向量正交性取结合式（6），得第72页，共148页，2023年，2月20日，星期一特征向量正交性进一步得正交性关系一个重要的关系第73页，共148页，2023年，2月20日，星期一刚度矩阵为半正定矩阵质量矩阵为半正定矩阵半正定矩阵的特征值问题特征向量特征值2自由度系统的例子第74页，共148页，2023年，2月20日，星期一Rayleigh商对给定的振型向量，可定义Rayleigh商Rayleigh商的上下界考虑式（12），可将式（15）表示为第75页，共148页，2023年，2月20日，星期一Rayleigh商将式（16）展开，得物理意义——假设模态相当于增加了约束，使得系统刚度提高第76页，共148页，2023年，2月20日，星期一误差估计设和是标准特征值问题的近似解，定义残差考虑下式得第77页，共148页，2023年，2月20日，星期一误差估计定义向量和矩阵的2范数对式（21）两边取2范数，整理得第78页，共148页，2023年，2月20日，星期一误差估计由于将式（24）代入式（23）得式（23）说明，只要残差的2范数很小，给定的就接近系统的某一阶特征值第79页，共148页，2023年，2月20日，星期一逆迭代法基本假定刚度矩阵非奇异迭代方程考虑到任意向量都可由特征向量表示为第80页，共148页，2023年，2月20日，星期一逆迭代法综合式（26）和式（27），得第81页，共148页，2023年，2月20日，星期一逆迭代法从式（30）中容易看出，当考虑Rayleigh商第82页，共148页，2023年，2月20日，星期一逆迭代法结合式（27）、（28）和（29）可将式（33）化为第83页，共148页，2023年，2月20日，星期一逆迭代法逆迭代法的一般步骤给定初始迭代向量；若已计算出前i-1阶模态，根据已得到的特征向量进行正交化处理；第84页，共148页，2023年，2月20日，星期一逆迭代法逆迭代法的一般步骤解方程；计算特征值第85页，共148页，2023年，2月20日，星期一逆迭代法几个要点正交化处理；归一化处理；当刚度矩阵奇异时的处理；系统包含近频时的处理移轴法第86页，共148页，2023年，2月20日，星期一带移轴的逆迭代法考虑特征值问题其特征值与原系统特征值之间满足如下关系由式（26）给出的迭代关系可得怎么定？？？第87页，共148页，2023年，2月20日，星期一Rayleigh商迭代法Rayleigh商迭代法——以Rayleigh商作为移轴量计算流程第88页，共148页，2023年，2月20日，星期一变换法变换法的基本思想是通过迭代的方法逐次构造变换矩阵，使矩阵和趋近对角形，其一般迭代思路如下式所示当时，应有第89页，共148页，2023年，2月20日，星期一变换法并且有两类应用广泛的变换法雅可比变换法Householder-QR法第90页，共148页，2023年，2月20日，星期一雅可比变换法适用范围——标准特征值问题，即考虑标准特征值问题取正交矩阵使得第91页，共148页，2023年，2月20日，星期一雅可比变换法雅可比法中选择的变换矩阵为旋转矩阵，在一次变换中使待变换矩阵中的一个非对角元变为零。例第92页，共148页，2023年，2月20日，星期一雅可比变换法将变换为零，即，由此得第93页，共148页，2023年，2月20日，星期一雅可比变换法定义阀值收敛条件第94页，共148页，2023年，2月20日，星期一雅可比变换法收敛性检验考虑特征值问题，其中由式（47）可计算得第95页，共148页，2023年，2月20日，星期一雅可比变换法经第一次循环，对非对角元（1，2）消零，得由式（52）可计算得第96页，共148页，2023年，2月20日，星期一雅可比变换法类似的，分别对非对角元（1，3）和（2，3）消零，得第97页，共148页，2023年，2月20日，星期一雅可比变换法雅可比变换法的一般步骤给定当前循环的阀值，一般取第m次迭代的阀值为10-2m；对上三角块的非对角元按式（48）计算，若计算结果大于给定的阀值，则按式（47）进行旋转变换；检查式（49）式是否满足，若满足则由迭代最终步骤计算出的矩阵给出特征值，计算特征向量，否则转向第1步，开始新一轮循环。第98页，共148页，2023年，2月20日，星期一广义雅可比变换法对雅可比变换矩阵进行修改，得转换矩阵第99页，共148页，2023年，2月20日，星期一广义雅可比变换法系数的选取应使刚度矩阵和质量矩阵在变换后其第i行j列的元素同时为零，即若，可取，否则第100页，共148页，2023年，2月20日，星期一广义雅可比变换法其中第101页，共148页，2023年，2月20日，星期一广义雅可比变换法例1，用雅可比法求解图示系统的固有频率第102页，共148页，2023年，2月20日，星期一广义雅可比变换法第103页，共148页，2023年，2月20日，星期一广义雅可比变换法精确解第104页，共148页，2023年，2月20日，星期一Householder-QR变换法Householder-QR变换法可大致分为三个步骤Householder变换；特征值的QR迭代计算；特征向量计算。第105页，共148页，2023年，2月20日，星期一Householder-QR变换法Householder变换的目的在于将矩阵变换为三对角阵，加速QR迭代收敛速度。QR迭代考虑到实对称矩阵可由正交矩阵和上三角矩阵分解为第106页，共148页，2023年，2月20日，星期一Householder-QR变换法由此定义迭代格式当时，，第107页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法若假设振型与系统的第i阶主振型相差一阶小量，则Rayleigh商与系统的第i阶特征值相差二阶小量。假设近似振型为其中里兹基向量里兹坐标第108页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法定义Rayleigh商取上式的驻值，得第109页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法由式（65）可以计算出新系统的特征值和特征向量满足则原系统的特征值和特征向量为第110页，共148页，2023年，2月20日，星期一给出的近似特征值前一半精度较高，因此，若需要求解前r阶模态，要先定义出2r阶假设振型；经变换后，系统规模大大降低；Rayleigh-Ritz法计算结果依赖于里兹基的选择；里兹基可通过静力问题求得。Rayleigh-Ritz法Rayleigh-Ritz法的特点第111页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法例2，用Rayleigh-Ritz法计算下述广义特征值问题的近似解。假设载荷为第112页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法由新系统特征值问题解得精确解第113页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法如果假设载荷为由新系统特征值问题解得第114页，共148页，2023年，2月20日，星期一Rayleigh-Ritz法进一步得第115页，共148页，2023年，2月20日，星期一子空间迭代法子空间迭代法假设r个初始向量同时进行迭代，计算前p个特征值和特征向量。子空间迭代法可以看为是逆迭代法的推广，算法步骤类似，其具体迭代步骤可表示为1、初始计算（1）形成刚度矩阵和质量矩阵；（2）按所给定的边界条件修正刚度矩阵；（3）三角分解，得第116页，共148页，2023年，2月20日，星期一子空间迭代法2、给定初始迭代向量矩阵3、对每次迭代，按下面步骤执行（1）赋值：（2）求解方程组：（3）赋值：（4）赋值：（5）赋值：（6）解广义特征值问题：第117页，共148页，2023年，2月20日，星期一子空间迭代法（7）检验特征值是否满足精度要求：若满足精度要求，则转到第4步，若不满足，则（8）赋值：（9）令：，并返回步骤（2）；4、赋值：5、输出结果：第118页，共148页，2023年，2月20日，星期一第四章系统动力学响应系统运动方程对一般的振动响应，考虑初值问题可建立系统运动方程如下求解方式振型叠加法直接积分法第119页，共148页，2023年，2月20日，星期一振型叠加法引入坐标变换其中广义位移由式（2）和式（1）得第120页，共148页，2023年，2月20日，星期一振型叠加法若阻尼矩阵可解耦，则方程可变换为利用Duhamel积分，得第121页，共148页，2023年，2月20日，星期一直接积分法直接积分法不对运动方程进行任何变换，直接对方程进行积分求解。其原理基于两个基本概念将求解域上任意时刻均需满足运动方程的要求更改为在一定条件下近似的满足运动方程，如，在若干离散的时间点上满足方程；在一定数目的时间间隔区域内，假设位移、速度和加速度的近似插值函数。误差来源截断误差——近似表达式中略去高阶小量；舍入误差——超过计算机字长的四舍五入。第122页，共148页，2023年，2月20日，星期一直接积分法解的稳定性无条件稳定——对任意初始条件，无论时间步长怎么选取，积分结果都不会无界增大；条件稳定——为保证积分结果的有界性，要求时间步长必须小于某个临界时间步长。考虑直接积分法递推格式在讨论稳定性时，可令，由式（7）得给定的初始条件递推矩阵第123页，共148页，2023年，2月20日，星期一直接积分法解的稳定性对递推矩阵进行谱分解，得递推矩阵为二阶非对称矩阵，且有重特征根，其约当标准形为式中为递推矩阵的约当标准形，其对角元为递推矩阵的特征值。当时，只有，才能使在时有界当时，只要，就能使在时有界第124页，共148页，2023年，2月20日，星期一直接积分法解的稳定性若递推矩阵为二阶非对称矩阵，没有重特征根，其约当标准形为易知，只要满足时，就有界递推矩阵的谱半径第125页，共148页，2023年，2月20日，星期一直接积分法解的稳定性稳定性准则若递推矩阵没有重特征根，要求其谱半径不大于1；若递推矩阵存在重特征根，要求其特征值的模小于1。数值阻尼（人工阻尼）——谱半径小于1.0时，随迭代增加，矩阵第126页，共148页，2023年，2月20日，星期一中心差分法将位移函数按Taylor级数展开，得由式（13）可得第127页，共148页，2023年，2月20日，星期一中心差分法将式（14）代入任意时刻的运动方程得其中第128页，共148页，2023年，2月20日，星期一中心差分法起始步的计算定义积分常数第129页，共148页，2023年，2月20日，星期一中心差分法计算步骤1、初始计算1）形成刚度矩阵和质量矩阵和阻尼矩阵；2）给定初始条件和，并计算；3）给定时间步长，计算积分常数；4）计算起步位移：；5）计算有效质量矩阵：；6）对有效质量矩阵进行三角分解：；第130页，共148页，2023年，2月20日，星期一中心差分法2、对每一个时间步1）计算时刻的有效载荷；2）计算时刻的位移：；3）如有需要，计算时刻的加速度和速度第131页，共148页，2023年，2月20日，星期一中心差分法中心差分法的数值稳定性对任意时刻的运动方程考虑其中心差分格式，有将式（21）写为递推矩阵格式有第132页，共148页，2023年，2月20日，星期一中心差分法其中对无阻尼系统，由式（23）计算出递推矩阵特征值为第133页，共148页，2023年，2月20日，星期一中心差分法根据收敛条件，可知从而得对有阻尼系统，中心差分法的临界时间步长为第134页，共148页，2023年，2月20日，星期一中心差分法由于对任意一个解耦后的运动方程来说均应满足式（27），因此，满足收敛条件的时间步长由系统中的最高阶模态对应的固有频率决定。容易看出，当满足式（25）时，谱半径为1，此时中心差分法没有数值阻尼。第135页，共148页，2023年，2月20日，星期一

Newmark法在位移和速度的Taylor级数中只保留一阶导数项，则可得若在式（28）中将导数的平均值代替时刻的导数值，则可将式（28）改写为第13