著者:佐々木 龍之介
今回は、数値解析アルゴリズムの一つである「ニュートン法」を紹介します。ニュートン法を利用すると、反復計算によって方程式の解を近似的に求められます。C言語を使ったサンプルプログラムを挙げながら、同手法について簡単に解説します。
シェルスクリプトマガジン Vol.80は以下のリンク先でご購入できます。
図3 「x2 – 1 = 0」という方程式の解を求めるCプログラムの例
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 |
#include <stdio.h> #include <math.h> #define e 0.000001 double f(double x) { return x*x - 1; } double df(double x) { return 2*x; } int main() { double x1 = 2.5, x2; while(1) { x2 = x1 - f(x1)/df(x1); if(fabs(f(x2)) < e) { break; } x1 = x2; } printf("%f\n", x2); } |
図4 「x2 – 1 = 0」という方程式の解を求めるCプログラムの改良版
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 |
#include <stdio.h> #include <math.h> #define e 0.000001 double f(double x) { return x*x - 1; } double df(double x) { return 2*x; } int main() { double x1 = 2.5, x2; while(1) { x2 = x1 - f(x1)/df(x1); if(fabs(x2 - x1) < e) { break; } x1 = x2; } printf("%f\n", x2); } |
図5 √2の近似解を求めるCプログラムの例
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 |
#include <stdio.h> #include <math.h> #define e 0.000001 double f(double x) { return x*x - 2; } double df(double x) { return 2*x; } double g(double x) { return x - f(x) / df(x); } int main() { double x1 = 2.0, x2; while(1){ x2 = g(x1); if(fabs(x2 - x1) < e) { break; } x1 = x2; printf("%f\n", x2); } printf("√2= %f\n", x2); } |