So-net無料ブログ作成
検索選択

今日のお休みソング、「創聖のアクエリオン」から『オムナマグニ』 [今日のアニソン]

今日のお休みソングは、「創聖のアクエリオン」から『オムナマグニ』です。



この歌詞は、ラテン語かな♪

「創聖のアクエリオン」といえば、やはり、この曲なので、これを埋め込まないわけにはいかないケロ。



nice!(0)  コメント(0)  トラックバック(0) 
共通テーマ:音楽

ニュートン法のプログラム [数値解析]

参考までに、C言語で書いたニュートン法のサンプルプログラムをを以下に示す。

#include <stdio.h>
#include <math.h>

double f(double x) {
    return tan(x)-1.;
}

double df(double x) {
    return 1/(cos(x)*cos(x));
}

double Newton(double x, double eps, int n) {
    double x0;
    int i;

    for (i = 0; i < n; i++) {
        x0=x;
        x=x-f(x)/df(x);
        if (fabs(x-x0) < eps ) return x;
    }
   
    return -1;
}

// f'(x)の計算で中心差分を用いたニュートン法
double Newton2(double x, double eps, int n) {
    double x0, h = 0.001;
    double df;
    int i;

    for (i = 0; i < n; i++) {
        x0=x;
        df=(f(x+h)-f(x-h))/(2.*h); // f'(x)の計算に中心差分を使う
        x=x-f(x)/df;
        if (fabs(x-x0) < eps ) return x;
    }
   
    return -1;
}   

main() {
    double x=0.0, eps =1.e-6;
    int n = 10;

//    f(x)=tan x  1=0を解く
    printf("x=%f\n",Newton(x,eps,n)); // ニュートン法
    x=0.0;
    printf("x=%f\n",Newton2(x,eps,n)); // f'(x)に中心差分を用いたニュートン法
}

収束判定として、|x−x₀|<10⁻⁶、または、計算の反復回数が10を越えたときを用いている。

解く方程式はtan x = 1で、0≦x<π/2とすると、この解は
 x=π/4≒0.785398163

反復回数6回で|x−x₀|<10⁻⁶を満たし、
 x=0.78598
という収束解に到達する。

上のプログラムでは、f'(x)の値を中心差分を用いて計算するニュートン法のプログラムもあわせて紹介している。
中心差分を用いると、微分係数f'(x)は

  

と近似できることを利用し、この値を実際の微分係数f'(x)の値のかわりに使用している。
つまり、中心差分を利用したニュートン法の場合、導関数f'(x)を与える必要がなく、f(x)だけをプログラムに記述すればよい。

このプログラムでは、あえてh=0.001と、hに粗い値を与えているが、微分係数と中心差分を用いて微分係数との誤差はおおよそh²=0.001²=10⁻⁶のオーダーなので、このように粗いhの値でも
 x=0.78598
という解が得られており、実用上、問題がないことがわかると思う。


ニュートン法 [数値解析]

ニュートン法

Newton-method-01.pngf(x)を微分可能な関数とする。

方程式

  

の解をx=αとする。

右の図のように適当な点x₁を選び、y=f(x)の点(x₁,y₁)における接線の方程式は

  

で、この接線とx軸との交点のx座標x₂は、上の式にy=0を代入することによって、

  

となる。

そして、同様ににおける接線の接線を引き、この接線とx軸との交点のx座標x₃を求めると、

  

となり、この操作を繰り返せば繰り返すほど、こうして求められたf(x)=0の近似値である

f(x)=0の解、x=αに近づいてゆくことが予想される。


これがニュートン法である。

漸化式の形で書けば、ニュートン法は次のようになる。

  


Newton-method-tab-01.pngx=√3の近似値をニュートン法を用いて求めることにする。

x=√3の両辺を2乗すると

  

これを

  

とおき、f(x)=x²−3=0x>0)とすれば、これはx=√3と同値。

f'(x)=2xだから、

  

計算開始のx₀=1として、表計算ソフトを使って計算したものは次の通り。

4、5回計算するだけで、x=√3≒1.732050808という近似値に到達している。


ニュートン法は前回の2分法よりも速く、しかも急速に収束することがわかると思う。


ただし、ニュートン法は、次の例のように、収束しないことがある。

  

計算の初期値としてx₀=1、または、x₀=2を取ると、

接線の方程式が、それぞれ、

  

となり、x=1x=2を交互に永遠に行き来する。

Newtonhou-02.png

この他にも、f'(x)=0になる点に差し掛かったとき、ゼロ割が発生するなど、危険な一面も有している。

こういうことは極まれにしか起きないけれど、運悪くこのような事態に遭遇することがある。

2分法と比較すると、ニュートン法は収束の速度は速いけれど、安定性に欠ける。


今日のアニソン、「花咲くいろは」から『ハナノイロ』 [今日のアニソン]

今日のアニソンは、「花咲くいろは」から『ハナノイロ』です。



ED曲も紹介するにゃ。



nice!(0)  コメント(0)  トラックバック(0) 
共通テーマ:音楽