数学研发网设为首页收藏本站
切换到窄版

 找回密码
 欢迎注册
数学研发论坛»论坛 › 【算法研发】 › 编程擂台 › 二进制32位整数快速平方根
1 ...234567891011... 15下一页
返回列表 发新帖
楼主: 无心人

[讨论] 二进制32位整数快速平方根

[复制链接]
liangbch
发表于 2009-2-13 09:38:58 | 显示全部楼层
我也利用54楼的思想写了一个函数。 这个函数包含2种做法，当 定义了 LOOP_UNROLL， 使用swich case 技巧做了循环展开，速度应该更快。当LOOP_UNROLL未定义，使用常规的循环算法，楼主可否测试一下这个版本的速度。 下面给出源代码。
  3. typedef unsigned long DWORD;
  4. typedef unsigned char BYTE;
  6. inline DWORD log2(DWORD n)
  7. {
  8. _asm
  9. {
  10. mov ecx,n
  11. bsr eax,ecx
  12. }
  13. }
  15. const unsigned char sqrtTab2[]=
  16. {
  17. 0, 1, 1, 1, 2, 2, 2, 2,
  18. 2, 3, 3, 3, 3, 3, 3, 3,
  19. 4, 4, 4, 4, 4, 4, 4, 4,
  20. 4, 5, 5, 5, 5, 5, 5, 5,
  21. 5, 5, 5, 5, 6, 6, 6, 6,
  22. 6, 6, 6, 6, 6, 6, 6, 6,
  23. 6, 7, 7, 7, 7, 7, 7, 7,
  24. 7, 7, 7, 7, 7, 7, 7, 7,
  25. 8, 8, 8, 8, 8, 8, 8, 8,
  26. 8, 8, 8, 8, 8, 8, 8, 8,
  27. 8, 9, 9, 9, 9, 9, 9, 9,
  28. 9, 9, 9, 9, 9, 9, 9, 9,
  29. 9, 9, 9, 9, 10, 10, 10, 10,
  30. 10, 10, 10, 10, 10, 10, 10, 10,
  31. 10, 10, 10, 10, 10, 10, 10, 10,
  32. 10, 11, 11, 11, 11, 11, 11, 11,
  33. 11, 11, 11, 11, 11, 11, 11, 11,
  34. 11, 11, 11, 11, 11, 11, 11, 11,
  35. 12, 12, 12, 12, 12, 12, 12, 12,
  36. 12, 12, 12, 12, 12, 12, 12, 12,
  37. 12, 12, 12, 12, 12, 12, 12, 12,
  38. 12, 13, 13, 13, 13, 13, 13, 13,
  39. 13, 13, 13, 13, 13, 13, 13, 13,
  40. 13, 13, 13, 13, 13, 13, 13, 13,
  41. 13, 13, 13, 13, 14, 14, 14, 14,
  42. 14, 14, 14, 14, 14, 14, 14, 14,
  43. 14, 14, 14, 14, 14, 14, 14, 14,
  44. 14, 14, 14, 14, 14, 14, 14, 14,
  45. 14, 15, 15, 15, 15, 15, 15, 15,
  46. 15, 15, 15, 15, 15, 15, 15, 15,
  47. 15, 15, 15, 15, 15, 15, 15, 15,
  48. 15, 15, 15, 15, 15, 15, 15, 15
  49. };
  51. const unsigned char sqrTab2[]=
  52. {
  53. 0, 1, 4, 9, 16, 25, 36, 49,
  54. 64, 81, 100,121,144,169,196,225
  55. };
  58. #define LOOP_UNROLL
  59. extern "C"
  60. DWORD __fastcall UintSqrt2(DWORD x)
  61. {
  62. DWORD bc,m1,m2,n,t,i;
  64. if (x<256)
  65. return sqrtTab2[x];
  67. bc=log2(x)/2;
  68. m2=(x>>(bc*2-6));
  69. n=sqrtTab2[m2];
  70. m1=sqrTab2[n];
  72. #ifdef LOOP_UNROLL
  73. switch(bc)
  74. {
  75. case 15:
  76. m2= (x>>22);
  77. n <<= 1; m1<<= 2;
  78. t = m1 + (n<<1)+1;
  79. if ( t<=m2) { m1=t; n++; }
  80. case 14:
  81. m2= (x>>20);
  82. n <<= 1; m1<<= 2;
  83. t = m1 + (n<<1)+1;
  84. if ( t<=m2) { m1=t; n++; }
  85. case 13:
  86. m2= (x>>18);
  87. n <<= 1; m1<<= 2;
  88. t = m1 + (n<<1)+1;
  89. if ( t<=m2) { m1=t; n++; }
  90. case 12:
  91. m2= (x>>16);
  92. n <<= 1; m1<<= 2;
  93. t = m1 + (n<<1)+1;
  94. if ( t<=m2) { m1=t; n++; }
  95. case 11:
  96. m2= (x>>14);
  97. n <<= 1; m1<<= 2;
  98. t = m1 + (n<<1)+1;
  99. if ( t<=m2) { m1=t; n++; }
  100. case 10:
  101. m2= (x>>12);
  102. n <<= 1; m1<<= 2;
  103. t = m1 + (n<<1)+1;
  104. if ( t<=m2) { m1=t; n++; }
  105. case 9:
  106. m2= (x>>10);
  107. n <<= 1; m1<<= 2;
  108. t = m1 + (n<<1)+1;
  109. if ( t<=m2) { m1=t; n++; }
  110. case 8:
  111. m2= (x>>8);
  112. n <<= 1; m1<<= 2;
  113. t = m1 + (n<<1)+1;
  114. if ( t<=m2) { m1=t; n++; }
  115. case 7:
  116. m2= (x>>6);
  117. n <<= 1; m1<<= 2;
  118. t = m1 + (n<<1)+1;
  119. if ( t<=m2) { m1=t; n++; }
  120. case 6:
  121. m2= (x>>4);
  122. n <<= 1; m1<<= 2;
  123. t = m1 + (n<<1)+1;
  124. if ( t<=m2) { m1=t; n++; }
  125. case 5:
  126. m2= (x>>2);
  127. n <<= 1; m1<<= 2;
  128. t = m1 + (n<<1)+1;
  129. if ( t<=m2) { m1=t; n++; }
  130. case 4:
  131. m2= x;
  132. n <<= 1; m1<<= 2;
  133. t = m1 + (n<<1)+1;
  134. if ( t<=m2) { m1=t; n++; }
  135. case 3:
  136. case 2:
  137. case 1:
  138. case 0:;
  139. }
  140. #else
  141. for (i=bc;i>=4;i--)
  142. {
  143. m2=(x>>(i*2-8));
  144. n <<= 1; m1 <<= 2;
  145. t= m1 + (n<<1)+1;
  146. if ( t<=m2)
  147. { m1=t; n++; }
  148. }
  149. #endif
  150. return n;
  151. }
复制代码
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 09:49:32 | 显示全部楼层
GxQ的代码 好像，并不快啊
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 09:51:26 | 显示全部楼层
GxQ代码的汇编版本 提高50%左右性能
  2. _declspec(naked)
  3. UINT32 __fastcall iSqrt_GxQ_asm( UINT32 n )
  4. {
  5. /*
  6. UINT32 t;
  7. UINT32 ret = 0;
  9. #define SQRT_CORE(x) \
  10. t = ret + (1UL<<(x)); \
  11. ret >>= 1; \
  12. if ( n >= t ) \
  13. { \
  14. n -= t; \
  15. ret += (1UL<<(x)); \
  16. }
  18. SQRT_CORE(30, 28, ..., ..., 0);
  19. */
  20. __asm
  21. {
  22. push ebx //t
  23. push esi
  24. push edi
  25. xor eax, eax //ret
  26. mov ebx, eax
  27. mov edx, 1
  28. shl edx, 30
  29. //30
  30. add ebx, edx
  31. shr eax, 1
  32. mov esi, ecx
  33. sub esi, ebx
  34. mov edi, eax
  35. add edi, edx
  36. cmp ecx, ebx
  37. cmovae ecx, esi
  38. cmovae eax, edi
  39. //28
  40. shr edx, 2
  41. mov ebx, eax
  42. add ebx, edx
  43. shr eax, 1
  44. mov esi, ecx
  45. sub esi, ebx
  46. mov edi, eax
  47. add edi, edx
  48. cmp ecx, ebx
  49. cmovae ecx, esi
  50. cmovae eax, edi
  51. //26
  52. shr edx, 2
  53. mov ebx, eax
  54. add ebx, edx
  55. shr eax, 1
  56. mov esi, ecx
  57. sub esi, ebx
  58. mov edi, eax
  59. add edi, edx
  60. cmp ecx, ebx
  61. cmovae ecx, esi
  62. cmovae eax, edi
  63. //24
  64. shr edx, 2
  65. mov ebx, eax
  66. add ebx, edx
  67. shr eax, 1
  68. mov esi, ecx
  69. sub esi, ebx
  70. mov edi, eax
  71. add edi, edx
  72. cmp ecx, ebx
  73. cmovae ecx, esi
  74. cmovae eax, edi
  75. //22
  76. shr edx, 2
  77. mov ebx, eax
  78. add ebx, edx
  79. shr eax, 1
  80. mov esi, ecx
  81. sub esi, ebx
  82. mov edi, eax
  83. add edi, edx
  84. cmp ecx, ebx
  85. cmovae ecx, esi
  86. cmovae eax, edi
  87. //20
  88. shr edx, 2
  89. mov ebx, eax
  90. add ebx, edx
  91. shr eax, 1
  92. mov esi, ecx
  93. sub esi, ebx
  94. mov edi, eax
  95. add edi, edx
  96. cmp ecx, ebx
  97. cmovae ecx, esi
  98. cmovae eax, edi
  99. //18
  100. shr edx, 2
  101. mov ebx, eax
  102. add ebx, edx
  103. shr eax, 1
  104. mov esi, ecx
  105. sub esi, ebx
  106. mov edi, eax
  107. add edi, edx
  108. cmp ecx, ebx
  109. cmovae ecx, esi
  110. cmovae eax, edi
  111. //16
  112. shr edx, 2
  113. mov ebx, eax
  114. add ebx, edx
  115. shr eax, 1
  116. mov esi, ecx
  117. sub esi, ebx
  118. mov edi, eax
  119. add edi, edx
  120. cmp ecx, ebx
  121. cmovae ecx, esi
  122. cmovae eax, edi
  123. //14
  124. shr edx, 2
  125. mov ebx, eax
  126. add ebx, edx
  127. shr eax, 1
  128. mov esi, ecx
  129. sub esi, ebx
  130. mov edi, eax
  131. add edi, edx
  132. cmp ecx, ebx
  133. cmovae ecx, esi
  134. cmovae eax, edi
  135. //12
  136. shr edx, 2
  137. mov ebx, eax
  138. add ebx, edx
  139. shr eax, 1
  140. mov esi, ecx
  141. sub esi, ebx
  142. mov edi, eax
  143. add edi, edx
  144. cmp ecx, ebx
  145. cmovae ecx, esi
  146. cmovae eax, edi
  147. //10
  148. shr edx, 2
  149. mov ebx, eax
  150. add ebx, edx
  151. shr eax, 1
  152. mov esi, ecx
  153. sub esi, ebx
  154. mov edi, eax
  155. add edi, edx
  156. cmp ecx, ebx
  157. cmovae ecx, esi
  158. cmovae eax, edi
  159. //8
  160. shr edx, 2
  161. mov ebx, eax
  162. add ebx, edx
  163. shr eax, 1
  164. mov esi, ecx
  165. sub esi, ebx
  166. mov edi, eax
  167. add edi, edx
  168. cmp ecx, ebx
  169. cmovae ecx, esi
  170. cmovae eax, edi
  171. //6
  172. shr edx, 2
  173. mov ebx, eax
  174. add ebx, edx
  175. shr eax, 1
  176. mov esi, ecx
  177. sub esi, ebx
  178. mov edi, eax
  179. add edi, edx
  180. cmp ecx, ebx
  181. cmovae ecx, esi
  182. cmovae eax, edi
  183. //4
  184. shr edx, 2
  185. mov ebx, eax
  186. add ebx, edx
  187. shr eax, 1
  188. mov esi, ecx
  189. sub esi, ebx
  190. mov edi, eax
  191. add edi, edx
  192. cmp ecx, ebx
  193. cmovae ecx, esi
  194. cmovae eax, edi
  195. //2
  196. shr edx, 2
  197. mov ebx, eax
  198. add ebx, edx
  199. shr eax, 1
  200. mov esi, ecx
  201. sub esi, ebx
  202. mov edi, eax
  203. add edi, edx
  204. cmp ecx, ebx
  205. cmovae ecx, esi
  206. cmovae eax, edi
  207. //0
  208. shr edx, 2
  209. mov ebx, eax
  210. add ebx, edx
  211. shr eax, 1
  212. mov esi, ecx
  213. sub esi, ebx
  214. mov edi, eax
  215. add edi, edx
  216. cmp ecx, ebx
  217. cmovae ecx, esi
  218. cmovae eax, edi
  219. //ret
  220. pop edi
  221. pop esi
  222. pop ebx
  223. ret
  225. }
  226. }
复制代码
FPU2 Version: 3589.069 ms iSqrt_16 Version: 26237.479 ms iSqrt_GxQ Version: 16068.268 ms iSqrt_GxQ_asm Version: 7884.029 ms
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 10:02:58 | 显示全部楼层
FPU2 Version: 3596.297 ms iSqrt_16 Version: 26235.439 ms iSqrt_GxQ Version: 16177.534 ms iSqrt_GxQ_asm Version: 7868.024 ms UintSqrt2 Version: 21584.053 ms ========================================= 急需各位机器的测试结果 我这里Core 2 XEON 1.6似乎还是浮点方式快的多
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 10:04:42 | 显示全部楼层
我的代码参考46 iSqrt_16在55 iSqrt_GxQ在56 iSqrtGxQ_在63 UintSqrt2在62
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 10:07:47 | 显示全部楼层
测试方法 1、从1循环到2^32-1 2、测试2100000000到2200000000
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2009-2-13 10:59:38 | 显示全部楼层
这个主题中，出现了各种不同版本开平方版本，我想对这些版本全部测一遍。 1.首先由作者提出需要测试的程序，在哪一楼，那个函数。 2.某个人负责将这些程序收集起来，做成一个测试程序，打包上传。 3.任何人就可以下载编译并测试了，测试完毕提交测试环境和测试结果。 2点说明： 1. 测试程序应包括正确性测试和性能测试。 2 我的3个整数版本中，运行速度和被开放数的大小有关。 所以我建议，在做性能测试时，增加一个case,测试n=1到65535时的速度。 我先开个头。 我总共编写了6个有效的程序，其中包含一个汇编版，我需要测试的程序有 1. 11楼的 C 整数版 2. 42楼的3个浮点汇编版 3. 61楼的C 整数非牛顿迭代版。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

无心人
 楼主| 发表于 2009-2-13 11:53:31 | 显示全部楼层
我需要测试46， 55， 63 这个测试很耗费时间 我的机器测试46三个版本是160秒左右 55版本是1200秒左右 63版本是330秒 56大概是670秒
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2009-2-13 14:23:39 | 显示全部楼层
代码整理中，请稍候。。。
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

liangbch
发表于 2009-2-13 14:56:27 | 显示全部楼层
全部函数和测试程序源码包。

Uint32sqrt.rar

8.21 KB, 下载次数: 19, 下载积分: 金币 -1 枚, 经验 1 点, 下载 1 次
毋因群疑而阻独见　　毋任己意而废人言
毋私小惠而伤大体　　毋借公论以快私情
回复 支持 反对

使用道具 举报

1 ...234567891011... 15下一页
返回列表 发新帖
高级模式
B Color Image Link Quote Code Smilies
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-11-22 00:50 , Processed in 0.025919 second(s), 17 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

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