#
B. 浮動小数点演算、その問題と制限

浮動小数点数は、計算機ハードウェアの中では、 基数を 2 とする (2進法の) 分数として表現されています。 例えば、小数

    0.125

は、 1/10 + 2/100 + 5/1000 という値を持ちますが、これと同様に、 2 進法の分数

    0.001

は 0/2 + 0/4 + 1/8 という値になります。これら二つの分数は同じ値を 持っていますが、ただ一つ、最初の分数は基数 10 で記述されており、 二番目の分数は基数 2 で記述されていることが違います。

残念なことに、ほとんどの小数は 2 進法の分数として正確に表わすことが できません。その結果、一般に、入力した 10 進の浮動小数点数は、 2 進法の浮動小数点数で近似された後、実際にマシンに記憶されます。

最初は基数 10 を使うと問題を簡単に理解できます。分数 1/3 を考えて みましょう。分数 1/3 は、基数 10 の分数として、以下のように近似する ことができます:

    0.3

さらに正確な近似は、

    0.33

です。さらに正確に近似すると、

    0.333

となり、以後同様です。何個桁数を増やして書こうが、結果は決して厳密な 1/3 にはなりません。しかし、少しづつ正確な近似にはなっていくでしょう。

同様に、基数を 2 とした表現で何桁使おうとも、10 進数の 0.1 は 基数を 2 とした分数で正確に表現することはできません。 基数 2 では、1/10 は循環小数 (repeating fraction)

    0.0001100110011001100110011001100110011001100110011...

となります。どこか有限の桁で止めると、近似値を得ることになります。 これこそが、以下のような事態:

    >>> 0.1
    0.10000000000000001

に出くわす理由です。

今日では、ほとんどのマシンでは、0.1 を Python のプロンプトから入力すると 上のような結果を目にします。そうならないかもしれませんが、これは ハードウェアが浮動小数点数を記憶するのに用いているビット数がマシンに よって異なり、Python は単にマシンに 2 進で記憶されている、真の 10 進の 値を近似した値を、されに 10 進で近似して出力するだけだからです。 ほとんどのマシンでは、Python が 0.1 を記憶するために 2 進近似した 真の値を 10 進で表すと、以下のような出力

    >>> 0.1
    0.1000000000000000055511151231257827021181583404541015625

になるでしょう！ Python プロンプトは、文字列表現を得るために 何に対しても `repr()` を使います。 浮動小数点数の場合、 `repr(float)` は真の 10 進値を 有効数字 17 桁で丸め、以下のような表示

    0.10000000000000001

を行います。

`repr(float)` が有効数字 17桁 の値を生成するのは、この値 が (ほとんどのマシン上で) 、全ての有限の浮動小数点数 x について `eval(repr(x)) == x` が成り立つのに十分で、 かつ有効数字 16 桁に丸めると成り立たないからです。

これは 2 進法の浮動小数点の性質です: Python のバグでも、ソースコードのバ グでもなく、浮動小数点演算を扱えるハードウェア上の、すべての言語で同じ 類の現象が発生します (ただし、言語によっては、デフォルトのモードや 全ての出力モードでその差を 表示しない かもしれません)。

Python の組み込みの `str()` 関数は有効数字 12 桁しか生成 しません。このため、この関数を代わりに使用したいと思うかもしれません。 この関数は `eval(str(x))` としたときに x を再現 しないことが多いですが、出力を目で見るには好ましいかもしれません:

    >>> print str(0.1)
    0.1

現実という考えからは、上の表示は錯覚であると気づくのは重要なことです: マシン内の値は厳密に 1/10 ではなく、単に真のマシン内の 表示される値 を丸めているだけなのです。

まだ驚くべきことがあります。例えば、以下

    >>> 0.1
    0.10000000000000001

を見て、`round()` 関数を使って桁を切り捨て、期待する 1 桁に したい誘惑にかられたとします。しかし、結果は依然同じ値です:

    >>> round(0.1, 1)
    0.10000000000000001

問題は、"0.1" を表すために記憶されている 2 進表現の浮動小数点数の値は、 すでに 1/10 に対する最良の近似になっており、値を再度丸めようとしても これ以上ましにはならないということです: すでに値は、`round()` で得られる値になっているというわけです。

もう一つの重要なことは、0.1 が正確に 1/10 ではないため、0.1 を 10 個加算すると厳密に 1.0 にはならないこともある、ということです:

    >>> sum = 0.0
    >>> for i in range(10):
    ...     sum += 0.1
    ...
    >>> sum
    0.99999999999999989

2 進の浮動小数点数に対する算術演算は、このような意外性をたくさん持って います。"0.1" に関する問題は、以下の "表現エラー" の章で詳細に説明します。 2 進法の浮動小数点演算にともなうその他のよく知られた意外な事象に関しては The Perils of Floating Point を参照してください。

究極的にいうと、``容易な答えはありません''。ですが、浮動小数点数の ことを過度に警戒しないでください！ Python の float 型操作における エラーは浮動小数点処理ハードウェアから受けついたものであり、 ほとんどのマシン上では一つの演算あたり高々 2**53 分の 1 です。 この誤差はほとんどの作業で相当以上のものですが、浮動小数点演算は 10 進の演算えはなく、浮動小数点の演算を新たに行うと、新たな 丸め誤差の影響を受けることを心にとどめておいてください。

異常なケースが存在する一方で、普段の浮動小数点演算の利用では、 単に最終的な結果の値を必要な 10 進の桁数に丸めて表示するのなら、 最終的には期待通りの結果を得ることになるでしょう。 こうした操作は普通 `str()` で事足りますし、よりきめ細かな 制御をしたければ、 Python の `%` 書式化演算子についての議論を 参照してください: `%g` 、`%f` 、および `%e` といった 書式化コードでは、浮動小数点数を表示用に丸めるための柔軟性のある、 簡単な手段を提供しています。

#
1 表現エラー

この章では、``0.1'' の例について詳細に説明し、このようなケースに 対してどのようにすれば正確な分析を自分で行えるかを示します。 ここでは、 2 進法表現の浮動小数点数についての基礎的な知識があるものとして 話を進めます。

表現エラー は、いくつかの (実際にはほとんどの) 10 進の小数が 2 進法 (基数 2 ) の分数として表現できないという事実に関係しています。 これは Python (あるいは Perl、 C、 C++、Japa、Fortran 、および その他多く) が期待通りの正確な 10 進数を表示できない主要な理由です:

    >>> 0.1
    0.10000000000000001

なぜこうなるのでしょうか？ 1/10 は 2 進法の分数で厳密に表現することが できません。今日 (2000年11月) のマシンは、ほとんどすべて IEEE-754 浮動小数点演算を使用しており、ほとんどすべてのプラットフォームでは Python の浮動小数点を IEEE-754 における "倍精度 (double precision)" に対応付けます。754 の double には 53 ビットの精度を持つ数が入るので、 計算機に入力を行おうとすると、可能な限り 0.1 を最も近い値の分数に変換 し、J/2**N の形式にしようと努力します。 J はちょうど 53 ビットの精度の整数です。

     1 / 10 ~= J / (2**N)

を書き直すと、

    J ~= 2**N / 10

となります。

J は厳密に 53 ビットの精度を持っている (`>= 2**52` だが `< 2**53` ) ことを思い出すと、 N として最適な値は 56 になります:

    >>> 2**52
    4503599627370496L
    >>> 2**53
    9007199254740992L
    >>> 2**56/10
    7205759403792793L

すなわち、56 は J をちょうど 53 ビットの精度のままに保つ N の唯一の値です。J の取りえる値はその商を丸めた ものです:

    >>> q, r = divmod(2**56, 10)
    >>> r
    6L

残りは 10 の半分以上なので、最良の近似は丸め値を一つ増やした (round up) ものになります:

    >>> q+1
    7205759403792794L

従って、754 倍精度における 1/10 の取りえる最良の近似は 2**56 以上 の値、もしくは

    7205759403792794 / 72057594037927936

となります。丸め値を 1 増やしたので、この値は実際には 1/10 より少し小 さいことに注意してください; 丸め値を 1 増やさない場合、商は 1/10 よりもわずかに小さくなります。しかし、どちらにしろ 厳密に 1/10 ではありません！

つまり、計算機は 1/10 を ``理解する'' ことは決してありません: 計算機が理解できるのは、上記のような厳密な分数であり、 754 の倍精度浮動小数点数で得られるもっともよい近似は:

    >>> .1 * 2**56
    7205759403792794.0

となります。

この分数に 10**30 を掛ければ、 有効数字 30 桁の十進数の (切り詰められた) 値を見ることができます:

    >>> 7205759403792794 * 10**30 / 2**56
    100000000000000005551115123125L

これは、計算機が記憶している正確な数値が、10 進数値 0.100000000000000005551115123125 にほぼ等しいということです。この値を 有効数字 17 桁で丸めると、Python が表示する値は 0.10000000000000001 に なります (もちろんこのような値になるのは、 IEEE 754 に適合していて、C ライブラリで可能な限り正確に値の入出力を行った場合だけです -- 読者の計算機ではそうではないかもしれません！)