Projective coordinates for prime curves
[u/mdw/catacomb] / exp.h
1 /* -*-c-*-
2 *
3 * $Id: exp.h,v 1.1.4.1 2004/03/20 00:13:31 mdw Exp $
4 *
5 * Generalized exponentiation
6 *
7 * (c) 2001 Straylight/Edgeware
8 */
9
10 /*----- Licensing notice --------------------------------------------------*
11 *
12 * This file is part of Catacomb.
13 *
14 * Catacomb is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU Library General Public License as
16 * published by the Free Software Foundation; either version 2 of the
17 * License, or (at your option) any later version.
18 *
19 * Catacomb is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU Library General Public License for more details.
23 *
24 * You should have received a copy of the GNU Library General Public
25 * License along with Catacomb; if not, write to the Free
26 * Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
27 * MA 02111-1307, USA.
28 */
29
30 /*----- Revision history --------------------------------------------------*
31 *
32 * $Log: exp.h,v $
33 * Revision 1.1.4.1 2004/03/20 00:13:31 mdw
34 * Projective coordinates for prime curves
35 *
36 * Revision 1.1 2001/06/16 13:00:59 mdw
37 * New generic exponentation code. Includes sliding-window simultaneous
38 * exponentiation.
39 *
40 */
41
42 #ifdef CATACOMB_EXP_H
43 # error "Multiple inclusion of <catacomb/exp.h>"
44 #endif
45
46 #define CATACOMB_EXP_H
47
48 #ifdef __cplusplus
49 extern "C" {
50 #endif
51
52 /*----- Header files ------------------------------------------------------*/
53
54 #include <stddef.h>
55
56 #include <mLib/alloc.h>
57
58 #ifndef CATACOMB_MP_H
59 # include "mp.h"
60 #endif
61
62 /*----- Data structures ---------------------------------------------------*/
63
64 typedef struct exp_simulscan {
65 mpw w;
66 size_t len;
67 const mpw *v;
68 } exp_simulscan;
69
70 typedef struct exp_simul {
71 unsigned b;
72 size_t o, n;
73 exp_simulscan *s;
74 } exp_simul;
75
76 /*----- Macros provided ---------------------------------------------------*/
77
78 /* --- Parameters --- */
79
80 #ifndef EXP_WINSZ /* Sliding window size */
81 # define EXP_WINSZ 4 /* Predefine if you need to */
82 #endif
83
84 /* --- These are determined from the window size --- */
85
86 #define EXP_TABSZ (1 << EXP_WINSZ)
87 #define EXP_THRESH (((MPW_BITS / EXP_WINSZ) << 2) + 1)
88
89 /* --- Required operations --- *
90 *
91 * The macros here are independent of the underlying group elements. You
92 * must provide the necessary group operations and other definitions. The
93 * group operation is assumed to be written multiplicatively.
94 *
95 * @EXP_TYPE@ The type of a group element, e.g., @mp *@.
96 *
97 * @EXP_COPY(d, x)@ Makes @d@ be a copy of @x@.
98 *
99 * @EXP_DROP(x)@ Discards the element @x@, reclaiming any
100 * memory it used.
101 *
102 * @EXP_MUL(a, x)@ Multiplies @a@ by @x@ (writing the result
103 * back to @a@).
104 *
105 * @EXP_FIX(x)@ Makes @x@ be a canonical representation of
106 * its value. All multiplications have the
107 * right argument canonical.
108 *
109 * @EXP_SQR(a)@ Multiplies @a@ by itself.
110 *
111 * @EXP_SETMUL(d, x, y)@ Sets @d@ to be the product of @x@ and @y@.
112 * The value @d@ has not been initialized.
113 *
114 * @EXP_SETSQR(d, x)@ Sets @d@ to be the square of @x@.
115 *
116 * Only @EXP_TYPE@, @EXP_MUL@ and @EXP_SQR@ are required for simple
117 * exponentation. Sliding window and simultaneous exponentation require all
118 * of the operations.
119 */
120
121 #ifndef EXP_TYPE
122 # error "EXP_TYPE not defined for <catacomb/exp.h>"
123 #endif
124
125 /* --- @EXP_SIMPLE@ --- *
126 *
127 * Arguments: @a@ = the result object, initially a multiplicative identity
128 * @g@ = the object to exponentiate
129 * @x@ = the exponent, as a multiprecision integer
130 *
131 * Use: Performs a simple left-to-right exponentiation. At the end
132 * of the code, the answer is left in @a@; @g@ and @x@ are
133 * unchanged.
134 */
135
136 #define EXP_SIMPLE(a, g, x) do { \
137 mpscan sc; \
138 unsigned sq = 0; \
139 \
140 /* --- Begin scanning --- */ \
141 \
142 mp_rscan(&sc, x); \
143 if (!MP_RSTEP(&sc)) \
144 goto exp_simple_exit; \
145 while (!MP_RBIT(&sc)) \
146 MP_RSTEP(&sc); \
147 \
148 /* --- Do the main body of the work --- */ \
149 \
150 EXP_FIX(g); \
151 for (;;) { \
152 EXP_MUL(a, g); \
153 sq = 0; \
154 for (;;) { \
155 if (!MP_RSTEP(&sc)) \
156 goto exp_simple_done; \
157 sq++; \
158 if (MP_RBIT(&sc)) \
159 break; \
160 } \
161 while (sq--) EXP_SQR(a); \
162 } \
163 \
164 /* --- Do a final round of squaring --- */ \
165 \
166 exp_simple_done: \
167 while (sq--) EXP_SQR(a); \
168 exp_simple_exit:; \
169 } while (0)
170
171 /* --- @EXP_WINDOW@ --- *
172 *
173 * Arguments: @a@ = the result object, initially a multiplicative identity
174 * @g@ = the object to exponentiate
175 * @x@ = the exponent, as a multiprecision integer
176 *
177 * Use: Performs a sliding-window exponentiation. At the end of the
178 * code, the answer is left in @a@; @g@ and @x@ are unchanged.
179 */
180
181 #define EXP_WINDOW(a, g, x) do { \
182 EXP_TYPE *v; \
183 EXP_TYPE g2; \
184 unsigned i, sq = 0; \
185 mpscan sc; \
186 \
187 /* --- Get going --- */ \
188 \
189 mp_rscan(&sc, x); \
190 if (!MP_RSTEP(&sc)) \
191 goto exp_window_exit; \
192 \
193 /* --- Do the precomputation --- */ \
194 \
195 EXP_FIX(g); \
196 EXP_SETSQR(g2, g); \
197 EXP_FIX(g2); \
198 v = xmalloc(EXP_TABSZ * sizeof(EXP_TYPE)); \
199 EXP_COPY(v[0], g); \
200 for (i = 1; i < EXP_TABSZ; i++) { \
201 EXP_SETMUL(v[i], v[i - 1], g2); \
202 EXP_FIX(v[i]); \
203 } \
204 EXP_DROP(g2); \
205 \
206 /* --- Skip top-end zero bits --- * \
207 * \
208 * If the initial step worked, there must be a set bit somewhere, so \
209 * keep stepping until I find it. \
210 */ \
211 \
212 while (!MP_RBIT(&sc)) \
213 MP_RSTEP(&sc); \
214 \
215 /* --- Now for the main work --- */ \
216 \
217 for (;;) { \
218 unsigned l = 1; \
219 unsigned z = 0; \
220 \
221 /* --- The next bit is set, so read a window index --- * \
222 * \
223 * Reset @i@ to zero and increment @sq@. Then, until either I read \
224 * @WINSZ@ bits or I run out of bits, scan in a bit: if it's clear, \
225 * bump the @z@ counter; if it's set, push a set bit into @i@, \
226 * shift it over by @z@ bits, bump @sq@ by @z + 1@ and clear @z@. \
227 * By the end of this palaver, @i@ is an index to the precomputed \
228 * value in @v@. \
229 */ \
230 \
231 i = 0; \
232 sq++; \
233 while (l < EXP_WINSZ && MP_RSTEP(&sc)) { \
234 l++; \
235 if (!MP_RBIT(&sc)) \
236 z++; \
237 else { \
238 i = ((i << 1) | 1) << z; \
239 sq += z + 1; \
240 z = 0; \
241 } \
242 } \
243 \
244 /* --- Do the squaring --- * \
245 * \
246 * Remember that @sq@ carries over from the zero-skipping stuff \
247 * below. \
248 */ \
249 \
250 while (sq--) EXP_SQR(a); \
251 \
252 /* --- Do the multiply --- */ \
253 \
254 EXP_MUL(a, v[i]); \
255 \
256 /* --- Now grind along through the rest of the bits --- */ \
257 \
258 sq = z; \
259 for (;;) { \
260 if (!MP_RSTEP(&sc)) \
261 goto exp_window_done; \
262 if (MP_RBIT(&sc)) \
263 break; \
264 sq++; \
265 } \
266 } \
267 \
268 /* --- Do a final round of squaring --- */ \
269 \
270 exp_window_done: \
271 while (sq--) EXP_SQR(a); \
272 for (i = 0; i < EXP_TABSZ; i++) \
273 EXP_DROP(v[i]); \
274 xfree(v); \
275 exp_window_exit:; \
276 } while (0)
277
278 /* --- @EXP_SIMUL@ --- *
279 *
280 * Arguments: @a@ = the result object, initially a multiplicative identity
281 * @f@ = pointer to a vector of base/exp pairs
282 * @n@ = the number of base/exp pairs
283 *
284 * Use: Performs a simultaneous sliding-window exponentiation. The
285 * @f@ table is an array of structures containing members @base@
286 * of type @EXP_TYPE@, and @exp@ of type @mp *@.
287 */
288
289 #define EXP_SIMUL(a, f, n) do { \
290 size_t i, j, jj, k; \
291 size_t vn = 1 << (EXP_WINSZ * n), m = (1 << n) - 1; \
292 EXP_TYPE *v = xmalloc(vn * sizeof(EXP_TYPE)); \
293 exp_simul e; \
294 unsigned sq = 0; \
295 \
296 /* --- Fill in the precomputed table --- */ \
297 \
298 j = 1; \
299 for (i = 0; i < n; i++) { \
300 EXP_COPY(v[j], f[n - 1 - i].base); \
301 EXP_FIX(v[j]); \
302 j <<= 1; \
303 } \
304 k = n * EXP_WINSZ; \
305 jj = 1; \
306 for (; i < k; i++) { \
307 EXP_SETSQR(v[j], v[jj]); \
308 EXP_FIX(v[j]); \
309 j <<= 1; jj <<= 1; \
310 } \
311 for (i = 1; i < vn; i <<= 1) { \
312 for (j = 1; j < i; j++) { \
313 EXP_SETMUL(v[j + i], v[j], v[i]); \
314 EXP_FIX(v[j + i]); \
315 } \
316 } \
317 \
318 /* --- Set up the bitscanners --- * \
319 * \
320 * Got to use custom scanners, to keep them all in sync. \
321 */ \
322 \
323 e.n = n; \
324 e.b = 0; \
325 e.s = xmalloc(n * sizeof(*e.s)); \
326 e.o = 0; \
327 for (i = 0; i < n; i++) { \
328 MP_SHRINK(f[i].exp); \
329 e.s[i].len = MP_LEN(f[i].exp); \
330 e.s[i].v = f[i].exp->v; \
331 if (e.s[i].len > e.o) \
332 e.o = e.s[i].len; \
333 } \
334 \
335 /* --- Skip as far as a nonzero column in the exponent matrix --- */ \
336 \
337 do { \
338 if (!e.o && !e.b) \
339 goto exp_simul_done; \
340 i = exp_simulnext(&e, 0); \
341 } while (!(i & m)); \
342 \
343 /* --- Now for the main work --- */ \
344 \
345 for (;;) { \
346 unsigned l = 1; \
347 unsigned z = 0; \
348 \
349 /* --- Just read a nonzero column, so read a window index --- * \
350 * \
351 * Clear high bits of @i@ and increment @sq@. Then, until either I \
352 * read @WINSZ@ columns or I run out, scan in a column and append \
353 * it to @i@. If it's zero, bump the @z@ counter; if it's nonzero, \
354 * bump @sq@ by @z + 1@ and clear @z@. By the end of this palaver, \
355 * @i@ is an index to the precomputed value in @v@, followed by \
356 * @n * z@ zero bits. \
357 */ \
358 \
359 sq++; \
360 while (l < EXP_WINSZ && (e.o || e.b)) { \
361 l++; \
362 i = exp_simulnext(&e, i); \
363 if (!(i & m)) \
364 z++; \
365 else { \
366 sq += z + 1; \
367 z = 0; \
368 } \
369 } \
370 \
371 /* --- Do the squaring --- * \
372 * \
373 * Remember that @sq@ carries over from the zero-skipping stuff \
374 * below. \
375 */ \
376 \
377 while (sq--) EXP_SQR(a); \
378 \
379 /* --- Do the multiply --- */ \
380 \
381 i >>= (z * n); \
382 EXP_MUL(a, v[i]); \
383 \
384 /* --- Now grind along through the rest of the bits --- */ \
385 \
386 sq = z; \
387 for (;;) { \
388 if (!e.o && !e.b) \
389 goto exp_simul_done; \
390 if ((i = exp_simulnext(&e, 0)) != 0) \
391 break; \
392 sq++; \
393 } \
394 } \
395 \
396 /* --- Do a final round of squaring --- */ \
397 \
398 exp_simul_done: \
399 while (sq--) EXP_SQR(a); \
400 for (i = 1; i < vn; i++) \
401 EXP_DROP(v[i]); \
402 xfree(v); \
403 } while (0)
404
405 /*----- Functions provided ------------------------------------------------*/
406
407 /* --- @exp_simulnext@ --- *
408 *
409 * Arguments: @exp_simul *e@ = pointer to state structure
410 * @size_t x@ = a current accumulator
411 *
412 * Returns: The next column of bits.
413 *
414 * Use: Scans the next column of bits for a simultaneous
415 * exponentiation.
416 */
417
418 extern size_t exp_simulnext(exp_simul */*e*/, size_t /*x*/);
419
420 /*----- That's all, folks -------------------------------------------------*/
421
422 #ifdef __cplusplus
423 }
424 #endif