Числа с плавающей точкой

Объем памяти, который отводится для хранения значений переменных, ограничен. На современных компьютерах целочисленная переменная может занимать 64 бита (8 байт), что достаточно для представления всех целых чисел до 10¹⁹. Этот диапазон значений достаточен для подавляющего большинства практических задач. С действительными числам ситуация принципиально другая: любой интервал содержит бесконечно много чисел, поэтому некоторые действительные (и рациональные) числа не могут быть представлены при конечном размере переменных. Это приводит к неизбежным ошибкам округления при вычислениях со значениями с плавующей точкой. Далее мы рассмотрим формат хранения чисел с плавающей точкой (стандарт IEEE-754), способы сравнения чисел на равенство и некоторые особенности вычислений с плавающей точкой.

Представление чисел

Нормализованным экспоненциальным представлением числа в системе счисления по основанию b является его запись в виде d₀.d₁d₂...d_k * b^e, где d_i суть цифры, причем d₀ > 0. Десятичное число 123.456 в нормализованном виде представляется как 1.23456 * 10².

Числом с плавающей точкой назвают действительное число, которое может быть представлено в указанном виде.

Число с плавающей точкой определяется тремя целыми числами: знаком (S), экспонентой (E) и мантиссой (M).

+---+-----------+---------------------+
| S |    exp    |     mantissa        |
+---+-----------+---------------------+
1бит    w бит           p бит  

Значение экспоненты может быть как положительным, так и отрицательным. Существует несколько способов кодирования знака: выделение отдельного бита, дополненное представление и смещение. При представлении чисел с плавающей точкой экспонента записывается со смещением. Для чисел одинарной точности экспонента содержит 8 бит и смещение равно 127. Это означает, что если последовательность бит в поле exp образует число k (если эта последовательность рассматривается как двоичное представлеие беззнакового целого числа), то значением экспоненты E является k-127.

Целые числа S, E и M опеределяют десятичное число (-1)^S * b^E * m.

Ошибки округления, распределение чисел

Сравнение значений на равенство

Ошибки округления могут приводить к ситуациям, когда (a + b) - a - b оказывается отличным от нуля. Например, при a=1.0, b=0.1 значение оказывается равным 2.2351742 10^-8 при вычислениях с однократной точностью, и 8.326672684688674 10^-17 при вычислении с двойной точностью. Поэтому сравнение числ с плавающей точкой на точное равенство с помощью оператора == не надежно.

Предположим, что мы хотим определить функцию int almostEqual(double a, double b); которая возвращает ненулевое значение, если числа a и b можно считать равными. Рассмотрим возможные варианты реализации такой функции.

Сравнение модуля разности с константой. Самым простым методом сравнения является следующий:

Плохо

#define EPSILON 1e-7
int almostEqual(double a, double b) {
    return fabs(a-b) < EPSILON;
}

Равными считаются числа, разность которых по модулю не превосходит заданной константы. В данном примере 10^-7.

Данный метод работает, если значения a и b близки к 1. Если абсолютные значения велики, то a и b могут быть соседними представимыми числами, но разность меджду ними будет больше, чем EPSILON.

Сравнение модуля разности с динамическим значением EPSILON. Второе решение состоит в динамическом изменении точности.

Плохо

#define EPSILON 1e-7
int almostEqual(double a, double b) {
    double largest = (fabs(a) > fabs(b) ? fabs(a) : fabs(b));
    return fabs(a-b) < (largest * EPSILON);
}

В данном случае при сравнении больших чисел значение "константы" будет увеличено. Этот метод лучше, но что будет, если a и b существенно меньше 1?

Сравнение модуля разности с динамическим значением EPSILON. Второе решение состоит в динамическом изменении точности.

Лучше

#define EPSILON 1e-7
#define MAX(a, b) (fabs(a) > fabs(b) ? fabs(a) : fabs(b))
#define MAX3(a, b, c) (MAX((a), MAX((b), (c))))
int almostEqual(double a, double b) {
    return fabs(a-b) < MAX3(a, b, 1.0) * EPSILON;
}

Этот метод достаточно надежен: при сравнении больших чисел использцется пропорционально увеличенное значение EPSILON, а при сравнении маленьких чисе – обычное значение.

Бесконечная арифметика

Стандарт IEEE-754 определяет специальные значения чисел с плавающей точкой: +Inf (положительная бесконечность), -Inf (отрицательная бесконечность) и NaN ("не число", Not a Number). Эти значения могут появиться в процессе выполнения арифметических операций. Например:

• a/0 = +Inf, если a > 0;
• a/0 = -Inf, если a < 0;
• 0/0 = Inf - Inf = Inf*0 = NaN.

Если переменные могут принимать специальные значения, то следующие два оператора не эквивалентны!

if (a > b) do_something(); else do_something_else(); 

if (b <= a) do_something_else(); else do_something(); 

Для проверки того, что значением переменной является NaN или бесконечность, можно использовать функции isnan и isinf соответсвенно, которые объявлены в заголовочном файле <math.h>.