1、有限元分析及应用有限元分析及应用 Finite Element Analysis and Application第七章 动力学问题的有限元法 在实际机械结构中，常作用于结构上的载荷是动在实际机械结构中，常作用于结构上的载荷是动载荷，即载荷随时间载荷，即载荷随时间t t相关，这时，结构上相应的位移，相关，这时，结构上相应的位移，应力和应变不仅随空间位置变化，还随时间应力和应变不仅随空间位置变化，还随时间t t而变化。而变化。结构动力学问题的有限元法的实质就是结构动力学问题的有限元法的实质就是将一个弹将一个弹性连续体的振动问题，离散为一个以有限个节点位移性连续体的振动问题，离散为一个以有限个节点位
2、移为广义坐标的多自由度系统的振动问题为广义坐标的多自由度系统的振动问题。其基本原理。其基本原理和分析方法类同静力学的有限元法，按杆梁、薄板等和分析方法类同静力学的有限元法，按杆梁、薄板等不同结构进行分析。不同的是，应用振动理论建立动不同结构进行分析。不同的是，应用振动理论建立动力学方程时，在单元分析中除需形成刚度矩阵外，还力学方程时，在单元分析中除需形成刚度矩阵外，还需形成质量矩阵，阻力需形成质量矩阵，阻力 矩阵；在整体分析中，不仅求矩阵；在整体分析中，不仅求动力响应，更多是求解特征值问题（结构振动的固有动力响应，更多是求解特征值问题（结构振动的固有频率及相应的振动型（或模态）频率及相应的振动
3、型（或模态）7-1 振动基本方程的建立 从静力学有限元法可知，有限元的基本思想是将弹性体离散成从静力学有限元法可知，有限元的基本思想是将弹性体离散成有限个单元，然后据各单元节点的位移协调和节点力平衡，建立整有限个单元，然后据各单元节点的位移协调和节点力平衡，建立整体刚度平衡方程：体刚度平衡方程：上式也同样适用于弹性体受动载荷作用的有限元分析，关于静上式也同样适用于弹性体受动载荷作用的有限元分析，关于静力问题和动力问题的区别，据达朗贝尔原理，动力学问题只要在外力问题和动力问题的区别，据达朗贝尔原理，动力学问题只要在外力在中计入惯性力后，便可按静力平衡处理。考虑到动力问题中的力在中计入惯性力后，便
4、可按静力平衡处理。考虑到动力问题中的载荷和位移均为时间的函数，上式可记为：载荷和位移均为时间的函数，上式可记为：由于动力载荷由于动力载荷 可为作用于弹性体上的动载荷可为作用于弹性体上的动载荷 ，也可，也可为弹性体的惯性力为弹性体的惯性力 ，也可为与速度相关的阻尼力，也可为与速度相关的阻尼力 ，即：，即：据惯性力定义表示为：据惯性力定义表示为：如阻尼力正比与速度，如阻尼力正比与速度，则动力学基本方程：则动力学基本方程：eKR()()eKtR t()()eTF tMt()R t()F t()TF t()cF t()()ecF tCt ()()()()eeeMtCtKtF t7-2 单元质量、阻尼、
5、刚阵计算 1 1、单元刚度阵、单元刚度阵 任取一个单元，单元节点位移为任取一个单元，单元节点位移为 ，节点速度和，节点速度和加速度为：加速度为：，则单元节点内任一点的位移，则单元节点内任一点的位移 NN为形函数，与时间为形函数，与时间t t无关，为无关，为X X、Y Y、Z Z的函数，它的函数，它与静力分析中一样；由于与静力分析中一样；由于NN与时间无关，则单元应与时间无关，则单元应变矩阵，应力矩阵仍与静力分析完全相同：变矩阵，应力矩阵仍与静力分析完全相同：则刚度矩阵同样与静力情况相同：则刚度矩阵同样与静力情况相同：()()etNt()()eett、()et TeVKBDB dV()()etB
6、t ()()()etDtDBt7-2 单元质量、阻尼、刚阵计算2 2、单元质量阵、单元质量阵设单元节点加速度为设单元节点加速度为 ，则单元内任一点的加速度：，则单元内任一点的加速度：设单元的质量密度为设单元的质量密度为 ，则单位体积中的惯性力为：，则单位体积中的惯性力为：负号表示惯性力与加速度相反。负号表示惯性力与加速度相反。显然，整个单元上惯性力即为上式的积分。显然，整个单元上惯性力即为上式的积分。如何将这个作用于如何将这个作用于单元上的惯性力移置到单元节点上单元上的惯性力移置到单元节点上，通常有两种方法：，通常有两种方法：1 1）虚功原理法）虚功原理法求得一致质量矩阵求得一致质量矩阵2 2
7、）直接分配法）直接分配法即按重心不变原则分配，求得集中质量矩即按重心不变原则分配，求得集中质量矩()et()()etNt()()Tp tt 7-2 单元质量、阻尼、刚阵计算 1 1）虚功原理法）虚功原理法 设单元中发生虚位移为设单元中发生虚位移为 则单元惯性力作的虚功为：则单元惯性力作的虚功为：单元节点上节点惯性力所作的功为：单元节点上节点惯性力所作的功为：将将 和和 代入可得代入可得这里这里MM为单元的一致质量矩阵。为单元的一致质量矩阵。显然，对于不同的单元，因形函显然，对于不同的单元，因形函数不同，则质量矩阵也是不同的。数不同，则质量矩阵也是不同的。如平面常应变三角形单元的一致如平面常应变
8、三角形单元的一致质量阵为：质量阵为：*()t*()et*()()etNt*()()TVUttdV *()()eeTWtF t()()etNt ()()()TeVeeeTNN dVtMFtt .50.250.250.50.250.25.50.2501.50.253.50.5eMA7-2 单元质量、阻尼、刚阵计算2 2）直接分配法）直接分配法将单元内分布质量按重心不变原将单元内分布质量按重心不变原则分配至单元节点上，所产生的则分配至单元节点上，所产生的质量矩阵是没有耦合项的对角矩质量矩阵是没有耦合项的对角矩阵。阵。如六自由度的平面三角形单元，如六自由度的平面三角形单元，单元总质量为单元总质量为W/
9、gW/g，则平均分配，则平均分配至三个节点上的质量所形成的质至三个节点上的质量所形成的质量阵为：量阵为：一般而言，一致质量较准确地反一般而言，一致质量较准确地反映了单元内质量分布的实际情况，映了单元内质量分布的实际情况，集中质量精度不如前者，但不存集中质量精度不如前者，但不存在耦合，使计算大大简化，是工在耦合，使计算大大简化，是工程中常用的方法。程中常用的方法。1000001000010001003101eWMg7-2 单元质量、阻尼、刚阵计算 3 3、单元阻尼阵、单元阻尼阵 单元阻尼力主要指结构阻尼力，它是由结构内部材单元阻尼力主要指结构阻尼力，它是由结构内部材料内摩擦引起的阻尼。设结构阻尼
10、系数为料内摩擦引起的阻尼。设结构阻尼系数为 ，则单位，则单位体积产生的阻尼力（即阻尼力密度）为：体积产生的阻尼力（即阻尼力密度）为：利用虚功原理同理可得：利用虚功原理同理可得：为了简化计算，工程中常将阻尼矩阵记为质量矩阵或刚为了简化计算，工程中常将阻尼矩阵记为质量矩阵或刚度矩阵的形式，如令阻尼系数度矩阵的形式，如令阻尼系数 则则()()cp tt ()()()eeceVeTNNFtdtCtV 是比例系数 eeCM7-2 单元质量、阻尼、刚阵计算 一旦单元刚阵、质量矩一旦单元刚阵、质量矩阵、阻尼矩阵求得，则阵、阻尼矩阵求得，则动力学方程中的整体刚动力学方程中的整体刚阵、质量阵等可类似静阵、质量阵
11、等可类似静力分析的刚度矩阵组装力分析的刚度矩阵组装得到：得到：111neneneKKMMCC7-3 固有频率和振型计算 计算结构的固有频率和振型是结构动力学分析的主要内计算结构的固有频率和振型是结构动力学分析的主要内容，也是分析结构动力响应和其它动力特性问题的基础。由容，也是分析结构动力响应和其它动力特性问题的基础。由于一般结构阻尼对结构的固有频率和振型影响极小，所以，于一般结构阻尼对结构的固有频率和振型影响极小，所以，求结构的固有频率和振型时，直接用无阻尼的自由振动方程求结构的固有频率和振型时，直接用无阻尼的自由振动方程求解。即求解。即因任意弹性体的自由振动都可分解为一系列的简谐振动的迭因任
12、意弹性体的自由振动都可分解为一系列的简谐振动的迭加：即结构上各节点位移为加：即结构上各节点位移为 为节点位移振幅向量（即振型），与时间为节点位移振幅向量（即振型），与时间t t无关的位移幅无关的位移幅值；值；为与该振型对应的频率。为与该振型对应的频率。()()0eeMtKt 0()ej tte 07-3 固有频率和振型计算1 1、特征方程、特征方程将节点位移代入动力方程，化简得广义特征值问题：将节点位移代入动力方程，化简得广义特征值问题：由于结构自由振动时，各个节点的振幅不可能全为零，则由于结构自由振动时，各个节点的振幅不可能全为零，则称为结构的称为结构的特征方程特征方程，即求结构的固有频率和
13、振型归结为特征值问，即求结构的固有频率和振型归结为特征值问题。设计结构的自由度为题。设计结构的自由度为n n，则特征方程为，则特征方程为 的的n n次代数方程，其次代数方程，其n n个根称为个根称为特征值特征值，记为，记为它们的平方根称为系统的固有频率，即它们的平方根称为系统的固有频率，即将这些固有频率从小到大依次排列为将这些固有频率从小到大依次排列为最低的频率最低的频率 称为基频，它是所有频率中最重要的一个。称为基频，它是所有频率中最重要的一个。20()0KM 20KM221,n 1,rrn 12n217-3 固有频率和振型计算 2 2、特征向量、特征向量 对应每个固有频率对应每个固有频率
14、，可有方程，可有方程 由此求得一组节点振幅不全为由此求得一组节点振幅不全为0 0的向量的向量 称称 为特征向量，也称为振型或模态向量。由为特征向量，也称为振型或模态向量。由于上述方程为齐次方程，显然解于上述方程为齐次方程，显然解 不唯一，也就是说：不唯一，也就是说：振型的形状是唯一的，但其振幅不是唯一的；振型的形状是唯一的，但其振幅不是唯一的；或或 一个特征值一个特征值 可对应有多个特征向量，但一个特征可对应有多个特征向量，但一个特征向量只对应一个特征值。向量只对应一个特征值。实际中，常选特征向量实际中，常选特征向量 使使 这个过程称之为正规化这个过程称之为正规化 利用正规化，可得利用正规化，
15、可得 2()0()0rKMr()()()()001020rrrrn()0 1,.,rrn()()001TrrM()0r()0r2r()()200TrrrK7-3 固有频率和振型计算3 3、特征向量的性质、特征向量的性质正交性：任意两个特征值对应的特征向量关于质量矩阵或刚度矩阵正正交性：任意两个特征值对应的特征向量关于质量矩阵或刚度矩阵正交。即设交。即设则有则有 （请同学们自己证明）（请同学们自己证明）若将所有的特征值若将所有的特征值 对应的特征向量对应的特征向量组装成特征向量矩阵，即组装成特征向量矩阵，即则对应所有的特征值问题则对应所有的特征值问题可简记为矩阵形式：可简记为矩阵形式：考虑到正规
16、化，考虑到正规化，可进一步记为：可进一步记为：2()2()00 rsrs()()000 TrsMrs 212nKM221,n(1)()00n()2()00 1,2,.rrrKMrn 212TnK TM(1)()00,n7-4 特征值问题的解法 结构固有频率和振型的计算归结为求结构固有频率和振型的计算归结为求的特征值和特征向量的特征值和特征向量 由于有限元法将结构离散为由于有限元法将结构离散为n n个自由度，个自由度，n n一般相当大，一般相当大，故故n n次特征方程的直接求解十分困难，常求其近似解，常用次特征方程的直接求解十分困难，常求其近似解，常用的求解方法有幂迭代法、逆迭代法、子空间迭代法
17、等的求解方法有幂迭代法、逆迭代法、子空间迭代法等1 1、幂迭代法、幂迭代法特点：用于计算最大（主）特征值十分有效。特点：用于计算最大（主）特征值十分有效。设设则则这里这里DD称为动力矩阵，也即一个变换矩阵，它可将任一特征称为动力矩阵，也即一个变换矩阵，它可将任一特征向量变换为一常数与其自身的乘积。向量变换为一常数与其自身的乘积。20()0KM21 1DKM 20000 KMD 广义特征值问题标准特征值问题7-4 特征值问题的解法 由于任两个特征值对应的特征向量是正交的，则由于任两个特征值对应的特征向量是正交的，则n n个特征向量可组成个特征向量可组成特征向量空间中的一个特征向量基，其特征向量空
18、间中的任一特征向量可特征向量空间中的一个特征向量基，其特征向量空间中的任一特征向量可表示为基向量的线性组合。即存在任一向量：表示为基向量的线性组合。即存在任一向量：设这个向量被设这个向量被DD变换后形成一新的特征向量为：变换后形成一新的特征向量为：类推，可得：类推，可得：(1)(2)()()010200011nnrnr ()()00()000211()()010111 rrrnrrDnnrrrrrrDD 11()0010111 ()npprrrppD7-4 特征值问题的解法 由于所有的特征值排列为：由于所有的特征值排列为：即即 存在存在 考虑到问题为齐次方程，特征向量前的系数考虑到问题为齐次方
19、程，特征向量前的系数 可以略去，则上式在可以略去，则上式在p p趋近无穷时，其第一项趋近无穷时，其第一项就趋近就趋近 实际计算，只需迭代有限次即可得精确解实际计算，只需迭代有限次即可得精确解12n111,()0 pprr 则当时1p(1)012n 11()(1)010011 ()pnpprrrp7-4 特征值问题的解法 幂法迭代格式幂法迭代格式 1 1、选初始特征向量、选初始特征向量 ，如单位向量，如单位向量 2 2、构造新特征向量，并归一化、构造新特征向量，并归一化 3 3、计算特征值近似值、计算特征值近似值 4 4、计算相邻两次迭代的特征值误差，、计算相邻两次迭代的特征值误差，检查是否收敛
20、检查是否收敛 若需计算二阶、三阶等特征值，则需构造新的动力矩阵若需计算二阶、三阶等特征值，则需构造新的动力矩阵 01 001000max1kkkikkkkDmm02110ikkik2211121kkk(1)(1)1100TrrrrrDDM7-4 特征值问题的解法 2 2、逆迭代法、逆迭代法 逆迭代法也称为反幂法，类似逆迭代法也称为反幂法，类似于幂法，特征值问题改写为：于幂法，特征值问题改写为：其具体迭代格式为：其具体迭代格式为：1 1）选初始向量）选初始向量 如单位向量如单位向量 2 2）计算中间向量）计算中间向量 3 3）求解线性方程组）求解线性方程组 4 4）归一化）归一化 5 5）计算特
21、征值近似值）计算特征值近似值 6 6）计算相邻两次迭代的特征值）计算相邻两次迭代的特征值误差，检查是否收敛误差，检查是否收敛 001KM 01 01kkYM 0kkKY 000max1kikkkkmm7-5 动力响应的计算 动力响应是指结构在外加激振力作用下，所产生的位移、速度、动力响应是指结构在外加激振力作用下，所产生的位移、速度、加速度。对于受迫振动，由于阻尼影响很大，不能忽略，即基本方加速度。对于受迫振动，由于阻尼影响很大，不能忽略，即基本方程为程为求解此方程通常有两种数值方法：求解此方程通常有两种数值方法：振型迭加法和逐次积分法振型迭加法和逐次积分法1 1、振型迭加法、振型迭加法振型迭
22、加法的基本思想是利用结构固有振型的正交性，把结构的复振型迭加法的基本思想是利用结构固有振型的正交性，把结构的复杂振动分解为一组相互独立的单自由度振动（即解耦），从而求得杂振动分解为一组相互独立的单自由度振动（即解耦），从而求得结构的位移响应。结构的位移响应。设结构无阻尼自由振动的各阶固有频率和相应的固有振型为：设结构无阻尼自由振动的各阶固有频率和相应的固有振型为：则结构任意时刻的受迫振动产生的位移可认为是则结构任意时刻的受迫振动产生的位移可认为是n n个固有振型为基的个固有振型为基的线性组合，即线性组合，即 为组合系数，是时间为组合系数，是时间t t的函数，也称为振形坐标。的函数，也称为振形坐
23、标。()()()()eeeMtCtKtF t(1)(2)()12000,nn(1)(2)()01020()()()()nnty ty ty t()1,iy tin7-5 动力响应的计算 上式可记为上式可记为 这里这里 代如动力学方程：代如动力学方程：左乘左乘 (1)(2)()00012 ()()()()nntYYy ty ty t ()MYCYKYF t ()TTTTMYCYKYF t T广义质量阵广义阻尼阵广义刚度阵广义激振力7-5 动力响应的计算 据正交性可知，这些广义矩阵均为对角矩阵，据正交性可知，这些广义矩阵均为对角矩阵，即表示方程各个变量之间是没有耦合项的，从而即表示方程各个变量之间
24、是没有耦合项的，从而动力方程转化为动力方程转化为n n个相互独立的单自由度振动的个相互独立的单自由度振动的动力方程，动力方程，分别求解这分别求解这n n个方程可求得个方程可求得 从而求得动力方程的位移解从而求得动力方程的位移解 进而可求得速度、加速度进而可求得速度、加速度 ()()()()()()000000()0()()()=()1,2,TTTrrrrrrrrrTrMy tCy tKy tF trn12(),(),()ny ty ty t(1)(2)()01020()()()()nnty ty ty t7-5 动力响应的计算2 2、逐次积分法、逐次积分法基本思想：将时间基本思想：将时间t t
25、离散为离散为n n个区间，并假设在一个个区间，并假设在一个 时间区间内，结时间区间内，结构的加速度响应为线性变化，由此，对加速度积分，可得速度和位移，构的加速度响应为线性变化，由此，对加速度积分，可得速度和位移，一旦所有区间计算完毕，则求出结构的动力响应。一旦所有区间计算完毕，则求出结构的动力响应。假设在假设在 至至t t的很小时间间隔内的很小时间间隔内 ，加速度线性变化：，加速度线性变化：对对 积分，并引入初始条件待定积分常数积分，并引入初始条件待定积分常数将将 代入代入t t时刻的动力方程时刻的动力方程并整理后即可逐步求解各时刻的加速度，然后求出各时刻的速度和位并整理后即可逐步求解各时刻的加速度，然后求出各时刻的速度和位移。移。ttt()()()()(0)ttttttt 2211221136ttttttttttttttttt (),(),()ttt ()()()()eeeMtCtKtF tt
