FC2ブログ

スポンサーサイト

上記の広告は1ヶ月以上更新のないブログに表示されています。
新しい記事を書く事で広告が消せます。
スポンサー広告 | --/--/-- --:--

2次ルンゲ・クッタ法(ホイン法・中点法)

お次は、2次ルンゲ・クッタ

導きたい最終的な式を\[
y_{n+1} = y_n + h(w_1 k_1 + w_2 k_2)
\tag{1}
\]とおく。ここで、
\[
k_1 = f(x_n, y_n)
\tag{2.1}
\]\[
k_2 = f(x_n + \alpha h, y_n + \beta k_1 h)
\tag{2.2}
\]とし、$w_1$、$w_2$、$\alpha$、$\beta$ は定数とする。

$k_2$ をテーラー展開。\[
k_2 = f + \alpha h f_x + \beta k_1 h f_y + \mathcal{O}(h^2)
\tag{3}
\]以下、単に、$f$ や $f_x$ などと書いてあるのはすべて、$(x_n, y_n)$ での値を表すものとする。

(2.1)と(3)を(1)に代入。\[
y_{n+1} = y_n + h(w_1 + w_2)f + h^2 w_2 (\alpha f_x + \beta f_y f) + \mathcal{O}(h^3)
\tag{4}
\]
一方、真の解 $y(x)$ のテーラー展開は、\[
y_{n+1} = y_n + h y'(x_n) + \frac{h^2}{2} y''(x_n) + \mathcal{O}(h^3)
\tag{5}
\]ここで、元の微分方程式から\[
y'(x_n) = f
\]\[
y''(x_n) = \frac{df}{dx} = f_x + f_y y' = f_x + f_y f
\]であるから、(5) は、\[
y_{n+1} = y_n + h f + \frac{h^2}{2} (f_x + f_y f) + \mathcal{O}(h^3)
\tag{6}
\]となる。 (4) と (6) を比較して、h の2次の範囲で両者が一致するためには、\[
w_1 + w_2 = 1
\tag{7.1}
\]\[
w_2 \alpha = \frac{1}{2}
\tag{7.2}
\]\[
w_2 \beta = \frac{1}{2}
\tag{7.3}
\]であればよい。
4つの未知数に対して、3つしか条件がないから、解は無数にある。

一例として、\[
w_1 = 0, \ w_2 = 1, \ \alpha = \beta = 1/2
\]とすると、中点法の式が得られる。
\[
k_1 = f(x_n, y_n) \\
k_2 = f(x_n + h/2, y_n + hk_1/2 ) \\
y_{n+1} = y_n + h k_2
\tag{8}
\]この方式では、中点における導関数の値を評価して、
それによって次の値を求めていることになる。

また、\[
w_1 = w_2 = 1/2, \ \alpha = \beta = 1
\]とすると、ホイン法の式が得られる。
\[
k_1 = f(x_n, y_n) \\
k_2 = f(x_n + h, y_n + k_1 h) \\
y_{n+1} = y_n + \frac{h}{2} ( k_1 + k_2 )
\tag{9}
\]この方式では、積分表示\[
y_{n+1} = y_n + \int_{x_n}^{x_{n+1}} f(x,y(x)) dx \tag{10}
\]において、積分を台形近似\[
\int_{x_n}^{x_{n+1}} f(x,y(x)) dx \simeq \frac{h}{2} \{ f(x_n,y_n) + f(x_{n+1}, y_{n+1}) \}
\]により求めていることになる。

この調子で、3次、4次と進んでいくわけですが、
このままいくと、計算がどんどん複雑になっていって、
とてもできる自信がありませんね・・・^^;
途中で挫折するかもです・・・


参考文献
[1] 川上一郎 「数値計算の基礎」 平成21年, http://www7.ocn.ne.jp/~kawa1/
[2] D. Young and R. T. Gregory, "A Survey of Numerical Mathematics Vol. I"
スポンサーサイト
ジャンル:[学問・文化・芸術]  テーマ:[数学
数値計算>微分方程式 | コメント(0) | 2014/10/27 22:39
コメント

管理者のみに表示

上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。