ラマヌジャンの円周率公式

今年も3月14日、3.14の日がやってきました。3.14といえば、もちろん円周率の近似値ですね！

円周率の近似値にちなんで、世界的には円周率の日（英語圏だとPIの日）と呼ばれているそうです。

毎年、この日にブログに書きたいと思っていた（できずにいた）話があります。それが次のラマヌジャンの円周率公式です。

$\displaystyle \frac{1}{\pi} = \frac{2\sqrt{2}}{9801} \sum_{n = 0}^{\infty} \frac{(4n)! (1103+26390n)}{(n!)^4 \cdot 396^{4n}} \tag{1}$

なんじゃこりゃと思うような公式ですね。

9801、1103、26390 といった謎めいた整数を複雑に絡み合わせた無限級数を計算すると、なんと円周率の逆数が出てくるというのです。

ラマヌジャンはインドが生んだ著名な数学者で、数学者の中でも群を抜く奇才として知られています。上の公式は、まさにラマヌジャンの奇才っぷりを詰め込んだような式になっていますね。

しかもこの公式、こんな見た目をしておきながらめちゃくちゃ収束が早いそうで、一時期は円周率を世界最高の精度で計算するプログラムに使われていたそうです。

ラマヌジャン自身はこの公式の証明を残しておらず、後の数学者によって証明は与えられています。しかし、彼独自の考え方はあったにせよ、こんな複雑な式を予想したというのはかえって不思議さが増してしまいますね。

こんな魅力的な公式なので、いつかはブログ記事にまとめたいと考えていました。とはいえ、いつか書きたいと思っていても、何もしなければ一生書けないとも思いました。

円周率の日にかこつけて、行けるところまでこの公式についての解説を試みたいと思います*1。

f:id:tsujimotter:20210314030658p:plain:w400

参考文献

先に参考文献を挙げておきます。今回は数学セミナーの2018年12月号に掲載された塚本真輝さんの記事「京都大学ガロア祭問題７」を参考にさせていただきました。

数学セミナー 2018年 12 月号 [雑誌]

発売日: 2018/11/12
メディア: 雑誌

まずは遊んでみよう

せっかくなので、ラマヌジャンの公式の精度の高さを味わってみましょう。

まず驚きなのは、無限級数を $n = 0$ まで足した結果です。 $n = 0$ までで打ち切ると、色々な項が消えて

$\displaystyle \frac{2\sqrt{2}\times 1103}{9801}$

となります。これをGoogle検索で計算してみましょう：

https://www.google.com/search?q=2*sqrt(2)*1103%2F9801

結果は

0.31830987844

となるのですが、この逆数をとると

3.14159273001

となり、なんと小数第6位まで一致しています！　凄まじいですね。まだ $n = 0$ までですよ。

Pythonで軽く計算してみると次のようになりました：

n = 0 まで: 3.1415927300133055
n = 1 まで: 3.1415926535897936
n = 2 まで: 3.141592653589793
n = 3 まで: 3.141592653589793

太字のところが円周率の正確な値と一致している部分です。たった $n = 2$ までを計算しただけで、表示される桁の分は計算しきってしまいました。精度が高すぎて、検証するためにはもっと高精度で計算できるプログラムを書かなければいけないようです。

STEP1: ルジャンドルの関係式

それでは、上のラマヌジャンの公式がいったいどうやって導かれるのか、考えていきたいと思います。

まず、円周率の公式を求めるにあたっては、円周率を表す既存の関係式からスタートするのが自然な発想かと思います。

今回用いるのは次のルジャンドルの関係式です。

$0 < k < 1$ 、 $k' = \sqrt{1-k^2}$ としたとき

$\displaystyle \frac{\pi}{2} = E(k) K(k') + E(k') K(k) - K(k) K(k') \tag{2}$

が成り立つ。

何はともあれ、左辺に円周率が入っているところが注目ポイントです。

$K(k), E(k)$ はそれぞれ第1種完全楕円積分、第2種完全楕円積分と呼ばれるものです。定義は次の通り：

$\displaystyle K(k) = \int_{0}^{1} \frac{dx}{\sqrt{(1-x^2)(1-k^2x^2)}} \tag{3}$

$\displaystyle E(k) = \int_{0}^{1} \sqrt{\frac{1-k^2x^2}{1-x^2}} dx \tag{4}$

このような積分は「楕円の孤の長さ」を計算する際に必要だったという経緯があります。

たとえば、円1/4の孤の長さは

$\displaystyle \frac{\pi}{2} = \int_{0}^{1} \frac{dx}{\sqrt{1-x^2}}$

のように表されます。これを楕円に変えると、上のような積分が現れるというわけですね。

なお、式 $(3)$ ・式 $(4)$ において $k = 0$ とすれば、円弧の積分と一致します。つまり、円弧の積分の拡張になっているわけです。

円弧の積分の積分範囲を $0$ から $x$ の範囲で動かすと $\sin$ の逆関数になります。一方、 $(3), (4)$ の積分において、積分範囲を同様に動かして得られる関数は、初等関数では表すことができません。そこで、新たに「楕円関数」という関数を導入することで、初等関数に留まらない新しい関数論の世界が広がっていくわけですね。

この辺が楕円積分が導入された背景です。語りたいことはたくさんありますが、先に進みましょう。

$K(k), E(k)$ についている $k$ は、モジュラスと呼ばれるパラメータです。この $k$ を動かすと、積分はさまざまな値をとるわけですね。この $k$ に対応するように

$k' = \sqrt{1-k^2}$

というパラメータを用意します。 $k'$ は $k$ に従属しているので、 $k$ が動くと、それに合わせて $k'$ も動きます。

この $k, k'$ を用いて式 $(2)$ の右辺のように楕円積分を計算すると、なんと円周率が出てくるというのですね。

しかも面白いことに、パラメータ $k$ によらず左辺が一定値になっているのです。これは面白いことですね。

証明の流れを簡単に追いたいと思います。関数 $EK' + E'K - KK'$ を $k$ で微分すると

$\displaystyle \frac{d}{dk}(EK' + E'K - KK') = 0$

となります。ただし、 $E = E(k), E' = E(k'), K = K(k), K' = K(k')$ と略記しています。 $\frac{dE'}{dk}, \frac{dK'}{dk}$ については、 $E', K'$ を $k$ の関数だと思って微分することになります。

簡単に書いていますが、この計算はやってみるとものすごく大変でした。

微分して $0$ になるということは、 $k$ を動かしても $EK' + E'K - KK'$ は定数として振る舞うということですね。

定数ということは、 $k$ に何を入れても同じ結果になるということです。ということは、一番計算しやすいパラメータを選んであげれば良いですね。実際、 $k\to 0$ とすると

$\displaystyle EK' + E'K - KK' \xrightarrow{k \to 0} \frac{\pi}{2}$

となることがわかり、ルジャンドルの関係式が得られるというわけです。

以上の計算過程は次の記事のp.19以降で読むことができます。なので、今回の記事では省略します。
http://www.oishi.info.waseda.ac.jp/~samukawa/GaussLeg2.pdf
ルジャンドルの関係式自体非常に面白い対象なので、また別途取り扱うときに詳細を紹介しましょう。

STEP2: 超幾何級数に帰着

ルジャンドルの関係式によって、円周率が楕円積分によって計算できることがわかりました。

楕円積分は積分の形になっていますが、一方のラマヌジャンの公式は級数の形です。したがって、楕円積分を級数の形に表す必要がありいます。どうすればよいでしょうか。

そんなときに使えるのが、超幾何級数です。超幾何級数は以前の記事でも紹介したことがありました：
tsujimotter.hatenablog.com

また会えたね。

超幾何級数は、以下の式で表される無限級数でした：

$\displaystyle F(a, b, c; x) = \sum_{n=0}^\infty \frac{(a)_n (b)_n}{n! (c)_n}x^n \tag{5}$

ここで $(a)_n$ はポッホハマー記号と呼ばれるもので

$(a)_n = a(a+1)(a+2) \cdots (a+n-1)$

と定義されるのでした。

このパラメータがたくさんあるお化けのような関数が超幾何級数です。とにかく $a, b, c$ と３つもパラメータがついているので、いろんな関数をシミュレートできます。

たとえば、対数関数や逆三角関数は

$\displaystyle F(1, 1, 2; -x) = \frac{\log(1+x)}{x}$

$\displaystyle F(\frac{1}{2}, \frac{1}{2}, \frac{3}{2}; x^2) = \frac{\arcsin(x)}{x}$

のように表すことができます。

同じ要領で、先ほどの楕円積分も次のように表すことができます：

$\displaystyle F(\frac{1}{2}, \frac{1}{2}, 1; k^2) = \frac{2K(k)}{\pi} \tag{6}$

$\displaystyle F(-\frac{1}{2}, \frac{1}{2}, 1; k^2) = \frac{2E(k)}{\pi} \tag{7}$

すごいですね！　最初この事実を知ったときは、超幾何関数はこんな風に使えたのか、と驚いた記憶があります。

ルジャンドルの関係式の両辺に $4/\pi^2$ をかけて、超幾何級数を代入すると次のようになります：

$\displaystyle \begin{align} \frac{2}{\pi} &= \frac{2E(k)}{\pi} \frac{2K(k')}{\pi} + \frac{2E(k')}{\pi} \frac{2K(k)}{\pi} - \frac{2K(k)}{\pi} \frac{2K(k')}{\pi} \\ &= F(-\frac{1}{2}, \frac{1}{2}, 1; k^2) F(\frac{1}{2}, \frac{1}{2}, 1; k'^2) + F(-\frac{1}{2}, \frac{1}{2}, 1; k'^2) F(\frac{1}{2}, \frac{1}{2}, 1; k^2) \\ & \;\;\;\;\;\;- F(\frac{1}{2}, \frac{1}{2}, 1; k^2) F(\frac{1}{2}, \frac{1}{2}, 1; k'^2) \end{align}$