第三节 Runge-Kutta方法 一、R一K方法的基本思想 二、N级R一K公式 三、小结
第三节 Runge-Kutta方法 一、R-K方法的基本思想 二、N级R-K公式 三、小结
一、R一K方法的基本思想 3.1 首先,从微分中值定理及方程(1.1)得 x)=)=y'5)=f5,》 x<g<xn 这里f(飞,(5)称为方程(1.1)的积分曲线(x)在区间 [xm,x+上的平均斜率。由此可知,只要对此平均斜率提 供一种算法,就可以得到一个相应的计算公式。下面, 我们来观察Euler格式和改进的Euler格式,将它们分别 写成
3.1 一、R-K方法的基本思想 首先,从微分中值定理及方程(1.1)得 ( ) ( , ( )) ( ) ( ) 1 y f y h y x y x n n = = + − n n+1 x x 这里 称为方程(1.1)的积分曲线 在区间 上的平均斜率。由此可知,只要对此平均斜率提 供一种算法,就可以得到一个相应的计算公式。下面, 我们来观察Euler格式和改进的Euler格式,将它们分别 写成 f (, y()) y(x) [ , ] n n+1 x x
ym-yn=f(xny) 前一式是用Xn点处的斜率f(xn,yn)来代替上面所说的平均 斜率的,后一式是用X,,X+1两点上的斜率的平均值来代 替平均斜率的。我们已经知道,uler格式是1阶方法,而 改进的Euler格式是2阶方法。 由此看来,如果在区间[xm,x]内多预报几个点的斜率值,然 后将它们加权平均斜率,就可以构造出更高阶的计算公式来
[ ( , ) ( , )] 2 1 ( , ) 1 1 1 1 + + + + = + − = − n n n n n n n n n n f x y f x y h y y f x y h y y 前一式是用 点处的斜率 来代替上面所说的平均 斜率的,后一式是用 , 两点上的斜率的平均值来代 替平均斜率的。我们已经知道,Euler格式是1阶方法,而 改进的Euler格式是2阶方法。 n x ( , ) n n f x y n x n+1 x 由此看来,如果在区间 内多预报几个点的斜率值,然 后将它们加权平均斜率,就可以构造出更高阶的计算公式来。 [ , ] n n+1 x x
二、N级R一K公式 Euler格式与改进的Euler格式可以改写成下面的形式 yn=yn +K K=hf(xnyn) y=+2K,+ K=hf(xn2y) K2=hf(xn+1>y+i) 显然,若在区间[xm,x+]内取N个不同的点,记积分曲 线y(x)在这N个点上的斜率分别为.K,K,,Kv于是可设 n+月=x,)+立C,K+0h) (3.1)
二、N级R-K公式 Euler格式与改进的Euler格式可以改写成下面的形式 = + = + ( , ) 1 n n n n K hf x y y y K = = = + + + + + ( , ) ( , ) 2 1 2 1 2 1 1 1 1 1 2 n n n n n n K hf x y K hf x y y y K K 显然,若在区间 内取 个不同的点,记积分曲 线 在这 个点上的斜率分别为, [ , ] n n+1 x x y(x) N K K KN , , , N 1 2 于是可设 ( ) ( ) ( ) 1 1 + + + = + + m N i n n i i y x h y x C K O h (3.1)
舍去误差项,便得到 +CK i+1 (3.2) K=hf(xnyn) K hf(xn +ahiyn+ biK,) i=2,3,…,W 这就是所谓的N级m阶的R一K公式。其中a,b,多是 待定系数,并且有a= ∑b,1=23,N (3.3 a,b,C,待定系数可用比较系数的方法求得
舍去误差项,便得到 = + + = = = + − = + + K hf x a h y b K i N K hf x y y y C K i j i n i i n i j j n n N i n n i i ( , ) 2,3, , ( , ) 1 1 1 1 1 (3.2) 这就是所谓的N级m阶的R-K公式。其中 多是 待定系数,并且有 i ij i a ,b ,c − = = = 1 1 , 2,3, , i j ai bij i N (3.3) ai ,bij ,ci 待定系数可用比较系数的方法求得
以N=2为例,说明待定系数的求法。当N=2时,由(3.1) 式有yxn+h)=y(xn)+C,K1+C2K2+O(h"m)(3.4) 将(3.4)式中的xn+h与K,K分别在x,处作Taylor展 开,有 ,+)=x)4wx.)+2,)+0h h K hf(y)=hy'( K2 hf[x +ah,y(x)+b2 hf(x y(x)] =hIf(xn,y(x))+hazf (x,y(x))+ hb2f(xn,xn》f(xn,y(xn》+O(h2) =hf +h'az (f:+ff)+O(h) hy'(x)+h'azy"(x)+O(h)
以 为例,说明待定系数的求法。当 时,由(3.1) 式有 (3.4) N = 2 ( ) ( ) ( ) 1 1 1 2 2 + + = + + + m n n y x h y x C K C K O h N = 2 将(3.4)式中的 与 分别在 处作Taylor展 开,有 y(x h) n + 1 2 K ,K n x ( ) ( ) ( ) ( ) ( ) ( , ( )) ( , ( )) ( ) [ ( , ( )) ( , ( )) [ , ( ) ( , ( ))] ( , ) ( ) ( ) ( ) 2 ( ) ( ) ( ) 3 2 2 3 2 2 2 2 1 2 2 2 2 1 1 3 2 hy x h a y x O h hf h a f f f O h hb f x y x f x y x O h h f x y x ha f x y x K hf x a h y x b hf x y x K hf x y hy x y x O h h y x h y x hy x n n x y n n y n n n n x n n n n n n n n n n n n n = + + = + + + + = + + = + + = = + = + + +
于是得到 C1+C2=1 1 22= 若取95,则可算得C=与4=-l,这时,由(32)式得 y4=.+2y)+j,+h.+x刀 称为修正的梯形公式。若取C,=1,则可算得C=0,ab=2 ,由(3.2)得 =x*x培} 称为修正的中矩形公式
于是得到 = + = 2 1 1 2 2 1 2 C a C C 若取 2 ,则可算得 ,这时,由(3.2)式得 1 C2 = , 1 2 1 C1 = a2 = b21 = [ ( , ) ( , ( , ))] 2 n 1 n n n n n n n f x y f x h y hf x y h y + = y + + + + 称为修正的梯形公式。若取 ,则可算得 ,由(3.2)得 1 C2 = 2 1 0, C1 = a2 = b21 = + = + + + ( , ) 2 , 2 n 1 n n n n n f x y h y h y y hf x 称为修正的中矩形公式
以上两个公式,都是在W=2及m=2的前提下构造出 来的。因此,它们都是2级2阶的R一K公式。注意, 上面求待定系数的方程组中,有一个自由参数,故2级 2阶的R一K公式有无穷多个。但是,在这些2级R一K 公式中,不可能存在高于2阶的方法。下面,我们给出 N级R一K公式可以达到的最高阶数: N, 当N=1,2,3,4时 1m= N-1, 当N=5,6,7时 <N-2, 当N≥8时
以上两个公式,都是在 及 的前提下构造出 来的。因此,它们都是2级2阶的R-K公式。注意, 上面求待定系数的方程组中,有一个自由参数,故2级 2阶的R-K公式有无穷多个。但是,在这些2级R-K 公式中,不可能存在高于2阶的方法。下面,我们给出 N级R-K公式可以达到的最高阶数: N = 2 m = 2 − − = = = 当 时 当 时 当 时 2, 8 1, 5,6,7 , 1,2,3,4 N N N N N N m
三、标准4级4阶R一K公式 y1=yn+(K1+2K2+2K3+K) K1=hf(xn2yn) K2=hf(xn+hyn+5K (3.5) 1 K,=hf(xn+5hyn+5K2》 2 Ka =hf(xn+h,yn +k3)
三、标准4级4阶R-K公式 = + + = + + = + + = + = + + + + ( , ) ) 2 1 , 2 1 ( ) 2 1 , 2 1 ( ( , ) ( 2 2 ) 6 1 4 3 3 2 2 1 1 1 1 2 3 4 K hf x h y K K hf x h y K K hf x h y K K hf x y y y K K K K n n n n n n n n n n (3.5)
例4用标准4级4阶R一K公式求解例1中给出的初值问 题,較02 解这里具体的计算公式如下: y=,+(K1+2K2+2K3+K) K=0.2ym K2=0.2ym+ K12(xn+0.1) 2 K3=0.2yn+ K22(xm+0.) 2 yn+ 2 K4=02yn+K- 2(xm+0.) n=0,l,…,5 y+K3
例4 用标准4级4阶R-K公式求解例1中给出的初值问 题,取 。 解 这里具体的计算公式如下: h = 0.2 = + + = + − + + = + − + + = + − = − + = + + + + 0,1, ,5 2( 0.1) 0.2 2 2( 0.1) 2 0.2 2 2( 0.1) 2 0.2 2 0.2 ( 2 2 ) 6 1 3 4 3 2 2 3 1 1 2 1 1 1 2 3 4 n y K x K y K K y K x K y K y K x K y y x K y y y K K K K n n n n n n n n n n n n n n