Мой ответ на главный вопрос: «Есть ли значение x с плавающей запятой, для которого x-x == 0 равно false?» это: по крайней мере, реализация с плавающей запятой на процессорах Intel делает арифметические потери NO в операциях "+" и "-", и поэтому вы не сможете найти x, для которого x-x == 0 равно false. То же самое верно для всех процессоров, которые поддерживают IEEE 754-2008 (см. Ссылки ниже).
Мой краткий ответ на другой ваш вопрос: если (xy == 0) точно так же безопасно, как если бы (x == y), значит, assert (xx == 0) в порядке, потому что арифметическое недопущение не будет быть произведенным в хх или (ху).
Причина в следующем. Число с плавающей запятой / двойное число будет храниться в памяти в виде мантиссы и двоичного показателя. В стандартном случае мантисса нормализована: это> = 0,5 и <1. В <code><float.h> вы можете найти некоторые константы из стандарта IEEE с плавающей запятой. Интересны сейчас для нас только следующие
#define DBL_MIN 2.2250738585072014e-308 /* min positive value */
#define DBL_MIN_10_EXP (-307) /* min decimal exponent */
#define DBL_MIN_EXP (-1021) /* min binary exponent */
Но не все знают, что вы можете иметь двойные числа меньше DBL_MIN. Если вы выполняете арифметические операции с числами в DBL_MIN, это число будет НЕ нормализовано, и поэтому вы будете работать с этими числами, как с целыми числами (только с мантиссой) без каких-либо "ошибок округления" .
Замечание : Лично я стараюсь не использовать слова "ошибки округления", поскольку в арифметических компьютерных операциях нет ошибок . Эти операции не совпадают с операциями +, -, * и / с теми же номерами компьютеров, что и плавающее число. Есть детерминированных операций над подмножеством чисел с плавающей запятой, которые могут быть сохранены в форме (мантисса, показатель степени) с четко определенным числом битов для каждого. Такое подмножество чисел с плавающей точкой мы можем назвать число с плавающей точкой компьютера . Таким образом, результат классической операции с плавающей запятой будет спроецирован обратно на набор чисел с плавающей запятой компьютера. Такая проекционная операция является детерминированной и имеет множество функций, например, если x1> = x2, то x1 * y> = x2 * y.
Извините за длинное замечание и вернемся к нашей теме.
Чтобы показать, что у нас есть, если мы оперируем числами, меньшими DBL_MIN, я написал небольшую программу на C:
#include <stdio.h>
#include <float.h>
#include <math.h>
void DumpDouble(double d)
{
unsigned char *b = (unsigned char *)&d;
int i;
for (i=1; i<=sizeof(d); i++) {
printf ("%02X", b[sizeof(d)-i]);
}
printf ("\n");
}
int main()
{
double x, m, y, z;
int exp;
printf ("DBL_MAX=%.16e\n", DBL_MAX);
printf ("DBL_MAX in binary form: ");
DumpDouble(DBL_MAX);
printf ("DBL_MIN=%.16e\n", DBL_MIN);
printf ("DBL_MIN in binary form: ");
DumpDouble(DBL_MIN);
// Breaks the floating point number x into its binary significand
// (a floating point value between 0.5(included) and 1.0(excluded))
// and an integral exponent for 2
x = DBL_MIN;
m = frexp (x, &exp);
printf ("DBL_MIN has mantissa=%.16e and exponent=%d\n", m, exp);
printf ("mantissa of DBL_MIN in binary form: ");
DumpDouble(m);
// ldexp() returns the resulting floating point value from
// multiplying x (the significand) by 2
// raised to the power of exp (the exponent).
x = ldexp (0.5, DBL_MIN_EXP); // -1021
printf ("the number (x) constructed from mantissa 0.5 and exponent=DBL_MIN_EXP (%d) in binary form: ", DBL_MIN_EXP);
DumpDouble(x);
y = ldexp (0.5000000000000001, DBL_MIN_EXP);
m = frexp (y, &exp);
printf ("the number (y) constructed from mantissa 0.5000000000000001 and exponent=DBL_MIN_EXP (%d) in binary form: ", DBL_MIN_EXP);
DumpDouble(y);
printf ("mantissa of this number saved as double will be displayed by printf(%%.16e) as %.16e and exponent=%d\n", m, exp);
y = ldexp ((1 + DBL_EPSILON)/2, DBL_MIN_EXP);
m = frexp (y, &exp);
printf ("the number (y) constructed from mantissa (1+DBL_EPSILON)/2 and exponent=DBL_MIN_EXP (%d) in binary form: ", DBL_MIN_EXP);
DumpDouble(y);
printf ("mantissa of this number saved as double will be displayed by printf(%%.16e) as %.16e and exponent=%d\n", m, exp);
z = y - x;
m = frexp (z, &exp);
printf ("z=y-x in binary form: ");
DumpDouble(z);
printf ("z will be displayed by printf(%%.16e) as %.16e\n", z);
printf ("z has mantissa=%.16e and exponent=%d\n", m, exp);
if (x == y)
printf ("\"if (x == y)\" say x == y\n");
else
printf ("\"if (x == y)\" say x != y\n");
if ((x-y) == 0)
printf ("\"if ((x-y) == 0)\" say \"(x-y) == 0\"\n");
else
printf ("\"if ((x-y) == 0)\" say \"(x-y) != 0\"\n");
}
Этот код выдает следующий вывод:
DBL_MAX=1.7976931348623157e+308
DBL_MAX in binary form: 7FEFFFFFFFFFFFFF
DBL_MIN=2.2250738585072014e-308
DBL_MIN in binary form: 0010000000000000
DBL_MIN has mantissa=5.0000000000000000e-001 and exponent=-1021
mantissa of DBL_MIN in binary form: 3FE0000000000000
the number (x) constructed from mantissa 0.5 and exponent=DBL_MIN_EXP (-1021) in binary form: 0010000000000000
the number (y) constructed from mantissa 0.5000000000000001 and exponent=DBL_MIN_EXP (-1021) in binary form: 0010000000000001
mantissa of this number saved as double will be displayed by printf(%.16e) as 5.0000000000000011e-001 and exponent=-1021
the number (y) constructed from mantissa (1+DBL_EPSILON)/2 and exponent=DBL_MIN_EXP (-1021) in binary form: 0010000000000001
mantissa of this number saved as double will be displayed by printf(%.16e) as 5.0000000000000011e-001 and exponent=-1021
z=y-x in binary form: 0000000000000001
z will be displayed by printf(%.16e) as 4.9406564584124654e-324
z has mantissa=5.0000000000000000e-001 and exponent=-1073
"if (x == y)" say x != y
"if ((x-y) == 0)" say "(x-y) != 0"
Итак, мы можем видеть, что если мы будем работать с числами, меньшими DBL_MIN, они не будут нормализованы (см. 0000000000000001
). Мы работаем с этими числами как с целыми числами и без каких-либо «ошибок». Таким образом, если мы присваиваем y=x
, то if (x-y == 0)
точно так же безопасен, как if (x == y)
, и assert(x-x == 0)
работает нормально. В этом примере z = 0,5 * 2 ^ (- 1073) = 1 * 2 ^ (- 1072). Это число действительно наименьшее число, которое мы можем сохранить в два раза. Вся арифметическая операция с числами, меньшими DBL_MIN, работает как целое число, умноженное на 2 ^ (- 1072).
Итак, у меня нет проблем с подвохом на компьютере с Windows 7 с процессором Intel. Если у кого-то есть другой процессор, было бы интересно сравнить наши результаты .
У кого-нибудь есть идеи, как можно получить арифметическое занижение с помощью - или + операций? Мои эксперименты выглядят так, что это невозможно.
РЕДАКТИРОВАНИЕ : Я немного изменил код для лучшей читаемости кода и сообщений.
ДОБАВЛЕННЫЕ ССЫЛКИ : Мои эксперименты показывают, что http://grouper.ieee.org/groups/754/faq.html#underflow абсолютно корректно на моем процессоре Intel Core 2. То, как это будет вычислено, дает без потерь в операциях "+" и "-" с плавающей запятой. Мои результаты не зависят от параметров строгого (/ fp: строгого) или точного (/ fp: точного) компилятора Microsoft Visual C (см. http://msdn.microsoft.com/en-us/library/e7s85ffb%28VS.80%29.aspx и http://msdn.microsoft.com/en-us/library/Aa289157)
Еще одна (вероятно, последняя) ссылка и мое последнее замечание : Я нашел хорошую ссылку http://en.wikipedia.org/wiki/Subnormal_numbers,, где описано то же самое, что я написал ранее. Включение денормальных чисел или денормализованных чисел (в настоящее время часто называемых субнормальными числами, например, в In IEEE 754-2008 ) следует следующему утверждению:
«Денормальные числа дают
гарантировать, что дополнение и
вычитание чисел с плавающей точкой
никогда не теряется; два рядом
Числа с плавающей точкой всегда имеют
представимая ненулевая разница.
Без постепенного снижения
вычитание a − b может опуститься и
производить ноль, хотя значения
не равны ».
Так что все мои результаты должны быть правильными на любом процессоре, который поддерживает IEEE 754-2008.