数值分析26_9.3Runge-kutta法
四阶龙格库塔法(Runge-Kutta)求解微分方程
四阶龙格库塔法（Runge-Kutta ）求解微分方程张晓颖（天津大学 材料学院 学号：1012208027）1 引言计算传热学中通常需要求解常微分方程。
这类问题的简单形式如下：{),(')(00y x f y y x y == （1）虽然求解常微分方程有各种各样的解析方法，但解析方法只能用来求解一些特殊类型的方程，实际问题中的多数微分方程需要采用数值解法求解。
初值问题（1）的数值解法有个基本特点，它们采取“步进式”，即求解过程顺着节点排序一步一步向前推进。
这类算法是要给出用已知信息y n 、 y n −i ……计算y n +1的递推公式即可。
2 龙格库塔法（Runge-Kutta ）介绍假设对于初值问题（1）有解 y = y (x ) ，用 Taylor 展开有：......)(!3)(!2)()()(321+'''+''+'+=+n n n n n x y h x y h x y h x y x y （2） 龙格库塔法（Runge-Kutta ）实质上是间接的使用 Taylor 级数法的一种方法。
对于差商hx y x y n n )()(1-+，根据微分中值定理，存在 0 < θ < 1 ，使得：)()()(1h x y hx y x y n n θ+'=-+ （3）于是对于 y = y (x ) ，可得到：))(,()()(1h x y h x hf x y x y n n n n θθ+++=+ （4）设))(,(*h x y h x f K n n θθ++=，为区间 [x n , x n +1 ] 上的平均斜率。
四阶龙格库塔格式中的*K 由下式计算得到：⎪⎪⎪⎪⎪⎩⎪⎪⎪⎪⎪⎨⎧++=++=++==++++=+),()2,2()2,2(),()22(6342312143211hK y h x f K K h y h x f K K h y h x f K y x f K K K K K h y y n n n n nn n n n n （5） 四阶龙格库塔法（Runge-Kutta ）的每一步需要四次计算函数值f ，其截断误差更低，计算的精度更高。
微分方程数值解的隐显式Runge-Kutta方法
的方法来求解初值问题。 事实上，如果一阶常微分方程初值问题：
，dy dx
=f（x，y），x∈ ， a，b
， ，
y（a）＝y0
（1.1）
的 解 和 它 的 右 端 函 数 f（x，y）是 充 分 光 滑 的 ，那 么 在 y（xn ）给 定 时 ， y
（x） 关 于 x 的 各 阶 导 数 （y′，y″… ） 在 点 xn 处 的 值 可 以 通 过 微 分 方 程
K）
） ）
1
＝f
xn ＋
1 2
－
姨3 6
h，yn
＋
1 4
hK1 ＋
1 4
－
姨3 6
hK2 ，
）
姨 姨 姨 姨 姨 姨 ）
K）
） ）
2
＝f
xn ＋
）
1 2
＋
姨3 6
h，yn ＋
1 4
＋
姨3 6
hK1
＋
1 4
hK2
．
实用对角隐式 Runge-Kutta 方法，s=2，p=3（2 级三阶）的方法为：
Runge-Kutta 的 常 用 格 式 一 般 为 s≤4。 Runge-Kutta 方 法 （1.2），
（1.3）是显式方法，类似可以构造隐式的 Runge-Kutta 方法，s 级的隐式
Runge-kutta 方法为：
）
s
Σ ）
y）
） n＋1
＝yn
＋h
bi Ki ，
）
i＝1
）
）
s
Σ ）
K）
）i
科技信息
○高校讲坛○
SCIENCE & TECHNOLOGY INFORMATION
龙格库塔方法
• 可以通过检查步长折半前后两次计
算结果的偏差
=
(h)
y2 n1
y(h) n1
变步长方
来判断所选取的步长是否合法适
（1）对于给定的精度，如果 ,则反复将步长折半 进行计算直至 为止，这时取折半以后的“ 新值” 作为结果；
四阶RungeKutta方法
这样，下一个值(yn+1)由现在的值(yn)加上时间间隔 (h)和一个估算的斜率的乘积 决定。该斜率是以下斜率的加权平均：
k1是时间段开始时的斜率； k2是时间段中点的斜率，通过欧拉法采用斜率 k1 来 决定 y在点 tn + h/2的值； k3也是中点的斜率，但是这次采用斜率 k2决定 y值； k4是时间段终点的斜率，其 y值用 k3 决定。 当四个斜率取平均时，中点的斜率有更大的权值：
75
0.0033
误差 0.45e-4 0.17e-4 0.15e-4 0.48e-4 0.25e-4 0.55e-4 0.14e-4
y(xn) 1.0954 1.1832 1.2649 1.3416 1.4142 1.4832 1.5492
改进Euler法一步需要计算两个函数值 (h=0.1) 四阶Runge-Kutta方法一步需要计算四 个函数值（h=0.2） 总计算量大致相当，但四阶RungeKutta方法精度更高
设从节点xn出发，先以h为步长求出一个近似值，记为yn(h)1 因为经典格式的局部截断误差为O(h5),因此有
y(xn+1)-y(nh)1 Ch5 其中C与y(5)(x)在[xn, xn1]内的值有关
将步长折半，即取
h 2
为步长从xn跨两步到xn+1，求得一个近似值yn( h21)
Runge-Kutta法
3½×RK· ¨
1.422
¸ ĽøEuler· ¨
1.42
1.418 0.2965 0.297 0.2975 0.298 0.2985 0.299 0.2995 0.3 0.3005 17
§ 7.2 Runge-Kutta法
龙格-库塔(Runge-Kutta)方法简称R-K法，是一种应用较广的 高精度的单步法。
所谓单步法就是在计算yi时只用到前一步信息yi-1的方法。
本节介绍R-K法的构造原理、常用公式。
1
一、Runge-Kutta方法的构造原理 对于常微分方程的初值问题
y f ( x , y ) y( a ) y0 a xb
另一方面
h2 y( xi 1 ) y( xi ) hy( xi ) y( xi ) o(h3 ) 2
9
式（6）的局部截断误差:
en (h) y ( xi 1 ) yi 1 1 h 1 (1 2 ) y( xi ) h2 22 y( xi ) o(h3 ) 2
1 2 1 1 2 2 2
有无穷多组解，从而可以得到许多具体的二阶R -K公式 如：
1 2 1，2 1
1 1 0， 2 1， 2 2
10
用类似的方法也可以构造三阶R-K公式
h yi yi 1 ( K1 4 K2 K3 ) 6 K1 f ( xi 1, yi 1 ) h h K 2 f ( xi 1 , yi 1 K1 ) 2 2 K3 f ( xi 1 h, yi 1 h(2 K2 K1 ))
因而方法(8)有4阶精度
12
例1. 使用高阶R-K方法计算初值问题
rungekutta方法求解常微分
rungekutta方法求解常微分嘿，朋友们！今天咱来聊聊这个神奇的 Runge-Kutta 方法求解常微分。
你说这常微分方程啊，就像是生活中的一个个小难题，有时候还真让人有点头疼呢！那这 Runge-Kutta 方法呢，就像是一把神奇的钥匙，能帮我们解开这些难题。
想象一下，常微分方程就像是一条弯弯曲曲的小路，我们要找到通过这条小路的最佳路径。
而 Runge-Kutta 方法呢，就像是一个经验丰富的向导，带着我们一步一步地往前走，每一步都走得稳稳当当。
它可不是随随便便就出手的哦！它会仔细地分析问题，然后给出精确的计算步骤。
就好像我们做菜一样，要精确地掌握各种调料的用量，才能做出美味的菜肴。
Runge-Kutta 方法里有好多不同的形式呢，就像不同口味的糖果，各有各的特点。
有时候我们要根据具体的问题来选择合适的那种形式，这可得有点眼力见儿呢！比如说，有些问题比较简单，那可能用简单一点的 Runge-Kutta 形式就够了。
但要是遇到复杂的大难题，那就得派出厉害一点的形式啦！这就好比对付小怪兽用普通武器就行，要是遇到大怪兽，那可得拿出厉害的法宝了。
在使用 Runge-Kutta 方法的时候，可不能马虎大意哦！每一个步骤都要认真对待，就像我们走路一样，一步一个脚印，踏踏实实地走。
要是有一步走错了，那可就麻烦啦，可能就找不到正确的答案了。
而且啊，这个方法还很灵活呢！它可以适应各种不同的情况，就像一个全能选手，不管遇到什么挑战都能应对自如。
那我们怎么才能用好这个神奇的 Runge-Kutta 方法呢？这就需要我们多练习啦！就像学骑自行车一样，刚开始可能会摔倒，但多练几次就熟练了。
我们要不断地去尝试，去探索，去发现它的奥秘。
总之呢，Runge-Kutta 方法是我们求解常微分方程的好帮手，只要我们用心去学，用心去用，就一定能在数学的海洋里畅游无阻！它就像是一盏明灯，照亮我们在常微分方程世界里前进的道路。
大家可别小瞧了它哦，赶紧去试试吧！。
Runge-Kutta方法
Runge-Kutta方法
由二元函数Taylor展开式:
k2 hf (xn , yn ) h2 fx(xn , yn ) hk1 f y(xn , yn ) O(h3) hy(xn ) h2 (fx(xn , yn ) f (xn , yn ) f y(xn , yn )) O(h3)
k1 hf (xn , yn )
k2
hf
( xnΒιβλιοθήκη 1 2h,yn
1 2
k1 )
这是二阶R - K方法.
(3)
取1
1 4
,
可得2
3 ,
4
2 ,又有算式 3
yn1 yn 14（k1 3k2 ) k1 hf (xn , yn )
k2
hf
( xn
2 3
h,
yn
2 3
k1 )
这也是二阶R - K方法。
k2
hf
( xn
(1 2
3 6
)h,
yn
(1 4
3 6
)k1
1 4
k2
)
数值计算方法
数值计算方法
Runge-Kutta方法
受改进的Euler方法启发，更一般算式可设为
yn1 yn 1k1 2k2
k1 hf (xn , yn )
(n 0,1,2,...)
k2 hf (xn h, yn k1)
适当选择参数1，2，，，使局部截断误差
Tn1 y(xn1) yn1 O(h3 ),这里仍假定yn y(xn )。
(h (h 1
2
1 2
(h
)
2
)
yn
(h)2 1 (h
4
)3
)
runge kutta法总结
四阶误差 0.0000 1.60E-7 3.40E-7 5.48E-7 7.69E-7 9.95E-7 1.20E-6 1.42E-6 1.68E-6 1.96E-6 2.32E-6
二阶误差 0.000000 2.38E-3 7.09E-3 1.44E-2 2.45E-2 3.78E-2 5.45E-2 7.51E-2 9.99E-2 1.29E-1 1.64E-1
此处即通过计算已知点的函数值(K1,K2)的线性组合代替高阶导数， 得到了较高的精度。
Runge-Kutta方法的推导
Runge-Kutta方法的一般形式：
r Yn 1 Yn h ci K i i 2,3,, r i 1 i 1 , K1 F (t n , Yn ) i ij i 1 j 1 K F (t h, Y h K ) n i n ij j i j 1
Runge-Kutta积分方法
dY F (t , Y )，一般的解法可以表示为：Yn 1 Yn h D(t n , Yn ). dt 其中D(t n , Yn )是动点在(t n , t n 1 )中的平均速度。设Yn是精确的， 对 h2 h3 Y (t n 1 ) Y (t n h) Y (t n ) hY (t n ) Y (t n ) Y (t n ) 2 3! 所以得到：
结果及比较
结果及比较
四阶显式Runge-Kutta方法
f2
f4
f3
f1
f
f 1 f1 2 f 2 2 f 3 f 4 6
xn
xn + h/2
xn + h
例 求解初值问题ODE : dy y x 2 ,
runge-kutta方法
runge-kutta方法Runge-Kutta方法是求解常微分方程组的一种常用数值方法，它是由德国数学家卡尔·龙格和马丁·库塔（Martin Kutta）发明的。
该方法是一种高阶精度的方法，可以求解一阶或高阶的常微分方程组。
$$ y_{n+1} = y_n + \sum_{i=1}^s b_i k_i $$$y_n$表示第$n$个步骤中的输出，$y_{n+1}$是下一个步骤的输出。
$k_i$表示第$i$个中间变量，$b_i$是权值，它们是根据所采用的方法计算得出的。
我们将在下文中详细介绍这些变量和权值的计算方法。
在数值计算中，我们通常采用自适应步长的方法，这需要我们给出一个误差容限值$\epsilon$，根据误差限定步长大小。
如果当前步长不满足误差要求，我们会将步长缩小再进行计算，或者选择更高阶的方法进行计算，以提高精度。
接下来我们将详细介绍如何使用Runge-Kutta方法求解常微分方程组。
1. Runge-Kutta方法的中间变量$k_i$的计算我们需要计算中间变量$k_i$的值。
它们的计算公式如下：$$ k_1 = hf(t_n, y_n) $$$$ k_2 = hf(t_n + c_2 h, y_n + a_{21}k_1) $$$$ ... $$$$ k_s = hf(t_n + c_s h, y_n + \sum_{j=1}^{s-1}a_{sj}k_j) $$$h$表示步长大小，$t_n$和$y_n$分别表示当前时间和状态，$c_i$和$a_{ij}$是与所采用的方法相关的常数。
$$ c_2 = \frac{1}{2} , a_{21} = \frac{1}{2}$$- 改进的Runge-Kutta方法（RK4）：根据具体选择的方法，我们可以计算出$k_i$的值。
$b_{ij}$也是与所采用的方法相关的常数，对于不同的方法而言，它们的计算公式不同。
根据计算出的中间变量$k_i$和权值$b_i$，我们可以进行Runge-Kutta方法的迭代，求解常微分方程组。
Runge-Kutta积分方法
微分方程
Runge-Kutta方法基于微分方程，用 于求解未知函数的导数或偏导数。
数值近似
通过在离散点上取值来近似连续函数 的积分，Runge-Kutta方法提供了一
种数值求解微分方程的方法。
线性插值
Runge-Kutta方法在离散点之间进行 线性插值，以估计函数在未知点的值。
方法的实现步骤
初始条件
为求解微分方程，需要提供初始条件，即函数在某一点的取值。
离散化
将微分方程的求解区间离散化，将连续的时间或空间划分为若干个离散点。
参数选择
根据微分方程的性质和精度要求，选择合适的Runge-Kutta方法参数。
迭代计算
根据离散点和参数，进行迭代计算，逐步逼近微分方程的解。
方法的收敛性与稳定性
收敛性
Runge-Kutta方法的收敛性是指 随着离散点数增加，数值解逐渐 接近真实解。
自适应步长Runge-Kutta方法
总结词
自适应步长Runge-Kutta方法是一种结合了 变步长策略的数值积分方法，通过自动调整 步长以适应解的性质，提高数值积分的精度 和稳定性。
详细描述
自适应步长Runge-Kutta方法在每一步计算 后，根据解的局部误差估计自动调整步长。 这种方法能够更好地处理解的局部变化，避 免因固定步长而产生的数值振荡或误差累积
改进方向
优化算法
通过改进Runge-Kutta方法的算法，减少计算量，提 高计算效率。
引入并行计算
利用并行计算技术，加快Runge-Kutta方法的计算速 度。
改进步长选择策略
通过改进步长选择策略，提高数值解的精度和稳定性。
05
Runge-Kutta方法的实例分 析
二阶偏微分方程组的runge kutta解法
二阶偏微分方程组的runge kutta解法二阶偏微分方程组的 Runge-Kutta 解法与一阶偏微分方程组的Runge-Kutta 解法类似，只需稍作扩展。
考虑一个二阶偏微分方程组：\[\frac{{\partial^2 u}}{{\partial t^2}} = f\left(t, u, \frac{{\partial u}}{{\partial t}}, \frac{{\partial^2 u}}{{\partial x^2}}\right),\]\[\frac{{\partial^2 v}}{{\partial t^2}} = g\left(t, v, \frac{{\partial v}}{{\partial t}}, \frac{{\partial^2 v}}{{\partial x^2}}\right).\]其中，\(t\) 为时间变量，\(x\) 为空间变量，\(u(t, x)\) 和 \(v(t,x)\) 为未知函数，\(f\) 和 \(g\) 为已知函数。
假设已知初始条件为：\[u(t_0, x) = u_0(x), \quad \frac{{\partial u}}{{\partial t}}(t_0, x)= u_1(x),\]\[v(t_0, x) = v_0(x), \quad \frac{{\partial v}}{{\partial t}}(t_0, x)= v_1(x).\]要求解该方程组，可以通过 Runge-Kutta 法进行离散化。
首先，将时间区间\([t_0, T]\) 离散化为\(N\) 个等距时间步长，每个时间步长为 \(\Delta t = \frac{{T-t_0}}{N}\)，对应时间节点为 \(t_n = t_0 + n \cdot \Delta t\)，其中 \(n = 0, 1, \ldots, N\)。
将空间区间离散化为 \(M\) 个等距空间步长，每个空间步长为\(\Delta x = \frac{{b-a}}{M}\)，对应空间节点为 \(x_j = a + j\cdot \Delta x\)，其中 \(j = 0, 1, \ldots, M\)，且给定边界条件为\(u(t, a) = u(t, b) = 0\) 和 \(v(t, a) = v(t, b) = 0\)。
四阶龙格库塔法(runge-kutta)求振动方程
四阶龙格库塔法(runge-kutta)求振动方程四阶龙格库塔法可以用于求解振动方程。
振动方程通常是一个二阶微分方程，可以转化为一个一阶微分方程组来使用四阶龙格库塔法求解。
假设我们要求解的振动方程为：m * d^2x/dt^2 + c * dx/dt + k * x = 0其中，m是质量，c是阻尼系数，k是弹性系数，x是位移，t 是时间。
首先，将振动方程转化为一阶微分方程组。
引入新的变量v来表示速度，则有：dx/dt = vdv/dt = (-c * v - k * x) / m接下来，使用四阶龙格库塔法求解该一阶微分方程组。
首先将时间区间分割成n个小区间，假设每个小区间的宽度为h。
则利用四阶龙格库塔法的公式可以得到如下更新规则：k1 = h * f(t, x, v)k2 = h * f(t + h/2, x + k1/2, v + l1/2)k3 = h * f(t + h/2, x + k2/2, v + l2/2)k4 = h * f(t + h, x + k3, v + l3)其中，f(t, x, v)表示的是微分方程组的右侧。
根据上述更新规则，可以得到新的位移和速度：x_new = x + (k1 + 2k2 + 2k3 + k4) / 6v_new = v + (l1 + 2l2 + 2l3 + l4) / 6重复以上步骤，可以一步一步地求解出整个时间区间内的位移和速度。
需要注意的是，初始条件x(0)和v(0)需要提供。
总结起来，使用四阶龙格库塔法求解振动方程的步骤如下：1. 将二阶微分方程转化为一阶微分方程组。
2. 设置初始条件x(0)和v(0)。
3. 将时间区间分割成n个小区间，计算每个小区间的宽度h。
4. 根据上述更新规则，逐步求解出整个时间区间内的位移和速度。
希望以上解答对您有所帮助！。
动力学问题的数值解法研究
动力学问题的数值解法研究引言动力学是研究物体在不同力的作用下的运动规律和相互作用关系的学科，它在物理学、工程学和生物学等领域具有广泛的应用。
在传统的理论框架下，解析解方法在求解动力学问题中扮演重要角色，但是在复杂的现实问题中，往往难以以解析方法求得准确的解。
此时，数值方法成为研究动力学问题的重要工具。
本文将探讨一些常见的动力学问题的数值解法及其应用。
一、常微分方程数值解法1. Euler方法Euler方法是最简单的一种常微分方程数值解法，它基于欧拉公式，将微分方程转化为差分方程。
具体而言，从初始条件出发，用差分逼近微分，然后迭代更新求解。
然而，Euler方法的精度较低，在解决一些复杂动力学问题时，往往无法满足精确性的要求。
2. Runge-Kutta方法Runge-Kutta方法是一类常用的常微分方程数值解法，其通过使用多步骤的迭代算法，提高了解的精确度。
相较于Euler方法，Runge-Kutta方法的误差更小，因此在实际问题中更受欢迎。
常见的二阶和四阶Runge-Kutta方法在动力学问题的数值求解中具有重要的应用。
二、偏微分方程数值解法偏微分方程是描述动力学问题中时空变化的数学模型，其数值解法相较于常微分方程更为复杂。
在动力学问题中，常见的偏微分方程数值解法包括有限差分法、有限元法和谱方法等。
下面将简要介绍其中两种方法。
1. 有限差分法有限差分法是一种将微分方程离散化的方法，它将求解区域划分为有限的网格，在格点上构造离散的差分方程，从而逼近原微分方程。
该方法通过选取合适的差分格式和网格剖分，能够得到模拟区域内各点的数值解。
有限差分法在求解热传导方程、波动方程和扩散方程等动力学问题中被广泛使用。
2. 有限元法有限元法是一种将求解域分割为许多小的子域，并在子域上构建局部基函数的方法。
通过在子域上选取适当的基函数，将原偏微分方程化简为子域上的代数方程组，再通过解代数方程组得到整个区域上的数值解。
runge定理
runge定理
Runge定理，又称Runge-Kutta定理，是数值分析中一种重要的定理，源于德国数学家K. Runge和G. Kutta的名字。
它是用于求解微分方程的一种重要方法。
该定理曾在20世纪初被流行，在数值计算中也有很重要的作用，它表明对于一阶精确的泰勒展开，任何具有正则步长的方法都不能实现更高精度。
Runge定理以Runge和Kutta两位在1901年发表的论文命名，它解决了流行的Euler鞍点问题。
Euler方法是一种用于数值解决微分方程的高效方法，但它有时会出现鞍点，即解会急剧变化。
在论文中，Runge和Kutta通过对方程进行拟合，建立出一个通用的流行的数值解析方法，解决了Euler鞍点问题，而这种方法就是现今流行的Runge定理。
Runge定理的主要思想是把微分方程的拟合拆分成两步，并利用它追踪解的流行。
首先，根据Euler鞍点问题，用某种方法求解解的初始值；其次，将方程分解为更小的步骤，用特定方法对其进行拟合，这样就可以得到需要的精确解。
Runge定理首先用于用于求解一阶微分方程组，但它通常也会被用于求解更高阶方程组，因此它在许多实际数值分析应用中有重要作用，包括科学研究和工程学。
它不仅能用于一般的拟合泰勒数展开，而且还可以作为对不同问题的拟合。
经典Runge-kutta方法数值分析实验报告
>> z=-v+u+1;
>> m=RungeK(z,0.1,0,1,1,[u v])
2 再画出方程准确解的图像，程序语句如下
在命令窗口输入
hold on
x=0:0.1:1;
y=zeros(11,1);
for i=1:11
y(i)=exp(-x(i))+x(i);
end
plot(x,y,'*b')
format long;
N=(b-a)/h;
y=zeros(N+1,1);
y(1)Байду номын сангаасy0;
g(1)=y0;
x=a:h:b;
for i=2:N+1
k(1)=subs(f,var,[x(i-1) y(i-1)]);
k(2)=subs(f,var,[x(i-1)+h/2 y(i-1)+h/2*k(1)]);
数值分析实验报告
姓名
学号
日期
实验项目
经典Runge-kutta方法
指导教师
一、上机实验的问题和要求（需求分析）：
实验要求
掌握求解一阶常微分方程初值问题经典Runge-kutta方法编程运算
实验内容：
用经典Runge-kutta方法在步长选为1/10时计算一阶常微分方程初值问题
（该问题的精确解为 ），
3 计算误差
在命令窗口输入
m-y
四、调试和运行程序过程中产生的问题及采取的措施：
开始画出的图像，忘记了用不同的颜色区分开来；经修改后，图像客观性增强；
五、运行输出结果及分析：
结果如下图
结果分析：
第4节 单步法 Runge-Kutta法
m
ki = f (t + hai , un + h∑ bij k j ) i = 1,L, m,
其中
i =1
m
tn = t0 + nh , ( n = 0, 1,L, ) , h为积分步长，un为u(tn)的
j =1
近似值， B=(bij)(i,j=1,…,m)称为Runge-Kutta矩阵，其中在 上面公式中未出现的元素定义为零，而a=(a1,a2,…,am)T 和 c=(c1,c2,…,cm)T 分别称为Runge-Kutta节点和Runge-Kutta 权。且这些系数要满足以下条件：
u ′ = f (t , u )
u (t0 ) = u0
的解充分光滑，将 u(t) 在t0处用Taylor公式展开：
h2 h p ( p) u(t ) = u(t 0 ) + hu′(t 0 ) + u′′(t 0 ) + L + u (t 0 ) + O h p+1 2! p!
其中
u(t0 ) = u0 ,
5.2 单步法的收敛性 将初值问题写成积分形式
u (t + h ) u (t ) = ∫
t +h t
f (τ , u (τ )) dτ
u (t0 ) = u0
(5.6)
如果有某已确定的函数F(t, u, h)（通过某种离散化），使初值 问题的任意解u(t) 满足 (5.7) u ( t + h ) u ( t ) = hφ ( t , u ( t ) , h ) + O ( h p +1 ) , 其中p≥1是使(5.7)成立的最大整数，则称算法
然后将导数和h=0.25代入（5.5）式，计算出：
线性多步法 - 龙格-库塔法
河北联合大学第9章常微分方程初值问题数值解法§9.3 龙格-库塔法 §9.4 线性多步法1. 什么是p 阶Runge-Kutta 法？写出经典的四阶Runge-Kutta 法。
答：考虑增加计算),(y x f 的次数，就有可能构造出精度较高的近似公式。
构造如下形式的计算公式：⎪⎪⎩⎪⎪⎨⎧=++==+=∑∑-==+),,3,2( ),( ),(11111p i K b h y h a x f K y x f K K c h y y i j j ij n i n i n n p i i i n n （1）其中i ij i c b a ,,都是参数，确定它们的原则是使近似公式在),(n n y x 处的Taylor 展开式与)(x y 在n x 处的Taylor 展开式的前面的项尽可能多地重合。
将式（1）称为p 阶显式Runge-Kutta 方法，简称R-K 方法。
取4=p ，通过更复杂的计算，式（1）成为⎪⎪⎪⎪⎩⎪⎪⎪⎪⎨⎧++=++=++==++++=+),( )2,2( )2,2( ),( )22(6342312143211hk y h x f K K h y h x f K K h y h x f K y x f K K K K K h y y n n n n n n n n n n 上式称为经典形式的四阶R-K 公式，是常用的四阶R-K 公式。
局部截断误差为5()O h 。
2. 如何导出线性多步法的公式？它与单步法有何区别？答：（1）利用Taylor 展开构造线性多步方法线性多步法(Linear-Multistep Method)的一般形式为∑∑-=-=-++=ri i n i r i i n i n f h y y 101βα （2）其中i i βα,均为常数，),(k k k y x f f =),,,1(r n n n k -+= 。
若01=-β，公式（2）为显式法；01≠-β，公式为隐式法。
龙格-库塔(Runge-Kutta)方法
——常微分方程及方程组的数值解 常微分方程及方程组的数值解
基本概念 只有一个自变量的微分方程为常微分方程 ODE, 否则称为偏微分方程 PDE。 。 方程中未知函数导数的最高阶数称为方程的阶
[1] [2] C. Runge, "Ueber die numerische Auflsung von Differentialgleichungen" Math. Ann. , 46 (1895) pp. 167–178 W. Kutta, "Beitrag zur naherungsweisen Integration von Differentialgleichungen" Z. Math. und Phys. , 46 (1901) pp. 435–453
2. 递推化 在具有唯一解的条件下， 在具有唯一解的条件下，通过步进法逐步计算出解 在一系列离散点上的值。从而得到原 的数值近似解。 在一系列离散点上的值。从而得到原ODE的数值近似解。 的数值近似解
初值问题解法 讨论一阶ODE（高阶可化为一阶ODEs） （高阶可化为一阶 讨论一阶 ） 的初值问题。初值问题可以一般地写成： 的初值问题。初值问题可以一般地写成：
欧拉（ 欧拉（Euler）方法 ）
Euler方法是求解初值问题的最简单方法， 方法是求解初值问题的最简单方法， 方法是求解初值问题的最简单方法 精度差。然而对理论分析很有用。 精度差。然而对理论分析很有用。Runge-Kutta 法是对Euler法的改进 法是对 法的改进
Euler方法： 方法： ym+1 = ym + hK1 + O(h ) K1 = f ( xm , ym )
K 1 = f ( xm , ym ) K i = f ( xm + hai , ym + h∑ bil K l )
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
和(1)式一致,即改进Euler公式,也称为二阶Runge-Kutta法
华长生制作
8
三、高阶Runge-Kutta方法
如果[ xn 1 , xn ]上增加一点x
1 n 1 2
x
n
1 2
h xn 1 2
且以y( x)在xn 1、x
1 n 2
和xn处的斜率K1、K2和K3的加权平均值
3
以下使用中值定理进行推导
为了同学们课后复习的方便,以下的内容将k写成n.
一、Runge-Kutta方法的导出 对于常微分方程的初值问题
y f ( x , y ) y( a ) y0 a xb
的解 y y(x), 在区间 xn 1 , xn ]上使用微分中值定理有 [ ,
K1 f ( xk 1 , yk 1 )
K2 f ( xk , yk 1 hK1 ) y0 y( x0 )
----------(1)
华长生制作
2
改进Euler法是由梯形公式和Euler公式复合而成
梯形公式具有2阶精度 同样可以证明,改进Euler法也具有2阶精度 形如(1)式的求解公式称为二阶Runge-Kutta法
如下图
y y(x)
则(4)式化为
yn yn 1 hf ( xn 1 , yn 1 ) ----(5)
K K
即Euler方法
由于
华长生制作
en (h) O(h2 )
xn 1
xn
x
Euler方法也称为一阶Runge-Kutta方法
6
2.如果以y( x)在xn 1和xn处的斜率K 1和K 2的算术平均值作为 y( x)在[ xn 1 , xn ]上的平均斜率
K1 y( xn 1 ) f ( xn 1 , yn 1 ) K2 y( xn ) [ xn , y( xn )] f ( xn , yn )
(由(5)式)
y
y y(x)
K2
K1
K
f ( xn , yn 1 hK1 )
K1 K2 K 2
xn 1
华长生制作
----(9)
因而三阶R-K方法(7)具有3阶精度 14
类似于(7)式,还可构造四阶(经典)Runge=Kutta方法
h yn yn 1 ( K 1 2 K 2 2 K 3 K 4 ) 6 K1 f ( xn 1 , yn 1 ) h h K 2 f ( xn 1 , yn 1 K 1 ) 2 2 h h K 3 f ( xn 1 , yn 1 K 2 ) 2 2 K4 f ( xn 1 h, yn 1 hK3 )
-----------(10)
y0 y( x0 )
en (h) O(h5 )
华长生制作
因而方法(10)有4阶精度
15
例1. 使用高阶R-K方法计算初值问题
y y 2 0 x 0.5 y(0) 1
取h 0.1.
解:
(1) 使用三阶R-K方法(7)
2 K1 y0 1
取
1 4 1 K K1 K2 K3 6 6 6
则(4)式化为
h yn yn 1 ( K 1 4 K 2 K 3 ) 6 K1 f ( xn 1 , yn 1 )
-----------(7) h h K 2 f ( xn 1 , yn 1 K 1 ) 2 2 K3 f ( xn 1 h, yn 1 h(2 K2 K1 ))
y( xn 1 ) hy( xn 1 ) h(2 K2 2 K1 ) f y( xn 1 , yn 1 )
h2 [ y( xn 1 ) ( 2 K 2 2 K1 ) f xy ( xn 1 , yn 1 ) 2 ( 4 K 2 4 K1 K 2 ) f y( xn 1 , yn 1 )] O(h3 )
(2) 如果使用四阶R-K方法(10)
n 1时
2 K1 y0 1
0 .1 K 2 ( y0 K 1 )2 1.1025 2
华长生制作 17
0.1 K 3 ( y0 K 2 )2 1.1133 2
K4 ( y0 0.1K3 )2 1.2351
y1 y0 0.1 ( K 1 2 K 2 2 K 3 K 4 ) 1.1111 6
y( xn ) y( xn 1 ) y(n 1 )(xn xn 1 )
n1 ( xn1 , xn )
即
华长生制作
y( xn ) y( xn1 ) hy(n1 )
----------(3)
4
引入记号
y( xn ) y( xn1 ) hK
K y( n 1 ) f [ n 1 , y( n 1 )] yn yn 1 hK
华长生制作
xn 1
xn
5
x
二、低阶Runge-Kutta方法
1.如果以y( x)在xn 1处的斜率作为 ( x)在[ xn 1 , xn ]上的平均斜率 y K
即
K y( xn 1 ) f [ xn 1 , y( xn 1 )] f ( xn 1 , yn 1 )
y
华长生制作 12
K3 f ( xn 1 h, yn 1 h(2 K2 K1 ))
f ( xn 1 , yn 1 ) hf x( xn 1 , yn 1 ) h(2 K2 K1 ) f y( xn 1 , yn 1 )
1 2 [ h f x( xn 1 , yn 1 ) h 2 ( 2 K 2 K 1 ) f xy ( xn 1 , yn 1 ) 2! h2 ( 2 K 2 K 1 )2 f y( xn 1 , yn 1 )] O( h 3 )
再作y y( x)在xn 1处的Tailor展开式
h2 h3 y( xn 1 h) y( xn 1 ) hy( xn 1 ) y( xn 1 ) y( xn 1 ) O( h 4 ) 2 3!
比较(7)(8)两式,可知
y( xn ) yn en (h) O(h4 )
令
K3 y( xn ) f ( xn 1 h, yn 1 h(2 K2 K1 ))
参照Simpson求解公式
h h h yn yn 1 [ f ( xn 1 , yn 1 ) 4 f ( xn 1 , y( xn 1 )) f ( xn , yn )] 6 2 2 10 华长生制作
华长生制作
xn 1
xn 1
h 2
xn
9
x
如果以xn 1处的斜率K1预测y( x
y( x
) yn 1 h K 1 1 n 2 2
1 n 2
)
前进Euler公式
h h 令 K 2 y( x 1 ) f ( xn 1 , yn 1 K 1 ) n 2 2 2 h 同样以 xn 1、xn 1 处的斜率 K 1、K 2预测 y( xn ) 2 y( xn ) yn 1 h(2 K2 K1 )
其余结果如下:
n 1.0000 2.0000 3.0000 4.0000 5.0000
华长生制作
xn k1 k2 k3 k4 yn 0.1000 1.0000 1.1025 1.1133 1.2351 1.1111 0.2000 1.2346 1.3756 1.3921 1.5633 1.2500 0.3000 1.5625 1.7639 1.7908 2.0423 1.4286 0.4000 2.0408 2.3428 2.3892 2.7805 1.6667 0.5000 2.7777 3.2600 3.3476 4.0057 2.0000
2 1 h2 h2 K1 h2 K1 [ f x( xn 1 , yn 1 ) f xy ( xn 1 , yn 1 ) f y( xn 1 , yn 1 )] O( h 3 ) 2! 4 4 4
h h2 y( xn 1 ) y( xn 1 ) y( xn 1 ) O( h 3 ) 2 8
xn
x
令
则(4)式化为
华长生制作 7
h yn yn 1 ( K 1 K 2 ) 2
K1 f ( xn 1 , yn 1 ) K2 f ( xn , yn 1 hK1 )
y0 y( x0 )
-----------(6)
en (h) O(h3 )
§ 9.2 龙格-库塔方法
k 1,2 ,3,
考虑改进Euler法
yk yk 1 hf ( xk 1 , yk 1 ) h yk yk 1 [ f ( xk 1 , yk 1 ) f ( xk , yk )] 2 h y k yk 1 ( K 1 K 2 ) 如果将其改成 2
y0 y( x0 )
(7)式称为三阶Runge-Kutta方法 华长生制作
11
分别作K1、K2、K3在xn 1处的Tailor展开式
K1 f ( xn 1 , yn 1 ) y( xn 1 )
K 2 f ( xn 1 h h , yn 1 K 1 ) 2 2 h hK 1 f ( xn 1 , yn 1 ) f x( xn 1 , yn 1 ) f y( xn 1 , yn 1 ) 2 2
RK.m
n 1时
0 .1 K 2 ( y0 K 1 )2 1.1025 2