找回密码
 欢迎注册
楼主: 无心人

[擂台] x86上128位二进制乘法最快速算法征解

[复制链接]
发表于 2008-3-26 07:33:49 | 显示全部楼层
原帖由 无心人 于 2008-3-25 21:50 发表
总感觉你这么安排指令会存在冲突


我在写汇编时,尽量让相同指令聚集,且尽量让相邻语句无关,即减少依赖性;
但确实有许多靠经验、凭感觉的成分在组织安排指令顺序。
算法是经过仔细推敲论证过的。

不太明白你说的“冲突”是指什么?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2008-3-26 08:00:23 | 显示全部楼层
比如连续4个乘
每个乘均是5延迟的
最多同时执行两条
那这4个就产生执行单元冲突和指令延时的耽误
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2008-3-26 08:07:21 | 显示全部楼层
ALU 指令讲究 U/V 通道的排布,
不知 SIMD 指令是否也有双通道问题?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2008-3-26 08:37:21 | 显示全部楼层
当然啊
你没仔细读你推荐的那些PDF么?
存在执行单元的限制啊
P4是5个执行单元,但并不是每个都能执行全部指令
每个都有固定的执行指令集合
能否并行是和执行单元数量有关的
指令在进入流水线解码就绪后是在一个地方集合等待
根据不同指令,要等不同执行单元空闲且操作数就绪后
才进入执行单元执行,最后在一个地方集合退出
每次执行单元的延迟数是不同的,虽然大部分都标1周期指令
但执行周期是不同的啊
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2008-3-26 08:39:34 | 显示全部楼层
一个函数如果指令数少
也许在解码后,可填满指令的cache
就是能保证执行前并不存在延时
但并不能保证执行冲突
1、执行单元限制
2、指令依赖
3、某些具体限制,比如AGI等
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2008-3-26 21:58:52 | 显示全部楼层
写了一个 使用 MMX 寄存器,SSE2 乘 和 加 指令的版本,速度大幅领先于 现有的所有版本。
  先贴上测试结果(PIV 2.6G 768M RAM),为了避免将函数调用放在后面容易占优势的嫌疑,将对该函数的调用放在前面。

Test function: UInt128x128To256_ANSI_C32(..) 10000000 times...
Elapsed time: 2035.286 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_76F(..) 10000000 times...
Elapsed time: 395.657 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_40F(..) 10000000 times...
Elapsed time: 1027.907 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_42F(..) 10000000 times...
Elapsed time: 762.286 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_54F(..) 10000000 times...
Elapsed time: 734.089 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_56F(..) 10000000 times...
Elapsed time: 607.176 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_58F(..) 10000000 times...
Elapsed time: 923.039 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

Test function: UInt128x128To256_SSE2_69F(..) 10000000 times...
Elapsed time: 529.028 ms
  92F84935 2DD0353C A7166445 8B52280D * ADBE65FA 33EC4DD6 F39E291C 6C1A072B
= 63BF184A DE964DAC 8CD871A5 3FB4A4E7   D031F0A2 D4767BD7 A7C40427 3337152F

再贴上所有源代码(注意:该版有bug,作者已在 96# 修正!--gxqcn)
  1. _declspec(naked) /* 注意:该版有bug,作者已在 96# 修正!*/
  2. void UInt128x128To256_SSE2_76F( UINT32 * const result,
  3.                                 const UINT32 * const left,
  4.                                 const UINT32 * const right )
  5. // 使用MMX寄存器, 但用SSE2指令
  6. // 指令的使用,乘法 使用 32bit * 32bit 的MMX指令,加法使用 64 bit + 64 bit 的SSE2指令

  7. // 16次乘法,分4趟, 操作数一律在mmx 寄存器中进行
  8. // 第1趟计算 left[] * right[0]
  9. // 第2趟计算 left[] * right[1]
  10. // 第3趟计算 left[] * right[2]
  11. // 第4趟计算 left[] * right[3]

  12. // 每趟的计算结果为 5 DWORD,最低位DWORD 直接存储到result, 前3趟的计算结果的高4个DWORD(中间结果)存储到一组寄存器

  13. // 寄存器的使用

  14. // 中间结果 放入 mm0,mm1,mm2,mm3

  15. // 掩码寄存器 mm7,  用来请64位寄存器的高32bit
  16. // 进位 寄存器 和 当前乘积寄存器 为 mm4 和 mm5,轮换使用
  17. // 乘数 寄存器 mm6
  18. // 可以看到8个mmx寄存器全部派上用场

  19. {
  20. #define LEFT_REG        eax
  21. #define RIGHT_REG        edx       
  22. #define RESULT_REG        ecx
  23. #define MASK_REG        mm7
  24. #define MUL_REG                mm6
  25.        
  26.         __asm
  27.     {
  28.         mov                    eax, 0xffffffff
  29.                 mov         RESULT_REG, dword ptr[esp + 04h]   ; result
  30.                 movd                MASK_REG, eax

  31.         mov         RIGHT_REG, dword ptr[esp + 0Ch]   ; right
  32.                 movq                MUL_REG,   [RIGHT_REG]
  33.                 mov         LEFT_REG,  dword ptr [esp + 08h]   ; left
  34.         
  35.                
  36.                 //--------第1次乘
  37.                 movd          mm4,        [LEFT_REG+0]
  38.                 pmuludq          mm4, MUL_REG
  39.                 movd          [RESULT_REG], mm4
  40.                 psrlq          mm4, 32                //carry

  41.                 //------第2次乘
  42.                 movd           mm5,        [LEFT_REG+4]
  43.                 pmuludq           mm5, MUL_REG
  44.                 paddq           mm5, mm4                //64bit + 32 bit carry
  45.                 movq           mm0, mm5
  46.                 psrlq           mm5, 32                        //carry
  47.                 pand           mm0, MASK_REG        //mm0,中间结果
  48.                
  49.                 //------第3次乘
  50.                 movd           mm4,        [LEFT_REG+8]
  51.                 pmuludq           mm4, MUL_REG
  52.                 paddq           mm4, mm5                        //64bit + 32 bit carry
  53.                 movq           mm1, mm4
  54.                 psrlq           mm4, 32                        //carry
  55.                 pand           mm1, MASK_REG        //mm1,中间结果
  56.                
  57.                 //------第4次乘
  58.                 movd           mm5,        [LEFT_REG+12]
  59.                 pmuludq           mm5, MUL_REG
  60.                 paddq           mm5, mm4                //64bit + 32 bit carry
  61.                 movq           mm2, mm5
  62.                 psrlq           mm5, 32                        //carry
  63.                 pand           mm2, MASK_REG        //mm2,中间结果
  64.                 movq           mm3, mm5                        //mm3,中间结果

  65.                 //----得到乘数
  66.                 psrlq          MUL_REG, 32                //得到right[1]

  67.                
  68.                 //---第二趟乘---------
  69.                
  70.                 //--------第5次乘
  71.                 movd          mm4,        [LEFT_REG+0]
  72.                 pmuludq          mm4, MUL_REG
  73.                 paddq          mm4, mm0                // 64bit + 32 bit 中间结果
  74.                 movd          [RESULT_REG+4], mm4
  75.                 psrlq          mm4, 32                //carry

  76.                 //------第6次乘
  77.                 movd           mm5,        [LEFT_REG+4]
  78.                 pmuludq           mm5, MUL_REG
  79.                 paddq           mm5, mm4                //64bit + 32 bit carry
  80.                 paddq           mm5, mm1                //64bit + 32 bit 中间结果
  81.                 movq           mm0, mm5
  82.                 psrlq           mm5, 32                        //carry
  83.                 pand           mm0, MASK_REG        //mm0,中间结果
  84.                
  85.                 //------第7次乘
  86.                 movd           mm4,        [LEFT_REG+8]
  87.                 pmuludq           mm4, MUL_REG
  88.                 paddq           mm4, mm5                        //64bit + 32 bit carry
  89.                 paddq           mm4, mm2                        //64bit + 32 bit 中间结果
  90.                 movq           mm1, mm4
  91.                 psrlq           mm4, 32                        //carry
  92.                 pand           mm1, MASK_REG        //mm1,中间结果
  93.                
  94.                 //------第8次乘
  95.                 movd           mm5,        [LEFT_REG+12]
  96.                 pmuludq           mm5, MUL_REG
  97.                 paddq           mm5, mm4                //64bit + 32 bit carry
  98.                 paddq           mm5, mm3                //64bit + 32 bit 中间结果
  99.                 movq           mm2, mm5
  100.                 psrlq           mm5, 32                        //carry
  101.                 pand           mm2, MASK_REG        //mm2,中间结果
  102.                 movq           mm3, mm5                        //mm3,中间结果

  103.                 //----得到乘数
  104.                 movq                MUL_REG,   [RIGHT_REG+8]
  105.                
  106.                 //--------第9次乘
  107.                 movd          mm4,        [LEFT_REG+0]
  108.                 pmuludq          mm4, MUL_REG
  109.                 paddq          mm4, mm0                // 64bit + 32 bit 中间结果
  110.                 movd          [RESULT_REG+8], mm4
  111.                 psrlq          mm4, 32                //carry

  112.                 //------第10次乘
  113.                 movd           mm5,        [LEFT_REG+4]
  114.                 pmuludq           mm5, MUL_REG
  115.                 paddq           mm5, mm4                //64bit + 32 bit carry
  116.                 paddq           mm5, mm1                //64bit + 32 bit 中间结果
  117.                 movq           mm0, mm5
  118.                 psrlq           mm5, 32                        //carry
  119.                 pand           mm0, MASK_REG        //mm0,中间结果
  120.                
  121.                 //------第11次乘
  122.                 movd           mm4,        [LEFT_REG+8]
  123.                 pmuludq           mm4, MUL_REG
  124.                 paddq           mm4, mm5                        //64bit + 32 bit carry
  125.                 paddq           mm4, mm2                        //64bit + 32 bit 中间结果
  126.                 movq           mm1, mm4
  127.                 psrlq           mm4, 32                        //carry
  128.                 pand           mm1, MASK_REG        //mm1,中间结果
  129.                
  130.                 //------第12次乘
  131.                 movd           mm5,        [LEFT_REG+12]
  132.                 pmuludq           mm5, MUL_REG
  133.                 paddq           mm5, mm4                //64bit + 32 bit carry
  134.                 paddq           mm5, mm3                //64bit + 32 bit 中间结果
  135.                 movq           mm2, mm5
  136.                 psrlq           mm5, 32                        //carry
  137.                 pand           mm2, MASK_REG        //mm2,中间结果
  138.                 movq           mm3, mm5                        //mm3,中间结果

  139.                
  140.                 //----得到乘数
  141.                 psrlq          MUL_REG, 32                //得到right[1]
  142.                
  143.                 //--------第13次乘
  144.                 movd          mm4,        [LEFT_REG+0]
  145.                 pmuludq          mm4, MUL_REG
  146.                 paddq          mm4, mm0                // 64bit + 32 bit 中间结果
  147.                 movd          [RESULT_REG+12], mm4
  148.                 psrlq          mm4, 32                //carry

  149.                 //------第14次乘
  150.                 movd           mm5,        [LEFT_REG+4]
  151.                 pmuludq           mm5, MUL_REG
  152.                 paddq           mm5, mm4                //64bit + 32 bit carry
  153.                 paddq           mm5, mm1                //64bit + 32 bit 中间结果
  154.                 movd           [RESULT_REG+16], mm5
  155.                 psrlq           mm5, 32                        //carry
  156.                
  157.                
  158.                 //------第15次乘
  159.                 movd           mm4,        [LEFT_REG+8]
  160.                 pmuludq           mm4, MUL_REG
  161.                 paddq           mm4, mm5                        //64bit + 32 bit carry
  162.                 paddq           mm4, mm2                        //64bit + 32 bit 中间结果
  163.                 movd           [RESULT_REG+20], mm4
  164.                 psrlq           mm4, 32                        //carry
  165.                
  166.                
  167.                 //------第16次乘
  168.                 movd           mm5,        [LEFT_REG+12]
  169.                 pmuludq           mm5, MUL_REG
  170.                 paddq           mm5, mm4                //64bit + 32 bit carry
  171.                 paddq           mm5, mm3                //64bit + 32 bit 中间结果
  172.                 movq           [RESULT_REG+24], mm5
  173.                 psrlq           mm5, 32                        //carry
  174.                
  175.                 //------保存最高DWORD
  176.                 movq           [RESULT_REG+28], mm5                        //mm3,中间结果

  177.                 //emms
  178.         ret
  179.     }
  180. }
复制代码
我的指令数为111条,指令数可能是所有版本中最小的。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2008-3-27 08:04:51 | 显示全部楼层
在我的 56# 中的指令数为 110 行(包含最后的那个 ret 指令),比ls的更少。

因为 mmx 寄存器与 fpu 寄存器物理上是共享的,
所以 ls 若将 emms 指令屏蔽掉,将导致程序后面无法正常切换到浮点单元,只不过测试程序通篇全整型,不存在这个问题罢了。

不知我这款 AMD CPU 是 MMX 优化得太好,还是对 SSE2 支持得太烂,测试结果如下:
UInt128x128To256_ANSI_C32(..): 640.631 ms
UInt128x128To256_SSE2_40F(..): 516.867 ms
UInt128x128To256_SSE2_42F(..): 810.261 ms
UInt128x128To256_SSE2_54F(..): 755.762 ms
UInt128x128To256_SSE2_56F(..): 646.375 ms
UInt128x128To256_SSE2_58F(..): 639.292 ms
UInt128x128To256_SSE2_69F(..): 633.678 ms
UInt128x128To256_SSE2_76F(..): 379.350 ms

在我这台机器上,76F 替代了 40F 坐上了头把交椅
现在,非常期待 Core 2 系列的测试结果,以了解CPU对指令集的优化程度和发展趋势。。。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2008-3-27 08:10:27 | 显示全部楼层
你厉害
谁帖Core2测试结果啊?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2008-3-27 09:34:31 | 显示全部楼层
P4 2.0测试结果
P4E P4D估计不会出现其他结果
验证了我说的MMX寄存器在P4上快的结论
现在差Core结果
40 1794
42 1907
54 1894
56  894
58  756
76  642
76带不带emms均超其他函数
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2008-3-27 10:27:33 | 显示全部楼层
准备再写一个完全的SSE2版的,主要思路:
1.一次SSE2乘法指令计算 2个 32bit 乘以 32bit 的积
2.累加使用SSE2 128bit 加指令,但高64bit不使用
3.和76楼一样,使用4个XMM寄存器保存中间结果
4.依然需要使用掩码寄存器清除高bit,但使用 掩码清除高位的指令的使用次数大大减少,只需要3次.
5.依然是8大XMM寄存器其上阵,最大限度利用系统资源.
6.预计指令数将更少,但不知速度能否超越76楼.
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

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

GMT+8, 2024-3-29 00:56 , Processed in 0.044419 second(s), 15 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

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