Ещё одна причуда Python, исследование её подноготной и попытка понять, почему так случается.
Недавно в сети X был популярен этот твит (см. скриншот), и я обратил внимание. Это очередной сюрприз в Python, связанный с характерными для него уникальными деталями реализации.
Известно, насколько приблизительна математика для чисел с плавающей точкой, по причине присущих ей пределов представления. В других языках int зачастую неявно повышаются до double, благодаря чему операции сравнения получаются непротиворечивыми. Но в Python этот процесс осложняется из-за того, что в языке присутствуют целые числа с неограниченной разрядностью, и их использование приводит к неожиданным результатам.
В Python целочисленное значение сравнивается с представлением числа с плавающей точкой (двойной точности). При этом может снижаться точность, что вызывает такие расхождения. В этой статье мы подробно разберём, как такие сравнения выполняются в CPython. Так у нас будет отличная возможность подробнее исследовать эти сложности.
Итак, данная статья состоит из следующих разделов:
Мы быстро вспомним, что представляет собой IEEE-754 — формат для чисел двойной точности. Именно он используется для представления в памяти чисел с плавающей точкой
Проанализируем представление трёх чисел в соответствии со стандартом IEEE-754
Изучим алгоритм CPython, предназначенный для сравнения целых чисел и чисел с плавающей точкой
Проанализируем три тестовых сценария в контексте алгоритма CPython
Освежим знания по формату IEEE-754
Внутрисистемно в CPython используется тип double
из C, и именно в нём хранится значение с плавающей точкой. Именно поэтому данный формат соответствует стандарту IEEE-754 для представления этих значений.
Быстро напомню, что сказано в стандарте о представлении значений с двойной точностью. Если вы всё это уже знаете, можете смело переходить к следующему разделу.
Значения с двойной точностью представляются каждое при помощи 64 разрядов, а именно:
1 разряд используется в качестве знакового бита
В 11 разрядах представлен порядок. При указании порядка используется смещение в 1023, т.e., фактическое значение порядка на 1023 меньше, чем значение, представленное с двойной точностью.
В последних 52 разрядах представлена мантисса.
Этот формат подкреплён обширной теорией, описывающей, как в данный формат преобразуются числа и каковы различные способы округления, применяемые в нём — на эту тему есть целая статья, которую вам стоит прочесть, если вы интересуетесь темой.
Здесь я освещу основы, но предупреждаю, что эта статья не является полным разбором формата IEEE-754. Я попробую объяснить, как числа преобразуются в формат IEEE-754 в обычном случае, когда рассматриваемое число больше единицы.
Преобразование числа с плавающей точкой в формат IEEE-754
Чтобы преобразовать число с плавающей точкой в его двоичное представление, соответствующее формату IEEE-754, нужно проделать следующие три шага. В качестве примера разберём эти шаги для числа 4.5.
Шаг 1: Преобразуем число в двоичный формат
Чтобы преобразовать значение с плавающей точкой из десятичной в двоичную систему, требуется отдельно преобразовать целочисленную часть и отдельно дробную часть, а затем в двоичном представлении записать их через двоичную точку (.).
4 в двоичном представлении будет 100, а 0.5 будет 1. Следовательно, можно представить 4.5 как 100.1.
Шаг 2: Нормализация двоичного представления
На следующем этапе нормализуем двоичное представление, полученное при предыдущем шаге. Нормализация сводится к тому, что у двоичного представления должен быть всего один разряд до двоичной точки.
В нашем примере двоичная точка расположена так: 100.1, то есть, у нас 3 разряда до двоичной точки. Чтобы нормализовать его, необходимо выполнить сдвиг вправо на 2 разряда. Таким образом, нормализованное представление примет вид 1.001 * 2^2.
Шаг 3: Извлечение мантиссы и порядка из нормализованного представления
Из нормализованного представления прямо следуют мантисса и порядок. В примере, который мы сейчас рассматриваем, нормализованное значение равно 1.001. Чтобы его получить, нам потребовалось сдвинуть исходное значение на два разряда вправо — следовательно, порядок равен 2.
В формате IEEE-754 мантисса имеет форму 1.M. В сущности, ведущая единица перед двоичной точкой не входит в двоичное представление, так как известно, что здесь всегда будет 1. Поэтому, не сохраняя её, экономим 1 разряд под нужные данные. Следовательно, мантисса здесь составит 001. Поскольку ширина мантиссы составит 52 разряда, остальные 49 разрядов здесь будут заполнены нулями.
Итоговое двоичное представление
Знаковый разряд: Число положительное, следовательно, знаковый разряд будет равен 0.
Порядок: порядок равен 2. Правда, при представлении числа в формате IEEE-754, к нему потребуется добавить сдвиг в 1023, что даёт нам 1025. Таким образом, двоичное представление порядка составит 10000000001.
Мантисса: мантисса равна 001000…000
Соответственно, двоичное представление числа 4.5 в соответствии с форматом IEEE-74 будет:
(0) (10000000001) (001000…000)
(Для наглядности я заключил разные компоненты числа в круглые скобки.)
Анализ представления IEEE-754 на трёх тестовых примерах
Давайте проверим, что мы выучили о формате IEEE-754 для чисел двойной точности и попробуем выяснить, как были бы представлены три числа из трёх тестовых примеров. При этом их будет особенно удобно проанализировать в контексте алгоритма CPython для сравнения чисел с плавающей точкой.
Представление числа 9007199254740992.0 в формате IEEE-754
Возьмём значение из первого тестового сценария: 9007199254740992.0.
Оказывается, что 9007199254740992 — это 2^53. Таким образом, в двоичном представлении это будет единица, за которой следуют 53 нуля: 10000000000000000000000000000000000000000000000000000.
Нормализованная форма в данном случае 1.0 * 2^53, мы получим её, сдвинув разряды на 53 позиции вправо и представим порядок как 53.
Также напомню, что ширина мантиссы составляет всего 52 разряда, а мы сдвигаем значение на 53 разряда вправо — и это значит, что ведущий знаковый бит в нашем исходном значении будет потерян. Но, поскольку здесь все разряды равны 0, ничего страшного не случится, и мы всё равно сможем представить 9007199254740992.0 с достаточной точностью.
Давайте обобщим, из чего состоит 9007199254740992.0:
Знаковый разряд: 0
Порядок: 53 + 1023 (сдвиг) = 1076 (10000110100 в двоичной системе)
Мантисса: 0000000000000000000000000000000000000000000000000000
Представление в стандарте IEEE-754: 0100001101000000000000000000000000000000000000000000000000000000
Представление 9007199254740993.0 в формате IEEE-754
Во втором тестовом сценарии рассмотрим значение 9007199254740993.0, на единицу выше, чем в предыдущем примере. Чтобы получить его двоичное представление, достаточно всего лишь добавить единицу к двоичному представлению предыдущего значения:
10000000000000000000000000000000000000000000000000001.0
Опять же, чтобы преобразовать это число в нормализованный вид, нужно выполнить сдвиг на 53 разряда вправо, в результате получится, что порядок равен 53.
Но ширина мантиссы составляет всего 52 разряда. При сдвиге нашего значения на 53 разряда наименьший значащий разряд (LSB), равный 1, будет потерян, и в мантиссе останутся одни нули. В результате получатся следующие компоненты:
Знаковый разряд: 0
Порядок: 53 + 1023 (сдвиг) = 1076 (10000110100 в двоичной системе)
Мантисса: 0000000000000000000000000000000000000000000000000000
Представление в стандарте IEEE-754: 0100001101000000000000000000000000000000000000000000000000000000
Обратите внимание, что представление в стандарте IEEE-754 для 9007199254740993.0 точно такое же, как для 9007199254740992.0. Таким образом, в памяти число 9007199254740993.0 на самом деле представлено как 9007199254740992.0.
Вот почему Python выдаёт для 9007199254740993 == 9007199254740993.0 результат False, ведь он трактует эту операцию как сравнение 9007199254740993 и 9007199254740992.0.
Об этой частности мы ещё вспомним, когда будем анализировать 2-й сценарий после того, как обсудим алгоритм CPython.
Представление 9007199254740994.0 в формате IEEE-754
В третьем и последнем тестовом сценарии разберём число 9007199254740994.0, которое на единицу больше предыдущего значения. Его двоичное представление можно получить, добавив 1 к двоичному представлению 9007199254740993.0. В результате получим:
10000000000000000000000000000000000000000000000000010.0
Здесь также 54 разряда до двоичной точки. Чтобы преобразовать это число в нормализованный вид, его нужно сдвинуть вправо на 53 разряда. Таким образом, порядок будет 53.
Обратите внимание: на этот раз второй наименьший значащий разряд равен 1. Следовательно, если сдвинуть это число вправо на 53 разряда, он превращается в наименьший значащий разряд мантиссы (ширина которой составляет 52 разряда).
Компоненты выглядят так:
Знаковый разряд: 0
Порядок: 53 + 1023 (сдвиг) = 1076 (10000110100 в двоичной системе)
Мантисса: 0000000000000000000000000000000000000000000000000001
Представление в стандарте IEEE-754: 0100001101000000000000000000000000000000000000000000000000000001
Представление числа 9007199254740994.0 в формате IEEE-754, чего не скажешь о 9007199254740993.0, достижимо без какой-либо потери точности. Следовательно, результат 9007199254740994 == 9007199254740994.0 на Python будет равен True.
Как в СPython реализуется сравнение чисел с плавающей точкой
Выше мы рассмотрели представление трёх чисел в формате IEEE-754 и упомянули, что в Python не выполняется неявное расширение целых чисел до чисел двойной точности. Именно этим и объясняются неожиданные результаты тестовых сценариев, которые обсуждались в сети X.
Но, если вы хотите подробнее разобраться в этой теме и узнать, как именно в Python сравниваются целые числа и числа с плавающей точкой — читайте далее.
В CPython есть функция, в которой числа с плавающей точкой сравниваются как объекты — это float_richcompare. Она определяется в файле floatobject.c.
Это довольно большая функция, так как в ней требуется обрабатывать огромное количество случаев. На скриншоте выше заметны комментарии, из которых следует, насколько же сложно сравнивать числа с плавающей точкой и целые числа. Я разделю эту функцию на более мелкие фрагменты и объясню их один за другим.
Часть 1: сравнение Float и Float
Python — это язык с динамической типизацией. Следовательно, когда интерпретатор вызывает эту функцию для обработки оператора ==, она «не представляет», каковы конкретные типы операндов. Функция должна определить конкретные типы и действовать соответственно. В первой части сравнивающей функции проверяется, относятся ли оба объекта к типу float (число с плавающей точкой), и в этом случае речь идёт о сравнении значений с двойной точностью, лежащих в их основе.
На следующем рисунке показан первый фрагмент функции:
Тут подсвечены и пронумерованы различные части. Давайте по порядку разберёмся, что здесь происходит:
В i и j содержатся базовые значения с двойной точностью для v и w, а в r содержится окончательный результат сравнения.
Когда эту функцию вызывает интерпретатор, уже известно, что аргумент v относится к типу double, для этого здесь и предусмотрено условие для проверки утверждения (assert condition). Причём, базовое значение двойной точности присваивается i.
Далее проверяется, является ли w также объектом float в Python. Если так, то можно просто присвоить базовое значение двойной точности переменной j и напрямую сравнить два этих значения.
Но, если w не является объектом float в Python, то нужно предусмотреть раннюю проверку, на тот случай, если в качестве значения v задано infinity. В таком случае переменной j присваивается значение 0.0, поскольку конкретное значение w непринципиально при сравнении его с бесконечностью.
Часть 2: Сравнение Float и Long
Если w в Python является не числом с плавающей точкой, а целым числом, то всё становится ещё интереснее. Алгоритм рассчитан на то, чтобы в целом ряде случаев избегать прямого сравнения. Рассмотрим все эти случаи по очереди.
В первом случае обрабатывается ситуация, когда v и w противоположны по знаку. В таком случае нет необходимости сравнивать конкретные значения v и w, именно так и устроен код в следующем листинге:
Давайте охарактеризуем каждую из пронумерованных частей этого листинга:
В этой строке мы убеждаемся, что w — это целочисленный объект Python.
vsign и wsign имеют тот же знак, что v и w соответственно.
Если два значения отличаются знаком, то их величины нет нужды сравнивать — ответ выводится по одному лишь знаку.
Часть 2.1: w — огромное целое число
Далее также рассмотрим простой случай, в котором число w настолько велико, что определённо окажется больше, чем любое значение двойной точности (кроме бесконечности, но случай с бесконечностью уже рассмотрен выше). Данная ситуация рассмотрена в следующем листинге:
Функция
_PyLong_NumBits
возвращает, сколько бит было использовано для представления целочисленного объекта w. Но, если целочисленный объект так велик, что в типе size_t не умещается нужное количество разрядов, то такая ситуация обозначается возвратом -1.Таким образом, если число w просто огромно, то оно определённо будет больше, чем любое значение двойной точности. В данном случае сравнение не составляет труда.
Часть 2.2: Когда w помещается в тип double
Следующий случай также прост. Если w умещается в тип double, то мы можем просто привести его к значению двойной точности, а затем сравнить два значения двойной точности. Код показан в следующем листинге:
Стоит подчеркнуть, что 48 разрядов используется в качестве критерия, по которому определяют, расширить ли до числа двойной точности целочисленное значение, лежащее в основе w. Не уверен, почему в данном случае было выбрано 48 разрядов, а не 54 (в последнем случае вышеописанного отказа не произошло бы) или какое-либо иное значение. Может быть, требовалась железная уверенность, что вероятность переполнения исключена на любой платформе? Не знаю.
Часть 2.3 Сравнение порядков
Большинство этих случаев специально подбиралось так, чтобы не приходилось фактически сравнивать значения int и float. Последний и окончательный выход, позволяющий избежать прямого сравнения чисел — сравнивать их порядки. Если порядки у них в нормализованной форме отличаются, то действительно, достаточно сравнить порядки. Если порядок числа больше, то и число больше (предполагается, что v и w положительны).
Чтобы извлечь порядок значения с двойной точностью в v, в коде CPython применяется функция frexp из стандартной библиотеки C. Функция frexp принимает значение двойной точности x и делит его на две части: нормализованное дробное значение f и порядок e, так что x = f * 2 ^ e. Диапазон нормализованной дробной части f равен [0.5, 1).
Например, рассмотрим значение 16.4. Его двоичное представление — 10000.0110011001100110011. Чтобы превратить его в нормализованную дробную часть в соответствии с определением frexp, нужно сдвинуть правый бит на 5 позиций, и в результате получится 0.100000110011001100110011 * 2^5. Следовательно, нормализованное дробное значение будет 0.512500, а порядок — 5.
Что касается w — это целое число, поэтому его порядок просто равен количеству разрядов в нём, то есть. nbits. Теперь давайте рассмотрим код для этого случая:
Давайте по отдельности разберём все пронумерованные элементы:
Сначала нужно убедиться, что значение v положительное.
Затем при помощи функции frexp число i разделяется на нормализованную часть и порядок
Далее проверяется, не является ли порядок v отрицательным или меньшим, чем порядок w. Если какое-то из этих условий выполняется, то v определённо меньше, чем w.
Напротив, если порядок v больше, чем порядок w, то v больше, чем w.
Часть 2.4: Порядок обоих чисел одинаков
Это последняя и наиболее сложная часть нашей функции сравнения. Здесь у нас есть только один выход — непосредственно сравнить два числа.
На данном этапе уже известно, что порядок v и w одинаков, и это значит, что в обоих этих числах для представления целочисленных значений используется одинаковое количество разрядов. Таким образом, чтобы сравнить v и w, достаточно сравнить целочисленное значение v с w.
Чтобы разделить v на целочисленную и дробную часть, в коде CPython используется функция modf из стандартной библиотеки C. Например, если функция modf получит на вход значение 3.14, то разделит его на целочисленную часть 3 и дробную 0.14.
Как правило, должно быть достаточно сравнить целочисленные значения v и w, но только не в случае, если эти значения равны. Например, если v=5.75 и w=5, то их целочисленные значения равны, но фактически v больше w. В таких случаях при сравнении требуется учитывать и дробную часть значения v.
Для этого в коде CPython делается простой фокус. Если дробная часть v ненулевая, то целочисленные значения v и w сдвигаются влево на 1, после чего ведущий знаковый бит целочисленного значения v устанавливается в 1.
Зачем нужна эта операция сдвига влево и установки ведущего знакового бита?
Сдвигая влево целочисленные значения v и w, мы просто умножаем их на 2. Если у v и w были одинаковые целочисленные значения, то от этого ничего не изменится. Если w было больше v, то в результате оно немного увеличится, но результат сравнения останется прежним.
Когда мы задаём ведущий знаковый бит целочисленного значения v, оно просто увеличивается на 1. Такой инкремент позволяет учесть то дробное значение, которое входило в состав v. В рассматриваемом случае у v и w были равные целочисленные значения. Поэтому после установки ведущего знакового бита целочисленное значение v будет на 1 больше, чем у w. Но если w было больше v, то такое увеличение целочисленного значения v на 1 ничего не изменит.
Например, если v=5.75 и w=5, то целочисленное значение у обоих равно 5. Сдвинув его влево на 1, получим 10. Если установить ведущий знаковый бит для целочисленного значения v, то станет равно 11, а w останется 10. Таким образом, мы получим верный результат, так как 11 > 10. С другой стороны, если v=5.75, а w=6, то, умножив их целочисленные значения на 2, получим 10 и 12 соответственно. Прибавив 1 к целочисленному значению v, получим 11. Но w всё равно останется больше, и мы получим верный результат.
Ниже приведён листинг, описывающий эту часть кода:
А теперь разбор:
fracpart и intpart содержат дробное и целочисленное значения v. vv и ww — это целочисленные объекты Python (с бесконечной точностью), в которых представлены целочисленные значения v и w.
Они вызывают функцию
modf
, чтобы извлечь целочисленную и дробную часть v.В этой части обрабатывается ситуация, когда у v ненулевая дробная часть. Как видите, vv и ww сдвигаются влево на 1, после чего ведущий знаковый бит vv устанавливается в 1.
Теперь, чтобы узнать возвращаемое значение, остаётся всего лишь сравнить два целочисленных значения vv и ww.
Резюмируем алгоритм сравнения чисел с плавающей точкой в Python
Давайте обобщим весь алгоритм, уже не обращаясь к коду:
Если и v, и w — это числа с плавающей точкой (объекты float), то Python просто сравнивает заложенные в них значения двойной точности.
Но, если объект w — это целое число, то:
Если числа противоположны по знаку, то достаточно сравнить их знаки.
Иначе, если w — это огромное целое число (в Python точность целых чисел не ограничена), то этап сравнения конкретных чисел также можно пропустить, поскольку w больше.
Если w умещается в 48 разрядов или менее, то Python преобразует w в число двойной точности, а затем напрямую сравнивает значения двойной точности чисел v и w.
Иначе, если порядки v и w в нормализованной дробной форме не равны, то для сравнения чисел достаточно сравнить их порядки.
Если ничего из вышеперечисленного не помогло, то сравниваем сами числа. Для этого разбиваем v на целочисленную и дробную часть, а затем сравниваем целочисленную часть v с w. (При этом мы учитываем дробную часть v).
Возвращаемся к нашей задаче
Опираясь на все наши знания о стандарте IEEE-754 и о том, как в CPython сравниваются числа с плавающей точкой, вернёмся к тестовым сценариям и попробуем порассуждать об их выводе.
Сценарий 1: 9007199254740992 == 9007199254740992.0
>>> 9007199254740992 == 9007199254740992.
True
>>>
Имеем v=9007199254740992.0 и w=9007199254740992.
Мы уже выяснили мантиссу и порядок 9007199254740992.0. Порядок этого числа — 53, а мантисса — 0000000000000000000000000000000000000000000000000000.
Кроме того, число 9007199254740992 — это 2^53, то есть, для представления этого числа в памяти необходимо 54 разряда (nbits=54)
Рассмотрим, какая часть алгоритма для сравнения чисел с плавающей точкой обрабатывает этот случай:
w — это целое число, поэтому часть 1 неприменима
Знак у v и w одинаков, поэтому часть 2 неприменима
w умещается в 54 разряда, т.e., оно не слишком велико, поэтому часть 2.1 неприменима
w требуется более 48 разрядов, поэтому часть 2.2 неприменима
Порядок как для v, так и для w в нормализованной дробной форме равен 54, т.e., порядки у них равные, поэтому часть 2.3 неприменима
Так мы приходим к последней части, 2.4. Давайте её разберём.
Нам требуется извлечь из v как целочисленную, так и дробную часть, прибавить единицу к целочисленной части, если дробная часть ненулевая, а затем сравнить целочисленное значение с w.
В данном случае у v=9007199254740992.0 целочисленная часть равна 9007199254740992, а дробная часть равна 0. Целочисленное значение w также равно 9007199254740992. Поэтому Python непосредственно сравнивает два целых числа и возвращает результат True.
Сценарий 2: 9007199254740993 == 9007199254740993.0
>>> 9007199254740993 == 9007199254740993.
False
>>>
Такого результата мы не ожидали. Давайте его проанализируем.
Имеем v=9007199254740993.0 и w=9007199254740993
Оба эти значения на 1 больше, чем значения, рассмотренные в предыдущем сценарии, т.е., v = w = 2^53 + 1.
В двоичном представлении 9007199254740993 это 100000000000000000000000000000000000000000000000000001.
Напомню, что, когда мы нашли представление, соответствующее стандарту IEEE-754, выяснилось, что это число невозможно представить точно, и его представление такое же, как у 9 007 199 254 740 992.0.
Поэтому в Python, фактически, сравниваются числа 9 007 199 254 740 993 и 9 007 199 254 740 992.0. Давайте посмотрим, что же происходит внутри алгоритма CPython.
Точно, как и в предыдущем сценарии, оказываемся в последней части алгоритма (2.4), где требуется извлечь целочисленную и дробную часть v, а затем сравнить целочисленную часть со значением w.
Внутрисистемно 9 007 199 254 740 993.0 фактически представлено как 9 007 199 254 740 992.0, поэтому целочисленная часть равна 9 007 199 254 740 992. При этом значение w равно 9 007 199 254 740 993. Поэтому при сравнении двух этих значений получим False.
В языке с неявным расширением численных типов w со значением 9 007 199 254 740 993 также было бы представлено как 9 007 199 254 740 992.0, и при сравнении мы получили бы результат True.
Сценарий 3: 9007199254740994 == 9007199254740994.0
>>> 9007199254740994 == 9007199254740994.
True
>>>
В отличие от предыдущего сценария, 9007199254740994.0 можно в точности представить в соответствии с форматом IEEE-754. Вот из чего оно состоит:
Знаковый разряд: 0
Порядок: 53 + 1023 (сдвиг) = 1076 (10 000 110 100 в двоичной системе)
Мантисса: 0000000000000000000000000000000000000000000000000001
Здесь мы также добираемся до последней части алгоритма, требующей сравнить целочисленные части v и w. В данном случае их значения одинаковы, поэтому, проверив их на равенство, получим результат True.
Заключение
Стандарт IEEE-754 и арифметика чисел с плавающей точкой по определению сложны, и сравнивать числа с плавающей точкой — нетривиальная задача. В некоторых языках, например, C и Java, реализуется неявное расширение типов, что позволяет преобразовывать целые числа в числа двойной точности, а затем побитово их сравнивать. Правда, и такой подход может приводить к неожиданным результатам из‑за потери точности.
Python в данном отношении уникален. В нём используются целые числа неограниченной разрядности, поэтому во многих случаях расширять типы попросту невозможно. Поэтому в Python для сравнения таких чисел применяется собственный специализированный алгоритм, и в этом алгоритме есть свои пограничные случаи.
Независимо от того, с каким языком вы работаете — C, Java, Python или каким‑нибудь другим, при сравнении чисел с плавающей точкой рекомендуется пользоваться библиотечными функциями, а не сравнивать числа напрямую — так можно обойти потенциальные подводные камни.
Честно говоря, до того, как писать эту статью, я никогда не задумывался, как именно в Python сравниваются целые числа и числа с плавающей точкой. Теперь, разобрав эту проблему, я чувствую, что стал увереннее ориентироваться в сложной математике чисел с плавающей точкой. Надеюсь, что и вы тоже!