イノベーション エンジニアブログ


株式会社イノベーションのエンジニアたちの技術系ブログです。ITトレンド・List Finderの開発をベースに、業務外での技術研究などもブログとして発信していってます!


このエントリーをはてなブックマークに追加

ニュートン法で学ぶプログラミングと数学の関連

はじめに

どうも、bigenです。
最近ガイアックスさんが主催されているクローズドの数学勉強会に参加させてもらっていて、
統計学のための数学入門30講 を勉強しています。

輪講形式で進めているのですが、やっぱ大学のとき数学やっててよかったなー、と思いました。
大学時代はあんなに嫌になった数学が、実際にビジネスの世界に出てみると今は楽しくて、
やっぱり何のために勉強するのかが分かっていると良いなぁと思います。

先週の勉強会でニュートン法というアルゴリズムが出てきたのですが、
プログラムを色々書くようになってからみると新しい気持ちで学び直せたので、紹介させてください。
すごいんですよ、ニュートン法

√2ってどうやって求めますか?

突然ですが、コンピュータで平方根ってどうやって求めますか?
例えば、√2 とか。

あっ、sqrt()はナシですよ!
今日はその裏側の話です。

まず、前提として√2に関して知ってることって普通はあんまりないですよね。

  • 1よりは大きそう

  • 2よりは小さそう

  • 2乗すると2になる

ぐらいでしょうか。
記憶力の良い人であれば「1.41421356…​.(人よ一夜に人見頃…​)ぐらいでしょ」
と暗記している人もいると思いますが、 これはどうやって求められたのでしょうか?

今回は小数8桁の精度で√2を求めることを考えてみます。

とにかく求めるだけであれば、片っ端から色んな答えを試してみるのも一手でしょう。
精度は小数8桁でいいので、

step1: 1.00000001はどうかな? → 1.00000001 * 1.00000001 < 2 → もう少し大きくてもいけそう
step2: 1.00000002はどうかな? → 1.00000002 * 1.00000002 < 2 → もう少し大きくてもいけそう
step3: 1.00000003はどうかな? → 1.00000003 * 1.00000003 < 2 → もう少し大きくてもいけそう
・
・
・
step41421356: 1.41421356はどうかな? → 1.41421356 * 1.41421356 < 2 → もう少し大きくてもいけそう
step41421357: 1.41421357はどうかな? → 1.41421357 * 1.41421357 > 2 → 大きくしすぎた!

という感じでいつかは答えにたどり着きます。 ソースコード等は省略しますが、4000万ステップの反復が必要になります。

二分探索法も発想としては難しくないですね。
片っ端から試さなくても、ざっくり半分に割って、どちらのグループに答えがあるか調べることで、
絞り込んでいこうというやつです。 答えは1~2の間にあるので、

step1: 1.5はどうかな? → 1.5 * 1.5 > 2 → 大きすぎた。答えは1~1.5の間にある!
step2:1.25はどうかな?→ 1.25 * 1.25 < 2 → 小さすぎた。答えは1.25~1.5の間にある!
.
.
.
step26:1.414213553071はどうかな? → 小さすぎた!答えは 1.414213553071〜1.414213579722の間にある!
step27:1.4142135605216はどうかな? → そろそろ計算終わり・・・

という感じです。
計算の結果、27ステップで答えが見つかりました。
実際に以下のコードを実行していただけると途中経過も見れます。

<?php
// 二分探索法で√2を求める
$upper = 2;
$lower = 1;
$next;
$step = 0;

while ( $upper - $lower > 0.00000001) {
	$step++;
	$mid = ($upper + $lower) / 2;

	if ( $mid*$mid < 2) {
		$lower = $mid;
	} else {
		$upper = $mid;
	}

	echo 'step' . $step . ': ' . $mid . PHP_EOL;
}

力任せ探索に比べて爆発的に早く求まりましたね!

当たり前ですが、力任せ探索は探索範囲が広ければ広いほど、反復回数が爆発的に増えます。
出発点が偶然答えに近い場合は早く終わりますが、最悪計算量は悲惨です。
数学的にはo(N) (Nはデータの数)と表現し、Nが100倍になったら計算量も100倍、Nが10000倍になったら計算量も10000倍という意味です。

二分探索法は一発アテにいくわけじゃないので、答えがどこにあっても、それなりの反復回数で解にたどり着けることが知られています。
数学的にはo(logN)になり、Nが100倍になっても計算量は10倍、Nが10000倍になっても計算量は100倍くらいで済みます。

3.ニュートン法

最後に、数学的な知識を使って解くニュートン法をご紹介します。
アルゴリズムは、超簡単です。
適当にスタートの値x[0]を与えて、次のステップでの値を

x[n+1] = (x[n] + 2/x[n])/2

として計算します。 これが答えに十分近ければ終了、近くなければ同じ計算を繰り返すだけです。

ソースコードは以下です。

<?php
// ニュートン法を使って√2を求めてみます
$step = 1;    // ステップ数
$before = 2;  // 適当な初期値
$after = next_val($before);  // 1ステップ目

echo 'step1: ' . $after . "\n";

while (abs($before - $after) > 0.00000001) {
	$step++;
	$before = $after;
	$after = next_val($before);

	echo 'step' . $step . ': ' . $after . "\n";
}

echo 'result: ' . $after . "\n";



function next_val($before) {
	return ($before + 2/$before) / 2;
}

結果はこちら

step1: 1.5
step2: 1.4166666666667
step3: 1.4142156862745
step4: 1.4142135623747
step5: 1.4142135623731
result: 1.4142135623731

たった5ステップで8桁精度の答えが見つかった・・・・すごい!

ニュートン法についての説明は Wikipedia を見ていただくのが確実ですが、 理系の方なら懐かしい「微分」の概念を使って計算します。
収束速度は数学的に |α-x[n+1]| = o(|α-x[n]|^2)と表現され、2次収束と呼ばれます。
答えと次のステップの誤差は、前のステップの誤差の二乗になるという意味で、前のステップで誤差0.01ぐらいだったら次のステップは誤差0.0001ぐらいになるよ、ということです。

もうまとめ

注目ポイントは、二分探索法とかは数学のこと知らなくても思い付けそうなんですが、 ニュートン法は全く無理ですよね〜っていうところ。
「次のステップを、前のステップの値と2を前のステップで割ったものの平均を取る?意味わからん!」
って普通はなるわけですよ。バビロニア人を除くと。
(実は上記の手法はバビロニアン法とも呼ばれていて、紀元前数世紀に既に発見されている)

微分をしっている人だけが思いつけるアルゴリズムで計算速度が何百倍も早くなるとすると、
やっぱ数学の勉強とエンジニアリングは切っても切り離せないなぁと思いました。

少し短いですが今日はこのへんで。