鱼C论坛

 找回密码
 立即注册
查看: 5295|回复: 14

[技术交流] 一个Sqrt函数引发的血案

[复制链接]
发表于 2011-11-25 18:41:50 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能^_^

您需要 登录 才可以下载或查看,没有账号?立即注册

x
好吧,我承认我标题党了,不过既然你来了,就认真看下去吧,保证你有收获。网上看到的
我们平时经常会有一些数据运算的操作,需要调用sqrt,exp,abs等函数,那么时候你有没有想过:这个些函数系统是如何实现的?就拿最常用的sqrt函数来说吧,系统怎么来实现这个经常调用的函数呢?
虽然有可能你平时没有想过这个问题,不过正所谓是“临阵磨枪,不快也光”,你“眉头一皱,计上心来”,这个不是太简单了嘛,用二分的方法,在一个区间中,每次拿中间数的平方来试验,如果大了,就再试左区间的中间数;如果小了,就再拿右区间的中间数来试。比如求sqrt(16)的结果,你先试(0+16)/2=8,8*8=64,64比16大,然后就向左移,试(0+8)/2=4,4*4=16刚好,你得到了正确的结果sqrt(16)=4。然后你三下五除二就把程序写出来了:
  1. float SqrtByBisection(float n) //用二分法
  2. {
  3.         if(n<0) //小于0的按照你需要的处理
  4.                 return n;
  5.         float mid,last;
  6.         float low,up;
  7.         low=0,up=n;
  8.         mid=(low+up)/2;
  9.         do
  10.         {
  11.                 if(mid*mid>n)
  12.                         up=mid;
  13.                 else
  14.                         low=mid;
  15.                 last=mid;
  16.                 mid=(up+low)/2;
  17.         }while(abs(mid-last) > eps);//精度控制
  18.         return mid;
  19. }
复制代码
然后看看和系统函数性能和精度的差别(其中时间单位不是秒也不是毫秒,而是CPU Tick,不管单位是什么,统一了就有可比性)

                               
登录/注册后可看大图
从图中可以看出,二分法和系统的方法结果上完全相同,但是性能上整整差了几百倍。为什么会有这么大的区别呢?难道系统有什么更好的办法?难道。。。。哦,对了,回忆下我们曾经的高数课,曾经老师教过我们“牛顿迭代法快速寻找平方根”,或者这种方法可以帮助我们,具体步骤如下:
  1. 求出根号a的近似值:首先随便猜一个近似值x,然后不断令x等于x和a/x的平均数,迭代个六七次后x的值就已经相当精确了。
复制代码
例如,我想求根号2等于多少。假如我猜测的结果为4,虽然错的离谱,但你可以看到使用牛顿迭代法后这个值很快就趋近于根号2了:(       4  + 2/4        ) / 2 = 2.25(     2.25 + 2/2.25     ) / 2 = 1.56944..( 1.56944..+ 2/1.56944..) / 2 = 1.42189..( 1.42189..+ 2/1.42189..) / 2 = 1.41423..…. file:///C:\Users\CaiXiaoyu\AppData\Roaming\Tencent\Users\470572797\QQ\WinTemp\RichOle\@2HIE]BM)GLBX~V}EM%N`VG.jpg 1.jpg 这种算法的原理很简单,我们仅仅是不断用(x,f(x))的切线来逼近方程x^2-a=0的根。根号a实际上就是x^2-a=0的一个正实根,这个函数的导数是2x。也就是说,函数上任一点(x,f(x))处的切线斜率是2x。那么,x-f(x)/(2x)就是一个比x更接近的近似值。代入 f(x)=x^2-a得到x-(x^2-a)/(2x),也就是(x+a/x)/2。相关的代码如下:
  1. float SqrtByNewton(float x)
  2. {
  3.         float val = x;//最终
  4.         float last;//保存上一个计算的值
  5.         do
  6.         {
  7.                 last = val;
  8.                 val =(val + x/val) / 2;
  9.         }while(abs(val-last) > eps);
  10.         return val;
  11. }
复制代码
然后我们再来看下性能测试:

                               
登录/注册后可看大图
哇塞,性能提高了很多,可是和系统函数相比,还是有这么大差距,这是为什么呀?想啊想啊,想了很久仍然百思不得其解。突然有一天,我在网上看到一个神奇的方法,于是就有了今天的这篇文章,废话不多说,看代码先:
  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 y0
  6.         x = *(float*)&i; // convert bits BACK to float
  7.         x = x*(1.5f-xhalf*x*x); // Newton step, repeating increases accuracy
  8.         x = x*(1.5f-xhalf*x*x); // Newton step, repeating increases accuracy
  9.         x = x*(1.5f-xhalf*x*x); // Newton step, repeating increases accuracy

  10.         return 1/x;
  11. }
复制代码
然后我们最后一次来看下性能测试:

                               
登录/注册后可看大图

这次真的是质变了,结果竟然比系统的还要好。。。哥真的是震惊了!!!哥吐血了!!!一个函数引发了血案!!!血案,血案。。。





                               
登录/注册后可看大图
该贴已经同步到 太子的微博
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2011-11-25 21:53:23 | 显示全部楼层
不错 看来这个方法很神奇,不过指针
我就知道我暂时看不懂就是了~~
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2011-11-27 16:10:11 | 显示全部楼层
make先,以后研究研究
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2011-11-28 11:03:56 | 显示全部楼层
好强的贴啊
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2011-11-28 19:48:16 | 显示全部楼层
问下 你是怎么测出时间的?
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-4-1 18:45:21 | 显示全部楼层
真是厉害,顶一个
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-4-1 20:17:09 | 显示全部楼层
picture no display
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-4-3 14:56:32 | 显示全部楼层
赞一个:lol
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-4-3 15:10:49 | 显示全部楼层
雷神之锤的那个 神秘的0x5f375a86

那是一个传奇
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-5-5 16:29:28 | 显示全部楼层
牛人啊!我只是路过打酱油的。
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-5-5 17:09:51 | 显示全部楼层
感恩无私的分享与奉献 :)
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
发表于 2013-10-27 12:17:10 | 显示全部楼层
现在虽然看不懂,不过也知道,这太强了。
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

发表于 2014-2-6 13:43:45 | 显示全部楼层
路过看一看= =!
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

发表于 2014-2-6 14:20:57 From FishC Mobile | 显示全部楼层
靠,后悔没高数没学了,这么牛叉
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

发表于 2014-2-6 15:14:18 | 显示全部楼层
楼主很聪明也很诚实。
想知道小甲鱼最近在做啥?请访问 -> ilovefishc.com
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

小黑屋|手机版|Archiver|鱼C工作室 ( 粤ICP备18085999号-1 | 粤公网安备 44051102000585号)

GMT+8, 2024-3-28 17:55

Powered by Discuz! X3.4

© 2001-2023 Discuz! Team.

快速回复 返回顶部 返回列表