あざらしとペンギンの問題

主に漫画、数値計算、幾何計算、TCS、一鰭旅、水族館、鰭脚類のことを書きます。

消したい過去が多すぎる

ということはいつも思っている。そんなこんなで、Twitter のすべてのツイートを消し去るアプリケーション rm -rf tweets を作ったという話である。本文は別に読まなくていい。

anhpms.org: rm -rf tweets

制作の動機はタイトルの通りである。自分のツイートを消す分にはそこまでやる必要はなかったのだが、Twitter のアプリケーション登録がいろいろめんどくさくなったので、多少なりとも需要はあるだろうと見て作ったまで。

ご存知の通り、現在の Twitter では過去のツイートを遡るのに制限があり、単純に API を叩くだけでは全てのツイートを消すことはできない。検索するとおそらくトップに出てくる某サイトでもやはりそのように書いてあった。そこで筆者は Twitter の設定のところから取得できるツイート履歴を利用することを考えた。

Twitter のアプリケーション認証を行うため、当然ながらサーバは HTTPS の保護された通信を使う。この程度はしないと Twitter アカウントを預けるのは無理である。

とりあえず手順はリンク先に書いてあるので、後はどうでもいい話をしよう。

Damnatio Memoriae する?

『ダムナティオ・メモリアエ(ラテン語:Damnatio Memoriae)』といえば、古代ローマにおいて国賊や反逆者と断じられた皇帝や政治家などを「歴史上存在しなかったことにする」処分である。訳して『記憶の破壊』などと言われることもある*1。対象とされた人物は、すべての文書、すべての偶像、その他その存在を記したすべての物体から跡形もなく削り取られる。

同様のことは古代ローマ以外でも各地で行われており、近代では赤い世界における事例が目立っているが、別に他ではあり得ないということではない。現代に至っても一部の世界では行われている可能性がある。

「過去改変」と言うと SF の話のように思われるが、逆に情報の拡散に制限があったからこそ可能だったとも言える。未来どころか現代でさえ一度インターネッツに拡散された情報は半永久的に残り続けると言っても過言ではない。もっとも、歴史が知るように、特定の人物に対して『記憶の破壊』が行われたことが史料から裏付けられていることは多く、やはりいつの時代どの世界においても、一度記憶に残ってしまった物事をキレイサッパリ抹殺する*2のは非常に困難なことであるようだ。

過去に『記憶の破壊』が行われた人物に対して、後から『記憶の復元』が行われたことも少なくない。というのも、『記憶の破壊』の対象となる人物 A は、それに値するだけの地位や名声を持った権力者であり、一方で『記憶の破壊』を実行する人物 B は概ね時の最高権力者である。絶えない政争によって B が失脚すれば、次に地位を得た C は「不当に貶められた」A を復権させることによって、B を始末することの正当性のプロパガンダに利用できると考えることもあるだろう。実際には私怨などの合理的でない要因が絡んだりするので、筋の通ったシナリオにならないことも多いのがヒトのアレなところではあるが。

ゲオルギー・ジューコフ*3という人物は旧ソ連の元帥にまでなった人物である。第二次世界大戦において対ドイツの最後の激戦である『ベルリン包囲戦』を指揮し、百万を超える死傷者とともにナチス・ドイツの総統ヒトラーを自殺に追い込んだといえば、母国においてどれだけの名声を得たか察することができよう。しかし、彼は二度も『記憶の抹消』に遭ってもいる。結局は復権されているが、あまり政治に関心を持たなかったらしい彼が『記憶の抹消』を受けた理由は、おそらくは大きすぎる名声と軍に対する絶大な影響力のためだと思われる。最終的に彼は膨大な回顧録を残したので、かえって記憶をしっかり留めることに成功したと言える。

影響力の大きな人物は自分がのし上がるために利用できるが、上に立った後は一転して最大の危険分子になり得るので排除するというのはよくある話で、あるいは、戦時の実力者は戦後の政治のステージに入ると爆弾を抱えたくない上層部から敬遠されるというのもよくある話で、さてこれは一体何の話だったか?

記録の消去

実行してはならないコマンド

次のコマンドを見たことがある人は手を挙げて、あざらしは鰭を挙げて:

# rm -rf /

root でこれを打ってッターン!すると一瞬にして HDD の全部のファイルが吹っ飛ぶことでお馴染みのコマンドである。筆者は何度かやったことがある。もちろんデータはちゃんとバックアップした上で、だが。実際には全部のファイルが消えるとまではいかないと思われるが、少なくとも再起不能にはなる。

件のアプリケーション rm -rf tweets の名前はもちろんこのコマンドに由来している。

筆者はノート PC に Arch Linux を入れて使っているので、OS のインストールは結構手作業でやったし、パーティションは何度か切り直している。そのためか、OS や HDD (SSD) のデータを吹き飛ばすことに対する抵抗が比較的少ない方だと思う。

Git の履歴の消去

私が実際に時々使っているコマンドに

$ rm -rf .git

がある。git は後からやり直しができるように差分を保存することに意義があるのに、それをリセットするのは宜しくないとは思ってはいるが、一度コミットしてしまった以上は復元できてしまうので、変なことをやってしまったと思ったときは『記憶の抹消』をしてしまいたいことがある。

もっとも、最初のコミットより後の変更ならば

$ git reset --soft abc123
$ git commit -m "Damnatio memoriae"

で途中のコミットを一纏めにしてその間をなかったことにできる。やや面倒だが途中のいくつかのコミットだけを一纏めにすることもできる。以前のコミットを push していないならそのまま push して問題は起こらないはず。

以前のコミットのひとつでも既に push してしまっている場合、普通に push するとリモートの記録と整合性が取れないので失敗するが、強制的に push することは可能:

$ git push origin +master

自分だけで作っている場合は構わないと思わないでもないが、さすがにお勧めはしない。

サービスの実装

ここで詳しくは述べないが、件のツイート全消去アプリケーションのフロントエンドは Scala 言語で Play! Framework 2.7.x を使って書いている。また、フロントエンド、バックエンドともアクターシステムおよび非同期処理のライブラリである Akka を使っている。元より Akka ベースの Play! Framework 2.7.x を使うことを前提として書いていたので、正式なリリースを待ってウェブサイトを公開することに決めた。

Play! Framework 2.7.x Documentation

Akka

Twitter 用のライブラリは自分で作った twiscala-twitter というものを使っている。今のところドキュメントは書いておらず、自分で使う以上のことは考えていない。当初から多数の bot アカウントで使うことを想定し、Akka を意識して設計したので、他の多数の Twitter ライブラリとはかなり異なった書き方をしていると思う。

そんなこんなで、細かい話はどうでもいいので、過去を消したくなったら rm -rf tweets を試してみてください。

*1:それぞれは英語の damn と memory の語源である。

*2:『抹殺』に kill の意味はない。

*3:Георгий Константинович Жуков: 1896年12月1日 - 1974年6月18日

平成最後のブログ更新になるかもしれない

今年で平成が終わるらしい。だから何だ、別に宇宙の時がそこまでで終わるわけではない、というところだが、そういえばそんな話もあった気がする。宇宙は置いといて日本という国で生活する以上は、瑣末なことと言うわけにはいかないのが問題だ。

元号が変わることについての感想を言えと問われたとして「面倒くさい」という他に答えようがない。今年に入って某手帳を更新したときに、期限は「平成31年」と書かれていたので、つまりは無期限と考えるべきだろうか?ということを考えたところで、結局はその時が来たら「XX元年」と書き替えられることになるだろう。

もし次の元号が『高志』*1だったとしたら書く度にいろいろ刺されることになるし、『齋藤』や『渡邊』では何が正解がわからないし、簡単な字でいかにもめでたい感じのする『金玉』ではいろいろアレだろう。現代風に『雲母』を「マザークラウド」と読ませるのもありかもしれない。

何はともあれ発表されないことには早くから書き慣れることもできず、改元後も「平成30年5月1日」などと書いてしまうこともあるだろう。という以前に現時点では、「平成30年5月1日」という日は存在しないが「平成30年5月1日」に改元される、という矛盾した書き方をせざるを得ないのである。

私は生まれこそ昭和だが、当時の記憶などは何一つというほどないので、眼鏡をかけたいかにも日本のおじさんらしい風貌のおじさんが『勝訴』ではなく『平成』と書かれた紙をカメラの前に晒す光景を、当人が首相になったときのニュースで見たくらいの印象しかない。なお、ここ数年とそれ以前の記憶は不思議な力によってかき消されているので、今書いたことも私の記憶違いだったかもしれない。

今日の時点で平成は残り3ヶ月となった。題名の通り、これが平成最後の更新になるかもしれない。と思って前回を見たら、約1年前に追加した以来だった。その記事については後で何度か書き直したりこそしているが、約1年も新しい記事がないという有様だ。別に残り3ヶ月でなくともこれが平成最後になる可能性は十分にあったと言えよう。

前回の記事はおそらく国公立大学の入学試験の時期に書かれたもので、私立の合格発表も概ね同時期であったと思われる。たとえ滑り止めにキープされるとしても金はぶん取らなければいけないので、発表を国公立の試験日あたりに行うのは合理的と言える。

そういえばセンター試験が近々廃止されると聞き及んでいる。その後のことは知らないし知る必要も多分ない。私がセンター試験に関して憶えていることでは、得点率が9割を超えたこと、苦手の英語がなぜか高得点だったこと、文系科目が全体的に良かったことがある。よく言われる文系(経済学を除く)にサインコサインタンジェントが不要かどうかは知らないが、理系で入学以来使ったことのない古文と漢文については、セ試国語の中では満点だったと思う。高校でも「古文、漢文、評論文は満点を取れ」と言われてきた憶えがある。

人が苦手だったことについて不要論を主張するのはよく見られることで、自己の正当化において邪魔になることを排除しようという意思が働くのだろうと勝手に思っているが、私においてのそれは紛れもなく英語、即ち English である。自分の Englishability の低さについてはかなりの自覚がある。

私が苦手とする教科では、体育、音楽、美術などは、苦手というより教科として存在することを憎んでいるくらいだった。今でも教科として理不尽に過ぎるので組織レベルでの改善が必要だとは思っている。しかしいずれもセ試にはない教科なので、あくまでセ試の範囲内では私には苦手とすることはなかった。英語もセ試に限れば判断ミスをすることがある程度だった。

話を戻して、私が英語について苦手だったからといって、日本において英語が不要であると主張することは無理がある。実を言えば、誰もが知るウェブ検索サービスの言語は English に設定しているし、その他の多くのサービスでも選択可能な言語に日本語があったとしても English を選択していることが多い。というのも、私が知りたいことについて日本語の情報があまりに少なすぎるか、あったとしてもお粗末に過ぎるものが多く、また、日本語訳の可読性が著しく低いことが多くかえって混乱する、といった事情があるからだ。

いわゆるオタク的な趣味*2に関するものを除けば、私が読んでいる言語の大半は英語である。英語で検索をかけて英語の文章を読まなければ、欲しい情報を得ることが概ね困難になる。だから私は英語を不要だとは言えないし、もっとまともに使えるようにならなければいけないと、割と最近になって思うようになった。

もう終わりに差し掛かって何だが、私は来年度、生物学系の専攻に転向することに決めた。というわけで、私は平成の次の元号では最初の入学となる(はず)。

そこでもやはり英語を使わないわけにはいかないし、専門用語のレベルは他の理工系諸分野と比較しても相当に高いように思う。どのみち入学試験では某公開テストのスコアが必要になるので避けられないことに違いはない。

というわけで、なんだかんだで志を高くして頑張る年にしよう。もし次の元号が『高志』になったとしても刺されないように。

それでは、次の更新では平成でなくなっているかもしれませんが、今年もよろしくお願いします。

*1:北陸地方は古くは『越(こし)』と呼ばれていたが、別の表記として『高志』があり、北陸地方では固有名詞に用いられることがある。

*2:ここでは趣味の範疇を指しているのであって、私自身がオタクであるとは言っていない。

ななついろ☆ドロップアウト

この時期にドロップアウトとは穏やかじゃないですが、筆者一押しの4コマ漫画である『スロウスタート』のアニメは遅れもなく始まり、私だけ取り残されました。私事しか書かないブログですが、今年もよろしくお願いします。(本当は1月中に公開する予定でしたが、諸事情により遅れました。)

さて、今回のタイトルの由来である『ななついろ☆ドロップス』とは、2006年に18禁 PC ゲームとして発売され、2007年にアニメ化された作品です。ディオメディアの初元請作品*1であり、筆者においては初めてくらいに観た UHF アニメでした。そんなわけで記念すべき作品のひとつにとして使いました。

今回取り上げるのは、計算機上で小数を扱う際に重要な問題となってくる『桁落ち』という現象です。

結論から言えば、『桁落ち』とは「近接した値同士で引き算をすると結果の不確かさが異常に大きくなる」ということです。かなり微妙な言い回しですが、ともかく「引き算には気をつけろ」といった意味だと今は思っておいて結構です。

小数と誤差

計算機上では普通は2進小数が使われるのですが、『桁落ち』および『情報落ち』は特に基数に依存するものではないので、我々にとって馴染み深い10進小数*2で説明します。

有効数字と誤差

現代において我々は整数でない数を小数で表し、特に無限小数となる場合は有限桁の小数で打ち切ることを日常的に行っています。例えば、

\displaystyle \sqrt{2} = 1.41421356\dots
\displaystyle \pi = 3.14159265\dots

のように。あるいは、打ち切る桁より下の位に丸め(四捨五入など)、切り上げ、切り下げのいずれかを行って

\displaystyle \sqrt{2} = 1.4142 \times 10^0
\displaystyle \pi = 3.1416 \times 10^0

のように表したりもします。これらの例は6桁目を四捨五入したものです。前にある小数は『仮数部』(mantissa)、10 の肩に乗っている数は『指数部』(exponent) と言います。仮数部は必ず 1 以上 10 未満の範囲に入っていなければなりません。この形式は自然科学や工学において不確さを含む値を表記するときと同じです。ただし、この等号は真の等号ではないことに注意が必要です。

上の例では5桁の数字に意味があるので、『有効数字』5桁の値と言います。

有効数字の桁数は先頭の0でない桁から数えます。例えば、

\displaystyle \frac{23}{4567} = 0.0050361287\dots = 5.0361 \times 10^{-3}

のように、元の小数の6桁目を丸めたのでなく、最初に0でない数字が出る桁から数えて6桁目を丸めた値を、有効数字5桁の値として採用します。なお、このような形式の小数を『浮動小数点数』と言います。計算機における小数の実装のほとんど(C言語の double など)は2進浮動小数点数です。

数に対して有効数字への丸めを行うと、当然ですが真の値に対して誤差が発生します。これは『丸め誤差』と呼ばれます。例えば、

\displaystyle \pi = 3.1416 \times 10^0

は真の値から 0.0000073464\dots だけ大きくなっています。あるいは、これを真の値から 0.00023384\dots \% だけ大きいと表現することもあります。前者を『絶対誤差』、後者を『相対誤差』と言います。それぞれの定義は次に回すとして、実数を有限の桁数で表す以上は、誤差は常について回ります。

計算を有限の桁数で行うことに起因する誤差は『打ち切り誤差』と呼ばれます。一般に、誤差には測定誤差などの元々の数値に含まれる誤差も含まれますが、それらは今回扱う問題の範囲外とします。

なお、切り下げ(非負の場合。負の場合は切り上げ)の場合を除いて、有効数字はその桁数の真の値の表示とは一致しないことがあります。上の例では5桁目の数字が異なっています。

このように「有効数字xx桁」という表現の意味するところは曖昧です。あくまで、精度を考える上で『不確かさ』を含む位以下を除いた数字を数えているだけであり、その桁数がどういう意味を持つかは、誤差の意味や計算上の丸めモードなどによって異なります。しかし、この記事においては以後も「有効数字xx桁」という表現を使うことにします。桁落ちを説明するにはそれで十分だからです。

絶対誤差と相対誤差

さて、少し前の文で『絶対誤差』と『相対誤差』という言葉が出てきました。例えば、

\displaystyle \sqrt{2} = 1.4142 \times 10^0
\displaystyle \pi = 3.1416 \times 10^0
\displaystyle \frac{23}{4567} = 5.0361 \times 10^{-3}

のそれぞれの右辺の値にはすべて丸め誤差が含まれています。このような誤差は大きく2種類に分けられるということです。

その2種類とはすなわち『絶対誤差』と『相対誤差』であり、次のように定義されます。ここで、x^* は真の値、 \tilde{x} は誤差を含む値とします。

絶対誤差:
\displaystyle \epsilon_a = \left| \tilde{x} - x^* \right|

相対誤差:
\displaystyle \epsilon_r = \frac{\epsilon_a}{x^*} = \frac{\left| \tilde{x} - x^* \right|}{x^*}

相対誤差に関する見やすい例では、真の値 100 に対する 100.003、それと真の値 100000 に対する 100003 を誤差の意味で同じと見なします。双方とも言葉では「真の値より 0.003% だけ大きい」と言うことができます。

これらの二つの誤差のうち、浮動小数点数では絶対誤差ではなく相対誤差の方が重要です。というのも、数の大きさが指数部に依存するためです。前文の例では、有効数字5桁として

\displaystyle 100.003 = (1.0000 + 0.00003) \times 10^{2}
\displaystyle 100003 = (1.0000 + 0.00003) \times 10^{4}

であり、両者は仮数部の誤差の意味で全く同じであることがわかります。ここで、1.0000 と 1 が異なる意味を持ち、両者を区別しなければならないことは、学校の理科の授業などで注意されたと思います。あくまで有効数字の意味の取り方ひとつとして、真の値 x^* と誤差を含む値 \tilde{x} との間で

\tilde{x} = 1.00000.99995 \le \tilde{x} < 1.00005
\tilde{x} = 10.5 \le \tilde{x} < 1.5

を意味すると仮定します。このとき、誤差を含む値の『不確かさ』を有効数字の桁数で表しているため、必然的に区別が必要なのです。ここでは有効数字5桁としたので、6桁目を四捨五入すると、

\displaystyle 100.003 = 1.0000 \times 10^{2}
\displaystyle 100003 = 1.0000 \times 10^{4}

となります。

この通り、浮動小数点数では誤差を相対誤差で見ることにより、仮数部に関するものに限定することができます。逆に、もし計算機が5桁までの値しか扱えないならば、仮数部の 0.00005 の不確かさを除去することができません。*3

このような計算機の有限性の制限のため生じる誤差の範囲を『マシンイプシロン』と呼んだりします。ただし、この言葉は絶対誤差の意味でも使われるため、曖昧な表現であることに注意が必要です。

ところで、今まで『誤差』と『不確かさ』という二つの言葉を使ってきましたが、この記事では両者を明確に区別して使うことにします。

『誤差』とは本節の最初に定義された『絶対誤差』または『相対誤差』のみを指すものとします。

対して、有効数字から外れるために誤差として扱われる部分の範囲を『不確かさ』と呼ぶことにします。有効数字という言葉が曖昧性を含むように、不確かさという言葉も曖昧性を含むことは注意しておきます。先にも断った通り、この曖昧性は本記事において重要ではありません。

計算による不確かさの拡大

さて、敢えて『不確かさ』という曖昧な言葉を使うことのいくつかの理由のひとつは、丸め誤差や打ち切り誤差を含む値同士で計算を行ったとき、一般に誤差の範囲が拡大するためです。例えば、有効数字5桁、6桁目を四捨五入するとして、小数を用いた計算では

\displaystyle \frac{1}{3} + \frac{1}{3} = 3.3333 \times 10^{-1} + 3.3333 \times 10^{-1} = 6.6666 \times 10^{-1}

となるところが、真の値を四捨五入すると

\displaystyle \frac{1}{3} + \frac{1}{3} = \frac{2}{3} = 6.6667 \times 10^{-1}

となることなどがあります。第一の計算値については

\displaystyle 6.6666\dots \times 10^{-1} < \left(\frac{2}{3} - 0.00005 \right) \times 10^{-1}

であるため、6桁目の四捨五入によって許容される誤差の範囲に収まっていないことが判ります。このように、有効数字の桁数はあくまで精度の目安であって、厳密な誤差の範囲を表すものではないことに注意が必要です。

ここまで、誤差を含む値を表すのに考えなく等号を使ってきましたが、当然のことながらその等号は正しくありません。また、計算による不確かさの拡大があるため、計算順序も重要な問題である。本題の『桁落ち』とは、最初に書いた通り、計算による不確かさの拡大を意味します。

そこで、計算の推移を \to、真の値から小数への変換を \Rightarrow 単に値が近いことを \approx、値が右に示される区間に入っていることを \in で表すことにします。上の例では、

\displaystyle \frac{1}{3} + \frac{1}{3} \Rightarrow 3.3333 \times 10^{-1} + 3.3333 \times 10^{-1} \to 6.6666 \times 10^{-1}
\displaystyle \frac{1}{3} + \frac{1}{3} \to \frac{2}{3} \Rightarrow 6.6667 \times 10^{-1}
\displaystyle \frac{2}{3} \approx 6.6667 \times 10^{-1}
\displaystyle \frac{2}{3} \in (6.6667 \pm 0.00005) \times 10^{-1}

のように表されます。

相対誤差の上界

ここでひとつ相対誤差について注意すべきことがあります。本来、相対誤差の定義では真の値から見た誤差を測るべきところが、これまでの例では誤差を含む値から見た真の値の存在する範囲を示しているのです。当然ですが両者は同じではありません。実を言えば我々は数値計算によって相対誤差それ自体を厳密に知ることはできません。というのも、我々が計算機上で扱う値は誤差を含む値の方だからです。

とは言っても、相対誤差の取りうる値の上界を見積もることは可能です。まず、誤差を含む値 \tilde{x}仮数部と指数部に分けて

\displaystyle \tilde{x} = m \times 10^e

とします。ここで 1 \le |m| < 10 かつ e は整数です。

このとき、相対誤差の上界は、\epsilon > 0 として

\displaystyle x^* \in (m \pm \epsilon) \times 10^e

のとき

\displaystyle \epsilon_r \le \frac{\epsilon \times 10^e}{(|m| - \epsilon) \times 10^e} = \frac{\epsilon}{|m| - \epsilon} < \frac{\epsilon}{1 - \epsilon}

となります。

ここで  |m| = 1 と置いて上界を取ったことはいささか乱暴に思われるでしょう。もし |m| = 9 ならば、誤差を9倍以上緩く見積もっていることになります。

この上界の緩さを狭めるには、2進小数を使うことが最良の方法となります。なぜなら、2進小数では |m| = 1 に限定されるためです。その場合、誤差の見積もりは高々2倍以内に収まります。計算機では2進小数が使われるのが普通ですが、機械的な実装だけでなく誤差の見積もりの意味でも、2進小数を使う方が有利なのです。

そのような事情はあるものの、この記事ではあくまでわかりやすさのために、我々にとって馴染み深い10進小数を使います。

10進で有効数字  k \ge 6*4(k + 1) 桁目を四捨五入するとした場合、\epsilon = 5 \times 10^{-(k + 1)} なので、相対誤差の上界を乱暴に計算すると

\displaystyle \epsilon_r < \frac{\epsilon}{1 - \epsilon} < 1.000006 \times 5 \times 10^{-k} =  5.00003 \times 10^{-(k + 1)}

となります。

桁落ち

さて、ここまでは準備で、本題はここからです。

最初に言ったように、『桁落ち』とは「相対的に近接した値同士で引き算をすると相対的な不確かさが異常に大きくなる」ということです。問題点をより明らかするために、『近接』と『不確かさ』の前に『相対的な』を付けたことに注意してください。

本節においても、小数は10進有効数字5桁で6桁目を四捨五入するものとします。このように有限個の数字で計算を行うことは、この問題に関して本質的です。

桁落ちの例

まず、次のような大きさの極端に異なる数の足し算を行ってみます。

\displaystyle \sqrt{2} + \frac{23}{4567} \Rightarrow (1.4142 + 0.0050361) \times 10^0 \to 1.41923 \times 10^0 \to 1.4192 \times 10^0

この計算により、誤差そのものは最初の誤差の上限 0.00005 より増えている可能性がありますが、この時点では相対的な不確かさの拡大は小さなものに収まっています。実際、

\displaystyle \sqrt{2} + \frac{23}{4567} = 1.4192496911228\dots \in (1.4192 \pm 0.00005) \times 10^0

なので、厳密な相対誤差の意味でも想定内に収まっていることが判ります。

問題は、先程求めた値 1.4192 \times 10^0 から \sqrt{2} = 1.4142 \times 10^0 を引くことによって発生します。注意すべきはこの時点で持っている値が 1.4192 \times 10^01.4142 \times 10^0 以外の何物でもないことです。あくまで数値のみで計算を行うため、記号としての \sqrt{2} などは忘却の彼方にあるというわけです。それでは実際に引いてみましょう。

\displaystyle 1.4192 \times 10^0 - 1.4142 \times 10^0 \to 0.0050 \times 10^0 = 5.0 \times 10^{-3}

はい。有効数字2桁になりました。5.0 \times 10^{-3}5.0000 \times 10^{-3} の違いを思い出してください。真の値を5桁に丸めたものは 5.0361 \times 10^{-3} なので、不確かさの拡大は厳密な相対誤差の拡大とも合致します。

改めて全体の計算を書くと、

\displaystyle \left( \sqrt{2} + \frac{23}{4567} \right) - \sqrt{2} \Rightarrow \left( (1.4142 + 0.0050361) - 1.4142 \right) \times 10^0
\displaystyle \to (1.4192 - 1.4142) \times 10^0 \to 5.0 \times 10^{-3}

となります。

以上のことをもって『桁落ち』が生じたと言います。

桁落ちはとにかく悪い

桁落ちの何が問題かというと、計算値の相対的な不確かさが他の四則演算に比べて異常に増えることです。

馴染み深い例は微分を数値的に計算するときに現れます。まず微分係数の定義式が

\displaystyle f'(x) = \lim_{h \to 0} \frac{f(x + h) - f(x)}{h}

であることは言うまでもありませんが、この式より(悪い)自明な計算方法として h を適当に小さな値に取って式を評価するということが考えられるでしょう。ここでも計算は全体を通して有効数字5桁で行われるとします。

以下、f(x) = \sqrt{x} として x = 2 ちょうどでの微分係数数値計算することを考えます。真の値は

\displaystyle f'(2) = \frac{1}{2 \sqrt {2}} = 0.35355339059 \in 3.5355339 \times 10^{-1}

です。

h は適当に小さな値として h = 0.001 ちょうどを取るとしましょう。このときの計算は

\displaystyle \frac{\sqrt{2 + 0.001} - \sqrt{2}}{0.001} \Rightarrow \frac{\sqrt{2.0000 + 0.0010} - \sqrt{2.0000}}{0.0010000}
\displaystyle \to (\sqrt{2.0010} - \sqrt{2.0000}) \times 1000.0 \to (1.4146 - 1.4142) \times 1000.0
\displaystyle \to 0.40000 = 4.0000 \times 10^{-1}

となり、大きく異なった値が出ることがわかります。

このように、桁落ちは不確かさが異常に大きくなることにより、どのような意味でも精度の悪い計算結果を生み出してしまうのです。

なお、よく知られた方法として、中心差分を用いた公式

\displaystyle f'(x) = \lim_{h \to 0} \frac{f(x + h) - f(x - h)}{2 h}

を用いて計算すると良い計算値が出せるというものがあります。これによれば、同じ有効数字5桁のもとで f'(2) は次のように計算されます:

\displaystyle \frac{\sqrt{2 + 0.001} - \sqrt{2 - 0.001}}{2 \times 0.001} \Rightarrow \frac{\sqrt{2.0000 + 0.0010000} - \sqrt{2.0000 - 0.0010000}}{2 \times 0.0010000}
\displaystyle \to (\sqrt{2.0010} - \sqrt{1.9990}) \times 500.00 \to (1.4146 - 1.4139) \times 500.00
\displaystyle \to 0.3500 = 3.5000 \times 10^{-1}

先程より少しマシな結果が得られていますね。結果が2桁まで合っていることが見て取れます。

それにはちゃんとした理由があります。次の図より直感的にはわかってもらえることを期待します。

f:id:azapen6:20151222103901p:plain

しかし、この設定でこれ以上良い結果を得ることはおそらく不可能でしょう。なぜなら、有効数字5桁では 2.000 + 0.0010000 \Rightarrow 2.0010 となるように、0.0010000 のうち2桁より下の情報が失われてしまうためです。それを次節で説明します。

桁落ちと情報落ち

『桁落ち』と関係する言葉として『情報落ち』というものがあります。これは、桁数が有限であるために元々の数が持っていた情報が計算によって失われることを指して言う言葉です。例えば、先に出した例:

\displaystyle \sqrt{2} + \frac{23}{4567} \Rightarrow (1.4142 + 0.0050361) \times 10^0 \to 1.41923 \times 10^0 \to 1.4192 \times 10^0

では後ろの桁の情報 361 が失われていることから『情報落ち』が発生しています。これは一般には『積み残し』と呼ばれ、桁違い*5に大きさの異なる数を足し引きしたときに発生するものです。

しかし、この段階では有効桁数はほとんど失われていないことに注意が必要です。つまり、『情報落ち』こそ発生していますが、『桁落ち』が発生したとまでは言えません。

『桁落ち』が発生するのは近接した値の引き算を行ったときする。上の結果から \sqrt{2} \in 1.4142 \times 10^0 を引いて

\displaystyle 1.4192 \times 10^0 - 1.4142 \times 10^0 \to 0.0050 \times 10^0 = 5.0 \times 10^{-3}

としたとき、初めて『桁落ち』が発生します。

『桁落ち』の原因は上位桁が等しい数の有効数字がキャンセルされてしまうことです。すなわち、データとして見れば『情報落ち』の範疇に含まれます。しかし、『桁落ち』は単なる情報落ちに留まらず、有効数字が著しく失われることを特に指して言う言葉です。有効数字がほとんど失われない『情報落ち』も避けるに越したことに違いはありませんが、『桁落ち』はそれよりずっと重要な問題であり、小さくない犠牲を払ってでも避けられるべきなのです。

桁落ちの回避

桁落ちを回避する方法は、とにもかくにも、その原因となる「近接した値同士の引き算」を行わないようにすることです。それには、計算する対象についての事前知識があると、有利に進めることができる場合があります。

完全に同じ値の引き算を記号的に除去する

前節のひとつ前例では、完全に同じとわかっている値を記号的に引いてしまってから、小数としての計算を行うことによって、桁落ちを回避できます。

\displaystyle \frac{1}{\left( \sqrt{2} + \frac{23}{4567} \right) - \sqrt{2}} \to \frac{1}{\left( \sqrt{2} - \sqrt{2} \right) + \frac{23}{4567}} \to \frac{4567}{23}
\Rightarrow \frac{4.5670 \times 10^3}{2.3000 \times 10^1} \to 1.9857 \times 10^2

馬鹿馬鹿しい例ですが、記号のまま処理できる部分を先に記号のまま処理しておくことは、桁落ちを回避する上で重要なことです。あるいは、数学的、物理的な知見からの不変性、対称性が利用できるかもしれません。

とにかく引き算を使わずに済ませる

桁落ちを回避するための最も一般的な方策は、足し算と引き算の選択肢では必ず足し算を選ぶこと、また、足し算と乗除算を組み合わせることで引き算を回避できるならばそれを行うことです。

最も馴染み深いのは二次方程式の解の公式でしょう。中学校で習うように、二次方程式

\displaystyle ax^2 + bx + c = 0\;(a \neq 0)

の二つの解は

\displaystyle x = \{ r, s \} = \frac{-b \pm \sqrt{b^2 - 4ac}}{2a}

という公式によって求められます。ここでは二つの解が実数である場合のみを扱うものとして、r \ge s とします。

以下において、コロン付き等号 := を用いて定義されている式は計算式であって、それを評価して得られる値とは区別されるものとします。

b > 0 のとき、r の直接の計算式

\displaystyle r_1 := \frac{-b + \sqrt{b^2 - 4ac}}{2a}

において、分子の計算

\displaystyle -b + \sqrt{b^2 - 4ac}

は引き算となります。ひとつの場合として b平方根の中の値に近い場合、すなわち b^2 \gg 4ac の場合に桁落ちが生じます。もうひとつの場合としては、b^2 \approx 4ac のときは平方根の中の引き算で桁落ちが発生します。ここでは前者について考えることにします。

とりあえず、一方の解 r を得るために r_1 を直接計算することは可能な限り避けなければなりません。

他方、s の計算

\displaystyle s_1 := \displaystyle x = \frac{-b - \sqrt{b^2 - 4ac}}{2a}

においては、分子は負の値同士の足し算となるため、桁落ちが生じないことが判ります。

よって、最初に s_1 だけを計算し、もう一方を解と係数の関係

\displaystyle rs = \frac{c}{a}

から

\displaystyle r_2 := \frac{c}{a s_1}

の計算によって求めることで、桁落ちを回避することが可能です。

例としてこの場合の問題が生じる二次方程式 A

\displaystyle x^2 + 10 x + 0.001 = 0

について、r_1r_2 の式を32ビット浮動小数点数(10進有効数字7桁程度)の範囲内で計算してみます。真の解の小数値は次の通り:

\displaystyle r = -1.0000100002000050\dots \times 10^{-4}, s = -9.9998999989999799\dots

コードと計算結果:

https://ideone.com/FFul9l

 r_1 = -1.001358 \times 10^{-4}, s_1 = -9.999900 \times 10^0
 r_2 = -1.000010 \times 10^{-4}

真の値と比較すると、r_1 はかなり精度が悪く、r_2 は10進有効数字7桁程度では十分な精度を持った値であることが判ります。s_1 について後ろの桁の値が真の値と異なっていますが、誤差の意味では精度から想定される範囲を逸脱するものではありません。

分母、分子の有理化

二次方程式の解の公式のような平方根の引き算を伴う計算では、分母または分子の有理化を行うことによって、桁落ちが回避できる場合があります。

\displaystyle r_1 := \frac{-b + \sqrt{b^2 - 4ac}}{2a}

の分子

\displaystyle -b + \sqrt{b^2 - 4ac}

を有理化するとは、分母と分子の両方に

\displaystyle -b - \sqrt{b^2 - 4ac}

を掛けることです。少し計算すると、分子は

\displaystyle b^2 - (b^2 - 4ac) = 4ac

分母は

\displaystyle -2a \left( b + \sqrt{b^2 - 4ac} \right)

となって、

\displaystyle r_3 := -\frac{2c}{b + \sqrt{b^2 - 4ac}}

という計算式が得られます。

再び前節の二次方程式 A

\displaystyle x^2 + 10 x + 0.001 = 0

について、r_1r_3 の式を32ビット浮動小数点数(10進有効数字7桁程度)の範囲内で計算してみます。

真の値: r = -1.0000100002000050\dots \times 10^{-4}

コードと計算結果:

https://ideone.com/FFul9l

 r_1 = -1.001358 \times 10^{-4}
 r_3 = -1.000010 \times 10^{-4}

期待通り r_3 は十分な精度の値になっていることが判ります。

ところで、 r_3 についてよく見ると、前節と同じ

\displaystyle s_1 := \frac{-b - \sqrt{b^2 - 4ac}}{2a}

を用いて

r_3 = \frac{2c}{2a s_1} = \frac{c}{a s_1}

と書けることが判ります。すなわち、r_3 は前節の r_2 をひとつの式にまとめたものと見なすことができます。計算手順の違い(コンパイラ依存なので特定はできない)はあるものの、桁落ちが生じないように改良された計算では大きな問題にはなっていません。

級数展開を使う

前の二次方程式の解の公式の問題を、級数展開を使って解決することも考えられます。三度、問題の引き算(平方根の外)

\displaystyle -b + \sqrt{b^2 - 4ac}

で桁落ちが発生するときというのは、b^2 \gg 4ac となるときです。少し変形すると、

\displaystyle -b + \sqrt{b^2 - 4ac} = -b + b \sqrt{1 - \frac{4ac}{b^2}} = b \left( \sqrt{1 - \frac{4ac}{b^2}} - 1 \right)

において

\displaystyle \frac{4ac}{b^2} \ll 1

となるときとも言えます。このとき、一般化二項定理により得られる級数展開

\displaystyle \sqrt{1 - \frac{4ac}{b^2}} = 1 - \frac{2ac}{b^2} - \frac{2 a^2 c^2}{b^4} \dots

は良い近似になっていると考えられます。式の中に引き算が出てきますが、近接した値でないことが保証されているので、桁落ちが問題になる場合ではありません。

というわけで、級数を上の3項で打ち切って代入すると、

\displaystyle -b + \sqrt{b^2 - 4ac} = -\frac{2ac}{b} - \frac{2 a^2 c^2}{b^3}

から

\displaystyle r_4 := -\frac{c}{b} \left( 1 + \frac{ac}{b^2} \right)

という計算式が得られます。この式において括弧の中は引き算になる可能性がありますが、仮定より近接した値でないことが保証されているので、それによる桁落ちは問題にはなりません。

というわけで、三度、問題の二次方程式 A

\displaystyle x^2 + 10 x + 0.001 = 0

について、r_1r_4 の式を32ビット浮動小数点数(10進有効数字7桁程度)の範囲内で計算してみます。

真の値:  r = -1.0000100002000050\dots \times 10^{-4}

コードと計算結果:

https://ideone.com/FFul9l

 r_1 = -1.001358 \times 10^{-4}
 r_4 = -1.000010 \times 10^{-4}

期待通り r_4 は十分な精度の値になっていることが判ります。計算式 r_4 の導出には近似を用いているので r_2, r_3 とは記号的にも異なっていますが、結果としては同程度の精度の値が出ていることは注目に値するでしょう。

加速法を使う

微分方程式の差分解法のような、微分係数の差分近似を含む計算を行うとき、桁落ちはほぼ宿命的に問題となります。この前の記事で微分係数を加速法によって高精度で求める方法を説明しています。

アクセルほど微分精度は加速していくよ! - あざらしとペンギンの問題

引き算には気をつけろ

桁落ちを回避する方法はここに書いた以外にも様々ありますが、どんな場合でも上手くいくような万能の処方箋はありません。とにもかくにも引き算に気をつけるという意識を持つことが第一です。可能ならば引き算を行わないようにすること、回避できない場合は可能な限り回数を減らす、また、ソートを行って順序を変えるなど、工夫して計算を行うことが求められます。

有限性の壁

計算機において有効数字、桁落ち、情報落ちといった概念が出てくる理由は、物質的な制約のために扱えるデータ量が限られることにあります。特には浮動小数点演算装置(FPU)のレジスタが(符号、指数を含めて)少ない桁数の分しかないということです。レジスタを大きくすれば桁数を増やすことはできますが、その分だけ数のデータを保持し、また、処理するためにより大きな、あるいは集積度の高い回路が必要となります。回路の素子は原子で構成されていて、宇宙には有限個の原子しかないので、無限の桁数を扱うことは当然ながら不可能です。それは、数のデータをチップの外、例えばハードディスクやクラウドストレージに入れたとしても同じです。

さらに言えば、データの読み出しはつまるところ測定なので、物理学的な不確かさが本質的に有効数字を限定しているとも考えられます。この不確かさは測定誤差ばかりでなく、量子力学的、統計力学的な不確かさ、つまり、この世界においてどうやっても克服できない不確かさが含まれます。桁数というものが名目上存在しないアナログ計算機でさえ、有効数字の束縛から逃れ得るものではないのです。

まとめ

桁落ちとは、浮動小数点数演算において近接した値の引き算を行うことによって有効数字が失われ、計算値の不確かさが異常に大きくなることでした。桁落ちは数値計算における重要な問題であり、可能な限り回避しなければならないものです。本記事ではいくつかの方法を示しました。前の記事に書いた加速法も桁落ちを回避する手段となり得ます。その他様々な方法がありますが、万能の処方箋はなく、問題を見て工夫すべき場合も少なくありません。本記事の教訓を一言で言えば、

引き算には気をつけろ

これは数値計算を行う上で常に意識しなければならないことです。

桁落ち、情報落ちは有限の世界で計算を行う以上は避けられない問題です。一見すると数値計算における細かく面倒な問題のようで、実は計算の本質に関わる問題であるというのが私の思うところです。

参考文献

1. 伊理 正夫、藤野 和建『数値計算の常識』

定評のある本で、初等的な数値計算において必ず考慮されるべき事柄について、様々な事例を題材に理論と実際の両面から考察を加えている。数値計算の『教科書』ではないが、予備知識は特に必要としない。

*1:当時の社名はスタジオ・バルセロナ

*2:他の記数法を用いる文化圏では、10進小数が馴染み深いとは言えないかもしれません。

*3:多倍長演算によって任意の桁数を扱うことが可能ですが、計算の速度は極端に落ちます。

*4:IEEE 754 の定める 32 ビット浮動小数点数は10進で有効数字7桁程度

*5:ここでは読んで字のままの意味