円周率を近似計算してみた with Python

Pythonロゴ

今まで全く気にしてこなかったのですが、円周率ってどうやって導き出すんだろうと思ったので、Pythonでプログラムしてみました。

色々な求め方があるんですね。。。
その中でなんとなくPythonで組めそうな感じがした以下の2つを試してみました。

グレゴリー級数

end = 1000
phi = 0
isTurnOver = True
for i in range(1, end, 2):
    if isTurnOver:
        g = 4 * 1 / i
    else:
        g = -1 * 4 * 1 / i
    phi += g

    isTurnOver = not isTurnOver

    # 最後の方だけ「パイ」を表示してみる
    if i > end*0.99:
        print('π=', i, phi)

π= 991 3.1395765266063047
π= 993 3.1436047239879765
π= 995 3.139584623485464
π= 997 3.143596659593789
π= 999 3.139592655589785

結果は・・・π= 3.139592655589785

モンテカルロ法

bR = .5 # 半径
sS = 0  # 四角のプロット数 
cS = 0  # 円のプロット数
count = 500

np.random.seed(0) # 乱数を固定
ps = np.random.uniform(-1*bR, bR, (count, 2))
x = ps[:, 0]
y = ps[:, 1]
r0 = np.sqrt(x**2 + y**2) # ランダムにプロットした点
for r in r0:
    if bR > r:
        cS += 1; # 円の中
    sS += 1;     # 四角形の中

print('π=', 4. * cS / sS)

結果は・・・π= 3.064

乱数を固定せずに何度かやれば、下3, 4桁位までは精度が出ることがあったので、誤差のせいか今回の500というプロット数は圧倒的に少ないということでしょうか・・・

意外と精度が出なかった

グレゴリー級数で1/1000まで、モンテカルロ法で500点で実験してみましたが意外と3.141592…とはならないんだなって思いました。バグがなければですが・・・
割り算した後足したり引いたりの所で、ごく僅かですが桁が落ちてるのが精度が出ない原因なのでしょうか・・・?
また時間があったら調べてみようかと思います。

実は

ファインマンポイントという円周率の下762桁目から「9」が6回続くらしいのですが、この「999999」が出るまでやってみたかったのですが、グレゴリー級数の方を10億までやってみたのですがまだ駄目でした。いつかチャレンジしてみます。

鹿児島県の出水市という所に住んでいまして、インターネット周辺で色々活動して行きたいと思ってるところです。 Webサイト作ったり、サーバ設定したり、プログラムしたりしている、釣りと木工好きなMacユーザです。 今はデータサイエンスに興味を持って競馬AI予想を頑張ってます。

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です

コメントする

このサイトはスパムを低減するために Akismet を使っています。コメントデータの処理方法の詳細はこちらをご覧ください