微分方程式の数値解析 #3
Heun法を使ってみる。
...snip double precision function nextx(t, x, h) double precision, intent(in) :: t, x, h double precision :: k1, k2 k1 = f(t, x) k2 = f(t + h, x + h * k1) nextx = x + h / 2 * (k1 + k2) end function nextx ...snip
計算してみると、
0.00000000 0.00000000 0.00000000 0.10000000 0.10033469 0.10033467 0.20000000 0.20271007 0.20271004 0.30000000 0.30933630 0.30933625 0.40000000 0.42279328 0.42279322 0.50000000 0.54630257 0.54630249 0.60000000 0.68413689 0.68413681 0.70000000 0.84228846 0.84228838 0.80000000 1.02963859 1.02963856 0.90000000 1.26015814 1.26015822 1.00000000 1.55740734 1.55740772 1.10000000 1.96475845 1.96475966 1.20000000 2.57214793 2.57215162 1.30000000 3.60208942 3.60210245 1.40000000 5.79781591 5.79788372 1.50000000 14.10022034 14.10141995
t=1.5で9e-5弱の相対誤差であり中点法と同程度の誤差を持つ。