04-1 振动数值仿真方法

1、燕山大学机械工程学院School of Mechanical Engineering, Yanshan University 振动分析的方法很多，数值仿真方法是进行振动振动分析的方法很多，数值仿真方法是进行振动分析的最直接的一类方法，它们可以应用于包括非线性分析的最直接的一类方法，它们可以应用于包括非线性振动在内的各种振动问题，这类方法是研究动态响应的振动在内的各种振动问题，这类方法是研究动态响应的有效手段之一。有效手段之一。第四章第四章 振动的仿真振动的仿真 从数学的观点来看，数值仿真方法是解微分方程从数学的观点来看，数值仿真方法是解微分方程边值问题和初值问题的逐步方法。在结构动力学响应计边

2、值问题和初值问题的逐步方法。在结构动力学响应计算方面，采用实用有效的数值仿真方法，可以对系统在算方面，采用实用有效的数值仿真方法，可以对系统在任意激励下的动态响应进行分析。任意激励下的动态响应进行分析。燕山大学机械工程学院School of Mechanical Engineering, Yanshan University在时间域内对响应的时间历程进行离散，把运动微分在时间域内对响应的时间历程进行离散，把运动微分方程分为各离散时刻的方程；方程分为各离散时刻的方程；将某时刻的速度和加速度用相邻时刻的各位移的线性将某时刻的速度和加速度用相邻时刻的各位移的线性组合表示，将系统的运动微分方程化为一个

3、由位移组成组合表示，将系统的运动微分方程化为一个由位移组成的某离散时刻的代数方程组；的某离散时刻的代数方程组；对耦合的系统运动微分方程进行逐步数值积分，从而对耦合的系统运动微分方程进行逐步数值积分，从而求出在一系列离散时刻上的响应值。求出在一系列离散时刻上的响应值。 数值仿真方法的特点数值仿真方法的特点 这种数值仿真方法称为逐步积分法这种数值仿真方法称为逐步积分法( (或直接积分法或直接积分法) )。 燕山大学机械工程学院School of Mechanical Engineering, Yanshan University中心差分法；中心差分法；侯博特侯博特(Houbolt)法；法；威尔逊威

4、尔逊(Wilson- )法；法；纽马克纽马克(Newmark- )法。法。 对于高频分量和低频分量混合的问题，采用无条件对于高频分量和低频分量混合的问题，采用无条件稳定的解法，可以提高计算效率。稳定的解法，可以提高计算效率。求解多自由度线性振动系统常用的方法有：求解多自由度线性振动系统常用的方法有：燕山大学机械工程学院School of Mechanical Engineering, Yanshan University 中心差分法是直接积分法的一种。中心差分法是直接积分法的一种。 它是将系统的运动微分方程在时间域内离散，化它是将系统的运动微分方程在时间域内离散，化成对时间的差分格式，然后根据

5、初始条件，利用逐步积成对时间的差分格式，然后根据初始条件，利用逐步积分求出在一系列离散时刻上的响应值。分求出在一系列离散时刻上的响应值。 4.14.1 中心差分法中心差分法离散系统的运动微分方程为离散系统的运动微分方程为 tRKxxCxM 式中式中M，C，K分别为系统的质量矩阵，阻尼矩阵和刚度分别为系统的质量矩阵，阻尼矩阵和刚度矩阵矩阵； , , ，x分别表示系统的加速度向量，速度向量分别表示系统的加速度向量，速度向量和位移向量；和位移向量；R(t)是外力向量。是外力向量。x x 燕山大学机械工程学院School of Mechanical Engineering, Yanshan Unive

6、rsity 在中心差分法中，按中心差分将速度和加速度向量在中心差分法中，按中心差分将速度和加速度向量离散化为离散化为ttttttxxx21nTt tttttttxxxx212 假定在假定在t=0时，位移、速度和加速度分别为已知时，位移、速度和加速度分别为已知的的 。求时间区间求时间区间0，T的解。的解。000,x xx的近似解。的近似解。目的：目的：确定时刻确定时刻Ttnttttttn,2, 0210, ,x x x 两式中，两式中，t t时刻的时刻的速度和加速度是速度和加速度是以相邻时刻的位以相邻时刻的位移表示的。移表示的。把时间全程把时间全程T划分为划分为n等份，即：等份，即：燕山大学机械

7、工程学院School of Mechanical Engineering, Yanshan University在在t时刻的动力方程为时刻的动力方程为 ttttRKxxCxM 式中式中ttttttttxCMxMKRR211222CMKtt211221212tttttttttttttxxxxxxxttRxK t燕山大学机械工程学院School of Mechanical Engineering, Yanshan University22211(3)2ttttttttRRKM xMC x211(2)2ttKMC(1)tttKxR求解方程式求解方程式(1)，可得，可得xt+ t。21212ttttt

8、ttttttttxxxxxxx02002xxxx ttt由式由式(3)可以看出，为求可以看出，为求xt+ t必须使用必须使用xt和和xt- t的值的值开始计算时，即开始计算时，即t=0时，要计算时，要计算x t的值，就需要已知的值，就需要已知的的x- t值，而值，而x- t是未知的。是未知的。需要一个起始技术，因而这种算法不是自起步的。需要一个起始技术，因而这种算法不是自起步的。由于由于 是已知的。根据是已知的。根据000,x xx燕山大学机械工程学院School of Mechanical Engineering, Yanshan University中心差分法的计算机实施格式中心差分法的计

9、算机实施格式22211(3)2ttttttttRRKM xMC x211(2)2ttKMC(1)tttKxRA. 初始计算初始计算 1. 形成质量矩阵形成质量矩阵M，阻尼矩阵，阻尼矩阵C和刚度矩阵和刚度矩阵K。 2. 给出初始值给出初始值 3. 选择时间步长选择时间步长t，t tcr，计算积分常数：，计算积分常数： 4. 计算计算 。 5. 形成有效刚度矩阵：形成有效刚度矩阵： 6. 对对 作三角分解：作三角分解：20120321,1 2,2,1atat aa aa0300 xxxx attCMK10aaTLDLK 000,x xxK燕山大学机械工程学院School of Mechanical

10、 Engineering, Yanshan UniversityB. 关于每一时间增量计算关于每一时间增量计算 1. 计算计算t时刻的有效载荷时刻的有效载荷 2. 计算计算tt时刻的位移时刻的位移 3. 如果需要，计算如果需要，计算t时刻的加速度和速度时刻的加速度和速度taaattttxCMxMKRR102ttRxLDLT tttattttxxxx20 ttatttxxx122211(3)2ttttttttRRKM xMC x211(2)2ttKMC(1)tttKxR燕山大学机械工程学院School of Mechanical Engineering, Yanshan University 中

11、心差分法是一种显式积分方法。中心差分法是一种显式积分方法。 使用中心差分法必须考虑积分的时间步长使用中心差分法必须考虑积分的时间步长t不不能大于临界值能大于临界值tcr，即，即式中式中Tn为离散系统的最小周期。为离散系统的最小周期。ncrTtt如果不满足上式，数值解将出现发散现象。如果不满足上式，数值解将出现发散现象。这种算法不是无条件稳定的。这种算法不是无条件稳定的。燕山大学机械工程学院School of Mechanical Engineering, Yanshan University 侯博特侯博特(Houbolt)法是法是Houbolt为研究飞机振动所为研究飞机振动所提出的方法。提出的

12、方法。 该方法以三级位移插值为基础的，通过四点的位该方法以三级位移插值为基础的，通过四点的位移建立三次式，用两个向后差分公式表示在时刻移建立三次式，用两个向后差分公式表示在时刻t+t的的速度和加速度，即速度和加速度，即4.24.2 侯博特法侯博特法222111189261254tttttttttt ttttttttttxxxxxxxxxx在在tt时刻的动力方程为时刻的动力方程为ttttttttRKxxCxM 燕山大学机械工程学院School of Mechanical Engineering, Yanshan University整理得关于整理得关于xt+ t的代数方程组的代数方程组式中式中1

13、ttttKxR（）MCKK22611tt222253431123ttttttttttt ttttRRMCxMCxMCx该方法不是自起步的，要用其它方法由该方法不是自起步的，要用其它方法由 起步，起步，例如可用中心差分法求出例如可用中心差分法求出x t和和x2 t后，才能使用后，才能使用Houbolt法的方程逐步求解。法的方程逐步求解。000 xxx， ，由上式可以看出，要计算由上式可以看出，要计算xt+ t时刻的解，必须使用前时刻的解，必须使用前三步的位移三步的位移xt, xt- t和和xt-2 t。 燕山大学机械工程学院School of Mechanical Engineering, Ya

14、nshan UniversityHoubolt法的计算机实施格式法的计算机实施格式A. 初始计算初始计算 1. 形成质量矩阵形成质量矩阵M，阻尼矩阵，阻尼矩阵C和刚度矩阵和刚度矩阵K。 2. 给出初始值给出初始值 x0， ， 。 3. 选择时间步长选择时间步长t，并计算积分常数：，并计算积分常数： ， ， ， ， ， ， 。 4.使用特殊的起始过程，计算使用特殊的起始过程，计算x t和和x2 t。 5. 形成有效刚度矩阵形成有效刚度矩阵 ： 6. 对对 作三角分解：作三角分解：0 x 0 x KKKCMK10aaTLDLK 202tata6111225tata33042aa235aa206aa

15、 937aa 燕山大学机械工程学院School of Mechanical Engineering, Yanshan UniversityB. 关于每一时间增量计算关于每一时间增量计算1. 计算计算t+t时刻的有效载荷时刻的有效载荷 2. 计算计算tt时刻的位移时刻的位移 3. 如果需要，计算如果需要，计算t+t时刻的加速度和速度时刻的加速度和速度24623572ttttttaaaaaa ttttttttRRMxxxCxxxttttaaaa26420tttttxxxxx ttttaaaa27531tttttxxxxx ttttRxLDLTHoubolt法和中心差分法的根本不同之处是刚度矩阵法和

16、中心差分法的根本不同之处是刚度矩阵K出现在方程出现在方程(1)的左端，因此的左端，因此Houbolt法是隐式积分格式，法是隐式积分格式，其舍入误差与步长其舍入误差与步长t的大小无关，所以的大小无关，所以Houbolt法是无法是无条件稳定的。条件稳定的。燕山大学机械工程学院School of Mechanical Engineering, Yanshan UniversityWilson- 法模型法模型4.3 威尔逊威尔逊- 法法 威尔逊威尔逊 (Wilson- )法法是假定在是假定在t, t+ t( 1)时间时间间隔内，加速度呈线性变化，间隔内，加速度呈线性变化，如图所示。令如图所示。令 为自

17、为自t时刻开始时刻开始的时间变量，适用于的时间变量，适用于0 t。根据线性加速度的假。根据线性加速度的假设，可得在此范围内的加速设，可得在此范围内的加速度为度为1tt ttttxxxx（）燕山大学机械工程学院School of Mechanical Engineering, Yanshan University若若 =t，由以上两式可得，由以上两式可得t +t瞬时的速度和位移瞬时的速度和位移 上式积分后得上式积分后得222tt tttttxxxxx（）321326tt ttttttxxxxxx（）22226ttttttt tttttttttxxxxxxxxxttttxxxx tt根据上式，将根

18、据上式，将t +t时刻的加速度和速度用位移表示时刻的加速度和速度用位移表示 。燕山大学机械工程学院School of Mechanical Engineering, Yanshan University22226ttttttt tttttttttxxxxxxxxx226624tttt tttttxxxxx（）tttttxxxxx 223tttt燕山大学机械工程学院School of Mechanical Engineering, Yanshan University在在t +t时刻的动力方程为时刻的动力方程为 式中式中ttttttttRKxxCxM ttttRRRRtt22662322tttt

19、tttt ttttttttttxxxxxxxxxx燕山大学机械工程学院School of Mechanical Engineering, Yanshan University整理得关于整理得关于xt+ t的线性方程组的线性方程组式中式中ttttRxKMCKK2263tt22663222tttttttt ttttttttRRRRMxxxCxxx求解上述代数方程组，可得求解上述代数方程组，可得xt+t+t t。燕山大学机械工程学院School of Mechanical Engineering, Yanshan University同样取同样取 = t，将式，将式(1)分别代入式分别代入式(2)和

20、式和式(3)，有，有 这样就完成了一步积分。这样就完成了一步积分。tttttxxxxx 3166222ttttttttxxxx ttt2tttttxxxxx 262tttt 求出求出t+t 瞬时的位移瞬时的位移xt+t后，代入式后，代入式(4)就可获就可获得得 。在式。在式(1)中取中取 = t，并将式，并将式(4)代入，有代入，有 ttx 1tt ttttxxxx（）226624tttt tttttxxxxx（）燕山大学机械工程学院School of Mechanical Engineering, Yanshan University 本方法的物理意义是：假定加速度在时刻本方法的物理意义是：

21、假定加速度在时刻tt+ t内为线性变化，首先计算内为线性变化，首先计算t, t t区间的近似解，但区间的近似解，但仅取其中前半部分仅取其中前半部分(到时刻到时刻tt)作为正式的近似解而作为正式的近似解而舍去后半部分舍去后半部分(时刻时刻tt以后的部分以后的部分)。这种巧妙的处。这种巧妙的处理并非出于物理的原因，而主要是数学计算技术的理由。理并非出于物理的原因，而主要是数学计算技术的理由。 在在Wilson- 中，只要中，只要 值取值取1.37以上，不管以上，不管t取怎取怎样的值都是稳定的样的值都是稳定的(即这种算法是无条件稳定的即这种算法是无条件稳定的)。实际实际上，上， 最好不要太大，否则精

22、度会下降最好不要太大，否则精度会下降(截断误差增加截断误差增加)。因此，因此，Wilson推荐的合理推荐的合理 值为值为1.4。 燕山大学机械工程学院School of Mechanical Engineering, Yanshan UniversityWilson- 法的计算机实施格式法的计算机实施格式A. 初始计算初始计算 1. 形成质量矩阵形成质量矩阵M，阻尼矩阵，阻尼矩阵C和刚度矩阵和刚度矩阵K。 2. 给出初始值给出初始值 x0， ， 。 3. 选择时间步长选择时间步长t，取，取 1.4，计算积分常数：，计算积分常数： ， ， ， ， ， ， ， 。 4. 形成有效刚度矩阵形成有效刚

23、度矩阵 ： 5. 对对 作三角分解：作三角分解：0 x 0 x KKKCMK10aaTLDLK 206tata31122aa 23ta04aa 25aa316a27ta628ta燕山大学机械工程学院School of Mechanical Engineering, Yanshan UniversityB. 关于每一时间增量计算关于每一时间增量计算 1. 计算计算t+t时刻的有效载荷时刻的有效载荷 2. 计算计算tt时刻的位移时刻的位移 3. 计算计算t+t时刻的加速度和速度和位移时刻的加速度和速度和位移 ttttRxLDLTttttRRRRttttttttxxxCxxxM 312022aaaa

24、tttttxxxxx 654aaattttttxxxx tta7tttttxxxxx 28ttat燕山大学机械工程学院School of Mechanical Engineering, Yanshan University Wilson- 是一种隐式积分方法，即每计算是一种隐式积分方法，即每计算一步，必须解一个线性代数方程组。一步，必须解一个线性代数方程组。 Wilson- 算法是自起步的，算法是自起步的，tt时刻的时刻的位移，速度和加速度都可由位移，速度和加速度都可由t时刻的变量表示，时刻的变量表示，不需要特别的起动技术。不需要特别的起动技术。燕山大学机械工程学院School of Mech

25、anical Engineering, Yanshan University 纽马克纽马克- (Newmark- )法同样也是假定在时间间隔法同样也是假定在时间间隔 t, tt内加速度呈线性变化，它的基本假定为内加速度呈线性变化，它的基本假定为 4.4 纽马克纽马克- 法法式中式中 和和 为按积分的精度和稳定性要求可以调整的参数。为按积分的精度和稳定性要求可以调整的参数。研究表明，当研究表明，当 1/2， 1/4(1/2+ )2时，时，Newmark- 法是法是无条件稳定的。无条件稳定的。1tttttttxxxx( 1)212(2)tt tttttttxxxxx根据以上两式，根据以上两式， 和

26、和 可用可用表示。表示。 ttx ttx , ,ttttxx x燕山大学机械工程学院School of Mechanical Engineering, Yanshan University1tttttttxxxx( 1)212(2)tt tttttttxxxxxtttttttxxxxx ttt121112(3)tttttttxxxxx tt211(4)Newmark- 法每步积分应满足法每步积分应满足t+t时刻的动力方程时刻的动力方程ttttttttRKxxCxM (5)将公式将公式(3)、(4)代入式代入式(5)，可得，可得关于关于xt+ t 的方程为：的方程为：燕山大学机械工程学院School of Mechanical Engineering, Yanshan University式中求解方程求解方程(6)就可得到就可得到xt+t ，然后根据式，然后根据式(3