インテジャーズ

INTEGERS

数、特に整数に関する記事。

arcsin^2のテイラー展開

この記事は日曜数学会のAdvent Calendarの17日目の記事です。

adventar.org

数日前、日曜数学会ミニが仙台で開催されました。この記事はそこでのライトニングトークの内容に基づいています。

和算の内容を含んでいますが、私は専門家ではなく十分な文献調査をしたわけでもないため、間違いや勘違いを多く含んでいると思われます。間違いが分かり次第その都度修正致しますが、ご了承いただけますと幸いです。

建部賢弘

建部賢弘(1664-1739)は和算家であり、『綴術算経』(1722年)に円周率の累遍増約術を用いた計算が書かれています。彼の術で小数第41位まで正確に求めることができます*1

直径1の円に内接する正2^k角形の周の長さをl_kとするとき、\lim_{k \to \infty}l_k=\piが成り立ちます。l_kを円周率の近似とすることは昔からなされていて、初等幾何学的に

\displaystyle l_{k+1}=\sqrt{2^{k+1}\left(2^k-\sqrt{4^k-l_k^2}\right)}

という漸化式の成立がわかるため、これを利用して円周率の近似計算をすることが原理的には可能です。ルドルフは正2^{62}角形の周の長さを計算して円周率を小数第35桁まで計算したそうです。

一方、関孝和はl_kまでのデータしか与えられなくとも、残りのl_{k'}(k' > k)の値を上手く近似的に補完することによって収束を加速させる増約術を開発しました(関の使ったデータと手法では円周率を小数第18位まで求められます)。

建部は関の増約術を繰り返し適用する累遍増約術を編み出しています。まず、l_1, \dots, l_{10}をあらかじめ計算しておきます。これを用いて

\displaystyle \varepsilon_k:=\frac{l_{k+1}-l_k}{l_k-l_{k-1}}

k=2, \dots, 9で計算します。その数値をみると、\varepsilon_kk \to \infty1/4に収束することが強く期待されます。

\begin{align}l_k&=l_9+(l_{10}-l_9)+(l_{11}-l_{10})+\cdots +(l_k-l_{k-1})\\
&=l_9+(l_{10}-l_9)(1+\varepsilon_{10}+\varepsilon_{10}\varepsilon_{11}+\cdots +\varepsilon_{10}\varepsilon_{11}\cdots\varepsilon_{k})
\end{align}


なので、\varepsilon_{10}以降を1/4と置き換えることによって

\displaystyle \pi ≒ l_9+(l_{10}-l_9)\sum_{i=0}^{\infty}\frac{1}{4^i}=\frac{1}{3}(4l_{10}-l_9)

と考えられます。これが増約術です。ここで終わらずに円周率の近似列

\displaystyle l_k^{(1)}:=\frac{1}{3}(4l_{k+1}-l_k) \quad (k=1, \dots, 9)

に対して同様の増約術を適用します。すると、\varepsilon_k^{(1)}1/4^2に収束することがみてとれ、

\displaystyle \pi ≒ l_8^{(2)}:=\frac{1}{15}(16l_{9}^{(1)}-l_8^{(1)})

を得ることができます。以下、同様に増約術を繰り返すことによって(累遍増約術)、円周率の近似

\displaystyle \pi ≒ l_1^{(9)}

に到達し、正しく計算すれば円周率の小数第41位まで正確に近似します([1])。

次に、円の弧長sを矢の長さhと直径dを用いて表すという問題を考えます。

f:id:integers:20181218180817p:plain:w300

ここで、弧ACBの長さをsとし、CDを矢と呼んで、その長さをhとしています。建部はd=10, h=10^{-5}のときに弧ACBを六十四等分して、定半背冪 \left(\frac{s}{2}\right)^2を累遍増約術で計算しました。その結果、

\displaystyle \left(\frac{s}{2}\right)^2=0.00010000003333335111112253969066667282347769479595875\cdots

を得ています。建部はこれが

\displaystyle \left(\frac{s}{2}\right)^2=1\times 10^{-4}+\frac{1}{3}\times 10^{-10}+\cdots

と級数展開されると見抜き、続く項の係数も零約術(=連分数展開)を用いることによって発見しました。例えば

\displaystyle \left(\frac{s}{2}\right)^2-\left(1\times 10^{-4}+\frac{1}{3}\times 10^{-10}\right)=a\times 10^{-16}

として(10^{-6}=h/d)、aを連分数展開して近似分数を求めると

\displaystyle \frac{1}{5}, \frac{1}{6}, \frac{2}{11}, \frac{3}{17}, \frac{5}{28}, \frac{8}{45}, \frac{34565}{194428}, \dots

となっており、8/45が係数であろうと見当をつけます。続けて

\displaystyle \left(\frac{s}{2}\right)^2-\left(1\times 10^{-4}+\frac{1}{3}\times 10^{-10}+\frac{8}{45}\times 10^{-16}\right)=b\times 10^{-22}

として、bを連分数展開して近似分数を求めると

\displaystyle \frac{1}{8}, \frac{1}{9}, \frac{3}{26}, \frac{4}{35}, \frac{40175}{351531}, \dots

となっており、4/35が係数であろうと考えます。このようにして、建部は

\displaystyle \left(\frac{s}{2}\right)^2=1\times 10^{-4}+\frac{1}{3}\times 10^{-10}+\frac{8}{45}\times 10^{-16}+\frac{4}{35}\times 10^{-22}+\frac{128}{1575}\times 10^{-28}+\frac{128}{2079}\times 10^{-34}+\cdots

なる展開を得ました*2。更に、一般には

\displaystyle s^2=4hd\left(1+\frac{1}{3}\cdot\frac{h}{d}+\frac{8}{45}\cdot\left(\frac{h}{d}\right)^2+\frac{4}{35}\cdot\left(\frac{h}{d}\right)^3+\frac{128}{1574}\cdot\left(\frac{h}{d}\right)^4+\frac{128}{2079}\cdot\left(\frac{h}{d}\right)^5+\cdots\right)

となっていることを理解します。最後に建部は綴術によって

\displaystyle \frac{1}{3}, \ \frac{8}{45}, \ \frac{4}{35}, \ \frac{128}{1574}, \ \frac{128}{2079}, \dots

\begin{align}&\frac{2^2}{3\cdot 4}, \quad \frac{2^2\cdot 4^2}{3\cdot 4\cdot 5 \cdot 6}, \quad \frac{2^2\cdot 4^2\cdot 6^2}{3\cdot 4\cdot 5 \cdot 6\cdot 7 \cdot 8}, \ \\
&\frac{2^2\cdot 4^2\cdot 6^2\cdot 8^2}{3\cdot 4\cdot 5 \cdot 6\cdot 7 \cdot 8\cdot 9 \cdot 10}, \quad \frac{2^2\cdot 4^2\cdot 6^2\cdot 8^2\cdot 10^2}{3\cdot 4\cdot 5 \cdot 6\cdot 7 \cdot 8\cdot 9 \cdot 10\cdot 11 \cdot 12}, \dots\end{align}

となっていることを見抜きます。そうして、建部は無限級数展開

\displaystyle s^2=4hd\left(1+\frac{2^2}{3\cdot 4}\cdot\frac{h}{d}+\frac{2^2\cdot 4^2}{3\cdot 4\cdot 5 \cdot 6}\cdot\left(\frac{h}{d}\right)^2+\frac{2^2\cdot 4^2\cdot 6^2}{3\cdot 4\cdot 5 \cdot 6\cdot 7 \cdot 8}\cdot\left(\frac{h}{d}\right)^3+\cdots\right)

に到達しました。

\angle COB=\thetaとすると弧長の公式より

\displaystyle \frac{s}{2} = \frac{d}{2}\theta

であり、一方

\displaystyle h=\frac{d}{2}-\frac{d}{2}\cos \theta= \sin^2\frac{\theta}{2}

なので、

\displaystyle s^2=4d^2\left(\arcsin\left(\sqrt{\frac{h}{d}}\right)\right)^2

なる関係があります。また、

\displaystyle \frac{2^2\cdot 4^2\cdots (2n-2)^2}{3\cdot 4 \cdot 5 \cdots (2n)}=\frac{4^{n-1}( (n-1)!)^2}{2^{-1}\cdot (2n)!}=\frac{2\cdot 4^{n-1}}{n^2\binom{2n}{n}}

なので、\sqrt{\frac{h}{d}}=xとおくことによって、建部の式は\arcsin^2のテイラー展開

\displaystyle 2(\arcsin x)^2 = \sum_{n=1}^{\infty}\frac{(2x)^{2n}}{n^2\binom{2n}{n}}

に他ならないことがわかりました。この発見はEulerより15年早いとのことです。x=1/2を代入すれば

\displaystyle \sum_{n=1}^{\infty}\frac{1}{n^2\binom{2n}{n}}=\frac{\pi^2}{18}

が得られます。

現代的な導出

|x| < 1とする。

\displaystyle \int_0^1\frac{\mathrm{d}t}{1-x^2+x^2t^2}=\frac{1}{x\sqrt{1-x^2}}\arctan\left(\frac{x}{\sqrt{1-x^2}}\right) = \frac{\arcsin x}{x\sqrt{1-x}^2}

なので、

\displaystyle \frac{\arcsin x}{\sqrt{1-x^2}}=\sum_{n=0}^{\infty}\left(\int_0^1(1-t^2)^n\mathrm{d}t\right)x^{2n+1}

を得る。Chu-Vandermondeの恒等式

\displaystyle \sum_{k=0}^n\frac{(-n)_k(b)_k}{(c)_kk!}=\frac{(c-b)_n}{(c)_n}

b=\frac{3}{2}, c=\frac{1}{2}を代入することによって

\displaystyle \sum_{k=0}^n\frac{(-1)^k\binom{n}{k}}{2k+1}=\frac{(2n)!!}{(2n+1)!!}

となるので*3

\displaystyle \int_0^1(1-t^2)^n\mathrm{d}t = \int_0^1\sum_{k=0}^n(-1)^k\binom{n}{k}t^{2k}\mathrm{d}t=\sum_{k=0}^n\frac{(-1)^k\binom{n}{k}}{2k+1}=\frac{(2n)!!}{(2n+1)!!}

と計算できる。従って、

\displaystyle \frac{\arcsin x}{\sqrt{1-x^2}}=\sum_{n=0}^{\infty}\frac{(2n)!!}{(2n+1)!!}x^{2n+1}

を得た。

\displaystyle \frac{\mathrm{d}(\arcsin x)^2}{\mathrm{d}x}=\frac{2\arcsin x}{\sqrt{1-x^2}}

なので、

\begin{align} 2(\arcsin x)^2 &= 4\int_0^x\frac{\arcsin t}{\sqrt{1-t^2}}\mathrm{d}t = 4\int_0^x\left(\sum_{n=0}^{\infty}\frac{(2n)!!}{(2n+1)!!}t^{2n+1}\right)\mathrm{d}t \\
&=2\sum_{n=1}^{\infty}\frac{(2n-2)!!}{n\cdot (2n-1)!!}x^{2n} = \sum_{m=1}^{\infty}\frac{(2x)^{2n}}{n^2\binom{2n}{n}}\end{align}

となり、所望のテイラー展開が示された。

参考文献

[1] 小川束, 円理の萌芽–建部賢弘の円周率計算–, 数理解析研究所講究録, 1019巻 (1997), 77–97.
[2] 村田全, 建部賢弘の数学とその思想, 科学図書館.
[3] 小松彦三郎, 綴術算経の異本と成立の順序, 数理解析研究所講究録, 1130巻 (2000) 229–244.

*1:建部が実際に記載した数値は第40位までが正しいらしいです。また、本に書かれている方法でl_k^2について計算を実行すると第37位までしか求められないらしいです。

*2:128/1575あたりから何番目の近似分数を採用するかが明確でなくなってきます。この辺については[3]で勉強できそうです。綴術算径には幾つかのヴァージョンがあってもっと数値を計算しているものもあるとのこと。

*3:二重階乗 - INTEGERS