πを求める公式

　　　　　　　　2002/12/20　後保範 ( Ushiro Yasunori )

1. 概要

　計算機によるπの計算は1949年Arctan公式で開始され、1983年以降は主に相加相乗法
になり、Arctan公式は使用されなくなた。その理由は、n桁のπの計算に要する演算量が
Arctan公式のO(n²)に対し相加相乗法はO(n(log₂n)²)と大幅に削減できるためである。
　しかし、DRM法 の発見によりarctan公式等の有利数級数や連分数も相加相乗法とほぼ
同等の演算量で計算可能なことが判明した。
　公式に対応するFORTRANとCのπ計算プログラムを順に掲載。本プログラムは内容を理解
してもらうことを目的とし、計算時間及び使用メモリ量の効率化は目的としていない。

2. 級数公式

2.1 Arctan公式

　arctan(1/p) = 1/p - 1/(3p³) + 1/(5p⁵) - 1/(7p⁷) + …… なる展開式で計算。
　多数あるが下記に有名な公式及び日本人の発見した公式を示す。
　(1) マチンの公式　：　　　　π/4 = 4arctan(1/5) - arctan(1/239)
　(2) ガウスの公式　：　　　　π/4 = 12arctan(1/18) + 8arctan(1/57) - 5arctan(1/239)
　(3) ｸﾘﾝｹﾞﾝｼｭﾃｨﾙﾅの公式：π/4 = 8arctan(1/10) - arctan(1/239) - 4arctan(1/515)
　(4) ｼｭﾃﾙﾏｰの公式　：　　　π/4 = 44arctan(1/57) + 7arctan(1/239) - 12arctan(1/682)
　　　　　　　　　　　　　　　　　　　 + 24arctan(1/12943)
　(5) 高野喜久雄の公式： π/4 = 12arctan(1/49) + 32arctan(1/57) - 5arctan(1/239)
　　　　　　　　　　　　　　　　　　　 + 12arctan(1/110443)

2.2 ラマヌジャン型公式

　(1) チュドノフスキーの公式
　　　1/π = 12/A・Σ_{n = 0}[Πⁿ_k=1 (2k-1)(6k-5)(6k-1)/(9k³B³)](-1)ⁿ(Cn+D)
　　　A = 640320^3/2　,　B = 106720　,　C = 545140134　,　D = 13591409
　(2) ラマヌジャンの公式
　　　1/π = 2^3/2/9801・Σ_{n = 0}[Πⁿ_k=1 (2k-1)(4k-3)(4k-1)/(2k³198⁴)](26390n+1103)
　(3) ラマヌジャンの無平方根公式
　　　1/π = 1/3528・Σ_n=0[Πⁿ_k=1 (2k-1)(4k-3)(4k-1)/(2k³3528²)](-1)ⁿ(21460n+1123)

3. 算術幾何平均(AGM, Arithmetic Geometric Mean)

3.1 ガウス・ルジャンドルの公式

　初期値：　a₀=1,b₀=2^-1/2,t₀=1/4 としn桁計算ではp=log₂nまで下記をk=1,2,……,pと反復計算
　　　a_k = (a_k-1+b_k-1)/2　, 　b_k = (a_k-1b_k-1)^1/2　, 　t_k = t_k-1-2^k-1(a_k-1-a_k)²
　
　すべての計算はn桁(+ガード桁)で行い、n桁のπは下記で求まる。
　　　π = (a_p+b_p)²/(4t_p)

3.2 ボールウエインの4次の公式

　初期値：　a₀=1,b₀=2^-1/4,t₀=(5-8^1/2)/8 としp=(log₂n)/2まで下記をk=1,2,……,pと反復計算
　　　a_k = (a_k-1+b_k-1)/2　, 　b_k = (a_k-1b_k-1(a_k-1²+ b_k-1²)/2)^1/4
　　　t_k = t_k-1+4^k{a_k⁴ -[(a_k²+b_k²)/2]²}
　
　すべての計算はn桁(+ガード桁)で行い、n桁のπは下記で求まる。
　　　π = a_p+1⁴/(1-2t_p)

4. 連分数公式

4.1 Γ関数系公式

　n桁計算でnと同程度の整数kを用いて下記のように計算する。
　　　π = { (2･4･6･……･(2k)) / (1･3･5･……･(2k-1)) }² / k ･ F(v)
　
　　　　　F(v) = 1 - 1 / [ 1 + v + (1･3) / { v + (3･5) / ( v + (5･7) / ( v + …… )) } ]
　　　　　　　,　v = 8k