机械系统动力学 第三章 机械系统运动微分方程的求解2

Runge-Kutta法是求解常微分方程应用最多的方法之一。 对于微分方程的定解问题，欧拉法求解，其截断误差 O(h2 ) 故具有1阶精度，改进欧拉法，由于预测了 结点的差商并 用 两个节点的差商的平均值来代替导数，可望达到2阶精
度。实际上，在区间[ xn , xn1 ]的等价积分形式为
xn 1
y(xn1) y(xn ) f (x, y(x))dx
改进的欧拉法以 Pn 和 Pn+1 两个节点的差商的平均值来 代替导数，由于 值为待求值，故计算 结点的差商采用 预测，其迭代公式
预测： 校正：
yn1 yn hf (xn , yn )
yn1
yn
h[f 2
(xn ,
yn )
f
( xn1, yn1)]
可以证明，欧拉法具有1阶精度，而改进的欧拉法具有2
其欧拉法的迭代公式为
x(t t) x(t) y(t)t
y(t
t)
y(t)
y(t)t
y(t) f (x, y,t)
x(t t) x(t) y(t)t
y(t
t)
y(t)
f
( x,
y, t )t
3-2机械系统的运动方程求解方法-数值法
3-2-1 欧拉法 改进的欧拉法的迭代公式为：
预报： xy((tt
yn1 yn h( xn , yn , h)
其中 称增量函数，可表示为
r
( xn , yn , h)= ciKi
K1
i 1
f (xn , yn )
i 1
Ki
f (xn
ih, yn
h ijK j )
j 1
i 2, , r
t) t)
x(t) y(t)
y(t)t f (x, y,
t)t
校正： xy((tt
t) t)
x(t) （y(t) y(t y(t) ( f (x, y,t)
t))t / 2 f (x(t t),
y(t
t),
t
t
))t
/
2
3-2机械系统的运动方程求解方法-数值法
3-2-2 Newmark- 法 Newmark- 法是线性加速度法之一。对于具有关于时间2 阶导数的单自由度机械系统运动微分方程式，其 x(t t) 的Talar展开式：
x(t
t)
x(t)
x(t)t
t 2 2!
x(t)
t 2 [ x(t +t )-x(t )]
x(t t)
x(t)
t [x(t+t)+x(t)] 2
式中 为调节公式特征的参数，一般取值范围为
0 1/2
3-2机械系统的运动方程求解方法-数值法
3-2-2 Newmark- 法 对于多自由度振动系统运动微分方程：
3-2机械系统的运动方程求解方法-数值法
3-2-1 欧拉法
对于具有关于时间2阶导数的单自由度机械系统运动微分 方程，形如
x f (x, x,t) x(0) x0, x(0) x0
可令 x y 将上式转化成1阶常微分方程组
y f (x, y,t) x y x(0) x0, y(0) y0 x0
2!
3! t
3-2机械系统的运动方程求解方法-数值法
3-2-2 Newmark- 法
x(t t) x(t) x(t)t x(t) t2 x (t) t3 o(t4)
2!
3!
线性加速度法的迭代公式 1
大致具有3阶精度，将上式的最后一项中
即为Newmark- 法。其迭代公式为
3!用
代替，
称为步长，一
般在计算时常取步长为定值,这时节点为
xn x0 nh, n 0,1, 2,
3-2机械系统的运动方程求解方法-数值法
3-2-1 欧拉法
初值问题3-2-1的数值解法的求解过程为：给出用已知
信息 yn , yn1, yn2
计算 yn1 的递推公式，从初
始条件出发，顺着节点排列的次序一步一步地向前推
xn
通过增加积分求积的结点数提高计算精度，故将右端 的积分表示为
xn 1
r
f (x, y(x))dx h ci f (xn ih,y(xn1 ih))
xn
i 1
一般来说，接点数越多，计算越准确
3-2机械系统的运动方程求解方法-数值法
3-2-3 Runge-Kutta法 仿照欧拉法的迭代公式，写成
x(t t) x(t) x(t)t x(t) t2 x (t) t3 o(t4)
2!
3!
上式中取前三项， 若认为加速度在区间[ t , t+t ]
为线性变化，则有
x (t)= x(t+t)-x(t) t
代入上式
x(t t) x(t) x(t)t t2 x(t) t3 x(t+t)-x(t)
MX (t) CX (t) KX (t) F(t)
t+t 时刻有关系式
MX (t+t) CX (t+t) KX (t+t) F (t+t)
代入式Newmark- 法迭代公式
MX (t+t) C[X (t) t X (t+t) X (t)] K{X (t) X (t)t t2 X (t)
进。即所谓“步进式”算法。
欧拉法以节点的差商代替导数值，构成的递推公式为：
yn1 yn xn1 xn
f (xn , yn )
即欧拉（Euler）公式： yn1 yn hf (xn 来自百度文库 yn )
3-2机械系统的运动方程求解方法-数值法
3-2-1 欧拉法
从图3-2-1（b）可以看出，由于欧 拉法是以差商代替导数，其误差较 大。为了提高计算精度，一种办法 是减小步长，但会导致累计误差增 大，当步长减小到一定程度后，计 算精度提高受限。另一种办法是改 进算法，如改进的欧拉法、RungeKutta法等。
2
2!
t2[X (t+t)-X (t)]} F(t+t)
整理移项：
X
(t +t )=[M+
C
t
(t)2 K]-1{F (t+t)
C[ X
(t)
t
X
(t)]
2
2
K[ X (t) tX (t) (1 )t2 X (t)]}
2
3-2机械系统的运动方程求解方法-数值法
3-2-3 Runge-Kutta法
3-2机械系统的运动方程求解方法-数值法
3-2-1 欧拉法
对于常微分方程的定解问题，形如
y f (x, y)
y(x0 )
y0
3-2-1
所谓数值解法, 就是寻求解 y(x) 在一系列离散节点
x1 x2 xn xn1 上的近似值 y1, y2 , , yn , yn1 。
相邻两个节点的间距 hn xn1 xn