円周率・その3:柳谷晃「円周率πの世界」
3月14日は、世界的に円周率の日である。
このように書くのは、去年おととしに次いで3度目だ。つまり、孫娘の2歳の誕生日を迎えることになる。去年の記事の冒頭で、次のように書いた。
まぁ、仕方あるまい。そのうちわかるだろう。
ということで、今年はブルーバックスの柳谷晃著「円周率πの世界」を読んだ。
2021年の7月に発行されたばかりなので、2019年の3月15日に発表されたGoogleの快挙もカバーされている。
Google、円周率計算31兆桁達成 世界記録更新 - ITmedia NEWS
ただし、執筆時点で追いついていなかったのか、2020年の1月にTimothy Mullican氏が 50兆桁計算したことについては触れられていない。さらに、今日ネットを検索したところによれば、2021年8月19日にスイスの Applied Science Universityが計算した 62 兆 8318 億 5307 万 1796 桁が現在の記録のようだ。計算に要した時間は 108日と9時間だったということである。2019年のGoogleとほぼ同じ計算時間で、2倍の桁数が計算されたことになるわけだが、当然、このことも含まれていない。
円周率計算で世界記録を大幅更新、スイス研究チームの高性能コンピューター - ZDNet Japan
第1章で、古代バビロニアやギリシャ・エジプト、中国やインドでの円周率の計算や意義から始まり、第2章で正多角形の周の長さから求めるアルキメデスの方法にふれ、第3章で、無限級数と微積分の話、2項定理に触れてから、グレゴリーの級数が紹介される。さらに第4章でオイラーの業績を紹介し、円周率が超越数であることのリンデマンの証明について触れたあとに、第5章でコンピュータによるπの計算の歴史が紹介される。
去年に、やはりブルーバックスの西岡久美子著「超越数とはなにか」と堀場芳数著「無理数の不思議」をそれぞれ読んだのだが、よい復習になった。
それとともに、こういうたぐいの教養書が難しいと思ったところは、やはり高校卒業から大学初等くらいの数学が身についていないと理解がむつかしいだろうな、という点だ。少し中途半端な感じがした。
とはいえ、興味に応じて、難しい理論がない1章と2章だけでも楽しめるしいろいろ発見があると思う。3章ー4章は数学アレルギーには少々難しいかもしれないが、それでもネットで調べながら人物名に親しむと、それなりに楽しめるだろうし、そこから深めることもできるだろう。5章の桁数の競争のあたりは面倒な理論の説明は少ないので面白く読めると思う。
わからないところについて、「そのうちにわかるだろう」と、いかにわからないままで放っておくか、という割り切りも、楽しみながら教養を身につけるうえで大事なことがある。
ただし、この5章の内容は、今回私が特に興味あったところだったので、もっともっとページ数をとって、理論や計算の技のアレコレをもっと詳しく解説してほしかったな、と個人的には思っている。
今回、Python (*1) を使ってそれぞれの級数を使ったπの計算を自分でして確かめながら読んでみた。
普通に何も考えずに浮動小数点で計算すると、もともと仮数部の桁数で制限があるから、それ以上の桁数を出すことはできない。そのうえ、たとえば、アルキメデスの方法だと非常に小さい数字を2乗して4.0から引く引き算が入るため、仮数の桁数の半分の桁数までしか計算が正しくできない。正24576角形を超えると、3.141592645 から精度が改善されなくなり、正393216角形では 3.141593669 となり、それ以上頂点を増やしていってもかえって真値から大きい方向に単調にはずれていき、正805306368角形で π = 6.0となり、その次の正805306368角形で π=0.0 となる。一辺の長さが小さくなりすぎて 0 になってしまったためである。
このへんは、数値計算を行うときに気を付けるべきところだろう。有効数字が何桁であるのか、桁落ちしていないか、あるいは2進数表記と10進数表記の綾で計算誤差が出てそれが累積していないのか、有効桁数を多くとりたいと思えば、ほんの初等の関数を扱うときも数値計算のアレコレを熟知していないとかえって変な結果がでてしまったりする。
有効数字について十分に理解し、連続な対象をデジタルで離散的に考えていることには常々注意が必要である。
5 年前くらいだったか、ある製品の設計シミュレーションツールを Python で自力でちまちま作っていたことがある。そのころは日本ではまだ Python はメジャーではなく、解説本もほとんどなくて少し苦労した。しばらく触ってなかったのだが、今回いい機会だと思い1か月ほど前からリハビリがてらパッケージをインストールして動かしてみている。環境が充実してだいぶん便利になったことを実感し感激している。
そして、Pythonの数値演算ライブラリをネットで改めて調べながら本書の式を実際に試してみて、いまさら知ったのが次の点であった。
(1) 整数型はメモリの許す限り多い桁数を扱える。
(2) 整数型の演算であれば、その桁数で演算が実行できる。
(3) 整数型の平方根が標準 math で用意されている。math.isquare()
(4) decimal ライブラリを使えば 10進数で有効桁数を明確にして計算できる。
(5) mpmath ライブラリを使うと、任意の桁数を扱うことができる。
(6) mpmath には π そのものがすでに用意されている。(100 mpmath one-liners for pi)
上のキャプチャ画面の左側の2番目の関数と3番目の関数はそれぞれ、ウオレスの式と、グレゴリーの式だが、(1) と (2) の点はすぐに発見したので、整数で計算がすべてできるようにインプリした。ウオレスの式は収束が悪いし、どんどん計算に必要な桁数が増えていき効率が悪すぎる。グレゴリーの式はそれに比べればよっぽどましだが、やはり収束は遅い。本書「円周率πの世界」にも書いてあるが、自分の手で実感できる。
そのほかの方法だと、巨大な桁数の整数の平方根が必要だし、普通に平方根をとると、浮動小数点で値が変えるしなぁ、どうしようかな、自分で作るしかないか、と思っていたのだが (3) をいまさら見つけた。もともと用意されていたのだった。あまりに迂闊だった。
「円周率πの世界」にははっきりと言及されていないが、ガウスの公式というのがあって、こちらは収束が速い。ネットに次の pdf があるのを見つけた。大学初等向けだろうか、とてもわかりやすいし、面白い。
円周率πを計算する –アルキメデス,和算,ガウスの方法– (2008年 上越教育大学 中川仁)
桁数の競争ではなくDIYの楽しみを Python で味わう目的ならば、整数型で計算ができるようにガウスの公式をインプリすればかなりいけそうだ。上記 (5), (6) も知ったので、mpmath で1万桁でも10万桁でもすぐに答えを出せるため、結果を比較することもできる。
本当は今日この記事を書くまでにやってしまおうと思っていたのだけど、整数の平方根を [ int(math.sqrt(x)) ;xは整数、桁数は任意 ] で求めらるというガセネタ(試すとわかるが、もっともらしいが誤った結果が得られる)にひっかかって時間をとられてしまってできなかった。
上述のアルキメデスの方法での計算もそうだったが、こんな簡単なことでも、デジタルの世界はいろいろ落とし穴があるものだ。
ガウスの公式で計算した結果やその計算のうんちくを披露しなければいけないだろうか。それとも他のもっと効率のよい公式について追及してみるべきだろうか。(*2)
いくら情熱と膨大なエネルギーを傾けて 62兆 を超える桁数の円周率を計算しても、その値は円周率に近いが円周率ではない。そして、実際には 3.14 でもほとんどの場合は十分だし、ほとんど人が暗記している。さらには、正しいはずのコンピュータもちょっとした不注意で大きく誤った結果を導く。
そんな計算を自分でしてみること自体、無意味だろうか。
いいではないか、デジタルな世の中を渡っていくには、割り切りと打ち切りと熱い心、そして寛容な精神が必要なのだ。
※タイトルの写真は、今年の正月に娘婿と孫娘と一緒に娘が帰って来たときに、孫娘がお絵描きしたものだ。いくつもの色のクレヨンを使い力強く好きなようにのびのびと線を描く姿を見るとほほえましい。
そして、右上のバイキンマンは私が描いたものだ。調子に乗って2枚ほど孫娘と楽しく遊んだ記録を残しておこう。
孫娘はこれ以来、私を尊敬のまなざしで見ている、。。。。ような気がする。それにしても、日本の子供の情操教育でアンパンマンが果たす役割は非常に高いと、感心しているところだ。また、大人が学ぶところだってあるはずだ。同様に、ブルーバックスは日本の少年少女の科学リテラシーを高めるし、大人が科学技術の分野で教養を身につけるにもよい。
■ 注記
(*1) Python
Webシステムでもよく使われているだけでなく、Excelファイルの操作の自動化など様々なツールの開発も可能、機械学習の分野でも広く使われているうえ、信頼できる科学技術計算ライブラリが豊富で、たいていのことができる。
(*2) アルゴリズムについて
Chudnovskyの公式とコンピュータでの計算については、次のような記事があった。
Chudnovsky の公式を用いた円周率の計算用メモ - Qiita
■ 関連 note 記事
この記事が気に入ったらサポートをしてみませんか?