浮点运算中如何处理过度精度？

Question

浮点运算中如何处理过度精度？

10

在我的数值模拟中，我有类似以下代码片段的代码

double x;
do {
  x = /* some computation */;
} while (x <= 0.0);
/* some algorithm that requires x to be (precisely) larger than 0 */

在某些编译器（例如gcc）和平台（例如linux、x87数学）上，可能会以高于双精度的精度（“超额精度”）计算x。（更新：这里所说的精度指的是精度和/或范围）。在这种情况下，即使下一次将x四舍五入为双精度后变成0，比较（x <= 0）仍有可能返回false。（而且不能保证x不会在任意时刻被四舍五入。）

有没有一种方法可以执行此比较：

是可移植的，
适用于得到内联的代码，
没有性能影响，并且
不排除一些任意范围（0，eps）？

我尝试使用（x < std::numeric_limits<double>::denorm_min()），但当使用SSE2数学时，这似乎会显着减慢循环。（我知道非规格化数可以减慢计算，但我没想到它们在移动和比较时会更慢。）

更新：一种替代方案是使用volatile强制将x写入内存，然后进行比较，例如：

`

} while (*((volatile double*)&x) <= 0.0);

然而，根据应用程序和编译器所应用的优化方式，这种解决方案可能会引入明显的开销。

更新： 任何公差问题的问题在于它非常主观，即它取决于特定的应用程序或上下文。我更倾向于在没有过度精度的情况下进行比较，这样我就不必做出任何额外的假设或将一些任意的 epsilon 引入到我的库函数的文档中。

- Stephan

有趣的问题，我第一次听说有人抱怨过度精确。 - dsimcha

当您在比较中添加显式转换，即 ((double)x) <= 0.0，会发生什么？ - Christoph

也许相关的链接：https://dev59.com/VnRC5IYBdhLWcg3wW_tk - Christoph

我认为C99标准规定将类型转换为其本身应该删除多余的精度。但是这并没有得到很好的实现... - user3458

但您最终将检查每个编译器支持的功能。 - Ismael

显示剩余3条评论

5个回答

2

在你的问题中，你提到使用volatile可以解决问题，但是会有很大的性能损失。那么在比较时只使用volatile变量，允许x保存在寄存器中，如何呢？

double x; /* might have excess precision */
volatile double x_dbl; /* guaranteed to be double precision */
do {
  x = /* some computation */;
  x_dbl = x;
} while (x_dbl <= 0.0);

您还应该检查是否可以通过显式使用long double并缓存最小的子规范值来加速比较，例如：

const long double dbl_denorm_min = static_cast<long double>(std::numeric_limits<double>::denorm_min());

然后进行比较。

x < dbl_denorm_min

我认为一个不错的编译器应该会自动完成这个操作，但谁也不能确定...

- Christoph

这有点取决于编译器和具体的循环。它似乎比与非规范化值进行比较要快，但在没有过度精度的平台上仍会引入明显的开销。 - Stephan

然而，在可能存在过度精度的平台上，我可能会使用一些宏/预处理器定义，该定义会扩展为 while (*((volatile double*)&x) <= 0.0);。 - Stephan

非规格化数的比较只会减缓SSE2数学运算速度，而不会影响x87数学运算，因为x87数学运算可以使用适当的长双精度浮点数。 - Stephan

@Stephan：我明白了。那么，我已经没有更多的想法了——除了使用特定于编译器的标志或降到汇编级别之外，我认为你无法做更多的事情了... - Christoph

1

我想知道你是否有正确的停止准则。听起来，x <= 0 是一个异常条件，而不是一个终止条件，而且终止条件更容易满足。也许应该在while循环内部设置一个break语句，当达到某个公差时停止迭代。例如，许多算法在两个连续迭代彼此足够接近时终止。

- John D. Cook

不，这种情况确实是一个终止条件（或者是终止条件的逻辑否定）。当然，我可以在稍后的时候检查是否发生了除以零的情况，但在这种情况下（以及许多类似情况），我真的需要我的函数不返回任何小于或等于0的值。 - Stephan

0

好的，GCC有一个标志，-fexcess-precision，会引起您讨论的问题。它还有一个标志，-ffloat-store，可以解决您讨论的问题。

“不要将浮点变量存储在寄存器中。这可以防止在诸如68000这样的机器上出现不必要的过度精度，因为浮点寄存器（68881）保留比double类型应该具有的更多的精度。”

我怀疑这种解决方案是否会对性能产生影响，但影响可能不是非常昂贵的。随意搜索建议它的成本约为20％。实际上，我认为没有一种既便携又没有性能影响的解决方案，因为强制芯片不使用超额精度通常会涉及一些非自由操作。但是，这可能是您想要的解决方案。

- Brian

在支持SSE2的基于x86架构的平台上，使用“-mfpmath = sse”应该可以解决问题，而且没有任何负面性能影响。然而，我认为指定非标准编译器选项并不是很可移植。 - Stephan

0

一定要将检查变为绝对值。它需要在零周围有一个上下限的 epsilon。

- duffymo

我不明白，使用(x < eps)作为do while循环的条件会拒绝负值，对吗？ - Stephan

问题不在于通常的浮点精度问题，而是原帖作者在循环条件中得到了“(x >= 0)”返回false，当“x == 0”在循环外为true。 - David Thornley

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- Christoph · Accepted Answer

如Arkadiy在评论中所述，显式转换((double)x) <= 0.0 应该可以使用-至少根据标准是这样的。

C99:TC3, 5.2.4.2.2 § 8：

除了赋值和强制类型转换(它们会删除所有额外的范围和精度)之外，带有浮点操作数和值的操作的值以及通常算术转换和浮点常量的值将计算为其范围和精度可能大于类型所需的格式。[...]

如果你在x86上使用GCC，你可以使用标志-mpc32、-mpc64和-mpc80来设置浮点操作的精度为单精度、双精度和扩展双精度。