Ste*_*sch 8 c++ floating-point precision
C++ has std::nextafter(),它返回给定浮点值f之后的下一个可表示值。在我的情况下,我想在尾数位中允许n位的斜率,因此 3 位的斜率需要在某个给定值f之后获得第 8 个下一个值。我可以打电话nextafter()八次,但是有没有更好的方法来处理这个问题?
对于大多数值,由于 IEEE 754 的布局,您可以通过将 FP 值强制转换为uint_64,添加容差(1<<3对于 3 位斜率),然后再强制转换回double。但是,这依赖于 IEEE 754 浮点数(一个很好的假设,但也不是坚如磐石)。
(对于背景,我想用它来提升光线-表面交点,由于 FP 不精确,这些交点偶尔位于表面内部。熟悉稳健浮点的人会明白为什么这epsilon是一个糟糕的解决方案。)
一种简单的方法可能是将某个值与下一个可表示浮点数之间的距离乘以 8,而不是调用 8 次std::nextafter
double advance_float(double x, int d)
{
double step = std::copysign((std::nextafter(x, x + d) - x) * d, d);
return x + step;
}
Run Code Online (Sandbox Code Playgroud)
这里有一些测试,但由您决定这是否适合您的用例。
编辑
正如史蒂夫·霍拉什( Steve Hollash)指出的那样,x可能是如此之大x + d == d。Daniel Jour建议利用frexp( 和ldexp),但在下面的尝试中,我将使用不同的方法来确定方向。
double advance_float(double x, int d)
{
const double to = std::copysign(std::numeric_limits<double>::infinity(), d);
const double next = std::nextafter(x, to);
return x + std::copysign(d * (next - x), d);
}
Run Code Online (Sandbox Code Playgroud)
请注意,它假设std::numeric_limits<double>::has_infinity == true, 否则必须使用::lowest()和。::max()
这些是一些结果
xd 上一个 x 下一个
-------------------------------------------------- ----------------------------------------------------
1 1 0x1.ffffffffffffffp-1 0x1p+0 0x1.0000000000001p+0
1 8 0x1.ffffffffffff8p-1 0x1p+0 0x1.0000000000008p+0
3.14159 8 0x1.921fb54442d1p+1 0x1.921fb54442d18p+1 0x1.921fb54442d2p+1
100.01 8 0x1.900a3d70a3d69p+6 0x1.900a3d70a3d71p+6 0x1.900a3d70a3d79p+6
-100.01 8 -0x1.900a3d70a3d79p+6 -0x1.900a3d70a3d71p+6 -0x1.900a3d70a3d69p+6
1e+67 8 0x1.7bd29d1c87a11p+222 0x1.7bd29d1c87a19p+222 0x1.7bd29d1c87a21p+222
1e-59 8 0x1.011c2eaabe7dp-196 0x1.011c2eaabe7d8p-196 0x1.011c2eaabe7ep-196
0 8 -0x0.0000000000008p-1022 0x0p+0 0x0.0000000000008p-1022
4.94066e-324 8 -0x0.0000000000007p-1022 0x0.0000000000001p-1022 0x0.0000000000009p-1022