第四讲 龙格-库塔方法.docVIP

龙格-库塔方法 3.2 Runge-Kutta法.2.1 显式 Runge-Kutta法的一般形式 上节已给出与初值问题 1.2.1 等价的积分形式 .2.1 只要对右端积分用不同的数值求积公式近似就可得到不同的求解初值问题 1.2.1 的数值方法，若用显式单步法 .2.2 当，即数值求积用左矩形公式，它就是Euler法 3.1.2 ，方法只有一阶，若取　　 3.2.3 就是改进Euler法，这时数值求积公式是梯形公式的一种近似，计算时要用二个右端函数f的值，但方法是二阶的若要得到更高阶的公式，则求积分时必须用更多的f值，根据数值积分公式，可将 3.2.1 右端积分表示为注意，右端f中还不能直接得到，需要像改进Euler法 3.1.11 一样，用前面已算得的f值表示为 3.2.3 ，一般情况可将 3.2.2 的表示为 3.2.4 其中 这里均为待定常数，公式 3.2.2 ， 3.2.4 称为r级的显式Runge-Kutta法，简称R-K方法它每步计算r个f值 即 ，而ki由前面 i-1 个已算出的表示，故公式是显式的例如当r 2时，公式可表示为 3.2.5 其中.改进Euler法 3.1.11 就是一个二级显式R-K方法参数取不同的值，可得到不同公式.2.2 二、三级显式R-K方法对r 2的显式R-K方法 3.2.5 ，要求选择参数，使公式的阶p尽量高，由局部截断误差定义 3.2.6 令，对 3.2.6 式在处按Taylor公式展开，由于将上述结果代入 3.2.6 得要使公式 3.2.5 具有的阶p 2,即，必须 3.2.7 即由此三式求的解不唯一.因r 2，，于是有解 3.2.8 它表明使 3.2.5 具有二阶的方法很多，只要都可得到二阶R-K方法若取，则，则得改进Euler法 3.1.11 ，若取，则得，此时 3.2.5 为 3.2.9 其中称为中点公式Euler法 3.1.11 及中点公式 3.2.9 是两个常用的二级R-K方法，注意二级R-K方法只能达到二阶，而不可能达到三阶因为r 2只有4个参数，要达到p 3则在 3.2.6 的展开式中要增加3项，即增加三个方程加上 3.2.7 的三个方程求4个待定参数是无解的当然r 2,p 2的R-K方法 3.2.5 当取其他数时，也可得到其他公式，但系数较复杂，一般不再给出对r 3的情形，要计算三个k值，即其中　将按二元函数在处按Taylor公式展开，然后代入局部截断误差表达式，可得可得三阶方法，其系数应满足方程  3.2.10 这是8个未知数6个方程的方程组，解也是不唯一的，通常.一种常见的三级三阶R-K方法是下面的三Kutta方法： 3.2.11 function y DELGKT3_kuta f, h,a,b,y0,varvec format long; N b-a /h; y zeros N+1,1 ; y 1 y0; x a:h:b; var findsym f ; for i 2:N+1 K1 Funval f,varvec,[x i-1 y i-1 ] ; K2 Funval f,varvec,[x i-1 +h/2 y i-1 +K1*h/2] ; K3 Funval f,varvec,[x i-1 +h y i-1 -h*K1+K2*2*h] ; y i y i-1 +h* K1+4*K2+K3 /6; %满足c1+c2+c3 1, 1/6 4/6 1/6 end format short; 3.2.3 四阶R-K方法及步长的自动选择利用二元函数Taylor展开式可以确定 3.2.4 中r 4,p 4的R-K方法， 其中，，而 共计13个参数待定，Taylor展开分析局部截断误差，使得精度达到四阶，即误差为。 于是，r 4，p 4的13个参数 c4不能为0 引出了多种方案和挑战，如： 参数优化使阶数增加到5阶，得到四阶五阶R-K方法，matlab中有程序ode45； 四级四阶R-K方法的步长自动选取； 结合新算法的应用算法构造； 适应于新的领域实现求解； …… 经典的四阶R-K方法是： 3.2.12 ，这里为 1/6 2/6 2/6 1/6 .它的局部截断误差，故p 4这是最常用的四阶R-K方法，数学库中都有用此方法求解初值问题的软件这种方法的优点是精度较高，缺点是每步要算4个右端函数值，计算量较大例.3 　用经典四阶R-K方法解例.1的初值问题，仍取h 0.1，计算到，并与改进Euler法、梯形法在处比较其误差大小.解 用四阶R-K方法公式 3.2.12 ，此处，于是当n 0时于是，按公式 3.2.12 可算出此方法误差： 改进Euler法误差： 梯形法误差： 可见四阶R-K方法的精度比二阶