logo标签列表

有没有一个浮点数x，使得x-x == 0为假？

floating-pointfloating-accuracyieee-754
21

在大多数情况下，我了解到浮点数比较测试应该在一定范围内实现（abs(x-y) < epsilon），但是自减操作是否意味着结果将为零？

// can the assertion be triggered?
float x = //?;
assert( x-x == 0 )

我猜nan/inf可能是特殊情况，但我更感兴趣的是针对简单值会发生什么。

编辑：

如果有人可以引用参考资料（IEEE浮点标准），我很乐意选择一个答案。

你已经接受了这个问题，但请也阅读一下我的回答https://dev59.com/o3E85IYBdhLWcg3wnU-p#2687323。它可以澄清（我希望如此）并解决你的问题。 - Oleg
6个回答
61
正如你所暗示的那样，inf - inf是NaN ，不等于零。同样，NaN - NaN也是NaN。然而，对于任何有限浮点数x，x-x == 0.0（取决于舍入模式，x-x的结果可能为负零，但是在浮点算术中，负零与0.0相等）。 编辑：这个问题很棘手，因为这是IEEE-754标准规则中涌现出来的性质。具体而言，它遵循第5条款中定义的操作必须正确舍入的要求。减法就是这样的运算（第5.4.1节 "算术运算"），x-x的正确舍入结果是适当符号的零（第6.3节，第3段）：
 

当两个相反符号的操作数之和（或两个相似符号的差）恰好为零时，在除roundTowardNegative外的所有舍入方向属性下，该和（或差）的符号都应为+0；在该属性下，精确零和（或差）的符号应为−0。

因此，x-x的结果必须是+/- 0，因此必须与0.0相等（第5.11节，第2段）：
 

比较应忽略零的符号。

Further Edit: 这并不意味着有错误的编译器不能导致该断言触发。您的问题含糊不清; 不存在任何有限浮点数x，使得 x - x == 0 为假。然而，这不是您发布的代码所检查的内容。它检查一个C语言风格中的某个表达式是否能够求值为非零值；特别地，在某些平台上，使用某些（考虑不周的）编译器优化，该表达式中的变量x的两个实例可能具有不同的值，导致断言失败（特别是如果x是某种计算结果，而不是常数可表示的值）。这是那些平台上数字模型中的错误，但是这并不意味着它不可能发生。
太棒了，正是我在寻找的。 - Andrew Walker
10
@Potatoswatter: 作为754标准的编辑，我花了几个月的时间，这有助于我知道在哪里查找这些信息。如果没有这方面的背景知识，我就不会知道该去哪里寻找。 Translated: @Potatoswatter: 在754标准的编辑角色上经历几个月的时间是非常有帮助的。如果没有那样的背景，我不会知道在哪里寻找这些内容。 - Stephen Canon
你能否评论一下我在 https://dev59.com/o3E85IYBdhLWcg3wnU-p#2687323 上的代码示例。谢谢。 - Oleg
1
当然，C和C++都不“需要”754。也许应该重新标记问题？ - MSalters
1
此外，如果x是全局的和易失性的，编译器可能会读取它两次，并且多线程可能会干扰。 - Robert Fraser
显示剩余2条评论
4
如果表示形式被转换（例如从x86的64位内存格式转换为80位内部寄存器格式），我认为在某些情况下断言可能会触发。
5
根据这个问题的措辞，这种情况可能是不可能的。但是，“x=a+b; assert(x-(a+b)==0)” 可能会引发它。 - Mark Ransom
我认为这是一个需要关注的关键问题 - 表达式 x - x 不太可能在实际代码中使用（你为什么要这样做呢？），但是减去（或比较）变量值与可能产生该值的表达式可能会发生，并且由于编译器可能处理中间值的精度而产生意外结果。请参见 https://dev59.com/ZUzSa4cB1Zd3GeqPpLwz，这是一个可能类似于真实世界中可能发生的示例。 - Michael Burr
3
我在主要问题上的回答：“是否存在一个浮点数x，使得x-x == 0不成立？”是：至少在英特尔处理器上的浮点实现在“+”和“-”运算中没有算术下溢，因此您将无法找到一个x，使得x-x == 0不成立。对于支持IEEE 754-2008的所有处理器也是如此（请参见下面的参考文献）。
我的简短回答另一个你的问题：如果(x-y == 0)就像(x == y)一样安全，所以assert(x-x == 0)是可以的，因为在x-x或(x-y)中不会产生算术下溢。
原因如下。浮点/双精度数将以尾数和二进制指数的形式保存在内存中。在标准情况下，尾数被规范化：它是>= 0.5且<1。在<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以下的数字进行算术运算，则该数字将不会被规范化，因此您将像处理整数一样（仅使用尾数进行操作），而没有任何“舍入误差”。 备注：我个人尝试不使用“舍入误差”一词，因为算术计算机操作中没有错误。这些操作与使用相同计算机数值的+、-、*和/操作不同。它们是浮点数子集上的确定性操作，可以以(mantissa，exponent)的形式保存，并为每个定义了位数。我们可以将这种浮点数子集称为计算机浮点数。因此，经典浮点运算的结果将被投影回计算机浮点数集。这样的投影操作是确定性的，并具有许多功能，例如如果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)。这个数字确实是我们可以保存在double中的最小数字。所有小于DBL_MIN的数字的算术运算都像乘以2 ^(-1072)的整数一样进行。
在我的Windows 7电脑上，配有Intel处理器，我没有遇到下溢问题。如果有人有另一种处理器，比较我们的结果会很有趣。
有人知道如何通过-或+操作产生算术下溢吗？我的实验似乎表明这是不可能的。
编辑：我稍微修改了代码，使代码和消息更易读。 添加链接: 我的实验表明，http://grouper.ieee.org/groups/754/faq.html#underflow 在我的英特尔 Core 2 CPU 上是绝对正确的。计算的方式在 "+" 和 "-" 浮点运算中产生了无下溢。我的结果与 Microsoft Visual C 编译器的 Strict (/fp:strict) 或 Precise (/fp:precise) 开关无关 (参见 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，其中描述了我之前所写的内容。包括非规格化数或非规范化数（现在通常称为子规范数，例如在IEEE 754-2008中），遵循以下声明：

“非规格化数确保浮点数的加法和减法永远不会下溢；两个相邻的浮点数总是有一个可表示的非零差异。没有逐渐下溢，减法a−b可以下溢并产生零，即使这些值不相等。”

因此，在支持IEEE 754-2008的任何处理器上，我的所有结果必须都是正确的。

3

是的，除了特殊情况，x-x总是等于0。但是x*(1/x)并不总是等于1 ;-)

2
他不是在问特殊情况吗？ - Frank Krueger
1
@Frank - 是的，但他忽略了ypnos所提到的两种特殊情况（inf和NaN）。 - Chris Lutz
3
自减操作应该总是得到0，除了一些特殊情况。当你在比较之前进行加、减、乘或除运算时，就会出现问题，其中指数和尾数被调整。当指数相同时，尾数相减，如果它们相等，那么结果就为0。更多信息请参考：http://grouper.ieee.org/groups/754/。
1
网页内容由stack overflow 提供, 点击上面的
可以查看英文原文，
原文链接