2003/12/02 先週の復習から 先週 オイラー法 今週 ルンゲクッタ法 常微分方程式の計算を行う場合に、テキストには、 \frac{du}{dx} = g(x,u) という形のものしか提示されていない。 しかし、二階の微分方程式 \frac{d^2u}{dx^2} + \alpha\frac{du}{dx} + \beta u = f(x,u) が与えられた時に、どうするか.. ? 簡単、新しい、変数を導入して、連立させればよい。 この例では、 \frac{du}{dx} = v とすれば、元の式は、 \frac{dv}{dx} + \alpha v + \beta u = f(x,u) より、 \frac{dv}{dx} = - \alpha v - \beta u + f(x,u) \frac{du}{dx} = v の連立で Okey 一般に、N 階の微分方程式は、N 本の 1 階の微分方程式の 連立方程式を解くことになる。 数値的には、各々の式にオイラー法を適用する # 物理現象では、余り、三次になることはないが、二次の # 方程式 ( 運動方程式 ) の連立になることはある。 == 初期値問題 始めから初期値が与えられていれば Okey 終値しか与えられていないこともある ( 境界値問題 ) 初期値問題の繰り返しで解く == 1 階の近似は、性質が良い ( オイラー法 ) 変なことが起きにくい 問題の性質が分ってない時に用いる しかし、精度が悪い 性質が分っている 精度が欲しい場合 高次の公式を利用する ルンゲ=クッタ ( 4 次の近似 ) 非線形性の強い 低精度で、Δx を小くする # 変化が激しいの、少しずつしか進めない (タフ、ロバスト) # 少ししか進めないなら、公式の精度上てもしょうがない 非線形性の弱い 高精度で、Δx を大きくする ( 速度を稼ぐ ) # 変化が少いの、沢山進める # 沢山進めても誤差がおきないように、高精度 # ex) 宇宙船の打ち上げ # 大気圈内 : 複雑 -> オイラーで # 宇宙: 単純 -> ルンゲクッタで オイラー n 次 -> n 回の f の計算が必要 ? ルンゲ=クッタでは.. 4 次 -> 4 回の f の計算で Okey 5 次 -> 5 回ではできない # 4 次が、精度と計算速度のトレードとして一番相応しい # 日本の Mr. ルンゲ=クッタ # 田中正次 先生 ( 日大工学を退官なされた.. ) == テイラー展開 f(x+\delta x) = f(x) + \frac{\delta x}{1!}f'(x) + .. + \frac{\delta x^n}{n!}f^{(n)}(x) + .. 純粋な線型 f'(x) が 0 ... これは、単純過ぎるので、扱わない 非線型が強い 高次の項目の係数が大きい \delta x の羃がかかっているので、 \delta x を小くすれば、無視できる == ルンゲ=クッタ 法 常微分方程式 \frac{dy}{dx} = f(x,y) を解く、 両辺を区間 [x_0 x_1] = [x_0 x_0+\delta x] で定積分をしてみる。 \int_{x_0}^[x_0+\delta] \frac{dy}{dx} dx = \int_{x_0}^[x_0+\delta] f(x,y) dx f(x_0+\delta x) - f (x_0) = \int_{x_0}^[x_0+\delta] f(x,y) dx より、 f(x_0+\delta x) = f (x_0) + \int_{x_0}^[x_0+\delta] f(x,y) dx つまり \int_{x_0}^[x_0+\delta] f(x,y) dx をどう近似するかが問題。 # これを、\delta x f(x_0,y_0) の一次式で近似すれば、 # オイラー法になる。 ルンゲ=クッタ 法 (ホインの方法 : 二次近似) では、 \int_{x_0}^[x_0+\delta] f(x,y) dx を数値積分と同じ手法で解く。 しかし、数値積分では、被積分関数の値が分っているので、それを利 用するが常微分方程式では、そもそも、その関数の値が解らない =>(誤差を含んでいるが..) 他の方法 ( 例えばオイラー法 ) で関数値の近似値を計算し、それを利用して数値積分を行う。 == ルンゲ=クッタ 法 の計算式 # text p.124-125 k_1 = h f(x_i,y_i) k_2 = h f(x_i+\frac{h}{2},y_i+\frac{k_1}{2}) k_3 = h f(x_i+\frac{h}{2},y_i+\frac{k_2}{2}) k_4 = h f(x_i+h,y_i+k_3) y_{i+1} = y_i + \frac{1}{6} ( k_1 + 2 k_2 + 2 k_3 + k_4 ) で計算する。 # 以下、原理 k_1 = h f(x_i,y_i) k_2 = h f(x_i+\alpha_1 h,y_i+\beta_1 k_1) k_3 = h f(x_i+\alpha_2 h,y_i+\beta_21 k_1+\beta_22 k_2) k_4 = h f(x_i+\alpha_3 h,y_i+\beta_31 k_1+\beta_32 k_2+\beta_33 k_3) y_{i+1} = y_i + \mu_1 k_1 + \mu_2 k_2 + \mu_3 k_3 + \mu_4 k_4 # 4 次の近似では、4 回の計算で近似できる # 5 次の近似では、5 回の計算では近似できない この \alpha_i, \beta_i を決めることにより、様々な公式が作れる # 一般的には \alpha_i < \alpha_j ( i < j ) とするの普通 \alpha, \beta のパラメータを変更して、 精度が高く 計算量が少なく という条件で、探せばよいが、これは解けない (任意パラメータが残 る..)。これに、更に、 「係数の簡単なもの」 という条件を追加したものの一つがルンゲ=クッタ # 古典的ルンゲ=クッタ では、 # \alpha_1 = \beta_1 = \frac{1}{2} # \alpha_2 = \frac{1}{2}, \beta_21 = 0, \beta_22 = \frac{1}{2} # \alpha_3 = 1, \beta_31 = 0, \beta_32 = 0, \beta_33 = 1 # で計算 ## \frac{1}{2} つまり、半分を取るのは、その時に最も精度が高い ## 場合が多いから.. # ルンゲ=クッタ=ジル # メモリを使わない方法 ## これ以外のものは、生き残らなかった.. ## 古典的ルンゲ=クッタ は美しい... == Text p.128 多段階法 # ルンゲ=クッタ法 ( オイラー法 ) は、単段階法 ルンゲ=クッタ法は、区間内の計算は、他と独立 # 区間毎に異なる方法を適用してもよい。 多段階法は、 新しい区間での計算の段階で、以前の計算結果 を利用する。 多段階法の例 ミルン法 ( アダムス法、アダムス=ノルトン法 アダムス=バッシュフォース法 ) # 予想子修正法 多段階法では、ルンゲ=クッタ法(単段階法) に比較して、 計算量を減らしている ( 過去の結果を利用しているので.. ) のに、それほど精度が落ちていない。 == [演習] 次の常微分方程式を、ルンゲ=クッタ法 ( Text. p.135 ) で解け [常微分方程式] \frac{du_1}{x} = -u_2 \frac{du_2}{x} = u_1 [初期値] u_1(0) = 1 u_2(0) = 0 # これの解は、u_1(x) = \cos{x}, u_2(x) = \sin{x} となる。 # 前回の、オイラー法の課題との違いは、二点 # 解法が変った # 1 次が 2 次元に変った # 問題を解く上での違いは、関数の値の計算方法と次数の違い # だけなので、次数の違いをマスターすれば、関数の値の計算の # 違いは、本質的でないので、これで常微分方程式は Okey