Mercurial > dropbear

comparison libtommath/bn_mp_toom_sqr.c @ 284:eed26cff980b

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
propagate from branch 'au.asn.ucc.matt.ltm.dropbear' (head 6c790cad5a7fa866ad062cb3a0c279f7ba788583) to branch 'au.asn.ucc.matt.dropbear' (head fff0894a0399405a9410ea1c6d118f342cf2aa64)
author Matt Johnston <matt@ucc.asn.au>
date Wed, 08 Mar 2006 13:23:49 +0000
parents
children 5ff8218bcee9
comparison
equal deleted inserted replaced
283:bd240aa12ba7 284:eed26cff980b
1 #include <tommath.h>
2 #ifdef BN_MP_TOOM_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 /* squaring using Toom-Cook 3-way algorithm */
19 int
20 mp_toom_sqr(mp_int *a, mp_int *b)
21 {
22 mp_int w0, w1, w2, w3, w4, tmp1, a0, a1, a2;
23 int res, B;
24
25 /* init temps */
26 if ((res = mp_init_multi(&w0, &w1, &w2, &w3, &w4, &a0, &a1, &a2, &tmp1, NULL)) != MP_OKAY) {
27 return res;
28 }
29
30 /* B */
31 B = a->used / 3;
32
33 /* a = a2 * B**2 + a1 * B + a0 */
34 if ((res = mp_mod_2d(a, DIGIT_BIT * B, &a0)) != MP_OKAY) {
35 goto ERR;
36 }
37
38 if ((res = mp_copy(a, &a1)) != MP_OKAY) {
39 goto ERR;
40 }
41 mp_rshd(&a1, B);
42 mp_mod_2d(&a1, DIGIT_BIT * B, &a1);
43
44 if ((res = mp_copy(a, &a2)) != MP_OKAY) {
45 goto ERR;
46 }
47 mp_rshd(&a2, B*2);
48
49 /* w0 = a0*a0 */
50 if ((res = mp_sqr(&a0, &w0)) != MP_OKAY) {
51 goto ERR;
52 }
53
54 /* w4 = a2 * a2 */
55 if ((res = mp_sqr(&a2, &w4)) != MP_OKAY) {
56 goto ERR;
57 }
58
59 /* w1 = (a2 + 2(a1 + 2a0))**2 */
60 if ((res = mp_mul_2(&a0, &tmp1)) != MP_OKAY) {
61 goto ERR;
62 }
63 if ((res = mp_add(&tmp1, &a1, &tmp1)) != MP_OKAY) {
64 goto ERR;
65 }
66 if ((res = mp_mul_2(&tmp1, &tmp1)) != MP_OKAY) {
67 goto ERR;
68 }
69 if ((res = mp_add(&tmp1, &a2, &tmp1)) != MP_OKAY) {
70 goto ERR;
71 }
72
73 if ((res = mp_sqr(&tmp1, &w1)) != MP_OKAY) {
74 goto ERR;
75 }
76
77 /* w3 = (a0 + 2(a1 + 2a2))**2 */
78 if ((res = mp_mul_2(&a2, &tmp1)) != MP_OKAY) {
79 goto ERR;
80 }
81 if ((res = mp_add(&tmp1, &a1, &tmp1)) != MP_OKAY) {
82 goto ERR;
83 }
84 if ((res = mp_mul_2(&tmp1, &tmp1)) != MP_OKAY) {
85 goto ERR;
86 }
87 if ((res = mp_add(&tmp1, &a0, &tmp1)) != MP_OKAY) {
88 goto ERR;
89 }
90
91 if ((res = mp_sqr(&tmp1, &w3)) != MP_OKAY) {
92 goto ERR;
93 }
94
95
96 /* w2 = (a2 + a1 + a0)**2 */
97 if ((res = mp_add(&a2, &a1, &tmp1)) != MP_OKAY) {
98 goto ERR;
99 }
100 if ((res = mp_add(&tmp1, &a0, &tmp1)) != MP_OKAY) {
101 goto ERR;
102 }
103 if ((res = mp_sqr(&tmp1, &w2)) != MP_OKAY) {
104 goto ERR;
105 }
106
107 /* now solve the matrix
108
109 0 0 0 0 1
110 1 2 4 8 16
111 1 1 1 1 1
112 16 8 4 2 1
113 1 0 0 0 0
114
115 using 12 subtractions, 4 shifts, 2 small divisions and 1 small multiplication.
116 */
117
118 /* r1 - r4 */
119 if ((res = mp_sub(&w1, &w4, &w1)) != MP_OKAY) {
120 goto ERR;
121 }
122 /* r3 - r0 */
123 if ((res = mp_sub(&w3, &w0, &w3)) != MP_OKAY) {
124 goto ERR;
125 }
126 /* r1/2 */
127 if ((res = mp_div_2(&w1, &w1)) != MP_OKAY) {
128 goto ERR;
129 }
130 /* r3/2 */
131 if ((res = mp_div_2(&w3, &w3)) != MP_OKAY) {
132 goto ERR;
133 }
134 /* r2 - r0 - r4 */
135 if ((res = mp_sub(&w2, &w0, &w2)) != MP_OKAY) {
136 goto ERR;
137 }
138 if ((res = mp_sub(&w2, &w4, &w2)) != MP_OKAY) {
139 goto ERR;
140 }
141 /* r1 - r2 */
142 if ((res = mp_sub(&w1, &w2, &w1)) != MP_OKAY) {
143 goto ERR;
144 }
145 /* r3 - r2 */
146 if ((res = mp_sub(&w3, &w2, &w3)) != MP_OKAY) {
147 goto ERR;
148 }
149 /* r1 - 8r0 */
150 if ((res = mp_mul_2d(&w0, 3, &tmp1)) != MP_OKAY) {
151 goto ERR;
152 }
153 if ((res = mp_sub(&w1, &tmp1, &w1)) != MP_OKAY) {
154 goto ERR;
155 }
156 /* r3 - 8r4 */
157 if ((res = mp_mul_2d(&w4, 3, &tmp1)) != MP_OKAY) {
158 goto ERR;
159 }
160 if ((res = mp_sub(&w3, &tmp1, &w3)) != MP_OKAY) {
161 goto ERR;
162 }
163 /* 3r2 - r1 - r3 */
164 if ((res = mp_mul_d(&w2, 3, &w2)) != MP_OKAY) {
165 goto ERR;
166 }
167 if ((res = mp_sub(&w2, &w1, &w2)) != MP_OKAY) {
168 goto ERR;
169 }
170 if ((res = mp_sub(&w2, &w3, &w2)) != MP_OKAY) {
171 goto ERR;
172 }
173 /* r1 - r2 */
174 if ((res = mp_sub(&w1, &w2, &w1)) != MP_OKAY) {
175 goto ERR;
176 }
177 /* r3 - r2 */
178 if ((res = mp_sub(&w3, &w2, &w3)) != MP_OKAY) {
179 goto ERR;
180 }
181 /* r1/3 */
182 if ((res = mp_div_3(&w1, &w1, NULL)) != MP_OKAY) {
183 goto ERR;
184 }
185 /* r3/3 */
186 if ((res = mp_div_3(&w3, &w3, NULL)) != MP_OKAY) {
187 goto ERR;
188 }
189
190 /* at this point shift W[n] by B*n */
191 if ((res = mp_lshd(&w1, 1*B)) != MP_OKAY) {
192 goto ERR;
193 }
194 if ((res = mp_lshd(&w2, 2*B)) != MP_OKAY) {
195 goto ERR;
196 }
197 if ((res = mp_lshd(&w3, 3*B)) != MP_OKAY) {
198 goto ERR;
199 }
200 if ((res = mp_lshd(&w4, 4*B)) != MP_OKAY) {
201 goto ERR;
202 }
203
204 if ((res = mp_add(&w0, &w1, b)) != MP_OKAY) {
205 goto ERR;
206 }
207 if ((res = mp_add(&w2, &w3, &tmp1)) != MP_OKAY) {
208 goto ERR;
209 }
210 if ((res = mp_add(&w4, &tmp1, &tmp1)) != MP_OKAY) {
211 goto ERR;
212 }
213 if ((res = mp_add(&tmp1, b, b)) != MP_OKAY) {
214 goto ERR;
215 }
216
217 ERR:
218 mp_clear_multi(&w0, &w1, &w2, &w3, &w4, &a0, &a1, &a2, &tmp1, NULL);
219 return res;
220 }
221
222 #endif