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
printmp(mpint * b,char * tag)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
timing(void)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
main(int argc,char ** argv)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