Python で円周率を求める

Python

Python で数値計算の一例として、円周率を求めるスクリプトについてまとめておきます。
以下の環境でスクリプトの動作確認をしています。

環境: Windows 10、Python 3.x

背景

先日、LaTex が使えるプラグインをサイトの WordPress に入れました。
サイト上で自由に数式が扱えるようになりましたので、何か計算をしたくなってきます。
そこで、一番、典型的な事例として、円周率を計算してみます。

計算式

円周率の計算は、いくつか方法が考えられます。
どのやり方でもよいのですが、ここでは以下の公式(グレゴリー級数)を使うことにします。

(式1)
$$ \frac{\pi}{4} = 1 – \frac{1}{3} + \frac{1}{5} – \frac{1}{7} + \frac{1}{9} – …$$
この式の証明は、末尾に概要を記載しておくことにします。
プログラミングのため、両辺に4を掛けて以下のように変形します。

(式2)
$$ \pi = \frac{4}{1} – \frac{4}{3} + \frac{4}{5} – \frac{4}{7} + \frac{4}{9} – …$$

サンプルコード

上記の(式2)をそのまま、プログラムにします。

pi1 = 0 
for i1 in range( 100000 ): 
    odd1 = 2.0*i1 + 1.0 
    v1 = ((-1.0)**i1) * 4.0/odd1
    pi1 = pi1 + v1 
    print( pi1 ) 

出力は、3.14… の数字がずらずらと出てきて、だんだん1つの値(真値)に漸近するだけですので省略します。
6行ほどでプログラムが書けることがわかります。

サンプルコードの説明

・ pi1 が円周率です。
・ (式2)より、各項を無限に足していけば、円周率が求められます。
・ range() 内の数値が級数の項の数です。(式2)からは、無限大の数値を入れることになりますが、所定の有効桁で求められればよいため、まずは 10000 としています。
・ 各項の分母は奇数となっていますので、最初に奇数 odd1 を作ります。
・ つぎに、(式2)の符号と分数の記載に基づき、各項の値を作ります。
各項を pi1 に足しこんでいけば、円周率を求めることができます。
・ range() 内の数値は、最初は少ない数値を入れてスクリプトを実行してみてください。パソコンの計算能力に応じて大きな値に修正し、計算精度を上げてください。
・ (式2)の右辺は、絶対値が次第に小さくなる符号の異なる数を足し続ける形となっています。
したがって、円周率の真の値は、プログラムの実行中に最後に示される2つの数の間に存在することになります。加えて、表示される数値の変動幅はかならず小さくなっていくことになります。
この結果、計算中に表示される数値は、かならず、真の値に近づいていき、確定した桁数が次第に増えていく動きとなります。

式の証明(概要)

まず、\(tan y\) を \( y \) で微分します。
$$ \frac{d}{dy}tan y = \frac{1}{cos^2 y }
= \frac{cos^2 y + sin^2 y}{cos^2 y } = 1+ tan^2 y $$
つぎに、上式で、\(tan y = x \) とおき、逆数を取ります。
$$ \frac{dy}{dx} = \frac{1}{1+ x^2} $$
ここで、\( y = tan^{-1} x \) であり、\( | x | < 1\) と仮定すると、右辺は等比数列の和の形となっているので、以下となります。
$$ \frac{d}{dx} tan^{-1} x= 1-x^2+x^4-x^6+ … $$
両辺が\(x \) で積分できるものと仮定(※)し、0 から \(x\) まで積分します。
$$ tan^{-1} x= x-\frac{1}{3} x^3+\frac{1}{5} x^5-\frac{1}{7} x^7+ … $$
\(x → 1\) とすると、 \(tan \frac{\pi}{4} = 1\) ですので、(式1)が得られます。
※ 途中で仮定を入れていますので、厳密な証明とはいえません。しかしながら、とくに問題なく円周率の計算はできます。

まとめ

Python での数値計算の事例として、円周率を計算できるようになりました。

円周率を求める程度であれば、6行ほどで書けることがわかります。
また、ざっと見たところ、LaTeX 形式で記載した数式もとくに違和感はなく、サイト上で数式が自由に扱えるようになりました。

関連リンク
・ WordPress で数式を扱う 【LaTeX】

タイトルとURLをコピーしました