1 #include <u.h> 2 #include <libc.h> 3 #include <mp.h> 4 #include <libsec.h> 5 6 char *sfactors[] = 7 { "3", "5", "17", "257", "641", "65537", "274177", "2424833", "6700417", "45592577", 8 "6487031809", "67280421310721", "1238926361552897", "59649589127497217", 9 "5704689200685129054721", "4659775785220018543264560743076778192897", 10 "7455602825647884208337395736200454918783366342657", 11 "93461639715357977769163558199606896584051237541638188580280321", 12 13 "741640062627530801524787141901937474059940781097519023905821316144415759504705008092818711693940737", 14 15 "130439874405488189727484768796509903946608530841611892186895295776832416251471863574140227977573104895898783928842923844831149032913798729088601617946094119449010595906710130531906171018354491609619193912488538116080712299672322806217820753127014424577" 16 }; 17 18 long start; 19 20 void 21 printmp(mpint *b, char *tag) 22 { 23 int n; 24 char *p; 25 26 print("%s (%d) ", tag, b->top); 27 p = mptoa(b, 10, nil, 0); 28 write(1, p, strlen(p)); 29 free(p); 30 print("\n"); 31 } 32 33 int 34 timing(void) 35 { 36 long now, span; 37 38 now = time(0); 39 span = now-start; 40 start = now; 41 42 return span; 43 } 44 45 int expdebug; 46 47 48 void 49 main(int argc, char **argv) 50 { 51 mpint *p, *k, *d, *b, *e, *x, *r; 52 int i; 53 54 start = time(0); 55 fmtinstall('B', mpconv); 56 mpsetminbits(2*Dbits); 57 58 x = mpnew(0); 59 e = mpnew(0); 60 r = mpnew(0); 61 p = mpnew(0); 62 63 // b = 2^32 64 b = mpcopy(mpone); 65 mpleft(b, 32, b); 66 67 // 2^29440 68 p = mpcopy(mpone); 69 mpleft(p, 29440, p); 70 // 2^27392 71 k = mpcopy(mpone); 72 mpleft(k, 27392, k); 73 // k = 2^29440 - 2^27392 74 mpsub(p, k, k); 75 // p = 2^29440 - 2^27392 + 1 76 mpadd(k, mpone, p); 77 78 // if(!probably_prime(p, 18)){ 79 // print("not a prime\n"); 80 // exits(0); 81 // } 82 // print("probably prime\n"); 83 84 mpright(k, 10, k); 85 printmp(k, "k ="); 86 87 expdebug = 1; 88 mpexp(b, k, p, x); 89 printmp(x, "x ="); 90 print("timing %d\n", timing()); 91 92 for(i = 0; i < nelem(sfactors); i++){ 93 d = strtomp(sfactors[i], nil, 10, nil); 94 // e = k/d 95 mpdiv(k, d, e, r); 96 printmp(r, "r ="); 97 98 // x = b^e mod p 99 mpexp(b, e, p, x); 100 printmp(x, "x ="); 101 print("timing %d\n", timing()); 102 } 103 } 104