2.4.2 Newton 法

ここでは非線形方程式を解くための代表的な方法である Newton 法を試して みましょう。 $f: \R\to\R$ という関数があった時に、$x$ についての方程式

$$f(x) = 0$$

を解くために、適当な初期値 $x_0$ を定めて漸化式

$$x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)}$$

で数列 $\{x_n\}$ を定め、その極限 $x=\dsp\lim_{n\to\infty}x_n$ を見出 して、それを解として採用する、というものです。

example4.c

/* example4.c -- Newton 法で正の数 d の平方根を求める */

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

int main()
{
    int n;
    double eps = 1.0e-15;
    double d, Xn, Xnp1;

    printf("Newton 法で正の数 d の平方根を求めます。\n");
    printf("d=");
    scanf("%lf", &d);
    if (d < 0.0) return 0;

    Xn = 1.0;

    for (n = 0; n < 100; n++) {
        Xnp1 = (Xn + d / Xn) / 2;
        printf("%20.15e\n", Xnp1);
        /* 次の繰り返しのため */
        Xn = Xnp1;
        /* 十分よい近似値が得られたかどうかチェックして、満足行く精度に
        なっていたら繰り返しを打ち切って END に飛ぶ */
        /* if (Xn == Xnp1) goto END; */
        /* if (fabs(Xn - Xnp1) < eps * fabs(Xn)) goto END; */
        if (fabs(Xn * Xn - d) / d < eps) break;
    }
    printf("Square of %20.15e = %20.15e\n", Xn,  Xn * Xn);
    return 0;
}

最後に求めた計算値を自乗してチェックするようになっています。このプログ ラムで $\sqrt{2}$ を計算すると次のようになります。

example4 の実行結果

oyabun%	./example4
Newton 法で正の数 d の平方根を求めます。
d=2
1.500000000000000e+00
1.416666666666667e+00
1.414215686274510e+00
1.414213562374690e+00
1.414213562373095e+00
Square of 1.414213562373095e+00 = 2.000000000000000e+00
oyabun%

上記のプログラム example4.c の for ループの部分は Xnp1 という変数を使わないで

      for (n = 0; n < 100; n++) {
        Xn = (Xn + d / Xn) / 2;
        printf("%20.15e\n", Xnp);
        if (fabs(Xn * Xn - d) / d < eps) goto END;
      }

のように書くことも出来ます (こちらの書き方が普通です)。