4次ルンゲ・クッタ法で常微分方程式を解く方法｜公式・アルゴリズム・C++実装例

ルンゲ・クッタ法（Runge-Kutta法）は、常微分方程式（ODE）を数値的に解くために広く使われている手法です。この方法では、x と y を含む導関数 dy/dx を用い、初期値 y(0)（x=0 のときの y の値）が必要になります。これらをもとに、任意の x に対する y の近似値を高精度に求めることができます。

常微分方程式を解くには、以下の公式に従います。

ここで h は区間の幅（刻み幅）を表します。

補足： これらの公式のうち、最初の2つの k1 と k2 だけを使えば、2次のルンゲ・クッタ法として ODE を解くことも可能です。

k1〜k4 の役割

4次ルンゲ・クッタ法では、1ステップあたり4つの傾きを計算し、それらを重み付き平均して y の増分を求めます。

• k1： 区間の始点における傾き
• k2： 区間の中点における傾き（k1 を利用して推定）
• k3： 区間の中点における傾き（k2 を利用して再推定）
• k4： 区間の終点における傾き（k3 を利用して推定）

最終的な増分は (k1 + 2k2 + 2k3 + k4)/6 として計算されます。この重み付けにより、4次ルンゲ・クッタ法は1ステップあたり O(h⁵)、全体で O(h⁴) の精度を持つことが知られており、オイラー法などよりもはるかに正確な近似解が得られます。

入力と出力

Input:
x0 と f(x0)：0 と 0
x の値 = 0.4
h の値 = 0.1
Output:
微分方程式の解：0.0213594

アルゴリズム

rungeKutta(x0, y0, x, h)

入力 − 初期値 x0 と y0、求めたい x の値、刻み幅 h。

出力 − 指定した x における y の値。

Begin
    iteration := (x − x0)/h
    y = y0
    for i := 1 to iteration, do
        k1 := h*f(x0, y)
        k2 := h*f((x0 + h/2), (y + k1/2))
        k3 := h*f((x0 + h/2), (y + k2/2))
        k4 := h*f((x0 + h), (y + k3))
        y := y + (1/6)*(k1 + 2k2 + 2k3 + k4)
        x0 := x0 + h
    done
    return y
End

C++による実装例

#include <iostream>
using namespace std;

double diffOfy(double x, double y) {
    return ((x*x)+(y*y)); // 関数 x^2 + y^2
}

double rk4thOrder(double x0, double y0, double x, double h) {
    int iteration = int((x - x0)/h);      // 反復回数を計算
    double k1, k2, k3, k4;
    double y = y0;       // 初期状態では y は f(x0)

    for(int i = 1; i<=iteration; i++) {
        k1 = h*diffOfy(x0, y);
        k2 = h*diffOfy((x0+h/2), (y+k1/2));
        k3 = h*diffOfy((x0+h/2), (y+k2/2));
        k4 = h*diffOfy((x0+h), (y+k3));

        y += double((1.0/6.0)*(k1+2*k2+2*k3+k4));       // Δy を用いて y を更新
        x0 += h;       // x0 を h ずつ進める
    }
    return y;       // f(x) の値
}

int main() {
    double x0, y0, x, h;
    cout << "x0 と f(x0) を入力："; cin >> x0 >> y0;
    cout << "x を入力："; cin >> x;
    cout << "h を入力："; cin >> h;
    cout << "微分方程式の解：" << rk4thOrder(x0, y0, x, h);
}

実行結果

x0 と f(x0) を入力：0 0
x を入力：0.4
h を入力：0.1
微分方程式の解：0.0213594

このように、4次ルンゲ・クッタ法は比較的少ない計算量で高い精度の近似解が得られるため、科学技術計算や物理シミュレーションなどの分野で標準的な数値解法として広く利用されています。