15. 浮點數運算：問題與限制

在計算機架構中，浮點數透過二進位小數表示。例如說，在十進位小數中：

0.125

可被分為 1/10 + 2/100 + 5/1000，同樣的道理，二進位小數 ：

0.001

可被分為 0/2 + 0/4 + 1/8。這兩個小數有相同的數值，而唯一真正的不同在於前者以十進位表示，後者以二進位表示。

不幸的是，大多數十進位小數無法精準地以二進位小數表示。一般的結果為，您輸入的十進位浮點數由實際存在計算機中的二進位浮點數近似。

在十進位中，這個問題更容易被理解。以分數 1/3 為例，您可以將其近似為十進位小數：

0.3

或者，更好的近似：

0.33

或者，更好的近似：

0.333

依此類推，不論你使用多少位數表示小數，最後的結果都無法精準地表示 1/3，但你還是能越來越精準地表示 1/3。

同樣的道理，不論你願意以多少位數表示二進位小數，十進位小數 0.1 都無法被二進位小數精準地表達。在二進位小數中， 1/10 會是一個無限循環小數：

0.0001100110011001100110011001100110011001100110011...

只要您停在任何有限的位數，您就只會得到近似值。而現在大多數的計算機中，浮點數是透過二進位分數近似的，其中分子從最高有效位元使開始用 53 個位元表示，分母則是以二為底的指數。在 1/10 的例子中，二進位分數為 3602879701896397 / 2 ** 55，而這樣的表示十分地接近，但不完全等同於 1/10 的真正數值。

由於數值顯示的方式，很多使用者並沒有發現數值是個近似值。Python 只會印出一個十進位近似值，其近似了儲存在計算機中的二進位近似值的十進位數值。在大多數的計算機中，如果 Python 真的會印出完整的十進位數值，其表示儲存在計算機中的 0.1 的二進位近似值，它將顯示為：

>>> 0.1
0.1000000000000000055511151231257827021181583404541015625

這比一般人感到有用的位數還多，所以 Python 將位數保持在可以接受的範圍，只顯示捨入後的數值：

>>> 1 / 10
0.1

一定要記住，雖然印出的數字看起來是精準的 1/10，但真正儲存的數值是能表示的二進位分數中，最接近精準數值的數。

有趣的是，有許多不同的十進位數，共用同一個最接近的二進位近似小數。例如說：數字 0.1 和 0.10000000000000001 和 0.1000000000000000055511151231257827021181583404541015625，都由 3602879701896397 / 2 ** 55 近似。由於這三個數值共用同一個近似值，任何一個數值都可以被顯示，同時保持 eval(repr(x)) == x。

歷史上，Python 的提示字元 (prompt) 與內建的 repr() 函式會選擇上段說明中有 17 個有效位元的數：0.10000000000000001。從 Python 3.1 版開始，Python（在大部分的系統上）可以選擇其中最短的數並簡單地顯示為 0.1。

注意，這是二進位浮點數理所當然的特性，並不是 Python 的錯誤 (bug)，更不是您程式碼的錯誤。只要有程式語言支持硬體的浮點數運算，您將會看到同樣的事情出現在其中（雖然某些程式語言預設不顯示差異，或者預設全部輸出）。

為求更優雅的輸出，您可能想要使用字串的格式化 (string formatting) 產生限定的有效位數：

>>> format(math.pi, '.12g')  # give 12 significant digits
'3.14159265359'

>>> format(math.pi, '.2f')   # give 2 digits after the point
'3.14'

>>> repr(math.pi)
'3.141592653589793'

要了解一件很重要的事，在真正意義上，浮點數的表示是一種幻覺：你基本上在捨入真正機器數值所展示的值。

這種幻覺可能會產生下一個幻覺。舉例來說，因為 0.1 不是真正的 1/10，把三個 0.1 的值相加，也不會產生精準的 0.3：

>>> .1 + .1 + .1 == .3
False

同時，因為 0.1 不能再更接近精準的 1/10，還有 0.3 不能再更接近精準的 3/10，預先用 round() 函式捨入並不會有幫助：

>>> round(.1, 1) + round(.1, 1) + round(.1, 1) == round(.3, 1)
False

雖然數字不會再更接近他們的精準數值，但 round() 函式可以對事後的捨入有所幫助，如此一來，不精確的數值就變得可以互相比較：

>>> round(.1 + .1 + .1, 10) == round(.3, 10)
True

二進位浮點數架構擁有很多這樣的驚喜。底下的「表示法錯誤」章節，詳細的解釋了「0.1」的問題。如果想要其他常見驚喜更完整的描述，可以參考 The Perils of Floating Point（浮點數的風險）。

As that says near the end, "there are no easy answers." Still, don't be unduly wary of floating-point! The errors in Python float operations are inherited from the floating-point hardware, and on most machines are on the order of no more than 1 part in 2**53 per operation. That's more than adequate for most tasks, but you do need to keep in mind that it's not decimal arithmetic and that every float operation can suffer a new rounding error.

While pathological cases do exist, for most casual use of floating-point arithmetic you'll see the result you expect in the end if you simply round the display of your final results to the number of decimal digits you expect. str() usually suffices, and for finer control see the str.format() method's format specifiers in Format String Syntax.

For use cases which require exact decimal representation, try using the decimal module which implements decimal arithmetic suitable for accounting applications and high-precision applications.

Another form of exact arithmetic is supported by the fractions module which implements arithmetic based on rational numbers (so the numbers like 1/3 can be represented exactly).

If you are a heavy user of floating point operations you should take a look at the Numerical Python package and many other packages for mathematical and statistical operations supplied by the SciPy project. See <https://scipy.org>.

Python provides tools that may help on those rare occasions when you really do want to know the exact value of a float. The float.as_integer_ratio() method expresses the value of a float as a fraction:

>>> x = 3.14159
>>> x.as_integer_ratio()
(3537115888337719, 1125899906842624)

Since the ratio is exact, it can be used to losslessly recreate the original value:

>>> x == 3537115888337719 / 1125899906842624
True

The float.hex() method expresses a float in hexadecimal (base 16), again giving the exact value stored by your computer:

>>> x.hex()
'0x1.921f9f01b866ep+1'

This precise hexadecimal representation can be used to reconstruct the float value exactly:

>>> x == float.fromhex('0x1.921f9f01b866ep+1')
True

Since the representation is exact, it is useful for reliably porting values across different versions of Python (platform independence) and exchanging data with other languages that support the same format (such as Java and C99).

Another helpful tool is the math.fsum() function which helps mitigate loss-of-precision during summation. It tracks "lost digits" as values are added onto a running total. That can make a difference in overall accuracy so that the errors do not accumulate to the point where they affect the final total:

>>> sum([0.1] * 10) == 1.0
False
>>> math.fsum([0.1] * 10) == 1.0
True

15.1. Representation Error

This section explains the "0.1" example in detail, and shows how you can perform an exact analysis of cases like this yourself. Basic familiarity with binary floating-point representation is assumed.

Representation error refers to the fact that some (most, actually) decimal fractions cannot be represented exactly as binary (base 2) fractions. This is the chief reason why Python (or Perl, C, C++, Java, Fortran, and many others) often won't display the exact decimal number you expect.

Why is that? 1/10 is not exactly representable as a binary fraction. Almost all machines today (November 2000) use IEEE-754 floating point arithmetic, and almost all platforms map Python floats to IEEE-754 "double precision". 754 doubles contain 53 bits of precision, so on input the computer strives to convert 0.1 to the closest fraction it can of the form J/2**N where J is an integer containing exactly 53 bits. Rewriting

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

as

J ~= 2**N / 10

and recalling that J has exactly 53 bits (is >= 2**52 but < 2**53), the best value for N is 56:

>>> 2**52 <=  2**56 // 10  < 2**53
True

That is, 56 is the only value for N that leaves J with exactly 53 bits. The best possible value for J is then that quotient rounded:

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

Since the remainder is more than half of 10, the best approximation is obtained by rounding up:

>>> q+1
7205759403792794

Therefore the best possible approximation to 1/10 in 754 double precision is:

7205759403792794 / 2 ** 56

Dividing both the numerator and denominator by two reduces the fraction to:

3602879701896397 / 2 ** 55

Note that since we rounded up, this is actually a little bit larger than 1/10; if we had not rounded up, the quotient would have been a little bit smaller than 1/10. But in no case can it be exactly 1/10!

So the computer never "sees" 1/10: what it sees is the exact fraction given above, the best 754 double approximation it can get:

>>> 0.1 * 2 ** 55
3602879701896397.0

If we multiply that fraction by 10**55, we can see the value out to 55 decimal digits:

>>> 3602879701896397 * 10 ** 55 // 2 ** 55
1000000000000000055511151231257827021181583404541015625

meaning that the exact number stored in the computer is equal to the decimal value 0.1000000000000000055511151231257827021181583404541015625. Instead of displaying the full decimal value, many languages (including older versions of Python), round the result to 17 significant digits:

>>> format(0.1, '.17f')
'0.10000000000000001'

The fractions and decimal modules make these calculations easy:

>>> from decimal import Decimal
>>> from fractions import Fraction

>>> Fraction.from_float(0.1)
Fraction(3602879701896397, 36028797018963968)

>>> (0.1).as_integer_ratio()
(3602879701896397, 36028797018963968)

>>> Decimal.from_float(0.1)
Decimal('0.1000000000000000055511151231257827021181583404541015625')

>>> format(Decimal.from_float(0.1), '.17')
'0.10000000000000001'