Renamed from `rsa-decrypt', since the name was no longer appropriate.
[u/mdw/catacomb] / mpx.h
1 /* -*-c-*-
2 *
3 * $Id: mpx.h,v 1.9 1999/12/22 15:49:07 mdw Exp $
4 *
5 * Low level multiprecision arithmetic
6 *
7 * (c) 1999 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: mpx.h,v $
33 * Revision 1.9 1999/12/22 15:49:07 mdw
34 * New function for division by a small integer.
35 *
36 * Revision 1.8 1999/12/11 10:57:43 mdw
37 * Karatsuba squaring algorithm.
38 *
39 * Revision 1.7 1999/12/11 01:51:28 mdw
40 * Change Karatsuba parameters slightly.
41 *
42 * Revision 1.6 1999/12/10 23:23:51 mdw
43 * Karatsuba-Ofman multiplication algorithm.
44 *
45 * Revision 1.5 1999/11/20 22:23:27 mdw
46 * Add function versions of some low-level macros with wider use.
47 *
48 * Revision 1.4 1999/11/17 18:04:43 mdw
49 * Add two's complement support. Fix a bug in MPX_UMLAN.
50 *
51 * Revision 1.3 1999/11/13 01:51:29 mdw
52 * Minor interface changes. Should be stable now.
53 *
54 * Revision 1.2 1999/11/11 17:47:55 mdw
55 * Minor changes for different `mptypes.h' format.
56 *
57 * Revision 1.1 1999/09/03 08:41:12 mdw
58 * Initial import.
59 *
60 */
61
62 #ifndef CATACOMB_MPX_H
63 #define CATACOMB_MPX_H
64
65 #ifdef __cplusplus
66 extern "C" {
67 #endif
68
69 /*----- The idea ----------------------------------------------------------*
70 *
71 * This file provides functions and macros which work on vectors of words as
72 * unsigned multiprecision integers. The interface works in terms of base
73 * and limit pointers (i.e., a pointer to the start of a vector, and a
74 * pointer just past its end) rather than base pointer and length, because
75 * that requires more arithmetic and state to work on.
76 *
77 * The interfaces are slightly bizarre in other ways. Try to use the
78 * higher-level functions where you can: they're rather better designed to
79 * actually be friendly and useful.
80 */
81
82 /*----- Header files ------------------------------------------------------*/
83
84 #include <string.h>
85
86 #ifndef CATACOMB_MPW_H
87 # include "mpw.h"
88 #endif
89
90 /*----- General manipulation ----------------------------------------------*/
91
92 /* --- @MPX_SHRINK@ --- *
93 *
94 * Arguments: @const mpw *v@ = pointer to vector of words
95 * @const mpw *vl@ = (updated) current limit of vector
96 *
97 * Use: Shrinks down the limit of a multiprecision integer vector.
98 */
99
100 #define MPX_SHRINK(v, vl) do { \
101 const mpw *_vv = (v), *_vvl = (vl); \
102 while (_vvl > _vv && !_vvl[-1]) \
103 _vvl--; \
104 (vl) = (mpw *)_vvl; \
105 } while (0)
106
107 /* --- @MPX_BITS@ --- *
108 *
109 * Arguments: @unsigned long b@ = result variable
110 * @const mpw *v@ = pointer to array of words
111 * @const mpw *vl@ = limit of vector (from @MPX_SHRINK@)
112 *
113 * Use: Calculates the number of bits in a multiprecision value.
114 */
115
116 #define MPX_BITS(b, v, vl) do { \
117 const mpw *_v = (v), *_vl = (vl); \
118 MPX_SHRINK(_v, _vl); \
119 if (_v == _vl) \
120 (b) = 0; \
121 else { \
122 unsigned long _b = MPW_BITS * (_vl - _v - 1) + 1; \
123 mpw _w = _vl[-1]; \
124 unsigned _k = MPW_BITS / 2; \
125 while (_k) { \
126 if (_w >> _k) { \
127 _w >>= _k; \
128 _b += _k; \
129 } \
130 _k >>= 1; \
131 } \
132 (b) = _b; \
133 } \
134 } while (0)
135
136 /* --- @MPX_OCTETS@ --- *
137 *
138 * Arguments: @size_t o@ = result variable
139 * @const mpw *v, *vl@ = pointer to array of words
140 *
141 * Use: Calculates the number of octets in a multiprecision value.
142 */
143
144 #define MPX_OCTETS(o, v, vl) do { \
145 const mpw *_v = (v), *_vl = (vl); \
146 MPX_SHRINK(_v, _vl); \
147 if (_v == _vl) \
148 (o) = 0; \
149 else { \
150 size_t _o = (MPW_BITS / 8) * (_vl - _v - 1); \
151 mpw _w = _vl[-1]; \
152 unsigned _k = MPW_BITS / 2; \
153 while (_k >= 8) { \
154 if (_w >> _k) { \
155 _w >>= _k; \
156 _o += _k >> 3; \
157 } \
158 _k >>= 1; \
159 } \
160 if (_w) \
161 _o++; \
162 (o) = _o; \
163 } \
164 } while (0)
165
166 /* --- @MPX_COPY@ --- *
167 *
168 * Arguments: @dv, dvl@ = destination vector base and limit
169 * @av, avl@ = source vector base and limit
170 *
171 * Use: Copies a multiprecision integer.
172 */
173
174 #define MPX_COPY(dv, dvl, av, avl) do { \
175 mpw *_dv = (dv), *_dvl = (dvl); \
176 size_t _dn = _dvl - _dv; \
177 const mpw *_av = (av), *_avl = (avl); \
178 size_t _an = _avl - _av; \
179 if (_av == _dv) { \
180 if (_dvl > _avl) \
181 memset(_dv, 0, MPWS(_dn - _an)); \
182 } else if (_an >= _dn) \
183 memmove(_dv, _av, MPWS(_dn)); \
184 else { \
185 memmove(_dv, _av, MPWS(_an)); \
186 memset(_dv + _an, 0, MPWS(_dn - _an)); \
187 } \
188 } while (0)
189
190 /* --- @MPX_ZERO@ --- *
191 *
192 * Arguments: @v, vl@ = base and limit of vector to clear
193 *
194 * Use: Zeroes the area between the two vector pointers.
195 */
196
197 #define MPX_ZERO(v, vl) do { \
198 mpw *_v = (v), *_vl = (vl); \
199 if (_v < _vl) \
200 memset(_v, 0, MPWS(_vl - _v)); \
201 } while (0)
202
203 /*----- Loading and storing -----------------------------------------------*/
204
205 /* --- @mpx_storel@ --- *
206 *
207 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
208 * @void *p@ = pointer to octet array
209 * @size_t sz@ = size of octet array
210 *
211 * Returns: ---
212 *
213 * Use: Stores an MP in an octet array, least significant octet
214 * first. High-end octets are silently discarded if there
215 * isn't enough space for them.
216 */
217
218 extern void mpx_storel(const mpw */*v*/, const mpw */*vl*/,
219 void */*p*/, size_t /*sz*/);
220
221 /* --- @mpx_loadl@ --- *
222 *
223 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
224 * @const void *p@ = pointer to octet array
225 * @size_t sz@ = size of octet array
226 *
227 * Returns: ---
228 *
229 * Use: Loads an MP in an octet array, least significant octet
230 * first. High-end octets are ignored if there isn't enough
231 * space for them.
232 */
233
234 extern void mpx_loadl(mpw */*v*/, mpw */*vl*/,
235 const void */*p*/, size_t /*sz*/);
236
237 /* --- @mpx_storeb@ --- *
238 *
239 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
240 * @void *p@ = pointer to octet array
241 * @size_t sz@ = size of octet array
242 *
243 * Returns: ---
244 *
245 * Use: Stores an MP in an octet array, most significant octet
246 * first. High-end octets are silently discarded if there
247 * isn't enough space for them.
248 */
249
250 extern void mpx_storeb(const mpw */*v*/, const mpw */*vl*/,
251 void */*p*/, size_t /*sz*/);
252
253 /* --- @mpx_loadb@ --- *
254 *
255 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
256 * @const void *p@ = pointer to octet array
257 * @size_t sz@ = size of octet array
258 *
259 * Returns: ---
260 *
261 * Use: Loads an MP in an octet array, most significant octet
262 * first. High-end octets are ignored if there isn't enough
263 * space for them.
264 */
265
266 extern void mpx_loadb(mpw */*v*/, mpw */*vl*/,
267 const void */*p*/, size_t /*sz*/);
268
269 /*----- Logical shifting --------------------------------------------------*/
270
271 /* --- @mpx_lsl@ --- *
272 *
273 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
274 * @const mpw *av, *avl@ = source vector base and limit
275 * @size_t n@ = number of bit positions to shift by
276 *
277 * Returns: ---
278 *
279 * Use: Performs a logical shift left operation on an integer.
280 */
281
282 extern void mpx_lsl(mpw */*dv*/, mpw */*dvl*/,
283 const mpw */*av*/, const mpw */*avl*/,
284 size_t /*n*/);
285
286 /* --- @mpx_lsr@ --- *
287 *
288 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
289 * @const mpw *av, *avl@ = source vector base and limit
290 * @size_t n@ = number of bit positions to shift by
291 *
292 * Returns: ---
293 *
294 * Use: Performs a logical shift right operation on an integer.
295 */
296
297 extern void mpx_lsr(mpw */*dv*/, mpw */*dvl*/,
298 const mpw */*av*/, const mpw */*avl*/,
299 size_t /*n*/);
300
301 /*----- Unsigned arithmetic -----------------------------------------------*/
302
303 /* --- @mpx_2c@ --- *
304 *
305 * Arguments: @mpw *dv, *dvl@ = destination vector
306 * @const mpw *v, *vl@ = source vector
307 *
308 * Returns: ---
309 *
310 * Use: Calculates the two's complement of @v@.
311 */
312
313 extern void mpx_2c(mpw */*dv*/, mpw */*dvl*/,
314 const mpw */*v*/, const mpw */*vl*/);
315
316 /* --- @mpx_ucmp@ --- *
317 *
318 * Arguments: @const mpw *av, *avl@ = first argument vector base and limit
319 * @const mpw *bv, *bvl@ = second argument vector base and limit
320 *
321 * Returns: Less than, equal to, or greater than zero depending on
322 * whether @a@ is less than, equal to or greater than @b@,
323 * respectively.
324 *
325 * Use: Performs an unsigned integer comparison.
326 */
327
328 #define MPX_UCMP(av, avl, op, dv, dvl) \
329 (mpx_ucmp((av), (avl), (dv), (dvl)) op 0)
330
331 extern int mpx_ucmp(const mpw */*av*/, const mpw */*avl*/,
332 const mpw */*bv*/, const mpw */*bvl*/);
333
334 /* --- @mpx_uadd@ --- *
335 *
336 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
337 * @const mpw *av, *avl@ = first addend vector base and limit
338 * @const mpw *bv, *bvl@ = second addend vector base and limit
339 *
340 * Returns: ---
341 *
342 * Use: Performs unsigned integer addition. If the result overflows
343 * the destination vector, high-order bits are discarded. This
344 * means that two's complement addition happens more or less for
345 * free, although that's more a side-effect than anything else.
346 * The result vector may be equal to either or both source
347 * vectors, but may not otherwise overlap them.
348 */
349
350 extern void mpx_uadd(mpw */*dv*/, mpw */*dvl*/,
351 const mpw */*av*/, const mpw */*avl*/,
352 const mpw */*bv*/, const mpw */*bvl*/);
353
354 /* --- @mpx_uaddn@ --- *
355 *
356 * Arguments: @mpw *dv, *dvl@ = source and destination base and limit
357 * @mpw n@ = other addend
358 *
359 * Returns: ---
360 *
361 * Use: Adds a small integer to a multiprecision number.
362 */
363
364 #define MPX_UADDN(dv, dvl, n) do { \
365 mpw *_ddv = (dv), *_ddvl = (dvl); \
366 mpw _c = (n); \
367 \
368 while (_c && _ddv < _ddvl) { \
369 mpd _x = (mpd)*_ddv + (mpd)_c; \
370 *_ddv++ = MPW(_x); \
371 _c = _x >> MPW_BITS; \
372 } \
373 } while (0)
374
375 extern void mpx_uaddn(mpw */*dv*/, mpw */*dvl*/, mpw /*n*/);
376
377 /* --- @mpx_usub@ --- *
378 *
379 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
380 * @const mpw *av, *avl@ = first argument vector base and limit
381 * @const mpw *bv, *bvl@ = second argument vector base and limit
382 *
383 * Returns: ---
384 *
385 * Use: Performs unsigned integer subtraction. If the result
386 * overflows the destination vector, high-order bits are
387 * discarded. This means that two's complement subtraction
388 * happens more or less for free, although that's more a side-
389 * effect than anything else. The result vector may be equal to
390 * either or both source vectors, but may not otherwise overlap
391 * them.
392 */
393
394 extern void mpx_usub(mpw */*dv*/, mpw */*dvl*/,
395 const mpw */*av*/, const mpw */*avl*/,
396 const mpw */*bv*/, const mpw */*bvl*/);
397
398 /* --- @mpx_usubn@ --- *
399 *
400 * Arguments: @mpw *dv, *dvl@ = source and destination base and limit
401 * @n@ = subtrahend
402 *
403 * Returns: ---
404 *
405 * Use: Subtracts a small integer from a multiprecision number.
406 */
407
408 #define MPX_USUBN(dv, dvl, n) do { \
409 mpw *_ddv = (dv), *_ddvl = (dvl); \
410 mpw _c = (n); \
411 \
412 while (_ddv < _ddvl) { \
413 mpd _x = (mpd)*_ddv - (mpd)_c; \
414 *_ddv++ = MPW(_x); \
415 if (_x >> MPW_BITS) \
416 _c = 1; \
417 else \
418 break; \
419 } \
420 } while (0)
421
422 extern void mpx_usubn(mpw */*dv*/, mpw */*dvl*/, mpw /*n*/);
423
424 /* --- @mpx_umul@ --- *
425 *
426 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
427 * @const mpw *av, *avl@ = multiplicand vector base and limit
428 * @const mpw *bv, *bvl@ = multiplier vector base and limit
429 *
430 * Returns: ---
431 *
432 * Use: Performs unsigned integer multiplication. If the result
433 * overflows the desination vector, high-order bits are
434 * discarded. The result vector may not overlap the argument
435 * vectors in any way.
436 */
437
438 extern void mpx_umul(mpw */*dv*/, mpw */*dvl*/,
439 const mpw */*av*/, const mpw */*avl*/,
440 const mpw */*bv*/, const mpw */*bvl*/);
441
442 /* --- @mpx_umuln@ --- *
443 *
444 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
445 * @const mpw *av, *avl@ = multiplicand vector base and limit
446 * @mpw m@ = multiplier
447 *
448 * Returns: ---
449 *
450 * Use: Multiplies a multiprecision integer by a single-word value.
451 * The destination and source may be equal. The destination
452 * is completely cleared after use.
453 */
454
455 #define MPX_UMULN(dv, dvl, av, avl, m) do { \
456 mpw *_dv = (dv), *_dvl = (dvl); \
457 const mpw *_av = (av), *_avl = (avl); \
458 mpw _c = 0; \
459 mpd _m = (m); \
460 \
461 while (_av < _avl) { \
462 mpd _x; \
463 if (_dv >= _dvl) \
464 break; \
465 _x = (mpd)_m * (mpd)*_av++ + _c; \
466 *_dv++ = MPW(_x); \
467 _c = _x >> MPW_BITS; \
468 } \
469 if (_dv < _dvl) { \
470 *_dv++ = MPW(_c); \
471 MPX_ZERO(_dv, _dvl); \
472 } \
473 } while (0)
474
475 extern void mpx_umuln(mpw */*dv*/, mpw */*dvl*/,
476 const mpw */*av*/, const mpw */*avl*/, mpw m);
477
478 /* --- @mpx_umlan@ --- *
479 *
480 * Arguments: @mpw *dv, *dvl@ = destination/accumulator base and limit
481 * @const mpw *av, *avl@ = multiplicand vector base and limit
482 * @mpw m@ = multiplier
483 *
484 * Returns: ---
485 *
486 * Use: Multiplies a multiprecision integer by a single-word value
487 * and adds the result to an accumulator.
488 */
489
490 #define MPX_UMLAN(dv, dvl, av, avl, m) do { \
491 mpw *_dv = (dv), *_dvl = (dvl); \
492 const mpw *_av = (av), *_avl = (avl); \
493 mpw _cc = 0; \
494 mpd _m = (m); \
495 \
496 while (_dv < _dvl && _av < _avl) { \
497 mpd _x; \
498 _x = (mpd)*_dv + (mpd)_m * (mpd)*_av++ + _cc; \
499 *_dv++ = MPW(_x); \
500 _cc = _x >> MPW_BITS; \
501 } \
502 MPX_UADDN(_dv, _dvl, _cc); \
503 } while (0)
504
505 extern void mpx_umlan(mpw */*dv*/, mpw */*dvl*/,
506 const mpw */*av*/, const mpw */*avl*/, mpw m);
507
508 /* --- @mpx_usqr@ --- *
509 *
510 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
511 * @const mpw *av, *av@ = source vector base and limit
512 *
513 * Returns: ---
514 *
515 * Use: Performs unsigned integer squaring. The result vector must
516 * not overlap the source vector in any way.
517 */
518
519 extern void mpx_usqr(mpw */*dv*/, mpw */*dvl*/,
520 const mpw */*av*/, const mpw */*avl*/);
521
522 /* --- @mpx_udiv@ --- *
523 *
524 * Arguments: @mpw *qv, *qvl@ = quotient vector base and limit
525 * @mpw *rv, *rvl@ = dividend/remainder vector base and limit
526 * @const mpw *dv, *dvl@ = divisor vector base and limit
527 * @mpw *sv, *svl@ = scratch workspace
528 *
529 * Returns: ---
530 *
531 * Use: Performs unsigned integer division. If the result overflows
532 * the quotient vector, high-order bits are discarded. (Clearly
533 * the remainder vector can't overflow.) The various vectors
534 * may not overlap in any way. Yes, I know it's a bit odd
535 * requiring the dividend to be in the result position but it
536 * does make some sense really. The remainder must have
537 * headroom for at least two extra words. The scratch space
538 * must be at least one word larger than the divisor.
539 */
540
541 extern void mpx_udiv(mpw */*qv*/, mpw */*qvl*/, mpw */*rv*/, mpw */*rvl*/,
542 const mpw */*dv*/, const mpw */*dvl*/,
543 mpw */*sv*/, mpw */*svl*/);
544
545 /* --- @mpx_udivn@ --- *
546 *
547 * Arguments: @mpw *qv, *qvl@ = storage for the quotient (may overlap
548 * dividend)
549 * @const mpw *rv, *rvl@ = dividend
550 * @mpw d@ = single-precision divisor
551 *
552 * Returns: Remainder after divison.
553 *
554 * Use: Performs a single-precision division operation.
555 */
556
557 extern mpw mpx_udivn(mpw */*qv*/, mpw */*qvl*/,
558 const mpw */*rv*/, const mpw */*rvl*/, mpw /*d*/);
559
560 /*----- Karatsuba multiplication algorithms -------------------------------*/
561
562 /* --- @KARATSUBA_CUTOFF@ --- *
563 *
564 * This is the limiting length for using Karatsuba algorithms. It's best to
565 * use the simpler classical multiplication method on numbers smaller than
566 * this.
567 */
568
569 #define KARATSUBA_CUTOFF 16
570
571 /* --- @KARATSUBA_SLOP@ --- *
572 *
573 * The extra number of words required as scratch space by the Karatsuba
574 * routines. This is a (generous) guess, since the actual amount of space
575 * required is proportional to the recursion depth.
576 */
577
578 #define KARATSUBA_SLOP 64
579
580 /* --- @mpx_kmul@ --- *
581 *
582 * Arguments: @mpw *dv, *dvl@ = pointer to destination buffer
583 * @const mpw *av, *avl@ = pointer to first argument
584 * @const mpw *bv, *bvl@ = pointer to second argument
585 * @mpw *sv, *svl@ = pointer to scratch workspace
586 *
587 * Returns: ---
588 *
589 * Use: Multiplies two multiprecision integers using Karatsuba's
590 * algorithm. This is rather faster than traditional long
591 * multiplication (e.g., @mpx_umul@) on large numbers, although
592 * more expensive on small ones.
593 *
594 * The destination and scratch buffers must be twice as large as
595 * the larger argument. The scratch space must be twice as
596 * large as the larger argument, plus the magic number
597 * @KARATSUBA_SLOP@.
598 */
599
600 extern void mpx_kmul(mpw */*dv*/, mpw */*dvl*/,
601 const mpw */*av*/, const mpw */*avl*/,
602 const mpw */*bv*/, const mpw */*bvl*/,
603 mpw */*sv*/, mpw */*svl*/);
604
605 /* --- @mpx_ksqr@ --- *
606 *
607 * Arguments: @mpw *dv, *dvl@ = pointer to destination buffer
608 * @const mpw *av, *avl@ = pointer to first argument
609 * @mpw *sv, *svl@ = pointer to scratch workspace
610 *
611 * Returns: ---
612 *
613 * Use: Squares a multiprecision integers using something similar to
614 * Karatsuba's multiplication algorithm. This is rather faster
615 * than traditional long multiplication (e.g., @mpx_umul@) on
616 * large numbers, although more expensive on small ones, and
617 * rather simpler than full-blown Karatsuba multiplication.
618 *
619 * The destination must be twice as large as the argument. The
620 * scratch space must be twice as large as the argument, plus
621 * the magic number @KARATSUBA_SLOP@.
622 */
623
624 extern void mpx_ksqr(mpw */*dv*/, mpw */*dvl*/,
625 const mpw */*av*/, const mpw */*avl*/,
626 mpw */*sv*/, mpw */*svl*/);
627
628 /*----- That's all, folks -------------------------------------------------*/
629
630 #ifdef __cplusplus
631 }
632 #endif
633
634 #endif