円周率を計算するPythonプログラムの書き方

πは重要な数である。とπについての計算や、πを使った角度の測定に使われる。 πは、無理数であるなど、興味深い性質を持っている。これは、繰り返しのパターンに一致しない桁が無限にあることを意味する。しかし、様々な方法でπを近似することができる。しかし、πを様々な方法で近似することができる。幸いなことに、コンピュータ・プログラムを書くのは難しくない。プログラミングの練習にもなるし、πという数についてより深く学ぶのにも良い方法だ。Pythonの基本的なプログラムを使ってπを計算する方法を学びましょう！

方法1

ニラカンタ級数を使う

1. ニラカンタ級数を理解しましょう。Nilakantha 級数は、π=3+42*3*4-44*5*6+46*7*8-48*9*10...{displaystyle \pi =3+{frac{4}{2*3*4}}-{frac{4}{4*5*6}}+{frac{4}{6*7*8}}-{frac{4}{8*9*10}}...}から始まります。
   となり、このパターンに従って続ける。つまり、書きたいアルゴリズムは以下のようになる：

   • 3を "答え "とし、数n=2{displaystyle n=2}から始める。
   • 4n∗(n+1)∗(n+2){displaystyle {frac {4}{n*(n+1)*(n+2)}}} を計算する。
   • その計算結果を答えに加減する。
   • 指定回数繰り返す。
   • 答えを返して表示する。
2. 新しいテキストファイルを作成する。好きなIDEを使ってもいいし、テキストエディタを使ってもいい。あなたのコンピュータがPythonのプログラムファイルとして認識できるように、ファイルに.pyという拡張子をつけてください。
3. decimal モジュールをインポートします。Pythonをこのモジュールや同様のライブラリなしで使うと、精度は17桁に制限されます。しかし、このモジュールを使えば、任意の精度で桁数を指定することができます。Pythonのデフォルト・ライブラリなので、別途インストールする必要はない。

   10進数から * をインポートする。
4. 小数の桁精度を設定します。何桁のπを計算したいかによって、どの程度の精度にするか決まります。例えば、100桁のπを計算する場合は、次の行を追加します：

   getContext().prec = 100
5. ニラカンタ級数の関数を定義する。Pythonの場合、この関数は以下のような構造になります：

   def nilakantha(reps): # 計算はここで行う return answer
6. answerは初期値3である。10進数ライブラリが提供する高精度が必要な数なので、必ず10進数にすること。この変数は後で足し算と引き算を交互に行うのに使う。

   def nilakantha(reps): answer = Decimal(3.0) op = 1 # 計算はここで行う return answer
7. forループを追加する。forループは変数nに2を初期設定する。そしてループの中で書かれていることを実行し、nの値を2ずつ増やし、上限値（2*reps+1）に達するまでこの処理を繰り返す。

   def nilakantha(reps): answer = Decimal(3.0) op = 1 for n in range(2, 2*reps+1, 2): # 計算はここで行う return answer
8. ニラカンタ級数の要素を計算し、答えに加える。分数の1つを10進数にするだけで十分で、Pythonはそれに応じて他の部分を変換します。式をプログラムするだけでなく、opを掛けます。

   • 最初のサイクルでは、opは1に設定されているので、掛け算しても何も起こりません。しかし、後で他の値に設定される。

   for n in range(2, 2*reps+1, 2): result += 4/Decimal(n*(n+1)*(n+2)*op)
9. opが1なら-1、-1なら1になる。負の数を足すことは、正の数を引くことと同じである。このように、このプログラムは足し算と引き算を交互に繰り返す。

   for n in range(2, 2*reps+1, 2): result += 4/Decimal(n*(n+1)*(n+2)*op) op *= -1
10. この関数のインターフェイスを書きましょう。おそらく、何回繰り返された級数を使用するかを入力する方法と、計算したπの近似値を表示する方法が必要でしょう。

    print("How many repetitions?") repetitions = int(input()) print(nilakantha(repetitions))

    • πの桁数が少ない場合は、計算結果と比較するために実際のπの始まりを表示することもできます。その場合は、次の行を追加する：
      print("3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679")
      (比較のためにπの桁数がもっと必要な場合は、インターネットからコピーできます)。
11. コードをチェックしてください。コード全体は次のようになるはずです（最後の行は省略できます）：

    from decimal import * getcontext().prec = 100 def nilakantha(reps): result = Decimal(3.0) op = 1 n = 2 for n in range(2, 2*reps+1, 2): result += 4/Decimal(n*(n+1)*(n+2)*op) op *= -1 return result print("How many repetitions?") repetitions = int(input()) print(nilakantha(repetitions)) print("3.1415926535897932384626433832795028841971693993751058209749445923078164062862089986280348253421170679")
12. プログラムを実行します。IDEの "Run "シンボルをクリックします。PythonのIDLEで、F5を押す。単純なテキストエディタで作業していた場合は、ファイルを保存し、Pythonで実行します。

    • 100のような少量の反復から始める。これでプログラムが動くかどうかがわかる。
    • πを何桁も求める場合は待つ覚悟が必要です。例えば、このシリーズを100万回反復すると、18桁のπが正しく得られ、約30秒かかります。

方法2

モンテカルロ法を使う

1. モンテカルロ法を理解しましょう。任意の長さの正方形を想像し、その中にその長さと同じ半径を持つ円の1/4を置く。プログラムは正方形の内側にランダムな点を生成し、それらが円の内側にもあるかどうかをチェックする。

   • たくさんの点がある場合、四分の一円の中の点の量を正方形の中の点の量で割ると、四分の一円の面積を正方形の面積で割るようなことになる。つまり
     AquartercircleAsquare=14πr2r2=14π{displaystyle {frac {A_{quartercircle}}{A_{square}}}={frac {{frac {1}{4}}pi r^{2}}{r^{2}}}={frac {1}{4}}pi }} である。
     でπを計算できる：
     4AquartercircleAsquare=π{displaystyle 4{frac {A_{quartercircle}}{A_{square}}=pi}}で計算できる。
   • 四分の一円の面積を計算するにはπが必要になるので、このプログラムは面積を直接使うことはできない。
   • これは効率的な方法ではない。例えばニラカンタ級数と同じ桁数のπを得るには、かなり長い時間を待たなければならない。しかし、想像しやすく視覚化しやすい方法である（その代償として、パフォーマンスはさらに落ちる）。
2. 必要なモジュールをインポートする。randomには乱数を生成する機能がある。mathには平方根など、点の距離を計算するのに必要な数学関数がある。このため動作は遅くなるが、メソッドを理解するのに役立ち、しばらく見ているだけでも面白い。πを速く計算したいのであれば、とにかく別の方法を選ぶべきである。

   import random import math import turtle
3. 何点計算するかをユーザーに尋ねる。これは以下のコードでできる：

   print("Insert number of points:") np = input() while not np.isdigit(): print("Insert number of points:") np = input() np = int(np)
4. タートルを高速化する。デフォルトでは、タートルはそれほど速くありません。タートルのスピードをfastestに設定することで、これを変更できます：

   turtle.speed("fastest")
5. 状況を描画する。長方形と1/4円がある座標系を描き、1/4円を描く。

   • まず、四角形の長さと四分の一円の半径をピクセル単位で格納する変数を定義する（これは同じ数なので、変数は1つでよい）。こうしておくと、四半円と正方形の大きさを変更するときに手間が省ける。
     length = 300 # ピクセル単位の円の半径と正方形の長さ
   • 次に、実際に座標軸と円を描く必要がある。このコードは長いが、これらの描画のためにタートルを動かしているだけだ。
     #draw y axis turtle.pensize(2) turtle.forward(length + 40) turtle.left(135) turtle.forward(20) turtle.back(20) turtle.left(90) turtle.forward(20) turtle.penup() turtle.home() turtle.pendown() #draw x axis turtle.left(90) turtle.forward(length + 40) turtle.左(135) turtle.forward(20) turtle.back(20) turtle.left(90) turtle.forward(20) turtle.penup() turtle.goto(0,length) turtle.left(45) turtle.left(180) turtle.pendown() #円の4分の1を描く turtle.pencolor("red") turtle.circle(length,-90)
6. 各ドットに必要な計算のループを作る。ループの前に、円の内側のドットの量（変数inside）を0に設定する。

   inside = 0 for i in range(0,np)：
7. ドットのランダムな位置を取得する。ドットのx位置とy位置の2つの乱数が必要である。計算を簡単にするために、前のステップでは1/4円の中心を(0,0)にしておいた。つまり、両方の数値が0から正方形の長さの間にある必要がある。このような数値をrandom.uniform()関数で求める：

   #ドットの位置を取得 x = random.randint(0,length) y = random.randint(0,length)
8. ドットが1/4円の中にあるかどうかをチェックします。点と中心の距離を求め、それが四半円の半径以下かどうかをチェックする必要がある。

   • 距離を計算するには、ピタゴラスの定理を使う必要がある。それは
     d=(x2−x1)2+(y2−y1)2{\displaystyle d={\sqrt {(x_{2}-x_{1})^{2}+(y_{2}-y_{1})^{2}}}}
     しかし、中心は(0,0)にあるので、x1とy1は両方とも0であり、無視できる。式はより簡単である：
     d=x22+y22{displaystyle d={sqrt {{x_{2}}^{2}+{y_{2}}^{2}}}}} です。
     Pythonコードでは（x2とy2は前のステップで得た座標です）：
     #中心からの距離を求める d = math.sqrt(x**2 + y**2)
   • 点が円の内側にある場合、円の内側にある点をカウントする変数を1増やす。より見やすくするために、円の内側にある点の色を赤に、円の外側にある点の色を青に設定する。
     if d <= length: inside += 1 turtle.pencolor("red") else: turtle.pencolor("blue")
9. 点を描く。これにはタートルを使う：

   #draw dot turtle.penup() turtle.goto(x,y) turtle.pendown() turtle.dot()
10. ループ終了後に結果を表示する。円の内側にいくつの点があったか、この計算がどのπの値を与えたかをユーザーに伝える：

    print("1/4円の内側:") print(inside) print("点の総量:") print(np) print("円周率はおよそ:") print((inside / np) * 4.0)
11. ユーザーが画面をクリックしたときだけ終了する。これはturtleモジュールのexitonclick()関数で行う。そうしないと、計算が終わったときに図面のあるウィンドウが閉じてしまい、ユーザーがそれを見る時間がなくなってしまうからです。行を追加します：

    turtle.exitonclick()
12. コードをチェックしてください。あなたのコード全体はこうなっているはずだ：

    import random import math import turtle print("Insert number of points:") np = input() while not np.isdigit(): print("Insert number of points:") np = input() np = int(np) turtle.speed("fastest") length = 300 # 円の半径と正方形の長さ（ピクセル） # y軸を描く turtle.pensize(2) turtle.forward(length + 40) turtle.left(135) turtle.forward(20) turtle.back(20) turtle.left(90) turtle.forward(20) turtle.penup() turtle.home() turtle.pendown() #x軸の描画 turtle.left(90) turtle.forward(length + 40) turtle.left(135) turtle.forward(20) turtle.back(20) turtle.left(90) turtle.forward(20) turtle.penup() turtle.goto(0,length) turtle.left(45) turtle.left(180) turtle.pendown() #円の4分の1を描く turtle.pencolor("red") turtle.circle(length,-90) inside = 0 for i in range(0,np): #点の位置を取得 x = random.uniform(0,length) y = random.uniform(0,length) #中心からの距離を求める d = math.sqrt(x**2 + y**2) if d <= length: inside += 1 turtle.pencolor("red") else: turtle.pencolor("blue") #点を描く turtle.penup() turtle.goto(x,y) turtle.pendown() turtle.dot() print("1/4円の内側:") print(内側) print("点の総量:") print(np) print("円周率はおよそ:") print((inside / np) * 4.0) turtle.exitonclick()
13. プログラムを実行します。IDEの "Run "シンボルをクリックします。PythonのIDLEでF5を押す。単純なテキストエディタで作業していた場合は、ファイルを保存し、Pythonで実行します。

    • 100個のような少量のドットから始める。これでプログラムが動くかどうかがわかる。
    • 非常に長く待つことを覚悟してください。1000点の計算でも約1分半かかり、πは数桁（1～2桁）しか得られません。10000点の計算には15分かかり、2～3桁のπが得られる。

この記事は、CC BY-NC-SAで公開されている " How to Write a Python Program to Calculate Pi " を改変して作成しました。特に断りのない限り、CC BY-NC-SAで利用可能です。