数値解析および プログラミング演習 [08 第 9 回目 ]
の解法 - 4. Ruge-Kua( ルンゲ クッタ 法 Ruge-Kua-Gill( ルンゲ クッタ ジル / ギル 法 5. 多段解法
解法の対象 常微分方程式 d( d 初期値条件 (, の変化に応じて変化する の値を求める. ( 0 ( 0 と 0 は,give 0 常微分方程式の初期値問題 と言う. 3
Ruge-Kua 法の導出 4. Ruge-Kua 法 ( まずは 次 目標 元の式 : d( d (, の情報から,D 後の を求める. 4
Ruge-Kua 法の導出 ( まずは 次 Euler 法の解を求める. D (, 4. Ruge-Kua 法 ~ Euler 法の解 D 5
Ruge-Kua 法の導出 ( まずは 次 Euler 法の解を求める. D (, b, ad と の間の点 4. Ruge-Kua 法 ( を考える. ~ Euler 法の解 D 6
Ruge-Kua 法の導出 ( まずは 次 Euler 法の解を求める. D (, 4. Ruge-Kua 法 と の間の点の傾きから解 を求める. (, D b ad ~ Euler 法の解 D 7
Ruge-Kua 法の導出 ( まずは 次 Euler 法の解を求める. D (, と の間での傾きから解 を求める. (, D b ad 3 と の増分から解を求める. w w 4. Ruge-Kua 法 Ruge-Kua 法の解 以降では, 適正な a,b,w,w を求める. wd, ~ Euler 法の解 ( w D 8
9 常微分方程式 ( b, ad を Tayor 展開すると, Ruge-Kua 法の導出 ( まずは 次 4. Ruge-Kua 法, ( D a b! ø ö ç ç è æ D D ø ö ç è æ D ( ( a b a b a b! ø ö ç ç è æ D ø ö ç è æ D ( a ab b a b D (! ø ö ç ç è æ D ø ö ç è æ D D 3 ( a ab b a b w w w w よって,, ( を と表記して, ( であるため, ( (, w w w w D D D a b
0 常微分方程式 Ruge-Kua 法の導出! 3 3! d d d d D D D ( まずは 次 一方,Taylor 展開より, d d d d ( ( d d ø ö ç è æ (- (- (-3 4. Ruge-Kua 法 (!! 3 3! D D D ø ö ç è æ (- と (-3 を (- に代入して, (-4
Ruge-Kua 法の導出 w ( まずは 次 ( と (-4 を比較.O(D 3 は無視. w b w a w w q とすると, w -q w q a b q 一意に決まらず, 無限の組み合わせが存在. 4. Ruge-Kua 法 ( (-4 æ ö D ç a è ø ( w w D w b! æ ç è ö ø D D 3! D 3 (!!
Ruge-Kua 法の導出 w w q とすると, w -q w q a b q ( まずは 次 ( と (-4 を比較.O(D 3 は無視. w b w a w 一意に決まらず, 無限の組み合わせが存在. 4. Ruge-Kua 法 a b であるため, 点 は, 直線 l 上. ~ ( b, ad D l
Ruge-Kua 法 (- と (-3 を (- に代入して,( と比較.O(D 3 は無視. q/ のとき, 修正オイラー法 ( 先週, 紹介 w / w / a b ( 次 4. Ruge-Kua 法 ~ l D 3
Ruge-Kua 法 q のとき, 修正オイラー法のバリエーション w 0 w a b / ( 次 4. Ruge-Kua 法 (- と (-3 を (- に代入して,( と比較.O(D 3 は無視. ~ l D 4
4. Ruge-Kua 法 Ruge-Kua 法 (4 次 同様の手法を用いて,Taylor 展開の 4 次までを使う. 係数の組合せは, 無限に存在. 5
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 3 始点 Aと点 B の中点 M の傾きから Aを始点に解 B 3 を求める. 3 D (, D 4. Ruge-Kua 法 4 点 B 3 の傾きから,Aを始点に解 B 4 を求める. D (, 3 D 4 5 以下を解とする. ( 3 6 4 B 4 B 3 B M B M A D 4 3 6
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 4. Ruge-Kua 法 B A D 7
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから, D (, D 4. Ruge-Kua 法 B A D M 8
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 4. Ruge-Kua 法 B B A D M 9
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 3 始点 A と点 B の中点 M の傾きから 4. Ruge-Kua 法 B A M D 0
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 4. Ruge-Kua 法 3 始点 Aと点 B の中点 M の傾きから Aを始点に解 B 3 を求める. 3 D (, D B 3 B 3 A M D
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 3 始点 Aと点 B の中点 M の傾きから Aを始点に解 B 3 を求める. 3 D (, D 4 点 B 3 の傾きから, 4. Ruge-Kua 法 B 3 A D
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 3 始点 Aと点 B の中点 M の傾きから Aを始点に解 B 3 を求める. 3 D (, D 4. Ruge-Kua 法 4 点 B 3 の傾きから,Aを始点に解 B 4 を求める. D (, 3 D 4 A D B 4 4 B 3 3
Ruge-Kua 法 (4 次 古典的な原型の場合 Euler 法の解 B を求める. D (, 始点 A と点 B の中点 M の傾きから,A を始点に解 B を求める. D (, D 3 始点 Aと点 B の中点 M の傾きから Aを始点に解 B 3 を求める. 3 D (, D 4. Ruge-Kua 法 4 点 B 3 の傾きから,Aを始点に解 B 4 を求める. D (, 3 D 4 5 以下を解とする. ( 3 6 4 B 4 4 B 3 B 3 M B M A D 4
4. Ruge-Kua 法 Ruge-Kua 法 (4 次 古典的な原型の場合式だけ抜き出すと. D (, D (, D 3 D (, D D (, 3 D 4 ( 3 6 4 数値積分法シンプソン則に類似 h S ( ( a 4 ( c ( b 6 5
Ruge-Kua-Gill 法 4. Ruge-Kua-Gill 法 計算機の記憶容量を少なくし, 情報落ちの誤差が少なくするように係数の組合せを工夫した方法. 初期の計算機では有効. 現在の計算機では利点にはならない. 6
7 常微分方程式式だけ抜き出すと., ( D, ( D D ø ö ç è æ - - D D, ( ( 3 ø ö ç è æ - D D, ( 0 3 4 4 3 6 ( ( 6 - Ruge-Kua-Gill 法 4. Ruge-Kua-Gill 法
Ruge-Kua-Gill 法 さらに工夫して, q 0 0 D (, D ( h, D 3 q q0 ( - q0-4. Ruge-Kua-Gill 法 h h - æ ö 3 D çh, D è ø q q 3( - ( - q - ( - h3 h ( ( 3 - q q 0 h ( - ( - q 8
Ruge-Kua-Gill 法 ( 続き ( h D, D 4 4. Ruge-Kua-Gill 法 q 3 q 3( ( 3 - q - ( h4 h3 ( 4 - q3 6 q4 q3 ( 4 - q3-4 h 4 3 9
5. 多段法 一段法 オイラー法, 修正オイラー法, ルンゲクッタ法多段法 に加えて, -, -, を用いて計算する方法. と - を用いて, を計算する方法 ( 段, 陽解法 について以下説明. 30
5. 多段法 一段法であるオイラー法では, ~ まで, (, 一定と考えている. * 真の解 解析解 D D 3
5. 多段法 多段法 ( 段の陽解法 では,, - における の補間関数 (Lagrage 補間法 を作成. 補間関数を用いて, ~ にて積分を実施. つまり, ò F( d ここで,F(: 補間関数. - F( ( -, - - - - - - (, - ( ( * 真の解 解析解 - - - 外挿積分 D D 3
33 常微分方程式 5. 多段法多段法 ( 段の陽解法 では, 式 (( より, 積分 ( に関する一次関数 を実施して,, (, ( 3 - - - D D
5. 多段法 このように,(,, -, -, の p 個の値に対して,p- 次の Lagrage 補間を行い, を求める方法を Admas 法と言う. を用いない場合 陽解法 を用いて, を, 収束計算により求める場合 陰解法 34
陽解法まとめ Euler 法陽解法 ( の情報で計算可能 最も簡潔なアルゴリズム微分方程式が安定していても, 刻み幅 Dを大きくすると, 安定した解が得られない ( 誤差が有界な範囲外に発散する ことがある. 累積誤差のオーダー :O(D. 修正 Euler 法 (Ruge-Kua 法 次 陽解法 ( の情報から計算が可能 簡便なアルゴリズム. 累積誤差のオーダー :O(D. Ruge-Kua 法 (4 次 陽解法 ( の情報から計算が可能 簡便なアルゴリズム. 累積誤差のオーダー :O(D 4. 35
演習について 演習問題 4の一つ目現状のコードは, 以下の挙動を示している. d 0. d d, d 0 とおくと, 方程式は d d, -0. d 意味は, バネ d 5 0 5 0 単振動 バネ定数 0. 質量 m -5-0 -5 0 0 40 60 / s 36
バネとダッシュポット d d c d d d, d 0 とおくと, 方程式は d d, -c - d d 意味は, バネとダッシュポット バネ定数 質量 m 粘性減衰定数 c 37
提出課題について d d 0 d d d, d 0 とおくと, 方程式は d d, - - 0 d d 意味は, バネとダッシュポット バネ定数 0 質量 m 粘性減衰定数 c ( 通常 c は定数, 今回は特別に, 時間とともに増加すると考える. 38