给定 x>0，sqrt(x) = 0 是否可能给出浮点错误？

Question

标题中的问题。

我有一段代码：

double ccss = c * c + s * s;
double sqrtCCSS = sqrt(ccss);
if (sqrtCCSS != 0)
{
    n = n1 / sqrtCCSS;
}

我只是想知道这是否安全：

double ccss = c * c + s * s;
if (ccss != 0)
{
    n = n1 / sqrt(ccss);
}

我的直觉告诉我是的，但浮点错误仍然有些神秘。

更新：

至少对于 Python 来说，到目前为止这似乎是不可能的：

double ccss = c * c + s * s;
double sqrtCCSS = sqrt(ccss);
if (sqrtCCSS != 0)
{
    n = n1 / sqrtCCSS;
}

如果 numpy 中可能的最小浮点数（假设 numpy 主要在 C 中，这似乎相关）不能导致零，则导致零的唯一可能数字必须是引入严重浮点错误的数字，但不是最小的数字浮点数。

更新：我更改了变量以满足评论

Answer 1

我只是想知道这是否安全：

double xxyy = x * x + y * y;
if (xxyy != 0)
{
    n = n1 / sqrt(xxyy);
}

它总是安全的，因为浮点数学不会陷入困境，并且您可以自由除以零甚至 Inf/NaN，除非您告诉编译器在这些情况下陷入困境

无论如何，如果您只是想检查分母是否为零，那么标题上的问题实际上与您在问题内容中讨论的内容不同。问题的答案

给定 x>0，有可能sqrt(x) 为零吗？

是，如果非正规数为零(DAZ) 和/或刷新为零(FTZ) 已打开。事实上，大多数非 x86 架构默认启用 DAZ/FTZ，或者甚至不支持非正规数并始终打开 DAZ/FTZ，因为对非正规数的操作非常慢（请参阅为什么将0.1f 更改为 0 会降低性能） 10 倍？）。在 x86 中，您还可以打开这些标志，并且某些编译器会默认打开它们

• 浮点运算中的清零行为
• 在代码级别禁用非规范浮动
• 设置 FTZ 和 DAZ 标志

例如下面的示例代码

f = 0x0Cp-1022;
_MM_SET_DENORMALS_ZERO_MODE(_MM_DENORMALS_ZERO_ON);
_MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON);
printf("%a %a\n", f, sqrt(f));

可能会打印，0x0.000000000000cp-1022 0x0p+0因为由于 DAZf被视为零sqrt

Godbolt 上的一些演示。以及ideone 上的 gcc 演示

上面的代码片段也可以打印，0x0p+0 0x0p+0因为由于 FTZ 的原因f，它被视为零printf，即使它包含非零非规范化值

如果您在分配之前设置了DAZ/FTZ，f如下所示

_MM_SET_DENORMALS_ZERO_MODE(_MM_DENORMALS_ZERO_ON);
_MM_SET_FLUSH_ZERO_MODE(_MM_FLUSH_ZERO_ON);
f = 0x0Cp-1022;
printf("%a %a\n", f, sqrt(f));

那么它也可能 print 0x0p+0 0x0p+0，但在这种情况下f可以被分配一个实际的零值，并且打印其位模式将导致所有零位

但是如果你问过这个

给定 x 和 y 不同时为零，是否有可能sqrt(x*x + y*y)为零

那么这个问题的答案是肯定的，即使 DAZ/FTZ 未启用，因为数字可能太小，以至于它们的乘积会下溢。尝试这样的事情，看看

double x = 0x1p-1021;
double y = 0x2p-1020; // or y = 0
double xxyy = x * x + y * y;
printf("%a %a\n", x, y);
printf("%a %a\n", xxyy, sqrt(xxyy));