Newton 法による方程式の近似解法

Newton 法とは, 次の図を根拠に近似解を計算する方法です. 函数 $y=f(x)$ が微分可能かつ下に凸であるとき, $f(x)=0$ の解$\alpha$を次の方法で求めます. まず, $\alpha$に近い値 $x_1$ を 1つとり, 関数 $y=f(x)$ のグラフ上の点 $(x_{1},f(x_1))$ における接線が，$x$軸と交わる点の $x$ 座標 を $x_2$ とします. このとき接線の方程式は,

$$y=f'(x_1)(x-x_1)+f(x_1)$$

となります. ここで $y=0$ とおいて, $x$ について解くと $x_2$ が求まります. すなわち

$$x_2 = x_1 -\frac{f(x_1)}{f'(x_1)}$$ (1)

次に, 点$(x_2,f(x_2))$における曲線の接線と $x$軸との交点の$x$座標を$x_3$とし, この操作を繰り返すと, 数列

$$x_1,x_2,x_3,\cdots ,x_n,\cdots$$ (2)

が漸化式

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

で定まります. この数列(2)は, 初期値 $x_1$ が $\alpha$ に十分近い時, 多くの場合急速に $\alpha$ に近づき, これを利用して近似解が効率よく求める事ができます. このようにして近似解を求める方法をニュートン法といいます.

ニュートン法では, 2分法と違って, 誤差の限界の正確な判定ができません. ここでは, 数列 (2) において, $\vert x_{n}-x_{n-1}\vert< \epsilon$ となったとき, $x_n$は誤差の限界が$\epsilon$の近似解であると考えます.

誤差の限界を $\epsilon$ として, 方程式 $f(x)=0$ の近似解を ニュートン法で求める手順は, 次のようになります.

1. 適当な初期値 $a$ を入力する.
2. $${b = a - \frac{f(a)}{f'(a)}}$$ を求める.
3. $\vert b-a\vert < \epsilon$を確かめ, 正しければ $b$ を近似解とする. 正しくなければ, $a$ に $b$ の値を入力して, 2に戻る.

次の例は, $x^3+x+1=0$ の唯一の実数解の近似値をニュートン法で 求めるプログラムです. 手抜きをして, 絶対誤差を評価したプログラムです.

/* File name 7-3.c */

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

#define EPSILON 1.0E-10

double next(double x);

main()
{
    double ans, ans1;

    ans=0.0;
    ans1=next(ans);
    while (fabs(ans-ans1) > EPSILON){
        ans1=ans;
        ans=next(ans1);
    }
    printf("%1.10f\n", ans);
}
double next(double x)
{
    return (x-(x*x*x+x+1)/(3*x*x+1));
}