使用 sqrt 和 floor 时的近似误差

Question

我必须计算一个方程的解，我知道 y < x *( sqrt(n) - 1 )，其中 x，y n 是整数。

我天真的方法是寻找 y 小于或等于 floor( x * ( sqrt( (float)n ) - 1 ) )。

我应该担心近似误差吗？
例如，如果我的表达式比整数m稍微大一点，我是否应该担心得到m-1 最后?
如何检测此类错误？

Answer 1

您绝对应该担心近似误差，但担心程度取决于您担心的 x 和 n 的取值范围关于。

在 IEEE 4 字节浮点表示中的计算将在 2^23 到 2^24 中大致有一部分的错误；对于 8 字节的表示（即 double），它大约是 2^52 到 2^53 的一部分。您可能会期望您需要使用 doubles 而不是 floats 来获得 32 位整数 x 和 的准确结果]n，即使是 double 也不足以用于 64 位整数。

例如，考虑代码：

template <typename F,typename V>
F approxub(V x,V n) {
    return std::floor(x*std::sqrt(F(n))-x);
}

uint64_t n=1000000002000000000ull; // (10^9 + 1)^2 - 1
uint64_t x=3;
uint64_t y=approxub<double>(x,n);

这给出了 y=3000000000 的值，但正确的值是 2999999999。

当 x 很大而 n 很小时更糟：64 位大整数在 IEEE doubles:

uint64_t n=9;
uint64_t x=5000000000000001111; // 5e18 + 1111
uint64_t y=approxlb<double>(x,n);

y 的正确值（抛开 n 何时是完全平方的问题——在这种情况下，真正的上限将少一）是 2 x = 10000000000000002222，即 1e19 + 2222。然而，计算出的 y 是 10000000000000004096.

避免浮点近似

假设您有一个函数 isqrt，它精确地计算了一个整数的平方根的整数部分。那么你可以说

y = isqrt(x*x*n) - x

并且假设乘积 x*x*n 适合您的整数类型，您将有一个精确的上限（如果 n 是完美的正方形。）编写 isqrt 函数的方法不止一种；这是基于 material at code codex:

的示例实现

template <typename V>
V isqrt(V v) {
    if (v<0) return 0;

    typedef typename std::make_unsigned<V>::type U;
    U u=v,r=0;

    constexpr int ubits=std::numeric_limits<U>::digits;
    U place=U(1)<<(2*((ubits-1)/2));

    while (place>u) place/=4;
    while (place) {
        if (u>=r+place) {
            u-=r+place;
            r+=2*place;
        }
        r/=2;
        place/=4;
    }
    return (V)r;
}

如果 x 太大怎么办？例如，如果我们最大的整数类型有 64 位，并且 x 大于 2^32。最直接的解决方案是进行二进制搜索，以 x r - x 和 x r 为边界，其中 r = [√n]是整数平方根.

使用 sqrt 和 floor 时的近似误差

Approximation error when using sqrt and floor

c++

numerical

approximation

避免浮点近似