实现一个小型的大数库

G-Spider

这两天在琢磨一个大数库，考虑了几个方面，一是语言，二是大数结构，三是基。 结合LibTomMath和本论坛的一些经验，特别是liangbch兄的一些测试，初步确定 大数结构为:

1. struct MP_INT
2. usd dd ? //实际使用块
3. alc dd ? //已分配块
4. sgn dd ? //符号0,-1
5. dat dd ? //数据块指针
6. ends

复制代码

alc 的值不可能小于1，dat总是指向一个预先分配有空间的数组。//初始值设为MP_PREC usd 的值不可能超过alc，并且必须大于等于0，其值暗示了数组dat中第(usd-1)不为0。数组dat中第usd个及 以上位置必须为0。 如果usd为0。那么sgn的值必须是0，这代表mp_int的值为0。 sgn 表示0或正数时其值为0，负数时其值为-1. 以2为基底，位在31~28 ，内测定义30位。 MP_MASK = 3FFFFFFFh ;base 30bit 28~31 MP_BIT = 30 打算纯汇编实现，(fasm 汇编语法)，x86结构平台支持，达到一定阶段后，转向X64.... 统一stdcall调用方式(x86很缺寄存器^_^)，对于返回的成功或失败由进位标志决定。 STC—Set Carry Flag 失败 (根源于申请内存错) CLC—Clear Carry Flag 成功. 稳定后，开源，以便大牛优化和更多的算法实现。 局部标签:

1. virtual at esi
2. a:
3. .usd dd ?
4. .alc dd ?
5. .sgn dd ?
6. .dat dd ?
7. end virtual
9. virtual at edi
10. b:
11. .usd dd ?
12. .alc dd ?
13. .sgn dd ?
14. .dat dd ?
15. end virtual
17. virtual at ebx
18. c:
19. .usd dd ?
20. .alc dd ?
21. .sgn dd ?
22. .dat dd ?
23. end virtual

复制代码

G-Spider

单位数乘法mp_mul_d，用于基的转换(各种基的格式化输入)，初始非优化可能有bug的雏形：

1. align 16
2. ;a=a*data ;data is positive
3. proc mp_mul_d _a,_data
4. push esi
5. push ebx
6. push edi
8. mov esi,[_a]
9. mov ecx,[a.usd]
10. jz .A02
11. and [_data],MP_MASK
12. cmp [_data],0
13. jnz .A00
14. stdcall mp_zero,esi
15. jmp .A02
16. .A00:
17. cmp [a.alc],ecx
18. ja .A01
19. lea eax,[ecx+1]
20. stdcall mp_grow,esi,eax
21. jc .exit
22. .A01:
23. xor ebx,ebx
24. mov esi,[a.dat]
25. xor edi,edi ;carry=0
26. @@:
27. mov eax,[esi+ebx*4]
28. mul [_data] ;edx:eax
29. add eax,edi
30. adc edx,edi
31. shld edx,eax,32-MP_BIT
32. and eax,MP_MASK
33. mov [esi+ebx*4],eax
34. mov edi,edx ;edi=carry
35. inc ebx
36. cmp ebx,ecx
37. jb @B
39. test edi,edi
40. jz .A02
41. mov [esi+ebx*4],edi
42. mov esi,[_a]
43. add [a.usd],1
44. .A02:
45. clc
46. .exit:
47. pop edi
48. pop ebx
49. pop esi
50. ret
51. endp

复制代码

wayne

高山仰止 望洋兴叹

gxqcn

楼主汇编功夫了得，若完全用汇编实现大数库，工作量非常巨大，且较难维护。 但本人非常支持楼主之创举。

G-Spider

在做前期工作中，分析LibTomMath的结构意义，usd有些强要示，不妥，修正：

1. struct MP_INT
2. usd dd ?
3. alc dd ?
4. sgn dd ?
5. dat dd ?
6. ends

复制代码

1.初始化后alc 的值不可能小于1，dat总是指向一个预先分配有空间的数组。初始值设为MP_PREC。 2.usd 的值不可能超过alc，并且必须大于等于0，其值暗示了数组dat中第(usd-1)不为0。 如果usd为0,这代表mp_int的值为0,此时sgn的值必须是0，数组dat中第usd个及以上位置不保证为0。 3.sgn 表示0或正数时其值为0，负数时其值为-1. MP_PREC = 8 ;byte 8 16 32 64 128 256 ... MP_MASK = 3FFFFFFFh ;base 30bit 28~31 MP_BIT = 30 相关标签：

1. virtual at esi
2. s:
3. .usd dd ?
4. .alc dd ?
5. .sgn dd ?
6. .dat dd ?
7. end virtual
9. virtual at edi
10. e:
11. .usd dd ?
12. .alc dd ?
13. .sgn dd ?
14. .dat dd ?
15. end virtual
17. virtual at eax
18. a:
19. .usd dd ?
20. .alc dd ?
21. .sgn dd ?
22. .dat dd ?
23. end virtual
25. virtual at ebx
26. b:
27. .usd dd ?
28. .alc dd ?
29. .sgn dd ?
30. .dat dd ?
31. end virtual
33. virtual at ecx
34. c:
35. .usd dd ?
36. .alc dd ?
37. .sgn dd ?
38. .dat dd ?
39. end virtual
41. virtual at edx
42. d:
43. .usd dd ?
44. .alc dd ?
45. .sgn dd ?
46. .dat dd ?
47. end virtual

复制代码

G-Spider

比如初始化一个mp_int 结构:

1. align 16
2. ;===================================================
3. proc mp_init _c
4. stdcall XMALLOC,MP_PREC*4
5. jc .exit
7. mov ecx,[_c]
8. mov [c.dat],eax
9. mov [c.alc],MP_PREC
10. mov [c.usd],0
11. mov [c.sgn],0
13. clc
14. ;---------------------
15. .exit:
16. ret
17. endp

复制代码

清除结构：

1. align 16
2. ;===================================================
3. proc mp_clear _a
4. mov eax,[_a]
5. cmp dword [a.dat],0
6. jz .exit
7. push eax
8. stdcall XFREE,[a.dat]
9. pop eax
10. xor ecx,ecx
11. mov [a.dat],ecx
12. mov [a.alc],ecx
13. mov [a.usd],ecx
14. mov [a.sgn],ecx
15. .exit:
16. ret
17. endp

复制代码

G-Spider

2楼mp_mul_d bug adc edx,edi 更正为adc edx,0 2~64进制的串读取: ========================================================================== ;fmp_read_radix(_a,_str,_radix) ;_radix -64 -2 ~2 64 取绝对值 若基不大于36时，字母大小写忽略; 串支持5个分隔符:空格，逗号，制表，回车，换行符; 串终止符二选一：零或字符'\$'; 支持负号'-'; ;_map="0123456789ABCDEFGHIJKLMNOPQRSTUVWSYZabcdefghijklmnopqrstuvwxyz+/",20h,',',09h,0ah,0dh ;最多10+26+26+2=64个有效字符+分隔号5个=69 ;预扫描终止符:0或'\$' ;输入串，以0结束或以\$结束: ;str[]="-12abc03"; ;str[]=" 1,2ab,c03"; ;str[]=" - 1, 2ab,c03\$afdf"; 预扫描字符' - 1, 2ab,c03\$' ;========================================================================= 为验证正确性，与tommath.lib得到的结果一致(base 30bit). mp_read_radix(&a,"-10000000000f000fefefadadefe000000f0000000Asdf00003fffffff",64); tommath不支持分隔符。 fmp_read_radix测试串szTest szTest db ' - 10000000, 000f000fefefadadefe000000f0000000 Asdf00,003 fffffff\$', 09h //stdcall fmp_read_radix,mp_t,szTest,64

1. 00000018 alc
2. 0000000c usd
3. ffffffff sgn
4. 001d3a70 dat
5. =================
6. 29a69a69 dat[0]
7. 00003a69
8. 369e9000
9. 0000000a
10. 00a40000
11. 00000000
12. 249e8a68
13. 29a29927
14. 00000a68
15. 00000029
16. 00000000
17. 00000040

复制代码

G-Spider

基于上面的结构，端午写了个64位的汇编版本，基2^62。api不再有返回错误，(只存在内存申请错，则终止)。内存申请可对齐到指定的字节，默认32字节。 简单做了个64位的static lib 作为与其它语言的调用测试(未严格测试正确性)，dll调用起来也容易与kernel32.dll无异:

1. ;==========================================================================================
2. ;mpz_read_radix(a,str,radix)
3. ;radix -64 -2~2 64
4. ;map="0123456789ABCDEFGHIJKLMNOPQRSTUVWSYZabcdefghijklmnopqrstuvwxyz+/",20h,',',09h,0ah,0dh
5. ;最多10+26+26+2=64个有效字符+分隔号5个=69
6. ;预扫描终止符:0或'\[code]#include <stdio.h>
7. typedef _int64 INT64;
8. typedef unsigned _int64 UINT64;
9. typedef unsigned char UCHAR;
11. typedef struct
12. {
13. UINT64 usd;
14. UINT64 alc;
15. INT64 sng;
16. UINT64 *dat;
17. }MP_INT;
19. void mpz_init(MP_INT *);
20. void mpz_read_radix(MP_INT *,UCHAR *str,UINT64 radix);
21. void mpz_clear(MP_INT *);
23. void main()
24. {
25. UCHAR *str0=" 7,ff ff,ffff,ffff,fff\$";
26. UCHAR *str1=" 1,2ab,c03";
27. UCHAR *str2=" - 10000000, 000f000fefefadadefe000000f0000000 Asdf00,003 fffffff";
28. MP_INT cT0,cT1,cT2;
29. UINT64 i=0;
31. mpz_init(&cT0);
32. mpz_init(&cT1);
33. mpz_init(&cT2);
35. mpz_read_radix(&cT0,str0,16);
36. mpz_read_radix(&cT1,str1,16);
37. mpz_read_radix(&cT2,str2,64);
38. printf("\ndat0=");
39. for(i=cT0.usd;i>0;i--)
40. printf("%016I64X ",cT0.dat[i-1]);
41. printf("\ndat1=");
42. if(cT1.sng < 0)
43. printf("-");
44. for(i=cT1.usd;i>0;i--)
45. printf("%016I64X ",cT1.dat[i-1]);
46. printf("\ndat2=");
47. if(cT2.sng < 0)
48. printf("-");
49. for(i=cT2.usd;i>0;i--)
50. printf("%016I64X ",cT2.dat[i-1]);
52. mpz_clear(&cT0);
53. mpz_clear(&cT1);
54. mpz_clear(&cT2);
56. system("pause");
57. }

复制代码

mpzlib.zip (1.26 KB, 下载次数: 7)

2013-6-13 23:32 上传

点击文件名下载附件
阅读权限: 5

;输入串，以0结束或以\$结束: ;str[]="-12abc03"; ;str[]=" 1,2ab,c03"; ;str[]=" - 1, 2ab,c03\$afdf"; 预扫描有效字符' - 1, 2ab,c03\ ;==========================================================================================[/code]

gogdizzy

G-Spider 发表于 2013-6-13 23:32
基于上面的结构，端午写了个64位的汇编版本，基2^62。api不再有返回错误，(只存在内存申请错，则终止)。内 ...

怎么第一个数字是7fff
输出是13fff

G-Spider

gogdizzy 发表于 2013-7-2 18:44
怎么第一个数字是7fff
输出是13fff

第一个数字0x7fffffffffffffff是2^4为基(逗号仅仅为了增加某些输入的可读性)。转化为2^62为基输出就是0x13fffffffffffffff(大数运算的内部表示)，最高位在bit63无法表示，进位到下一个limb中，如果是2^64为基，则为0x7fffffffffffffff。