ルンゲ＝クッタ法のリスト

数値解析 > 常微分方程式の数値解法 > 線型多段法 > ルンゲ＝クッタ法 > ルンゲ＝クッタ法のリスト

ルンゲ＝クッタ法は、以下の形の常微分方程式の初期値問題の解を数値で近似計算する方法である。

y'=f(t,y),\;y(t_{0})=y_{0}

一般的に、ルンゲ＝クッタ法は以下の形で与えられる。

y_{n+1}=y_{n}+h\sum _{i=1}^{s}b_{i}k_{i}\,

ただし、

k_{i}=f\left(t_{n}+c_{i}h,y_{n}+h\sum _{j=1}^{i-1}a_{ij}k_{j}\right),

以下のリストで記述するすべての計算方法は、それに対応するブッチャー配列で与えられる。ある一つの方法に対する係数をブッチャー配列で以下の形で表す。

{\begin{array}{c|cccc}c_{1}&a_{11}&a_{12}&\dots &a_{1s}\\c_{2}&a_{21}&a_{22}&\dots &a_{2s}\\\vdots &\vdots &\vdots &\ddots &\vdots \\c_{s}&a_{s1}&a_{s2}&\dots &a_{ss}\\\hline &b_{1}&b_{2}&\dots &b_{s}\\\end{array}}

また、陽的ルンゲ＝クッタ法に対応するルンゲ＝クッタ行列は狭義の下三角行列であるので上三角成分の表記は省略される。

陽的ルンゲ＝クッタ法[編集]

（前進）オイラー法[編集]

オイラー法は1次の方法である。安定性と精度が低いため、オイラー法は入門の例でしか使われない（実際は大規模な計算では今でも使われることがある）。

{\begin{array}{c|c}0&\\\hline &1\\\end{array}}

陽的中点法[編集]

（陽的）中点法は2段2次の方法である（陰的中点法も参照）。

{\begin{array}{c|cc}0&&\\1/2&1/2&\\\hline &0&1\\\end{array}}

ホイン法[編集]

ホイン法（英語版）も2段2次の方法である (陽的台形公式としても知られている)。

{\begin{array}{c|cc}0&&\\1&1&\\\hline &1/2&1/2\\\end{array}}

Ralston法[編集]

Ralston法は2段2次の方法のうちで局所誤差の上界が最小のものである^[1]。

{\begin{array}{c|cc}0&&\\2/3&2/3&\\\hline &1/4&3/4\\\end{array}}

一般的な2段2次の方法[編集]

{\begin{array}{c|ccc}0&&\\x&x&\\\hline &1-{\frac {1}{2x}}&{\frac {1}{2x}}\\\end{array}}

クッタの3次の方法[編集]

{\begin{array}{c|ccc}0&&&\\1/2&1/2&&\\1&-1&2&\\\hline &1/6&2/3&1/6\\\end{array}}

^[2]

古典的ルンゲ＝クッタ法[編集]

{\begin{array}{c|cccc}0&&&&\\1/2&1/2&&&\\1/2&0&1/2&&\\1&0&0&1&\\\hline &1/6&1/3&1/3&1/6\\\end{array}}

クッタの3/8公式[編集]

この方法は、上記の古典的方法と同じ論文で提出されたが^[3]、古典的方法に比べるとあまり用いられていない。

{\begin{array}{c|cccc}0&&&&\\1/3&1/3&&&\\2/3&-1/3&1&&\\1&1&-1&1&\\\hline &1/8&3/8&3/8&1/8\\\end{array}}

埋め込み型ルンゲ＝クッタ法[編集]

詳細は「ルンゲ＝クッタ法#埋め込み型ルンゲ＝クッタ法」を参照

埋め込み型の方法はルンゲ＝クッタ法の局所誤差を推定するために開発された方法である。それらの方法は誤差を制御するために刻み幅を調整する。

埋め込み型方法に対応するブッチャー配列は以下のように与えられる。

{\begin{array}{c|cccc}c_{1}&a_{11}&a_{12}&\dots &a_{1s}\\c_{2}&a_{21}&a_{22}&\dots &a_{2s}\\\vdots &\vdots &\vdots &\ddots &\vdots \\c_{s}&a_{s1}&a_{s2}&\dots &a_{ss}\\\hline &b_{1}&b_{2}&\dots &b_{s}\\&b_{1}^{*}&b_{2}^{*}&\dots &b_{s}^{*}\\\end{array}}

ここで、上側の段の係数 $b i$ は $p$ 次陽的方法に対応するものであり、下側の段の係数 $b * i$ は $p -1$ 次陽的方法に対応するものである。

ホイン・オイラー法[編集]

この方法は、2次のホイン法と1次のオイラー法を組み合わせる方法であり、もっとも単純な埋め込み型方法である。

{\begin{array}{c|cc}0&\\1&1\\\hline &1/2&1/2\\&1&0\end{array}}

フェールベルグ RK1(2)[編集]

フェールベルグ法は3段で、次数2と1の方法を用いる^[4]。

{\begin{array}{c|ccc}0&\\1/2&1/2\\1&1/256&255/256\\\hline &1/512&255/256&1/512\\&1/256&255/256&0\end{array}}

Bogacki–Shampine法[編集]

Bogacki–Shampine法（英語版）は4段で、次数3と2の方法を用いる^[5]。MATLABのコマンド ode23 はこの方法の実装である^[6]。

{\begin{array}{c|ccc}0&\\1/2&1/2\\3/4&0&3/4\\1&2/9&1/3&4/9\\\hline &2/9&1/3&4/9&0\\&7/24&1/4&1/3&1/8\end{array}}

ルンゲ＝クッタ＝フェールベルグ法[編集]

ルンゲ＝クッタ＝フェールベルグ法は5段で、次数5と4の方法を用いる^[4]。

{\begin{array}{c|cccccc}0&\\1/4&1/4\\3/8&3/32&9/32\\12/13&1932/2197&-7200/2197&7296/2197\\1&439/216&-8&3680/513&-845/4104\\\hline &16/135&0&6656/12825&28561/56430&-9/50&2/55\\&25/216&0&1408/2565&2197/4104&-1/5&0\end{array}}

Cash-Karp法[編集]

Cash-Karp法（英語版）はフェールベルグの最初のアイディアを変型した方法である^[7]。フェールベルグの方法と同じく5段で、次数5と4の方法を用いる。

{\begin{array}{c|cccccc}0&\\1/5&1/5\\3/10&3/40&9/40\\3/5&3/10&-9/10&6/5\\1&-11/54&5/2&-70/27&35/27\\\hline &37/378&0&250/621&125/594&0&512/1771\\&2825/27648&0&18575/48384&13525/55296&277/14336&1/4\end{array}}

ドルマン＝プリンス法[編集]

ドルマン＝プリンス法は6段で、次数5と4の方法を用いる^[8]。MATLABのコマンド ode45 はこの方法を実装したものである^[6]。

{\begin{array}{c|ccccccc}0&\\1/5&1/5\\3/10&3/40&9/40\\4/5&44/45&-56/15&32/9\\8/9&19372/6561&-25360/2187&64448/6561&-212/729\\1&9017/3168&-355/33&46732/5247&49/176&-5103/18656\\\hline &35/384&0&500/1113&125/192&-2187/6784&11/84&0\\&5179/57600&0&7571/16695&393/640&-92097/339200&187/2100&1/40\end{array}}

陰的ルンゲ＝クッタ[編集]

上述の埋め込み型方法と同じく、ブッチャー配列に二つの方法が含まれている場合、表で重ねて書かれた下側の段の係数の方法が誤差をコントロールするためのものとなる。

後退オイラー法（陰的オイラー法）[編集]

後退オイラー法は1次の方法である。この方法は、偏微分方程式である線型拡散方程式の時間方向の離散化に用いた場合には無条件に安定で非振動的な方法である。

{\begin{array}{c|c}1&1\\\hline &1\\\end{array}}

陰的中点法[編集]

陰的中点法は2次方法である。選点法であり、以下のガウス・ルジャンドル法（英語版）の最も簡単な場合である。

{\begin{array}{c|c}1/2&1/2\\\hline &1\end{array}}

ガウス・ルジャンドル法[編集]

これらの方法はガウス求積法に基づいた方法であり、高い次数を持つ（ $s$ 段ガウス・ルジャンドル法の次数は $2 s$ である）。

4次の方法は以下のブッチャー配列で与えられる^[9]。

{\begin{array}{c|cc}{\frac {1}{2}}-{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}&{\frac {1}{4}}-{\frac {\sqrt {3}}{6}}\\{\frac {1}{2}}+{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}+{\frac {\sqrt {3}}{6}}&{\frac {1}{4}}\\\hline &{\frac {1}{2}}&{\frac {1}{2}}\\&{\frac {1}{2}}+{\frac {1}{2}}{\sqrt {3}}&{\frac {1}{2}}-{\frac {1}{2}}{\sqrt {3}}\\\end{array}}

さらに6次の方法に対応する配列は以下で与えられる^[9]。

{\begin{array}{c|ccc}{\frac {1}{2}}-{\frac {\sqrt {15}}{10}}&{\frac {5}{36}}&{\frac {2}{9}}-{\frac {\sqrt {15}}{15}}&{\frac {5}{36}}-{\frac {\sqrt {15}}{30}}\\{\frac {1}{2}}&{\frac {5}{36}}+{\frac {\sqrt {15}}{24}}&{\frac {2}{9}}&{\frac {5}{36}}-{\frac {\sqrt {15}}{24}}\\{\frac {1}{2}}+{\frac {\sqrt {15}}{10}}&{\frac {5}{36}}+{\frac {\sqrt {15}}{30}}&{\frac {2}{9}}+{\frac {\sqrt {15}}{15}}&{\frac {5}{36}}\\\hline &{\frac {5}{18}}&{\frac {4}{9}}&{\frac {5}{18}}\\&-{\frac {5}{6}}&{\frac {8}{3}}&-{\frac {5}{6}}\end{array}}

Lobatto法[編集]

Lobatto法^[10] は主に IIIA、 IIIB と IIIC（古典的な文献によって、記号 I と II は二種類のRadau法にのみ使われる）と呼ばれる三種類の方法を指している。方法の名称は Rehuel Lobatto にちなむ。それらの方法はすべて陰的であり、次数 $2 s -2$ を持ち、係数に対し条件 $c 1 = 0$ と $c s = 1$ を満たす。

Lobatto IIIA法[編集]

Lobatto IIIA法はコロケーション法である。

2次の方法は陰的台形公式として知られ^[11]、以下の配列で与えられる。

{\begin{array}{c|cc}0&0&0\\1&1/2&1/2\\\hline &1/2&1/2\\\end{array}}

さらに4次の方法は以下の配列で与えられる^[12]。

{\begin{array}{c|ccc}0&0&0&0\\1/2&5/24&1/3&-1/24\\1&1/6&2/3&1/6\\\hline &1/6&2/3&1/6\\\end{array}}

これらの方法はどれもA-安定であるが、L-安定やB-安定ではない。

Lobatto IIIB法[編集]

Lobatto IIIB法はコロケーション法ではないけど、非連続的コロケーション法として見ることができる。

2次の方法は以下の配列で与えられる^[13]。

{\begin{array}{c|cc}0&1/2&0\\1&1/2&0\\\hline &1/2&1/2\\\end{array}}

さらに4次の方法は以下の配列であたえられる^[12]。

{\begin{array}{c|ccc}0&1/6&-1/6&0\\1/2&1/6&1/3&0\\1&1/6&5/6&0\\\hline &1/6&2/3&1/6\\\end{array}}

これらの方法はどれもA-安定であるが、L-安定やB-安定ではない。

Lobatto IIIC法[編集]

Lobatto IIIC法も非連続的コロケーション法である。

2次の方法は以下の配列で与えられる^[11]。

{\begin{array}{c|cc}0&1/2&-1/2\\1&1/2&1/2\\\hline &1/2&1/2\\&1&0\\\end{array}}

さらに4次の方法は以下の配列で与えられる^[12]。

{\begin{array}{c|ccc}0&1/6&-1/3&1/6\\1/2&1/6&5/12&-1/12\\1&1/6&2/3&1/6\\\hline &1/6&2/3&1/6\\\end{array}}

これらの方法は、すべてL-安定であり、さらに代数的安定（よってB-安定）でもある。そのため、硬い方程式に対する適切な方法である。

Lobatto IIIC*法[編集]

Lobatto IIIC*法^[13] は文献によって、Lobatto III法^[14]、ブッチャーのLobatto法やLobatto IIIC法としても知られている。

2次の方法は上述の陽的台形公式にあたり（よって陰的方法ではない）、以下の配列で与えられる^[11]。

{\begin{array}{c|cc}0&0&0\\1&1&0\\\hline &1/2&1/2\\\end{array}}

さらに4次の方法は以下の配列で与えられる^[12]。

{\begin{array}{c|ccc}0&0&0&0\\1/2&1/4&1/4&0\\1&0&1&0\\\hline &1/6&2/3&1/6\\\end{array}}

これらの方法は、どれもA-安定でも、L-安定でも、B-安定でもない。

一般化Lobatto法[編集]

上述のLobatto法の係数はすべて一意に定められるため、方法の線型結合も考えられる^[13]。一般的に、3つの実数パラメータ $(\alpha _{A},\alpha _{B},\alpha _{C})$ からなるLobatto係数を以下のようにする。

a_{i,j}(\alpha _{A},\alpha _{B},\alpha _{C})=\alpha _{A}a_{i,j}^{A}+\alpha _{B}a_{i,j}^{B}+\alpha _{C}a_{i,j}^{C}+\alpha _{C*}a_{i,j}^{C*}

但し、

\alpha _{C*}=1-\alpha _{A}-\alpha _{B}-\alpha _{C}

ここで、 $a A i,j$ はLobatto IIIA法に対するルンゲ＝クッタ行列である。

$α A = 2$ 、 $α B = 2$ 、 $α C = -1$ のとき、対応する方法はLobatto IIID法であり、Lobatto IIINW法とも呼ばれる。

2次の方法は以下の配列であたえられる。

{\begin{array}{c|cc}0&1/2&1/2\\1&-1/2&1/2\\\hline &1/2&1/2\\\end{array}}

さらに4次の方法は以下の配列であたえられる。

{\begin{array}{c|ccc}0&1/6&0&-1/6\\1/2&1/12&5/12&0\\1&1/2&1/3&1/6\\\hline &1/6&2/3&1/6\\\end{array}}

これらの方法もすべてL-安定であり、さらに代数的安定（よってB-安定）である。

Radau法[編集]

Radau法^[10] は次数 $2 s -1$ を持ち、A-安定である。しかし、方法に対する計算コストが高い以上、方法の次数が落ちるという恐れもある。

Radau IA法[編集]

Radau IA法の係数 $c i$ は方程式

P_{s}(2x-1)+P_{s-1}(2x-1)=0

の解である。ここで、 $P s$ は $s$ 次ルジャンドル多項式である。

3次の方法は以下の配列で与えられる^[15]。

{\begin{array}{c|cc}0&1/4&-1/4\\2/3&1/4&5/12\\\hline &1/4&3/4\\\end{array}}

さらに5次の方法は以下の配列で与えられる^[15]。

{\begin{array}{c|ccc}0&{\frac {1}{9}}&{\frac {-1-{\sqrt {6}}}{18}}&{\frac {-1+{\sqrt {6}}}{18}}\\{\frac {3}{5}}-{\frac {\sqrt {6}}{10}}&{\frac {1}{9}}&{\frac {11}{45}}+{\frac {7{\sqrt {6}}}{360}}&{\frac {11}{45}}-{\frac {43{\sqrt {6}}}{360}}\\{\frac {3}{5}}+{\frac {\sqrt {6}}{10}}&{\frac {1}{9}}&{\frac {11}{45}}+{\frac {43{\sqrt {6}}}{360}}&{\frac {11}{45}}-{\frac {7{\sqrt {6}}}{360}}\\\hline &{\frac {1}{9}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}\\\end{array}}

Radau IIA法[編集]

Radau IIA法の係数 $c i$ は方程式

P_{s}(2x-1)-P_{s-1}(2x-1)=0

の解である。

3次の方法は以下の配列で与えられる^[15]。

{\begin{array}{c|cc}1/3&5/12&-1/12\\1&3/4&1/4\\\hline &3/4&1/4\\\end{array}}

さらに5次の方法は以下の配列で与えられる^[11]。

{\begin{array}{c|ccc}{\frac {2}{5}}-{\frac {\sqrt {6}}{10}}&{\frac {11}{45}}-{\frac {7{\sqrt {6}}}{360}}&{\frac {37}{225}}-{\frac {169{\sqrt {6}}}{1800}}&-{\frac {2}{225}}+{\frac {\sqrt {6}}{75}}\\{\frac {2}{5}}+{\frac {\sqrt {6}}{10}}&{\frac {37}{225}}+{\frac {169{\sqrt {6}}}{1800}}&{\frac {11}{45}}+{\frac {7{\sqrt {6}}}{360}}&-{\frac {2}{225}}-{\frac {\sqrt {6}}{75}}\\1&{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {1}{9}}\\\hline &{\frac {4}{9}}-{\frac {\sqrt {6}}{36}}&{\frac {4}{9}}+{\frac {\sqrt {6}}{36}}&{\frac {1}{9}}\\\end{array}}

脚注[編集]

^ Ralston 1962
^ Iserles 2008, p. 40
^ Hairer, Nørsett & Wanner (1993, p. 138) refer to Kutta (1901).
^ ^a ^b Fehlberg 1970
^ Bogacki & Shampine 1989
^ ^a ^b Moler 2014
^ Cash & Karp 1990
^ Dormand & Prince 1980
^ ^a ^b Iserles 2008, p. 47
^ ^a ^b Butcher (2008) 小節344にLobattoとRadau求積について詳述がある。
^ ^a ^b ^c ^d Butcher 2008, p. 226
^ ^a ^b ^c ^d Butcher 2008, p. 227
^ ^a ^b ^c Jay
^ Butcher 2008
^ ^a ^b ^c Butcher 2008, p. 225

参考文献[編集]

Bogacki, Przemyslaw; Shampine, Lawrence F. (1989), “A 3(2) pair of Runge–Kutta formulas”, Applied Mathematics Letters 2 (4): 321–325, doi:10.1016/0893-9659(89)90079-7, ISSN 0893-9659 .
Butcher, John C. (2008), Numerical Methods for Ordinary Differential Equations, New York: John Wiley & Sons, ISBN 978-0-470-72335-7 .
Cash, J. R.; Karp, Alan H. (1990). “A variable order Runge-Kutta method for initial value problems with rapidly varying right-hand sides”. ACM Transactions on Mathematical Software (ACM) 16 (3): 201-222. doi:10.1145/79505.79507. .
Dormand, J. R.; Prince, P. J. (1980), “A family of embedded Runge-Kutta formulae”, Journal of Computational and Applied Mathematics 6 (1): 19–26, doi:10.1016/0771-050X(80)90013-3 .
Fehlberg, Erwin (1970). “Klassische Runge-Kutta-Formeln vierter und niedrigerer Ordnung mit Schrittweiten-Kontrolle und ihre Anwendung auf Wärmeleitungsprobleme”. Computing (Springer-Verlag) 6 (1): 61-67. doi:10.1007/BF02241732. .
Hairer, Ernst; Nørsett, Syvert Paul; Wanner, Gerhard (1993), Solving ordinary differential equations I: Nonstiff problems, Berlin, New York: Springer-Verlag, ISBN 978-3-540-56670-0 .
Iserles, Arieh (2008), A First Course in the Numerical Analysis of Differential Equations (Second Edition), Cambridge University Press, ISBN 978-0-521-73490-5 .
Jay, Laurent O.. “Lobatto Methods”. 2016年12月31日閲覧。
Kutta, Martin Wilhelm (1901), “Beitrag zur näherungsweisen Integration totaler Differentialgleichungen”, Zeitschrift für Mathematik und Physik 46: 435–453 .
Moler, Cleve (2014年). “Ordinary Differential Equation Solvers ODE23 and ODE45”. 2016年12月31日閲覧。
Ralston, Anthony (1962). “Runge-Kutta Methods with Minimum Error Bounds”. Mathematics of Computation 16 (80): 431-437. doi:10.2307/2003133. .

[1] Ralston 1962

[2] Iserles 2008, p. 40

[3] Hairer, Nørsett & Wanner (1993, p. 138) refer to Kutta (1901).

[#1-4] Fehlberg 1970

[5] Bogacki & Shampine 1989

[Moler-6] Moler 2014

[7] Cash & Karp 1990

[8] Dormand & Prince 1980

[Iserles200847-9] Iserles 2008, p. 47

[Butcher2008LR-10] Butcher (2008) 小節344にLobattoとRadau求積について詳述がある。

[Butcher2008226-11] Butcher 2008, p. 226

[Butcher2008227-12] Butcher 2008, p. 227

[Jay-13] Jay

[14] Butcher 2008

[Butcher2008225-15] Butcher 2008, p. 225

[1]

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9]

[10]

[11]

[12]

[13]

[14]

[15]