c++ - 是否只允许在一个表达式中使用 float 缩写？

Question

这个 C/C++ 简化测试用例:

int x = ~~~;
const double s = 1.0 / x;
const double r = 1.0 - x * s;

assert(r >= 0); // fail

在数值上不稳定，并命中断言。原因是最后的计算可以用FMA完成，将r带入负数。

Clang 默认启用 FMA ( since version 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

这是 IEEE754/FP_CONTRACT 的保证行为，还是这是在玩火，人们应该找到一种数值更稳定的方法？我找不到任何迹象表明 fp 收缩仅意味着“局部”发生(在 one 表达式内)，并且像上面这样的简单拆分足以防止它们。

(当然，在适当的时候，也可以考虑用数值更稳定的算法替换算法。或者在 [0.0, 1.0] 范围内添加一个 clamp，但这感觉很老套。)

Answer 1

C++ 标准允许以额外的范围和精度计算浮点表达式，因为 C++ 2020 草案 N4849 7.1 [expr.pre] 6 说:

The values of the floating-point operands and the results of floating-point expressions may be represented in greater precision and range than that required by the type; the types are not changed thereby.

然而，注释 51 告诉我们:

The cast and assignment operators must still perform their specific conversions as described in 7.6.1.3, 7.6.3, 7.6.1.8and 7.6.19.

这句话的意思是赋值或强制转换必须将值转换为标称类型。因此，如果使用了额外的范围或精度，则在执行对 double 的赋值时，必须将该结果转换为实际的 double 值。 (我希望，为此目的，赋值包括定义中的初始化。)

所以 1.0 - x * s 可以使用融合乘加，但是 const double tmp = x * s; const double r = 1.0 - tmp; 必须计算 x * s 的 double 结果，然后从 1.0 中减去该 double 结果.

请注意，这并不排除 const double tmp = x * s; 使用额外的精度来计算 x * s 然后再次舍入以获得 双 结果。在极少数情况下，这会产生双舍入错误，结果与通过舍入 x•s 的实数算术结果略有不同直接到一个double。这在实践中不太可能发生； C++ 实现没有理由以更高的精度计算 x * s，然后将其四舍五入为 double。

另请注意，C 和 C++ 实现不一定符合 C 或 C++ 标准。