InvSqrt

来源:互联网 发布:阿里云 centos 中文 编辑:程序博客网 时间:2024/06/07 03:01

(转载自飞羽的博客http://blog.csdn.net/zhengjiaxiang135/article/details/6673639)

先解释下InvSqrt函数吧!InvSqrt(value)函数相当于1.0/sqrt(value),所以你大概应该明白了它是什么意思了吧!由于计算机图形学程序里经常要求一个面或点的单位法线,就不可避免地要用到1.0/sqrt(length)这样的式子。你可能会说:1.0/sqrt(value)不是挺好的么?这还有什么好拿出来的?!

       但是!下面的这个InvSqrt函数比传统的1.0/sqrt(value)平均要快4倍,迭代一次就可以得到相应的结果,而且精度很robust。

下面是源代码:

[cpp] view plaincopy
  1. float invSqrt(float x)  
  2. {  
  3.     float xhalf = 0.5f * x;  
  4.     int i = *(int*)&x;          // get bits for floating value  
  5.     i =  0x5f375a86 - (i>>1);    // gives initial guess  
  6.     x = *(float*)&i;            // convert bits back to float  
  7.     x = x * (1.5f - xhalf*x*x); // Newton step  
  8.     return x;  
  9. }  

         乍看一下上面的代码,你或许还没明白是怎么回事。事实上,说白了,这就是个牛顿迭代法的应用。考虑函数,牛顿迭代法求f(y)=0的迭代式是。所以代入f(y)和f ' (y),就很自然地就明白了
x = x * (1.5f - xhalf*x*x); // Newton step
        这一步的由来。所以这里的关键在于如何选择一个好的初始值来进行迭代才能以尽可能少的次数就达到足够的精度呢?代码中的整形数 i 转换成其机器码对应的浮点数就是一个很好的初始值,事实上它非常地接近x的平方根分之一,因此这也是为什么这里只需迭代一次就可以得到很好的结果的原因。
        这个代码牛逼也正好牛逼在选择了一个好的数据:0x5f375a86-(i>>1),这里的i 是浮点数x对应的机器码转换成的整型数据。然后i>>1就是相当于i/2,统统左移一位,再用一个常量来减i>>2,得到另一个整型数据,这个整形数据对应的浮点数就是一个满足要求的初始值了。

        问题的关键在于,你怎么知道要选择0x5f375a86这么好的一个数据呢?对此,可以后面的参考论文[InvSqrt.pdf],由于里面的数学知识太多,在这个页面里的排版不太方面,所以就不贴出来了~~论文中对于IEEE754浮点数据的格式进行了分析,并将浮点部分与指数部分分开讨论,对指数的奇偶作了一番探讨,最后将R-(i>>1)中的R范围给大致地确定出来,并分析不同取值的结果。大家如果不喜欢看这些头疼的数学,了解一下这个代码到底是怎么回事就行了,没有必

0 0