130 likes | 273 Views
微分方程式の数値解法 (具体例その1: オイラー (Euler) の方法 ) 近似誤差を Order( h 2 ) とすると ( h 2 以上の高次の h の関係式を無視すると ) , F(x, y, h ) = h f(x, y) であり, y(x + h ) = y(x) + Φ(x, y, h ) y(x + h ) ≒ y(x) + F(x, y, h ) Y(x + h ) ≒ Y(x) + h f(x, y) と近似可能.. そこで, h 2 以上の高次の h の関係式を無視する 微分方程式の
E N D
微分方程式の数値解法 (具体例その1:オイラー(Euler)の方法) 近似誤差を Order(h2) とすると (h2以上の高次のhの関係式を無視すると), F(x, y, h) = h f(x, y) であり, y(x + h) = y(x) + Φ(x, y, h) y(x + h) ≒ y(x) + F(x, y, h) Y(x + h) ≒ Y(x) + h f(x, y) と近似可能.
そこで,h2以上の高次のhの関係式を無視する 微分方程式の • 数値解析的解法の手法をアルゴリズムで記述すると, • アルゴリズム(オイラーの方法)は, • x = a = x0, x = b = xnの区間情報より, きざみ h =(b-a)/n を決定 • x = x0, y(x) = y(x0) = y0を基に Y1 = y0 + h f(x0, y0) より値を決定 • xi=xi-1 + h より,きざみを1つ前進 • Yi+1 = Yi + h f(xi , Yi) より,近似値 Ynを順次求める • と表現可能.
解の計算精度を向上させるため,近似関数の精度を修正して,近似誤差を Order(h3) とすると(すなわち,h3以上の高次のhの関係式を無視する方法を考えると), 近似関数 F(x, y, h) は F(x, y, h) = h { f(x, y) + (h /2!) df(x, y)/dx } と表現可能. 何故かという疑問には正確に答え難いものの,ここで, 近似関数 F(x, y, h) を F(x, y, h) = h { α f(x, y) + β f(x + γh, y + γh f(x, y)) } なる漸化式で表現できるようα,βおよびγを定める
テーラー展開を利用して,未定係数法などにより,上記のα, βおよび γを求めると,これら間には,次のような関係式が 成立する.すなわち, • α + β = 1 • β γ = 1/2 • なる方程式が得られる.2つの方程式では3つの未知数, α , βおよびγを一意には決定できないので, パラメータ λ を導入して, β = λ とし,パラメータ λを介して, α,βおよびγを • α = 1 - λ β = λ γ = 1 / (2λ) • と表現可能.
そこで,近似関数 F(x, y, h) は y(x + h) = y(x) + h { α f(x, y) + β f(x + γ h, y + γ h f(x, y) ) } = y(x) + h { (1 - λ) f(x, y) + λ f(x + h/(2λ), y + h/(2λ) f(x, y) ) } なる漸化式で表現可能
微分方程式の数値解法 (具体例その2:修正オイラー(Modified Euler)の方法) λ = 1 とすると,修正オイラーの方法が導出. この時,近似関数 F(x, y, h)は, y(x + h) = y(x) + h f(x + h/2, y + h/2f(x, y) ) と表現可能.
アルゴリズム(修正オイラー(Modified Euler)の方法) は, • x = x0の時,Y0 = y(x0) = y0を初期値とする • μ0 = h f(x0, Y0)ν0 = h f( x0 + h/2, Y0 + μ0/2 ) • 従って,Y1 = Y0 + ν0 • xi=xi-1 + h より,きざみを1つ前進 • μi-1 = h f(xi-1, Yi-1)νi-1 = h f( xi-1 + h/2, Yi-1 + μi-1/2 ) • Yi = Yi-1 + νi-1より,近似値 Ynを順次求める • と記述できる.
微分方程式の数値解法 (具体例その3:ホイン(Heun)の方法) λ = 1/2 とすると, ホイン の方法が導出. この時,近似関数 F(x, y, h)は, y(x + h) = y(x) + h { 1/2 f(x, y) + 1/2 f(x + h, y + h f(x, y) ) } と記述.
アルゴリズム(ホインの方法) は, • x = x0の時,Y0 = y(x0) = y0を初期値とする • ξ0 = h f(x0, Y0) η0 = h f( x0 + h, Y0 + ξ0 ) • 従って,Y1 = Y0 + 1/2 (ξ0 + η0) • xi = xi-1 + h より,きざみを1つ前進 • ξi-1 = h f(xi-1, Yi-1)ηi-1 = h f( xi-1 + h, Yi-1 + ξi-1 ) • Yi = Yi-1 + 1/2 (ξi-1 + ηi-1) より,近似値 Ynを順次求める • と記述.
微分方程式の数値解法 • (具体例その4:ルンゲ・クッタ(Runge-Kutta)の方法) • 資料を見て確認してもらうため,ここではアルゴリズムの基に なる漸化式のみを掲載. • ルンゲ・クッタ(Runge-Kutta)の方法は微分方程式の数値解析 的解法を 統一的に扱ったもので,現在最も一般的に利用. • これまで述べてきたように,近似関数のもつ誤差をどのよう に 減少されるかを検討し,原理的には,Order(h5)以上の誤 差を 無視して近似関数を求める手法を採用.
これから説明するルンゲ・クッタ法では,近似関数 Y(x + h) が Y(x + h) = Y(x) + h/6 ( k1 + 2 k2 + 2 k3 + k4 ) 但し, k1 = f(x, Y(x)) k2 = f(x + h/2, Y(x) + h/2 k1) k3 = f(x + h/2, Y(x) + h/2 k2) k4 = f(x + h, Y(x) + hk3) なる近似関数で与えられる. これを一般にルンゲ・クッタ法と呼び,次ページに示すような漸化式で微分方程式の解を 計算することになる.
アルゴリズム(ルンゲ・クッタの方法) は, • x = x0の時,Y0 = y(x0) = y0を初期値とする • k1 = f(x0, Y0)k2 = f(x0 + h/2, Y0 + h/2 k1)k3 = f(x0 + h/2, Y0 + h/2 k2)k4 = f(x0 + h, Y0 + h k3) • 従って,Y1 = Y0 + h/6 ( k1 + 2 k2 + 2 k3 + k4 ) • xi = xi-1 + h より,きざみを1つ前進 • 以下,同様にして,Yiを漸化式より,順次求めていく. • と記述される.
この方法は,次数を明示して,4次のルンゲ・クッタ法と呼ばれることがある.この方法は,次数を明示して,4次のルンゲ・クッタ法と呼ばれることがある. • 同様に, • オイラーの方法: 1次のルンゲ・クッタ法 • 修正オイラーの方法,あるいは,ホインの方法: 2次のルンゲ・ クッタ法 • と呼ばれることがある. 一応,これで,数値解析としては,一通りのメニューを終了することになりました.