15. Арифметика с плавающей точкой: проблемы и ограничения

Числа с плавающей точкой представлены в компьютерном оборудовании в виде дробей с основанием 2 (двоичных). Например, десятичная дробь 0.625 имеет значение 6/10 + 2/100 + 5/1000, и таким же образом двоичная дробь 0.101 имеет значение 1/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 и встроенная функция repr() выбирали вариант с 17 значащими цифрами, 0.10000000000000001. Начиная с Python 3.1 на большинстве систем теперь выбирается самый короткий из них, и просто отображается 0.1.

Обратите внимание, что это заложено в самой природе двоичных чисел с плавающей точкой — это не ошибка в Python или в вашем коде. То же самое вы увидите во всех языках, использующих арифметику с плавающей точкой аппаратного уровня (хотя некоторые языки могут не отображать разницу по умолчанию или во всех режимах вывода).

Для более красивого вывода вы можете использовать форматирование строк для получения ограниченного количества значащих цифр:

>>> format(math.pi, '.12g')  # получить 12 значащих цифр
'3.14159265359'

>>> format(math.pi, '.2f')   # получить 2 цифры после точки
'3.14'

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

Важно понимать, что в действительности это иллюзия: вы просто округляете отображаемое истинное значение, хранящегося в машине.

Одна иллюзия может порождать другую. Например, поскольку 0.1 не является точно 1/10, сумма трёх значений 0.1 может не дать точно 0.3:

>>> 0.1 + 0.1 + 0.1 == 0.3
False

Кроме того, поскольку 0.1 не может стать ближе к точному значению 1/10, а 0.3 — к точному значению 3/10, предварительное округление с помощью функции round() не может помочь:

>>> round(0.1, 1) + round(0.1, 1) + round(0.1, 1) == round(0.3, 1)
False

Хотя числа невозможно сильнее приблизить к их предполагаемым точным значениям, функция math.isclose() может быть полезна для сравнения неточных значений:

>>> math.isclose(0.1 + 0.1 + 0.1, 0.3)
True

В качестве альтернативы для сравнения грубых приближений можно использовать функцию round():

>>> round(math.pi, ndigits=2) == round(22 / 7, ndigits=2)
True

Двоичная арифметика с плавающей точкой таит в себе много подобных сюрпризов. Проблема с «0.1» подробно описана ниже, в разделе «Ошибка представления». См. Примеры проблем с плавающей точкой для наглядного обзора работы двоичной арифметики с плавающей точкой и проблем, наиболее часто встречающихся на практике. Также см. Опасности плавающей точки для более полного описания других типичных сюрпризов.

Как сказано почти в самом конце: «простых решений нет». Тем не менее, не стоит чрезмерно опасаться чисел с плавающей точкой! Ошибки при операциях с float в Python наследуются от аппаратного обеспечения с плавающей точкой и на большинстве машин составляют не более 1 части из 2**53 за операцию. Этого более чем достаточно для большинства задач, но нужно помнить, что это не десятичная арифметика и каждая операция с float может привести к новой ошибке округления.

Хотя патологические случаи действительно встречаются, при обычном использовании арифметики с плавающей точкой вы, как правило, получите ожидаемый результат, если просто округлите его до нужного количества десятичных цифр. Функции str() обычно достаточно, а для более точного контроля используйте спецификаторы формата метода str.format() из Format string syntax.

Для случаев, когда требуется точное десятичное представление, попробуйте использовать модуль decimal, реализующий десятичную арифметику, подходящую для бухгалтерских задач и приложений с высокой точностью.

Другая форма точной арифметики поддерживается модулем fractions, который реализует арифметику, основанную на рациональных числах (поэтому числа типа 1/3 могут быть представлены точно).

Если вы активно пользуетесь вычислениями с плавающей точкой, вам следует взглянуть на пакет NumPy и многие другие пакеты для математических и статистических операций, поставляемые проектом SciPy. См. <https://scipy.org>.

Python предоставляет инструменты, которые могут помочь в тех редких случаях, когда вы действительно хотите узнать точное значение числа с плавающей точкой. Метод float.as_integer_ratio() выражает значение float в виде дроби:

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

Поскольку это отношение точное, его можно использовать для воссоздания исходного значения без потерь:

>>> x == 3537115888337719 / 1125899906842624
True

Метод float.hex() выражает число с плавающей точкой в шестнадцатеричном виде (по основанию 16), снова предоставляя точное значение, сохраненное на вашем компьютере:

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

Это точное шестнадцатеричное представление можно использовать для восстановления значения float в точности:

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

Поскольку представление является точным, оно полезно для надежного переноса значений между разными версиями Python (независимость от платформы) и обмена данными с другими языками, поддерживающими тот же формат (например, Java и C99).

Ещё одним полезным инструментом является функция sum(), которая помогает уменьшить потерю точности во время суммирования. Она использует расширенную точность для промежуточных шагов округления по мере добавления значений к сумме. Это может повлиять на общую точность, предотвращая накопление ошибок до уровня, когда они влияют на итоговый результат:

>>> 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 + 0.1 == 1.0
False
>>> sum([0.1] * 10) == 1.0
True

Функция math.fsum() идёт дальше и отслеживает все «потерянные цифры» при добавлении значений к сумме, так что результат подвергается только одному округлению. Это медленнее, чем sum(), но даёт большую точность в редких случаях, когда слагаемые большой величины в основном взаимно компенсируют друг друга, оставляя итоговую сумму близкой к нулю:

>>> arr = [-0.10430216751806065, -266310978.67179024, 143401161448607.16,
...        -143401161400469.7, 266262841.31058735, -0.003244936839808227]
>>> float(sum(map(Fraction, arr)))   # Точное суммирование с одним округлением
8.042173697819788e-13
>>> math.fsum(arr)                   # Одно округление
8.042173697819788e-13
>>> sum(arr)                         # Многократное округление с расширенной точностью
8.042178034628478e-13
>>> total = 0.0
>>> for x in arr:
...     total += x                   # Многократное округление со стандартной точностью
...
>>> total                            # Прямое сложение не даёт корректный результат!
-0.0051575902860057365

15.1. Ошибка представления

В этом разделе подробно объясняется пример «0.1» и показывается, как можно самостоятельно выполнить точный анализ подобных случаев. Предполагается базовое знакомство с двоичным представлением чисел с плавающей точкой.

Ошибка представления относится к тому факту, что некоторые (на самом деле большинство) десятичные дроби не могут быть представлены точно как двоичные дроби (по основанию 2). Это основная причина, по которой Python (или Perl, C, C++, Java, Fortran и многие другие) часто не отображает точное десятичное число, которое вы ожидаете.

Почему так происходит? Число 1/10 нельзя точно представить как двоичную дробь. По крайней мере, с 2000 года почти все машины используют двоичную арифметику с плавающей точкой стандарта IEEE 754, и почти все платформы сопоставляют числа float в Python со значениями IEEE 754 binary64 «двойной точности». Эти значения имеют точность 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 <=  2**56 // 10  < 2**53
True

То есть 56 — единственное значение для N, при котором в J остаётся ровно 53 бита. Наилучшее возможное значение для J — это округленное частное:

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

Поскольку остаток больше половины от 10, наилучшее приближение получается округлением в большую сторону:

>>> q+1
7205759403792794

Поэтому наилучшее приближение к 1/10 в двойной точности IEEE 754:

7205759403792794 / 2 ** 56

Деление числителя и знаменателя на два уменьшает дробь до:

3602879701896397 / 2 ** 55

Обратите внимание: поскольку мы округлили в большую сторону, на самом деле это немного больше, чем 1/10; если бы мы не округляли в большую сторону, частное было бы немного меньше 1/10. Но ни в коем случае оно не может быть ровно 1/10!

Таким образом, компьютер никогда не «видит» 1/10: он работает с точной дробью, указанной выше — лучшим приближение двойной точности IEEE 754, которое он может получить:

>>> 0.1 * 2 ** 55
3602879701896397.0

Если мы умножим эту дробь на 10**55, мы увидим значение из 55 десятичных цифр:

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

это означает, что точное число, хранящееся в компьютере, равно десятичному значению 0,1000000000000000055511151231257827021181583404541015625. Вместо отображения полного десятичного значения во многих языках (включая старые версии Python) результат округляется до 17 значащих цифр:

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

Модули fractions и decimal упрощают эти вычисления:

>>> 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'