pep*_*ppe 2 c++ floating-point precision
这个 C/C++ 简化测试用例:
int x = ~~~;
const double s = 1.0 / x;
const double r = 1.0 - x * s;
assert(r >= 0); // fail
Run Code Online (Sandbox Code Playgroud)
数值不稳定,并且符合断言。原因是最后的计算可以用FMA来完成,这就带来了r负面影响。
Clang 默认启用了 FMA(从版本 14 开始),因此它导致了一些有趣的回归。这是一个运行版本:https://godbolt.org/z/avvnnKo5E
有趣的是,如果将最后一个计算一分为二,则不会发出 FMA 并且结果始终为非负:
int x = ~~~;
const double s = 1.0 / x;
const double tmp = x * s;
const double r = 1.0 - tmp;
assert(r >= 0); // OK
Run Code Online (Sandbox Code Playgroud)
这是 IEEE754 / FP_CONTRACT 的保证行为,还是这是在玩火,应该找到一种数值更稳定的方法?我找不到任何迹象表明 fp 收缩仅意味着“本地”发生(在一个表达式内),并且像上面这样的简单分割足以防止它们。
(当然,在适当的时候,人们也可以考虑用一种数值更稳定的算法来替换该算法。或者在 [0.0, 1.0] 范围内添加一个限制,但这感觉有点hackish。)
C++ 标准允许以额外的范围和精度计算浮点表达式,因为 C++ 2020 草案 N4849 7.1 [expr.pre] 6 表示:
\n\n\n浮点操作数的值和浮点表达式的结果可以用比类型要求更高的精度和范围来表示;类型不会因此改变。
\n
然而,注释 51 告诉我们:
\n\n\n强制转换和赋值运算符仍必须执行其特定转换,如 7.6.1.3、7.6.3、7.6.1.8\n和 7.6.19 中所述。
\n
其含义是赋值或强制转换必须将值转换为名义类型。double因此,如果使用了额外的范围或精度,则在分配给 a 时,必须将该结果转换为实际值double因此,如果使用了额外的范围或精度,则在执行(我希望,为此目的,赋值包括定义中的初始化。)
因此1.0 - x * s可以使用融合乘加,但const double tmp = x * s; const double r = 1.0 - tmp;必须计算double的结果x * s,然后double从 1.0 中减去该结果。
请注意,这并不排除const double tmp = x * s;使用额外的精度进行计算x * s,然后再次舍入以获得double结果。在极少数情况下,这可能会产生双舍入错误,其结果与将x\xe2\x80\xa2的实数算术结果s直接舍入为double. 这在实践中不太可能发生;C++ 实现没有理由x * s以额外的精度进行计算,然后将其舍入为double。
另请注意,C 和 C++ 实现不一定符合 C 或 C++ 标准。
\n| 归档时间: |
|
| 查看次数: |
161 次 |
| 最近记录: |