二进制32位整数快速平方根

仙剑魔

我抄过来改了下，只和sqrt()比较了下正确性

3. #ifdef USING_ASM
4. inline DWORD __fastcall Log2_Floor(DWORD x)
5. {
6. __asm
7. {
8. bsr eax,x
9. jnz RETURN
10. mov eax,0ffffffffh
11. RETURN:
12. }
13. }
14. #else
15. inline DWORD Log2_Floor(DWORD x)
16. {
17. DWORD l=0;
18. if(x&0xffff0000){l|=16;x>>=16;}
19. if(x&0x0000ff00){l|= 8;x>>= 8;}
20. if(x&0x000000f0){l|= 4;x>>= 4;}
21. if(x&0x0000000c){l|= 2;x>>= 2;}
22. if(x&0x00000002){l|= 1;/*x>>= 1;*/}
23. return l;
24. }
25. #endif
27. DWORD __fastcall Sqrt(DWORD n)
28. {
29. static WORD sqrtTab[]=
30. {
31. 0, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3,
32. 4, 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 5, 5,
33. 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6,
34. 6, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7,
35. 256,258,260,262,264,266,268,270,272,274,276,278,280,282,284,286,
36. 288,288,290,292,294,296,298,300,300,302,304,306,308,310,310,312,
37. 314,316,318,320,320,322,324,326,326,328,330,332,334,334,336,338,
38. 340,340,342,344,346,346,348,350,352,352,354,356,356,358,360,362,
39. 362,364,366,366,368,370,370,372,374,374,376,378,378,380,382,384,
40. 384,386,386,388,390,390,392,394,394,396,398,398,400,402,402,404,
41. 406,406,408,408,410,412,412,414,416,416,418,418,420,422,422,424,
42. 424,426,428,428,430,430,432,434,434,436,436,438,438,440,442,442,
43. 444,444,446,448,448,450,450,452,452,454,454,456,458,458,460,460,
44. 462,462,464,464,466,468,468,470,470,472,472,474,474,476,476,478,
45. 480,480,482,482,484,484,486,486,488,488,490,490,492,492,494,494,
46. 496,496,498,498,500,500,502,502,504,504,506,506,508,508,510,512
48. };
49. static const DWORD Table2[]=
50. {
51. 0<<10,1<<10,4<<10,9<<10,16<<10,25<<10,36<<10,49<<10,
52. };
54. if(n<64) return sqrtTab[n];
56. DWORD r;
57. if(n<65536)
58. {
59. r=sqrtTab[n>>10];
60. n-=Table2[r];
61. r<<=10;
62. DWORD c=0x100|r;r>>=1;
63. if(n>=c){n-=c;r|=0x100;}
64. c=0x40|r;r>>=1;
65. if(n>=c){n-=c;r|=0x40;}
66. c=0x10|r;r>>=1;
67. if(n>=c){n-=c;r|=0x10;}
68. c=0x04|r;r>>=1;
69. if(n>=c){n-=c;r|=0x04;}
70. c=0x01|r;r>>=1;
71. if(n>=c) r|=0x01;
72. return r;
73. }
75. switch(Log2_Floor(n)&0xfe)
76. {
77. #define SQRT_ADJUST {r=(r+n/r)>>1;if(r*r>n) --r;break;}
78. case 16:r=sqrtTab[n>>10] ;SQRT_ADJUST;
79. case 18:r=sqrtTab[n>>12]<<1;SQRT_ADJUST;
80. case 20:r=sqrtTab[n>>14]<<2;SQRT_ADJUST;
81. case 22:r=sqrtTab[n>>16]<<3;SQRT_ADJUST;
82. case 24:r=sqrtTab[n>>18]<<4;SQRT_ADJUST;
83. case 26:r=sqrtTab[n>>20]<<5;SQRT_ADJUST;
84. case 28:r=sqrtTab[n>>22]<<6;SQRT_ADJUST;
85. case 30:r=sqrtTab[n>>24]<<7;SQRT_ADJUST;
86. #undef SQRT_ADJUST
87. }
88. return r;
89. }

复制代码

sprime

多谢分享！

G-Spider

很多时间会跑来回顾下这里的经典算法，#137 应该是适合游戏中小浮点数的近似开方，然后也看到一段类似又有意思的话： "人们很早就在Quake3源代码中发现了如下的C代码，它可以快速的求1/sqrt(x)，在3D图形向量计算方面应用很广。 float InvSqrt(float x){ float xhalf=0.5f*x; long i=*(long*)&x; i=0x5f3759df - (i>>1); x=*(float *)&i; x=x*(1.5f-xhalf*x*x); return x; } Beyond3D.com的Ryszard Sommefeldt一直在想到底是哪个家伙写了这些神奇的代码？2003年Chris Lomont还写了一篇文章（PDF)对这些代码进行了分析。毫无疑问写出这些代码的人绝对是天才。" 是John Carmack？Michael Abrash？John Carmack在邮件回复中明确表示不是他，也不是Michael，可能是Terje Matheson。于是侦探Ryszard又向Terje Mathisen寻求答案。Terje说他写过类似的高效代码，但上面的不是。他猜测可能在MIT的旧文档中可以找到。Ryszard的求证之路显然无止境。