ネムネコ法、２次のルンゲクッタ法に勝つ！！（差分の基礎）

ネムネコ法、２次のルンゲクッタ法に勝つ！！（差分の基礎）

 

ネムネコが新たに提案した、差分を用いた１階常微分方程式の初期値問題に対する数値解法は、実際、どれくらい使い物になるのか、次の非線形常微分方程式の初期値問題に使ってみた。

　　

 

ネムネコが提案した方法は、上の微分方程式を

　　

と差分方程式に置き換えて解く方法。

この方法をネムネコ法と、ひとまず、仮称することにする。

（２）式の右辺第１項のは、前進差分法（陽解法）

　　

に由来するものであり、第２項のは後退差分（陰解法）

　　

に由来するもので、（２）式の打ち切り誤差はO(h³）で２次のルンゲ・クッタ法や修正オイラー法の誤差と同程度と考えられる。

 

ネムネコ法は、陽解法である修正オイラー法や２次ルンゲ・クッタ法とは異なり、の値を求めるためには超越方程式（２）を解く必要がある陰解法である。

（２）を（数値的に）解く必要があるので、修正オイラー法や２次のルンゲ・クッタ法のように表計算ソフトなどを使って簡単に微分方程式（１）の数値的な近似解を求められないという欠点がある。

なのですが、自分でプログラムを作ってこの問題を解くという場合には、これは必ずしも欠点にならない。というのは、（２）はニュートン法などを用いることによって簡単に数値的に解くことができるからだ。

 

そこで、プログラムを作り、解いてみた。

そして、Δx=0.1としてx=0〜1まで計算してみた結果は次のとおり。

 

 

厳密解とよく合っているじゃないのよ。我ながら驚きの結果！！

グラフにしてしまえば、厳密解との差がわからないほどよく合っているからね〜。

 

 

しかも、この計算結果は、Δx=0.1として２次のルンゲ・クッタ法を用いて解いた計算結果よりも精度が高いんだケロよ。

CASIOさんの高精度計算サイトの計算結果はコチラ。

 

 

スゴイにゃ、（仮称）ネムネコ法！！

そして、どうも、（仮称）ネムネコ法は、線形常微分方程式よりも非線形常微分方程式で強みを発揮するみたいだね。

 

この計算に使用したプログラムは次のとおり。

 

! ネムネコ法(^^ゞ
parameter (n=10)
real x(0:n) , y(0:n)

eps = 1.e-6
x = 0. ; y=0. ! 初期化
x(0)=0.; y(0)=0. ! 初期条件
dx =0.1 ! 分割の幅

do i = 1, n
    x(i)= i*dx
    yo = y(i-1);
   
    ! ニュートン法で非線形方程式解く
    do k=1, 10
!        yn = yo - (yo - dx*sin(yo)-y(i-1)-dx)/(1-dx*cos(yo))
        yn= yo- (yo-0.5*dx*(sin(yo)+sin(y(i-1)))-y(i-1)-dx)/(1-0.5*dx*cos(yo))
        err = abs(yn-yo)
        if (abs(yn-yo).lt.eps*abs(yo)) exit ! 収束判定
        yo=yn
    end do
   
    y(i)= yn ! 計算結果をセット
end do

! 計算結果の出力
write(*,*) '    x       y(calc)    y(exact)'
do i=0, n
    write(*,100) x(i), y(i) , 2.*atan(x(i)/(2-x(i)))
end do

100 format(f10.7,1x,f10.7,1x,f10.7)

end

 

ちなみに、（１）の厳密解は、

　　

 

などと求めることができる。

なお、

　　

 

画像元：YouTubeの上の動画

さらに、ネムネコ極悪モード！！

画像元：YouTubeの下の動画

わかったら、ほれほれ！！