优化 sqrt(n) - sqrt(n-1)

Question

这是我每秒调用多次的函数：

static inline double calculate_scale(double n) { //n may be int or double
    return sqrt(n) - sqrt(n-1);
}

循环调用如下：

for(double i = 0; i < x; i++) {
    double scale = calculate_scale(i);
    ...
}

而且速度太慢了。优化此功能以获得尽可能准确的输出的最佳方法是什么？

参数n：从1开始，实际不限，主要用于1-10范围内的小数。它是整数（整数），但它可能既是 int 也可能是 double，具体取决于哪个表现更好。

Answer 1

你可以尝试用下面的近似值代替

sqrt(n) - sqrt(n-1) == 
(sqrt(n) - sqrt(n-1)) * (sqrt(n) + sqrt(n-1)) / (sqrt(n) + sqrt(n-1)) ==
(n - (n + 1)) / (sqrt(n) + sqrt(n-1)) ==
1 / (sqrt(n) + sqrt(n-1))

对于足够大的 n，最后一个方程非常接近 1 / (2 * sqrt(n))。所以你只需要调用 sqrt 一次。还值得注意的是，即使没有近似值，就较大的 n.

的相对误差而言，最后一个表达式在数值上更稳定

Answer 2

你说 n 主要是一个小于 10 的数字。你可以对小于 10 的数字使用预先计算的 table，甚至更多，因为它很便宜，然后回退到实际计算如果数字更大。

代码看起来像这样：

static inline double calculate_scale(double n) { //n may be int or double
    if (n <= 10.0 && n == floor(n)) {
        return precomputed[(int) n]
    }
    return sqrt(n) - sqrt(n-1);
}

Answer 3

首先感谢大家的建议。我做了一些研究并发现了一些有趣的实现和事实。

1。在循环中或使用预计算 table

（感谢@Ulysse BN）您可以通过简单地保存以前的 sqrt(n) 值来优化循环。以下示例演示了此优化用于设置预计算 table.

    /**
     * Init variables
     *      i       counter
     *      x       number of cycles (size of table)
     *      sqrtI1  previous square root = sqrt(i-1)
     *      ptr     Pointer for next value
     */
    double i, x = sizeof(precomputed_table) / sizeof(double);
    double sqrtI1 = 0;

    double* ptr = (double*) precomputed_table;

    /**
     * Optimized calculation
     * In short:
     *      scale = sqrt(i) - sqrt(i-1)
     */
    for(i = 1; i <= x; i++) {
        double sqrtI = sqrt(i);
        double scale = sqrtI - sqrtI1; 
        *ptr++ = scale;
        sqrtI1 = sqrtI;
    }

使用预先计算的 table 是可能是最快的方法，但它的缺点可能是它的大小有限。

static inline double calculate_scale(int n) {
    return precomputed_table[n-1];
}

2。使用平方根反比对 BIG 数进行近似

所需的反（倒数）平方根函数rsqrt

此方法对大数字的结果最准确。小数字有错误：

1    2     3      10       100     1000
0.29 0.006 0.0016 0.000056 1.58e-7 4.95e-10

上面是我用来计算结果的JS代码：

function sqrt(x) { return Math.sqrt(x); } function d(x) { return (sqrt(x)-sqrt(x-1))-(0.5/sqrt(x-0.5));} console.log(d(1), d(2), d(3), d(10), d(100), d(1000));

您还可以在单图中看到与 two-sqrt 版本相比的准确性：https://www.google.com/search?q=(sqrt(x)-sqrt(x-1))-(0.5%2Fsqrt(x-0.5))

用法：

static inline double calculate_scale(double n) {
    //Same as: 0.5 / sqrt(n-0.5)
    //but lot faster
    return 0.5 * rsqrt(n-0.5);
}

在一些较旧的 cpu 上（硬件平方根较慢或没有硬件平方根），使用来自 Quake 的 floats 和 Fast inverse square root 可能会更快：

static inline float calculate_scale(float n) {
    return 0.5 * Q_rsqrt(n-0.5);
}

float Q_rsqrt( float number )
{
    long i;
    float x2, y;
    const float threehalfs = 1.5F;

    x2 = number * 0.5F;
    y  = number;
    i  = * ( long * ) &y;                       // evil floating point bit level hacking
    i  = 0x5f3759df - ( i >> 1 );               // what the fuck? 
    y  = * ( float * ) &i;
    y  = y * ( threehalfs - ( x2 * y * y ) );   // 1st iteration
//  y  = y * ( threehalfs - ( x2 * y * y ) );   // 2nd iteration, this can be removed

    return y;
}

有关实施的详细信息，请参阅 https://en.wikipedia.org/wiki/Fast_inverse_square_root and http://www.lomont.org/Math/Papers/2003/InvSqrt.pdf . Not recommended to use on modern cpus with hardware reciprocal square root。

不总是解：0.5 / sqrt(n-0.5)

请注意，在某些处理器上（例如 ARM Cortex A9, Intel Core2）除法几乎与硬件平方根的时间相同，所以最好使用具有 2 个平方根的原始函数 sqrt(n) - sqrt(n-1) 或者如果存在 0.5 * rsqrt(n-0.5)，则乘以倒数平方根。

3。将预计算 table 与回退

结合使用

此方法是前两种解决方案之间的良好折衷。它兼具良好的准确性和性能。

static inline double calculate_scale(double n) {
    if(n <= sizeof_precomputed_table) {
        int nIndex = (int) n;
        return precomputed_table[nIndex-1];
    }
    //Multiply + Inverse Square root
    return 0.5 * rsqrt(n-0.5);
    //OR
    return sqrt(n) - sqrt(n-1);
}

在我的例子中，我需要非常准确的数字，所以我预先计算的 table 大小是 2048。

欢迎任何反馈。

优化 sqrt(n) - sqrt(n-1)

Optimizing sqrt(n) - sqrt(n-1)

c

optimization

mathematical-optimization

1。在循环中或使用预计算 table

2。使用平方根反比对 BIG 数进行近似

3。将预计算 table 与回退