AZNPRUNGEKUTTAMETHOD：　　数値計算の説明　ルンゲ・クッタ法



Return to Top　数値計算の説明　オイラー法　ルンゲ・クッタ法　差分法　図形表現のページ　参考文献　
Last Updated　2002-6/01


数値計算の説明　ルンゲ・クッタ法

数値計算を行う方法には、
　　　　オイラー法、ルンゲ・クッタ法、差分法、有限要素法などがあり、計算対象とする
　　　　現象を表現する方程式の種類や、必要とする精度等に応じて使い分けます。

常微分方程式は、オイラー法やルンゲ・クッタ法で数値解を求めることができます。
常微分方程式には、

　　　  １階微分方程式　　　
　　         dy/dx=2(x+1)y

　　　  連立１階非線形微分方程式　　　
    　　    dx/dt=ax-bxy
        　　dy/dt=-cy+dxy

　　　  ２階微分方程式
　　        Md2x/dt2+pdx/dt+qx=0

　　　  連立２階微分方程式
　　　　　　M1d2x/dt2+hdx/dt+k1x-k2(y-x)=Asinλt
        　　M2d2y/dt2           +k2(y-x)=0　　　

があります。


-------------------------------------------------------------------------------

ルンゲ・クッタ法：

導関数　dy/dx=f(x,y)　が与えられているとき、関数　ｙを求めることを目的とします。

ルンゲ・クッタ法の基本的な、大きい手順は、オイラー法と同じですが、
オイラー法に較べて、傾きの求め方を調整・工夫します。

すなわち、例えば、y3の値は、位置x2におけるy2の値が前回計算で求められてあり、
位置x2における傾きは与えられた導関数dy/dxから 調整した 傾き k を算出することとするので、
刻み幅　h　を調整した傾き k に乗算することにより、
y2をこれに加えて、y3の値とする。以下同様の手順を限界まで繰り返す。

[ y3の値 = 位置x2におけるy2の値 + 位置x2における調整した傾き k * 刻み幅h ]　

この手順において、x0におけるy0の値を初期値として与える必要がある。
初期値を与えれば、全てのxに対してyの値が算出可能である。
なお、初期条件の個数としては、１階の微分方程式では１つ。２階の微分方程式では２つ必要です。

計算手順：

１．x1=x0+h, y1=y0+k0・h
２．x2=x1+h, y2=y1+k1・h
３．x3=x2+h, y3=y2+k2・h
　　・・・
Ｎ．xn+1=xn+h, yn+1=yn+kn・h

　　　ただし、調整する傾き kn は下式にて求める。
　　　kn=(k1+2k2+2k3+k4)/6
　　　ここで、

　　　k1=f(xn+0/1,yn+000/1)
　　　k2=f(xn+h/2,yn+kn1/2)
　　　k3=f(xn+h/2,yn+kn2/2)
　　　k4=f(xn+h/1,yn+kn3/1)
　　　の計算を付加的に行う必要がある。
　　　
必要回数にて、計算を終了する。

算出された(x1,y1),(x2,y2),(x3,y3),(x4,y4),,,(xi,yi),,,(xn,yn)を
平面上にプロットすると図形が描かれる。
こうして解としての関数 y　が図形上に求められたことになる。

以上がルンゲ・クッタ法です。考え方や、プログラミングはやや複雑ですが、
計算精度を必要とする場合に使われます。


-------------------------------------------------------------------------------


参考文献

//数値計算

//鷹尾「数値計算のはなし」日科技連
//登坂「偏微分方程式の数値シミュレーション」東京大学出版会
//水本「ＦＯＲＴＲＡＮによる数値計算法入門」近代科学社
//山崎「偏微分方程式の数値解法入門」森北出版

AZURE  Return to Top