若干快速算法总结

1.快速开方

1.1 二分法

假定中值为最终解,定义下限为0，上限为x，然后求中值；然后比较中值的平方和x的大小，并根据大小修改下限或者上限；重新计算中值，开始新的循环，直到前后两次中值的距离小于给定的精度为止。需要注意的一点是，如果x小于1，将上限置为1。缺点：收敛过慢。

float SqrtByBisection(float n)
{
    float low,up,mid,last; 
    low=0,up=(n<1?1:n); 
    mid=(low+up)/2; 
    do
    {
        if(mid*mid>n)
            up=mid; 
        else 
            low=mid;
        last=mid;
        mid=(up+low)/2; 
    }while(fabsf(mid-last) > eps);

    return mid; 
}

1.2 牛顿迭代法

也就是将中值替换为切线方程的零根作为最终解。存在初始值的问题，但迭代次数更加少。

比如现在已知a，对a开方，则可以等价为求下示与x轴的交点：

                                                                                        

从上图可以得到每一次更新的收敛公式可以为：

                                                   

代码如下：

float SqrtByNewton(float x)
{
    float val=(((*(int *)&x)&0xff7fffff)>>1)+(64<<23);;//初始值
    float last;
    do
    {
        last = val;
        val =(val + x/val) / 2;
    }while(fabsf(val-last) > eps);
    return val;
}

1.3卡马克算法

本质上是在求取开方的倒数。只需要单次牛顿迭代，能够达到相当高的精度。

                                                                                 

同样得到迭代公式：

                                        

代码如下：

//版本1：3次迭代算法
float SqrtByCarmack( float number )
{
    int i;
    float x2, y;
    const float threehalfs = 1.5F;

    x2 = number * 0.5F;
    y  = number;
    i  = * ( int * ) &y;     
    i  = 0x5f375a86 - ( i >> 1 ); 
    y  = * ( float * ) &i;
    y  = y * ( threehalfs - ( x2 * y * y ) ); 
    y  = y * ( threehalfs - ( x2 * y * y ) );   
    y  = y * ( threehalfs - ( x2 * y * y ) ); 
    return number*y;
}

//原始版本
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;
}

                                                                              

Reference：

[1]快速浮点开方运算