移動式のブログ

ガジェット、アニメ、プログラミング、考えたことその他色々・・・特にこれといったテーマはないカオスなブログです。

円周率求めるアルゴリズム

ガウスルジャンドルアルゴリズム

idoushiki.hatenablog.com

以前、円周率を暗記するゲーム「PIMEMORIZER」を作成した。
円周率って無限に続くけど、どうやって円周率を求めるのだろう?どんなふうに求められているのかなーと思った。
「円周率 アルゴリズム」で検索したら「ガウスルジャンドルアルゴリズム」というものがあったので、Pythonでプログラム組んで計算してみた。


ガウス=ルジャンドルのアルゴリズム - Wikipedia



f:id:idoushiki:20170504202350p:plain
まず、これが初期値




f:id:idoushiki:20170504202411p:plain
そして、この4つの式を何度も計算する



f:id:idoushiki:20170504202537p:plain
四つの式を何度も反復して計算して得た答えをこの式に当てはめて、円周率を求める


これら式をPythonで書いてみる。

Python3で組んだコード

適当に19回ぐらい繰り返して計算させることにした。
桁数は20桁分まで見れるように"{0:.20f}".format(pi)

実行結果

2.91421356237309492343
3.14057925052216857509
3.14159264621354283875
3.14159265358979400418この行から答え同じ↓
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418
3.14159265358979400418

計算するごとに、3.14159....とより正確な円周率に近づいていることがわかる。
しかし、四回目から、同じ計算結果しか出なくなってしまった。

円周率πの計算:ガウス=ルジャンドルのアルゴリズム
このサイトは、ガウスの円周率求める式を計算したようだ。

このサイトによると、64ビット浮動小数点で計算しているので、14桁分まで求めるのが限界だそうだ。
ウイキペディアにも、
倍精度浮動小数点数 - Wikipedia
指数部が全て 0(ゼロ、および非正規化数)の場合を除き、通常の数(正規化された数)の精度は、53ビット相当、十進に直すと約16桁であると書いてある。
「3」と「.」も合わせれば、確かに16桁分になる。円周率も、3.14159265358979までは、あっている(円周率表)。


正確に求めるには、もっと高度なプログラムを組まないといけないようだな。