Subject: ルンゲ=クッタ法の刻み幅調節
Date: Tue, 16 Jul 2002 14:22:46 +0900 (JST)
To: staff
From: OSHIRO Naoki <oshiro>

大城です．

うちの研究室では数値シミュレーションでルンゲ=クッタ法（あるいは場合に
よってはオイラー法）を使っていると思いますが，どちらにせよ算出精度が刻
み幅に依存してしまいますし，近似計算が収束しなくなる場合もあります．そ
こで値によって刻み幅の調節が可能だというルンゲ=クッタ=フェールベルグ法
というものを試してみました．

  参考：ＣとJavaで学ぶ数値シミュレーション：峰村吉泰，森北出版，1999.

y'=y で y=exp(x) の場合の結果です．他にオイラー法などとの比較も出して
います．オイラー法が３つ，４次のルンゲ=クッタがあります．・=フェールベ
ルグ法の場合は４次の場合に比べ少し複雑になるのですが，４次と５次の近似
解を一度に求められるのでその両者の比較で誤差の評価を行なって，これを元
に刻み値を調節（半分，２倍に）するのだそうです．刻みが小さくなるので演
算回数は増えますが，当然その分の精度が向上します．また，変化のない箇所
では刻みが大きくなるので不要な演算が減るはずです．また，５次の値がある
ので刻み制御なしでも通常のルンゲ=クッタ法よりもいい結果が出ています．

示した結果は exp() を近似したものなので，まだそれほど複雑なものではあ
りません（近似誤差が出にくい）が，制御などが入ってくる場合にはここで示
したようなより精度の高い方法を検討する必要があるでしょう．

パラメータ詳細・実行ファイルは ~oshiro/rungekutta.tcl にあります．

   x      y      解析解との誤差
Euler formula:
  0.1 1.10000000 5.17e-03
  0.2 1.21000000 1.14e-02
  0.3 1.33100000 1.89e-02
  0.4 1.46410000 2.77e-02
  0.5 1.61051000 3.82e-02

improved Euler formula:
  0.1 1.10500000 1.71e-04
  0.2 1.22102500 3.78e-04
  0.3 1.34923262 6.26e-04
  0.4 1.49090205 9.23e-04
  0.5 1.64744677 1.27e-03

modified Euler formula:
  0.1 1.10500000 1.71e-04
  0.2 1.22102500 3.78e-04
  0.3 1.34923262 6.26e-04
  0.4 1.49090205 9.23e-04
  0.5 1.64744677 1.27e-03

Runge-Kutta method:
  0.1 1.10517083 8.47e-08
  0.2 1.22140257 1.87e-07
  0.3 1.34985850 3.11e-07
  0.4 1.49182424 4.58e-07
  0.5 1.64872064 6.32e-07

   x      y      解析解との誤差（近似解同士の誤差）
Runge-Kutta-Fehlberg method:
  0.1 1.10517092 9.26e-10(1.23e-08)
  0.2 1.22140276 2.05e-09(1.36e-08)
  0.3 1.34985880 3.40e-09(1.51e-08)
  0.4 1.49182469 5.00e-09(1.67e-08)
  0.5 1.64872126 6.91e-09(1.84e-08)
    ..
  9.5 13359.72576387 1.07e-03(1.49e-04)
  9.6 14764.78037533 1.19e-03(1.65e-04)
  9.7 16317.60586885 1.33e-03(1.82e-04)
  9.8 18033.74344378 1.48e-03(2.01e-04)
  9.9 19930.36878138 1.66e-03(2.22e-04)
 10.0 22026.46394523 1.85e-03(2.46e-04)

   x      y      解析解との誤差（近似解同士の誤差） 刻み切り替え回数 次の値までの演算回数
Runge-Kutta-Fehlberg method (h control):
  0.1 1.10517092 5.65e-12(0.00e+00) 10 18
  0.2 1.22140276 1.70e-13(0.00e+00) 7 15
  0.3 1.34985881 1.40e-11(0.00e+00) 7 15
  0.4 1.49182470 8.73e-12(0.00e+00) 7 15
  0.5 1.64872127 1.99e-11(1.00e-11) 7 15
    ..
  9.5 13359.72682870 9.62e-07(0.00e+00) 17 49
  9.6 14764.78156427 1.28e-06(0.00e+00) 11 43
  9.7 16317.60719640 1.62e-06(0.00e+00) 9 41
  9.8 18033.74492625 1.53e-06(0.00e+00) 19 55
  9.9 19930.37043651 1.73e-06(0.00e+00) 23 61
 10.0 22026.46579295 1.91e-06(0.00e+00) 22 62

---