Крайние случаи чисел с плавающей точкой

Предпосылки: Числа с плавающей точкой.

Числа с плавающей точкой | Кодирование текста

IEEE 754 даёт связное объяснение для «обычных» чисел: три поля, смещение, скрытый бит, округление до ближайшего представимого значения. Но обычная арифметика редко остаётся в рамках обычных чисел. Массив физических измерений содержит ноль, деление может упереться в него, корень — потребовать отрицательный аргумент, цикл — просуммировать миллион значений, а веб-API — передать число через JSON. На каждом из этих поворотов формат либо отвечает отдельным правилом, либо подводит тихой ошибкой.

Округление отвечает, что сохраняется при каждой операции. Специальные значения (NaN, ±∞, ±0) появляются на границах диапазона и в неопределённых ситуациях. Денормализованные числа описывают поведение у самого нуля. Катастрофическая потеря значимости и накопление ошибки показывают, почему невинные вычисления ломаются на длинной дистанции. Из всего этого складываются практические правила, позволяющие не попадаться в ловушки формата.

Округление

Когда результат вычисления не укладывается точно в мантиссу, его нужно округлить. IEEE 754 определяет пять режимов, по умолчанию используется round to nearest, ties to even — округление к ближайшему, при равном расстоянии к чётному. Если результат ровно посередине между двумя представимыми числами, выбирается то, у которого младший бит мантиссы = 0.

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

Помимо выбора нейтрального режима, стандарт фиксирует и саму точность. IEEE 754 гарантирует: результат любой базовой операции (сложение, вычитание, умножение, деление, квадратный корень) — это точный математический результат, округлённый по выбранному режиму. То есть операция выполняется так, как если бы у процессора была бесконечная точность, и только потом результат округляется до мантиссы. Ошибка каждой отдельной операции — не больше половины шага (0.5 ulp, unit in the last place — единица в последнем разряде мантиссы).

Остальные режимы (toward zero, toward +∞, toward −∞, toward nearest ties away from zero) нужны для специальных задач — например, интервальной арифметики, где важна гарантия «результат не превысил верхнюю границу». В обычном коде их переключают редко.

Специальные значения: ±0, ±∞, NaN

Округление решает проблему конечной точности для «обычных» чисел. Но что делать, когда результат выходит за пределы диапазона — например, при делении на ноль? Или когда операция математически не определена — корень из отрицательного числа? Для таких случаев стандарт резервирует специальные комбинации битов.

Два значения экспоненты зарезервированы: все нули (00000000) и все единицы (11111111 для float32). Они кодируют случаи, которые не вписываются в обычную схему.

Ноль. Экспонента = 0, мантисса = 0. Формально 1.0 × 2⁻¹²⁷ не равно нулю, поэтому ноль закодирован как исключение. Бит знака при этом работает: существуют +0 (все биты нули) и −0 (знаковый бит = 1, остальное нули). По стандарту +0 = −0 при сравнении, но знак сохраняется — он имеет значение при делении: 1.0 / (+0) = +∞, 1.0 / (−0) = −∞. Знаковый ноль полезен в геометрии и комплексных функциях, где направление подхода к нулю определяет знак результата.

Бесконечность. Экспонента = все единицы, мантисса = 0. Результат операций, выходящих за диапазон: 1.0 / 0.0 = +∞, log(0) = −∞, 10³⁰⁰ × 10³⁰⁰ = +∞. Бесконечность участвует в дальнейших вычислениях по арифметическим правилам: ∞ + 1 = ∞, ∞ × (−1) = −∞, 1 / ∞ = 0. Это не математическая бесконечность, а сигнальное значение, которое можно распознать (isinf(x)) и обработать.

NaN (Not a Number — не число). Экспонента = все единицы, мантисса отлична от нуля. Результат неопределённых операций: 0.0 / 0.0, √(−1), ∞ − ∞, 0 × ∞. Зачем отдельное значение, а не немедленная ошибка? Чтобы длинная цепочка вычислений не прерывалась на середине. NaN «пропитывает» результат: любая арифметическая операция с NaN даёт NaN. Программа доходит до конца, и программист видит NaN в итоге — сигнал, что где-то раньше произошла неопределённая операция.

NaN обладает уникальным свойством: NaN не равен ничему, включая самого себя. Выражение NaN == NaN возвращает ложь. Это единственное значение, для которого x != x — истина, и этим свойством пользуются для проверки: если x != x, значит x — это NaN. Почему такой странный выбор? Потому что NaN — это «ответа нет». Если спросить «равен ли ответ-которого-нет самому себе», единственный непротиворечивый ответ — нет.

Специальные значения float32:
 
Экспонента   Мантисса    Значение
00000000     = 0         ±0 (знак из бита S)
00000000     != 0        денормализованные числа (см. ниже)
11111111     = 0         ±∞
11111111     != 0        NaN
01..11111110 любая       обычные (нормализованные) числа

Денормализованные числа: плавный переход к нулю

Наименьшее нормализованное число в float32 — 1.0 × 2⁻¹²⁶ ≈ 1.18 × 10⁻³⁸. Следующее за ним в сторону нуля — это ноль? Тогда между нулём и этим числом была бы пустота: значение вроде 0.5 × 10⁻³⁸ невозможно представить. Хуже того, разность двух близких маленьких чисел может дать ноль, хотя числа не равны — вычитание потеряет информацию.

Денормализованные числа (denormalized, subnormals) заполняют этот промежуток. У них экспонента = 0 и мантисса отлична от нуля. Отличие от нормализованных: скрытый бит не подставляется. Вместо 1.мантисса число интерпретируется как 0.мантисса × 2⁻¹²⁶.

Нормализованные:     1.мантисса × 2^(экспонента − 127)
                     скрытый бит = 1
 
Денормализованные:   0.мантисса × 2^(−126)
                     скрытый бит = 0, экспонента фиксирована

Наименьшее денормализованное float32 — 2⁻¹⁴⁹ ≈ 1.4 × 10⁻⁴⁵. Числа плавно уменьшаются от наименьшего нормализованного к нулю, теряя по одному биту точности на каждом шаге. Это называется постепенное исчезновение (gradual underflow): вместо резкого прыжка «число → ноль» точность деградирует плавно, и свойство «если a ≠ b, то a − b ≠ 0» сохраняется.

Без денормализованных чисел шкала выглядит так:

... smallest_normal ... 0    <-- пустота между 0 и smallest_normal

С денормализованными:

... smallest_normal ... denorm ... denorm ... denorm ... 0
                        менее точные, но заполняют пробел

Цена платится временем. Многие процессоры обрабатывают денормализованные числа отдельным, медленным путём — операция над такими операндами может быть в десятки раз медленнее обычной. В аудио-DSP и других горячих циклах это проявляется как внезапная просадка производительности, когда затухающий сигнал уходит в область денормализованных амплитуд. Стандартный обходной путь — флаг FTZ (flush to zero) или DAZ (denormals are zero): любой денормализованный результат или операнд обнуляется. Это быстрее, но ломает свойство плавного исчезновения, и применять его можно только там, где потеря крошечных значений безопасна.

Катастрофическая потеря значимости

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

Пример в десятичной системе с 7 значащими цифрами:

a = 1.000001?????   (цифры после 7-й неизвестны)
b = 1.000000?????
 
a − b = 0.000001????

Семь значащих цифр превратились в одну. Результат почти полностью состоит из хвоста ошибок округления. В двоичной системе — то же самое, но с 24 или 53 значащими битами.

Это называется catastrophic cancellation — катастрофическая потеря значимости. Название точное: одна операция вычитания способна уничтожить все значащие цифры, которые с трудом набирались в предыдущих шагах.

Классический пример — формула корней квадратного уравнения:

x = (−b ± √(b² − 4ac)) / 2a

Когда сильно больше 4ac, подкоренное выражение ≈ , и √(b² − 4ac)|b|. Один из двух корней получается как разность двух почти равных чисел — и теряет точность. Численный трюк: вычислять сначала тот корень, где знаки складываются (−b − √... при b > 0), а второй восстанавливать через теорему Виета: x₁ × x₂ = c/a. Вторая формула даёт ту же математическую величину без опасного вычитания.

Ассоциативность сложения тоже нарушается из-за похожего механизма. Если a = 10²⁰, b = −10²⁰, c = 1.0, то (a + b) + c = 0 + 1 = 1.0, но a + (c + b) может дать не 1.0 — при сложении a + c число 1.0 теряется на фоне 10²⁰ (оно меньше одного шага мантиссы при таком масштабе). a + b + c ≠ c + b + a для float — порядок операций имеет значение.

Накопление ошибок и Kahan summation

Одиночная ошибка вычитания — частный случай, его ловят переписыванием формулы. Но даже без катастрофического вычитания ошибки могут копиться просто от количества операций. Каждое округление добавляет не более половины ulp, и на одной операции это ничтожно. На миллионе операций — уже заметно.

Наивная сумма миллиона чисел порядка 1.0 в float64 может отклониться от правильного ответа на несколько единиц. Причина: когда накопленная сумма выросла до, скажем, 10⁶, прибавление очередного 1.0 округляется, теряя часть значения в хвосте. Чем длиннее сумма, тем больше накопленная погрешность — в худшем случае линейно от количества слагаемых, O(N).

Алгоритм Кэхана (Kahan summation) превращает линейное накопление в константное. Идея: завести отдельную переменную c для «потерянных на округлении» битов и компенсировать их на каждом шаге.

sum = 0.0
c = 0.0  # компенсация потерь
for x in values:
    y = x - c          # x, исправленный на прошлую потерю
    t = sum + y        # новая сумма (теряет часть y)
    c = (t - sum) - y  # восстанавливаем, что именно потерялось
    sum = t
return sum

Три лишние операции на шаг — плата за то, что суммарная ошибка становится O(1) независимо от длины суммы. Научные библиотеки редко используют именно Кэхана, но часто применяют родственные схемы с пониженной ошибкой. NumPy, например, для np.sum использует partial pairwise summation — при axis=None гарантированно, при сведении по конкретной оси зависит от оси (см. документацию). Это одна из причин, почему np.sum на большом массиве обычно точнее ручного Python-цикла с +=.

Практические рекомендации

Не сравнивайте float через ==. Две математически одинаковые формулы дают разные округления. Правильный шаблон — относительная погрешность: |a − b| ≤ ε × max(|a|, |b|). Абсолютный ε работает только если масштаб чисел заранее известен и постоянен.

Проверяйте NaN через x != x или штатную функцию. x == NaN всегда false по определению NaN. В Python — math.isnan(x), в JavaScript — Number.isNaN(x), в C99 — isnan(x).

Деньги храните в целых. Сумма 0.10 + 0.20 в float даёт 0.30000000000000004. Умноженная на миллион транзакций, эта прибавка превращается в расхождение баланса. Стандартные решения: целочисленные копейки/центы; тип DECIMAL/NUMERIC в базе данных (хранит число в десятичной записи с фиксированной шкалой, арифметика над ним точна для всех дробей с конечной десятичной записью); специализированные money-библиотеки в языке (BigDecimal в Java, decimal в Python).

JSON и целочисленные ID опасны после 2⁵³. Спецификация JSON не указывает тип чисел, но большинство парсеров (JavaScript в первую очередь) читают их в float64. Мантисса float64 — 53 бита, и все целые от 0 до 2⁵³ = 9 007 199 254 740 992 представляются точно. За этой границей шаг между соседними представимыми числами становится больше 1, и, например, 2⁵³ + 1 округляется обратно до 2⁵³. Идентификаторы баз данных, Twitter Snowflake, Discord snowflakes (все 64-битные) в JSON безопасны только если передаются строками. Twitter столкнулся с этим в 2011, когда клиенты начали терять младшие биты ID твитов — с тех пор API отдаёт id_str рядом с id.

Осторожно с катастрофической потерей. Если формула содержит разность близких чисел (особенно с большими компонентами), проверьте альтернативную запись: теорема Виета, формулы половинного угла, разложение в ряд около нуля, log1p/expm1 вместо log(1+x)/exp(x)−1. Большинство стандартных библиотек содержат «численно устойчивые» варианты ровно для таких случаев.

Для больших сумм — Kahan или pairwise. В научных вычислениях и финансовой агрегации не суммируйте миллион значений простым циклом +=. Пользуйтесь библиотечным sum (NumPy, Julia и аналоги), который реализует устойчивую схему — обычно pairwise summation, опционально Kahan.


Числовая арифметика разобрана: целые, дробные и крайние случаи IEEE 754 укладываются в байты. Но экран показывает буквы, файлы содержат имена, сеть передаёт сообщения — а договор между байтом и символом пока не обсуждался. Какой именно договор превращает последовательность битов в A или é, и почему этот договор сам по себе — источник множества тонких багов?

Sources


Числа с плавающей точкой | Кодирование текста