|
282
|
1 #include <tommath.h> |
|
|
2 #ifdef BN_MP_KARATSUBA_SQR_C |
|
|
3 /* LibTomMath, multiple-precision integer library -- Tom St Denis |
|
|
4 * |
|
|
5 * LibTomMath is a library that provides multiple-precision |
|
|
6 * integer arithmetic as well as number theoretic functionality. |
|
|
7 * |
|
|
8 * The library was designed directly after the MPI library by |
|
|
9 * Michael Fromberger but has been written from scratch with |
|
|
10 * additional optimizations in place. |
|
|
11 * |
|
|
12 * The library is free for all purposes without any express |
|
|
13 * guarantee it works. |
|
|
14 * |
|
|
15 * Tom St Denis, [email protected], http://math.libtomcrypt.org |
|
|
16 */ |
|
|
17 |
|
|
18 /* Karatsuba squaring, computes b = a*a using three |
|
|
19 * half size squarings |
|
|
20 * |
|
|
21 * See comments of karatsuba_mul for details. It |
|
|
22 * is essentially the same algorithm but merely |
|
|
23 * tuned to perform recursive squarings. |
|
|
24 */ |
|
|
25 int mp_karatsuba_sqr (mp_int * a, mp_int * b) |
|
|
26 { |
|
|
27 mp_int x0, x1, t1, t2, x0x0, x1x1; |
|
|
28 int B, err; |
|
|
29 |
|
|
30 err = MP_MEM; |
|
|
31 |
|
|
32 /* min # of digits */ |
|
|
33 B = a->used; |
|
|
34 |
|
|
35 /* now divide in two */ |
|
|
36 B = B >> 1; |
|
|
37 |
|
|
38 /* init copy all the temps */ |
|
|
39 if (mp_init_size (&x0, B) != MP_OKAY) |
|
|
40 goto ERR; |
|
|
41 if (mp_init_size (&x1, a->used - B) != MP_OKAY) |
|
|
42 goto X0; |
|
|
43 |
|
|
44 /* init temps */ |
|
|
45 if (mp_init_size (&t1, a->used * 2) != MP_OKAY) |
|
|
46 goto X1; |
|
|
47 if (mp_init_size (&t2, a->used * 2) != MP_OKAY) |
|
|
48 goto T1; |
|
|
49 if (mp_init_size (&x0x0, B * 2) != MP_OKAY) |
|
|
50 goto T2; |
|
|
51 if (mp_init_size (&x1x1, (a->used - B) * 2) != MP_OKAY) |
|
|
52 goto X0X0; |
|
|
53 |
|
|
54 { |
|
|
55 register int x; |
|
|
56 register mp_digit *dst, *src; |
|
|
57 |
|
|
58 src = a->dp; |
|
|
59 |
|
|
60 /* now shift the digits */ |
|
|
61 dst = x0.dp; |
|
|
62 for (x = 0; x < B; x++) { |
|
|
63 *dst++ = *src++; |
|
|
64 } |
|
|
65 |
|
|
66 dst = x1.dp; |
|
|
67 for (x = B; x < a->used; x++) { |
|
|
68 *dst++ = *src++; |
|
|
69 } |
|
|
70 } |
|
|
71 |
|
|
72 x0.used = B; |
|
|
73 x1.used = a->used - B; |
|
|
74 |
|
|
75 mp_clamp (&x0); |
|
|
76 |
|
|
77 /* now calc the products x0*x0 and x1*x1 */ |
|
|
78 if (mp_sqr (&x0, &x0x0) != MP_OKAY) |
|
|
79 goto X1X1; /* x0x0 = x0*x0 */ |
|
|
80 if (mp_sqr (&x1, &x1x1) != MP_OKAY) |
|
|
81 goto X1X1; /* x1x1 = x1*x1 */ |
|
|
82 |
|
|
83 /* now calc (x1-x0)**2 */ |
|
|
84 if (mp_sub (&x1, &x0, &t1) != MP_OKAY) |
|
|
85 goto X1X1; /* t1 = x1 - x0 */ |
|
|
86 if (mp_sqr (&t1, &t1) != MP_OKAY) |
|
|
87 goto X1X1; /* t1 = (x1 - x0) * (x1 - x0) */ |
|
|
88 |
|
|
89 /* add x0y0 */ |
|
|
90 if (s_mp_add (&x0x0, &x1x1, &t2) != MP_OKAY) |
|
|
91 goto X1X1; /* t2 = x0x0 + x1x1 */ |
|
|
92 if (mp_sub (&t2, &t1, &t1) != MP_OKAY) |
|
|
93 goto X1X1; /* t1 = x0x0 + x1x1 - (x1-x0)*(x1-x0) */ |
|
|
94 |
|
|
95 /* shift by B */ |
|
|
96 if (mp_lshd (&t1, B) != MP_OKAY) |
|
|
97 goto X1X1; /* t1 = (x0x0 + x1x1 - (x1-x0)*(x1-x0))<<B */ |
|
|
98 if (mp_lshd (&x1x1, B * 2) != MP_OKAY) |
|
|
99 goto X1X1; /* x1x1 = x1x1 << 2*B */ |
|
|
100 |
|
|
101 if (mp_add (&x0x0, &t1, &t1) != MP_OKAY) |
|
|
102 goto X1X1; /* t1 = x0x0 + t1 */ |
|
|
103 if (mp_add (&t1, &x1x1, b) != MP_OKAY) |
|
|
104 goto X1X1; /* t1 = x0x0 + t1 + x1x1 */ |
|
|
105 |
|
|
106 err = MP_OKAY; |
|
|
107 |
|
|
108 X1X1:mp_clear (&x1x1); |
|
|
109 X0X0:mp_clear (&x0x0); |
|
|
110 T2:mp_clear (&t2); |
|
|
111 T1:mp_clear (&t1); |
|
|
112 X1:mp_clear (&x1); |
|
|
113 X0:mp_clear (&x0); |
|
|
114 ERR: |
|
|
115 return err; |
|
|
116 } |
|
|
117 #endif |
