71092fd8c18b1462c2afa813d5ac9f584bd2592c
3 * $Id: mp.h,v 1.12 2001/06/16 12:57:43 mdw Exp $
5 * Simple multiprecision arithmetic
7 * (c) 1999 Straylight/Edgeware
10 /*----- Licensing notice --------------------------------------------------*
12 * This file is part of Catacomb.
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.
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.
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,
30 /*----- Revision history --------------------------------------------------*
33 * Revision 1.12 2001/06/16 12:57:43 mdw
34 * Move the @mpmont_factor@ structure and rename it now that it's used for
35 * Barrett simultaneous exponentiation too.
37 * Revision 1.11 2001/04/03 19:36:05 mdw
38 * Add some simple bitwise operations so that Perl can use them.
40 * Revision 1.10 2000/10/08 12:03:16 mdw
41 * Provide @mp_eq@ and @MP_EQ@ for rapidly testing equality of two
44 * Revision 1.9 2000/07/29 17:03:31 mdw
45 * Add support for left-to-right bitscanning, for use in modular
48 * Revision 1.8 2000/06/22 19:02:01 mdw
51 * Revision 1.7 2000/06/17 11:45:09 mdw
52 * Major memory management overhaul. Added arena support. Use the secure
53 * arena for secret integers. Replace and improve the MP management macros
54 * (e.g., replace MP_MODIFY by MP_DEST).
56 * Revision 1.6 1999/12/10 23:19:46 mdw
57 * Minor bugfixes. New interface for suggested destinations.
59 * Revision 1.5 1999/11/22 20:50:37 mdw
60 * Add support for computing Jacobi symbols.
62 * Revision 1.4 1999/11/21 22:13:02 mdw
63 * Add mp version of MPX_BITS.
65 * Revision 1.3 1999/11/19 13:19:14 mdw
66 * Fix const annotation.
68 * Revision 1.2 1999/11/17 18:02:16 mdw
69 * New multiprecision integer arithmetic suite.
80 /*----- Header files ------------------------------------------------------*/
87 #ifndef CATACOMB_MPW_H
91 #ifndef CATACOMB_ARENA_H
95 #ifndef CATACOMB_MPARENA_H
99 #ifndef CATACOMB_MPX_H
103 /*----- Data structures ---------------------------------------------------*/
105 /* --- A multiprecision integer --- */
108 mpw
*v
, *vl
; /* Vector of digits, current limit */
109 size_t sz
; /* Size of digit buffer in words */
110 mparena
*a
; /* Arena for buffer allocation */
111 unsigned f
; /* Flags (see below) */
112 unsigned ref
; /* Reference counter */
115 #define MP_NEG 1u /* Negative (signed magnitude) */
116 #define MP_BURN 2u /* Secret (viral flag) */
117 #define MP_CONST 4u /* Uses strange memory allocation */
118 #define MP_UNDEF 8u /* Contains nothing interesting */
119 #define MP_DESTROYED 16u /* Has been destroyed */
121 /* --- A factor for simultaneous exponentation --- *
123 * Used by the Montgomery and Barrett exponentiators.
126 typedef struct mp_expfactor
{
131 /*----- Useful constants --------------------------------------------------*/
133 extern mp mp_const
[];
135 #define MP_ZERO (&mp_const[0])
136 #define MP_ONE (&mp_const[1])
137 #define MP_TWO (&mp_const[2])
138 #define MP_THREE (&mp_const[3])
139 #define MP_FOUR (&mp_const[4])
140 #define MP_FIVE (&mp_const[5])
141 #define MP_TEN (&mp_const[6])
142 #define MP_256 (&mp_const[7])
143 #define MP_MONE (&mp_const[8])
145 #define MP_NEW ((mp *)0)
146 #define MP_NEWSEC (&mp_const[9])
148 /*----- Trivial macros ----------------------------------------------------*/
150 /* --- @MP_LEN@ --- *
152 * Arguments: @mp *m@ = pointer to a multiprecision integer
154 * Returns: Length of the integer, in words.
157 #define MP_LEN(m) ((m)->vl - ((m)->v))
159 /*----- Memory management and reference counting --------------------------*/
161 /* --- @mp_new@ --- *
163 * Arguments: @size_t sz@ = size of vector required
164 * @unsigned f@ = flags to set
166 * Returns: Pointer to a new MP structure.
168 * Use: Allocates a new multiprecision integer. The data space is
169 * allocated from either the standard global or secret arena,
170 * depending on the initial flags requested.
173 extern mp
*mp_new(size_t /*sz*/, unsigned /*f*/);
175 /* --- @mp_create@ --- *
177 * Arguments: @size_t sz@ = size of vector required
179 * Returns: Pointer to pristine new MP structure with enough memory
182 * Use: Creates a new multiprecision integer with indeterminate
183 * contents. The integer has a single reference.
186 extern mp
*mp_create(size_t /*sz*/);
188 /* --- @mp_createsecure@ --- *
190 * Arguments: @size_t sz@ = size of vector required
192 * Returns: Pointer to pristine new MP structure with enough memory
195 * Use: Creates a new multiprecision integer with indeterminate
196 * contents. The integer has a single reference. The integer's
197 * data space is allocated from the secure arena. Its burn flag
201 extern mp
*mp_createsecure(size_t /*sz*/);
203 /* --- @mp_build@ --- *
205 * Arguments: @mp *m@ = pointer to an MP block to fill in
206 * @mpw *v@ = pointer to a word array
207 * @mpw *vl@ = pointer just past end of array
211 * Use: Creates a multiprecision integer representing some smallish
212 * number. You must provide storage for the number and dispose
213 * of it when you've finished with it. The number is marked as
214 * constant while it exists.
217 extern void mp_build(mp */
*m*/
, mpw */
*v*/
, mpw */
*vl*/
);
219 /* --- @mp_destroy@ --- *
221 * Arguments: @mp *m@ = pointer to a multiprecision integer
225 * Use: Destroys a multiprecision integer. The reference count isn't
226 * checked. Don't use this function if you don't know what
227 * you're doing: use @mp_drop@ instead.
230 extern void mp_destroy(mp */
*m*/
);
232 /* --- @mp_copy@ --- *
234 * Arguments: @mp *m@ = pointer to a multiprecision integer
236 * Returns: A copy of the given multiprecision integer.
238 * Use: Copies the given integer. In fact you just get another
239 * reference to the same old one again.
242 extern mp
*mp_copy(mp */
*m*/
);
244 #define MP_COPY(m) ((m)->ref++, (m))
246 /* --- @mp_drop@ --- *
248 * Arguments: @mp *m@ = pointer to a multiprecision integer
252 * Use: Drops a reference to an integer which isn't wanted any more.
253 * If there are no more references, the integer is destroyed.
256 extern void mp_drop(mp */
*m*/
);
258 #define MP_DROP(m) do { \
261 if (_mm->ref == 0 && !(_mm->f & MP_CONST)) \
265 /* --- @mp_split@ --- *
267 * Arguments: @mp *m@ = pointer to a multiprecision integer
269 * Returns: A reference to the same integer, possibly with a different
272 * Use: Splits off a modifiable version of the integer referred to.
275 extern mp
*mp_split(mp */
*m*/
);
277 #define MP_SPLIT(m) do { \
279 if ((_m->f & MP_CONST) || _m->ref > 1) { \
280 size_t _len = MP_LEN(_m); \
281 mp *_mm = mp_new(_len, _m->f); \
282 if (!(_m->f & MP_UNDEF)) \
283 memcpy(_mm->v, _m->v, MPWS(_len)); \
290 /* --- @mp_resize@ --- *
292 * Arguments: @mp *m@ = pointer to a multiprecision integer
293 * @size_t sz@ = new size
297 * Use: Resizes the vector containing the integer's digits. The new
298 * size must be at least as large as the current integer's
299 * length. This isn't really intended for client use.
302 extern void mp_resize(mp */
*m*/
, size_t /*sz*/);
304 #define MP_RESIZE(m, ssz) do { \
306 size_t _sz = (ssz); \
307 mparena *_a = (_m->f & MP_BURN) ? MPARENA_SECURE : MPARENA_GLOBAL; \
309 size_t _len = MP_LEN(_m); \
310 assert(((void)"can't make size less than length", _sz >= _len)); \
311 _v = mpalloc(_a, _sz); \
312 if (!(_m->f & MP_UNDEF)) \
313 memcpy(_v, _m->v, MPWS(_len)); \
314 if (_m->f & MP_BURN) \
315 memset(_m->v, 0, MPWS(_m->sz)); \
316 mpfree(_m->a, _m->v); \
319 _m->vl = _v + _len; \
322 /* --- @mp_ensure@ --- *
324 * Arguments: @mp *m@ = pointer to a multiprecision integer
325 * @size_t sz@ = required size
329 * Use: Ensures that the integer has enough space for @sz@ digits.
330 * The value is not changed.
333 extern void mp_ensure(mp */
*m*/
, size_t /*sz*/);
335 #define MP_ENSURE(m, ssz) do { \
337 size_t _ssz = (ssz); \
338 size_t _len = MP_LEN(_m); \
339 if (_ssz >= _len) { \
341 mp_resize(_m, _ssz); \
342 if (!(_m->f & MP_UNDEF) && _ssz > _len) \
343 memset(_m->vl, 0, MPWS(_ssz - _len)); \
344 _m->vl = _m->v + _ssz; \
348 /* --- @mp_dest@ --- *
350 * Arguments: @mp *m@ = a suggested destination integer
351 * @size_t sz@ = size required for result, in digits
352 * @unsigned f@ = various flags
354 * Returns: A pointer to an appropriate destination.
356 * Use: Converts a suggested destination into a real destination with
357 * the required properties. If the real destination is @d@,
358 * then the following properties will hold:
360 * * @d@ will have exactly one reference.
362 * * If @m@ is not @MP_NEW@, then the contents of @m@ will not
363 * change, unless @f@ has the @MP_UNDEF@ flag set.
365 * * If @m@ is not @MP_NEW@, then he reference count of @m@ on
366 * entry is equal to the sum of the counts of @d@ and @m@ on
369 * * The size of @d@ will be at least @sz@.
371 * * If @f@ has the @MP_BURN@ flag set, then @d@ will be
372 * allocated from @MPARENA_SECURE@.
374 * Understanding this function is crucial to using Catacomb's
375 * multiprecision integer library effectively.
378 extern mp
*mp_dest(mp */
*m*/
, size_t /*sz*/, unsigned /*f*/);
380 #define MP_DEST(m, ssz, f) do { \
382 size_t _ssz = (ssz); \
384 _m = mp_dest(_m, _ssz, _f); \
388 /*----- Size manipulation -------------------------------------------------*/
390 /* --- @mp_shrink@ --- *
392 * Arguments: @mp *m@ = pointer to a multiprecision integer
396 * Use: Reduces the recorded length of an integer. This doesn't
397 * reduce the amount of memory used, although it can improve
398 * performance a bit. To reduce memory, use @mp_minimize@
399 * instead. This can't change the value of an integer, and is
400 * therefore safe to use even when there are multiple
404 extern void mp_shrink(mp */
*m*/
);
406 #define MP_SHRINK(m) do { \
408 MPX_SHRINK(_mm->v, _mm->vl); \
413 /* --- @mp_minimize@ --- *
415 * Arguments: @mp *m@ = pointer to a multiprecision integer
419 * Use: Reduces the amount of memory an integer uses. It's best to
420 * do this to numbers which aren't going to change in the
424 extern void mp_minimize(mp */
*m*/
);
426 /*----- Bit scanning ------------------------------------------------------*/
428 #ifndef CATACOMB_MPSCAN_H
432 /* --- @mp_scan@ --- *
434 * Arguments: @mpscan *sc@ = pointer to bitscanner block
435 * @const mp *m@ = pointer to a multiprecision integer
439 * Use: Initializes a bitscanner on a multiprecision integer.
442 extern void mp_scan(mpscan */
*sc*/
, const mp */
*m*/
);
444 #define MP_SCAN(sc, m) do { \
445 const mp *_mm = (m); \
446 mpscan *_sc = (sc); \
447 MPSCAN_INITX(_sc, _mm->v, _mm->vl); \
450 /* --- @mp_rscan@ --- *
452 * Arguments: @mpscan *sc@ = pointer to bitscanner block
453 * @const mp *m@ = pointer to a multiprecision integer
457 * Use: Initializes a reverse bitscanner on a multiprecision
461 extern void mp_rscan(mpscan */
*sc*/
, const mp */
*m*/
);
463 #define MP_RSCAN(sc, m) do { \
464 const mp *_mm = (m); \
465 mpscan *_sc = (sc); \
466 MPSCAN_RINITX(_sc, _mm->v, _mm->vl); \
469 /* --- Other bitscanning aliases --- */
471 #define mp_step mpscan_step
472 #define mp_bit mpscan_bit
473 #define mp_rstep mpscan_rstep
474 #define mp_rbit mpscan_rbit
476 #define MP_STEP MPSCAN_STEP
477 #define MP_BIT MPSCAN_BIT
478 #define MP_RSTEP MPSCAN_RSTEP
479 #define MP_RBIT MPSCAN_RBIT
481 /*----- Loading and storing -----------------------------------------------*/
483 /* --- @mp_octets@ --- *
485 * Arguments: @const mp *m@ = a multiprecision integer
487 * Returns: The number of octets required to represent @m@.
489 * Use: Calculates the external storage required for a multiprecision
493 extern size_t mp_octets(const mp */
*m*/
);
495 /* --- @mp_bits@ --- *
497 * Arguments: @const mp *m@ = a multiprecision integer
499 * Returns: The number of bits required to represent @m@.
501 * Use: Calculates the external storage required for a multiprecision
505 extern unsigned long mp_bits(const mp */
*m*/
);
507 /* --- @mp_loadl@ --- *
509 * Arguments: @mp *d@ = destination
510 * @const void *pv@ = pointer to source data
511 * @size_t sz@ = size of the source data
513 * Returns: Resulting multiprecision number.
515 * Use: Loads a multiprecision number from an array of octets. The
516 * first byte in the array is the least significant. More
517 * formally, if the bytes are %$b_0, b_1, \ldots, b_{n-1}$%
518 * then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
521 extern mp
*mp_loadl(mp */
*d*/
, const void */
*pv*/
, size_t /*sz*/);
523 /* --- @mp_storel@ --- *
525 * Arguments: @const mp *m@ = source
526 * @void *pv@ = pointer to output array
527 * @size_t sz@ = size of the output array
531 * Use: Stores a multiprecision number in an array of octets. The
532 * first byte in the array is the least significant. If the
533 * array is too small to represent the number, high-order bits
534 * are truncated; if the array is too large, high order bytes
535 * are filled with zeros. More formally, if the number is
536 * %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
537 * then the array is %$b_0, b_1, \ldots, b_{n-1}$%.
540 extern void mp_storel(const mp */
*m*/
, void */
*pv*/
, size_t /*sz*/);
542 /* --- @mp_loadb@ --- *
544 * Arguments: @mp *d@ = destination
545 * @const void *pv@ = pointer to source data
546 * @size_t sz@ = size of the source data
548 * Returns: Resulting multiprecision number.
550 * Use: Loads a multiprecision number from an array of octets. The
551 * last byte in the array is the least significant. More
552 * formally, if the bytes are %$b_{n-1}, b_{n-2}, \ldots, b_0$%
553 * then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
556 extern mp
*mp_loadb(mp */
*d*/
, const void */
*pv*/
, size_t /*sz*/);
558 /* --- @mp_storeb@ --- *
560 * Arguments: @const mp *m@ = source
561 * @void *pv@ = pointer to output array
562 * @size_t sz@ = size of the output array
566 * Use: Stores a multiprecision number in an array of octets. The
567 * last byte in the array is the least significant. If the
568 * array is too small to represent the number, high-order bits
569 * are truncated; if the array is too large, high order bytes
570 * are filled with zeros. More formally, if the number is
571 * %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
572 * then the array is %$b_{n-1}, b_{n-2}, \ldots, b_0$%.
575 extern void mp_storeb(const mp */
*m*/
, void */
*pv*/
, size_t /*sz*/);
577 /*----- Simple arithmetic -------------------------------------------------*/
581 * Arguments: @mp *d@ = destination
584 * Returns: Result, @a@ converted to two's complement notation.
587 extern mp
*mp_2c(mp */
*d*/
, mp */
*a*/
);
591 * Arguments: @mp *d@ = destination
594 * Returns: Result, @a@ converted to the native signed-magnitude
598 extern mp
*mp_sm(mp */
*d*/
, mp */
*a*/
);
600 /* --- @mp_lsl@ --- *
602 * Arguments: @mp *d@ = destination
604 * @size_t n@ = number of bits to move
606 * Returns: Result, @a@ shifted left by @n@.
609 extern mp
*mp_lsl(mp */
*d*/
, mp */
*a*/
, size_t /*n*/);
611 /* --- @mp_lsr@ --- *
613 * Arguments: @mp *d@ = destination
615 * @size_t n@ = number of bits to move
617 * Returns: Result, @a@ shifted left by @n@.
620 extern mp
*mp_lsr(mp */
*d*/
, mp */
*a*/
, size_t /*n*/);
624 * Arguments: @const mp *a, *b@ = two numbers
626 * Returns: Nonzero if the numbers are equal.
629 extern int mp_eq(const mp */
*a*/
, const mp */
*b*/
);
631 #define MP_EQ(a, b) \
632 ((((a)->f ^ (b)->f) & MP_NEG) == 0 && \
633 mpx_ueq((a)->v, (a)->vl, (b)->v, (b)->vl))
635 /* --- @mp_cmp@ --- *
637 * Arguments: @const mp *a, *b@ = two numbers
639 * Returns: Less than, equal to or greater than zero, according to
640 * whether @a@ is less than, equal to or greater than @b@.
643 extern int mp_cmp(const mp */
*a*/
, const mp */
*b*/
);
645 #define MP_CMP(a, op, b) (mp_cmp((a), (b)) op 0)
647 /* --- @mpx_and@, @mpx_or@, @mpx_xor@, @mpx_not@ --- *
649 * Arguments: @mp *d@ = destination
650 * @mp *a, *b@ = sources
652 * Returns: The result of the obvious bitwise operation.
655 extern mp
*mp_and(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
656 extern mp
*mp_or(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
657 extern mp
*mp_xor(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
658 extern mp
*mp_not(mp */
*d*/
, mp */
*a*/
);
660 /* --- @mp_add@ --- *
662 * Arguments: @mp *d@ = destination
663 * @mp *a, *b@ = sources
665 * Returns: Result, @a@ added to @b@.
668 extern mp
*mp_add(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
670 /* --- @mp_sub@ --- *
672 * Arguments: @mp *d@ = destination
673 * @mp *a, *b@ = sources
675 * Returns: Result, @b@ subtracted from @a@.
678 extern mp
*mp_sub(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
680 /* --- @mp_mul@ --- *
682 * Arguments: @mp *d@ = destination
683 * @mp *a, *b@ = sources
685 * Returns: Result, @a@ multiplied by @b@.
688 extern mp
*mp_mul(mp */
*d*/
, mp */
*a*/
, mp */
*b*/
);
690 /* --- @mp_sqr@ --- *
692 * Arguments: @mp *d@ = destination
695 * Returns: Result, @a@ squared.
698 extern mp
*mp_sqr(mp */
*d*/
, mp */
*a*/
);
700 /* --- @mp_div@ --- *
702 * Arguments: @mp **qq, **rr@ = destination, quotient and remainder
703 * @mp *a, *b@ = sources
705 * Use: Calculates the quotient and remainder when @a@ is divided by
709 extern void mp_div(mp
**/
*qq*/
, mp
**/
*rr*/
, mp */
*a*/
, mp */
*b*/
);
711 /* --- @mp_odd@ --- *
713 * Arguments: @mp *d@ = pointer to destination integer
714 * @mp *m@ = pointer to source integer
715 * @size_t *s@ = where to store the power of 2
717 * Returns: An odd integer integer %$t$% such that %$m = 2^s t$%.
719 * Use: Computes a power of two and an odd integer which, when
720 * multiplied, give a specified result. This sort of thing is
721 * useful in number theory quite often.
724 extern mp
*mp_odd(mp */
*d*/
, mp */
*m*/
, size_t */
*s*/
);
726 /*----- More advanced algorithms ------------------------------------------*/
728 /* --- @mp_sqrt@ --- *
730 * Arguments: @mp *d@ = pointer to destination integer
731 * @mp *a@ = (nonnegative) integer to take square root of
733 * Returns: The largest integer %$x$% such that %$x^2 \le a$%.
735 * Use: Computes integer square roots.
737 * The current implementation isn't very good: it uses the
738 * Newton-Raphson method to find an approximation to %$a$%. If
739 * there's any demand for a better version, I'll write one.
742 extern mp
*mp_sqrt(mp */
*d*/
, mp */
*a*/
);
744 /* --- @mp_gcd@ --- *
746 * Arguments: @mp **gcd, **xx, **yy@ = where to write the results
747 * @mp *a, *b@ = sources (must be nonzero)
751 * Use: Calculates @gcd(a, b)@, and two numbers @x@ and @y@ such that
752 * @ax + by = gcd(a, b)@. This is useful for computing modular
753 * inverses. Neither @a@ nor @b@ may be zero.
756 extern void mp_gcd(mp
**/
*gcd*/
, mp
**/
*xx*/
, mp
**/
*yy*/
,
757 mp */
*a*/
, mp */
*b*/
);
759 /* --- @mp_jacobi@ --- *
761 * Arguments: @mp *a@ = an integer less than @n@
762 * @mp *n@ = an odd integer
764 * Returns: @-1@, @0@ or @1@ -- the Jacobi symbol %$J(a, n)$%.
766 * Use: Computes the Jacobi symbol. If @n@ is prime, this is the
767 * Legendre symbol and is equal to 1 if and only if @a@ is a
768 * quadratic residue mod @n@. The result is zero if and only if
769 * @a@ and @n@ have a common factor greater than one.
772 extern int mp_jacobi(mp */
*a*/
, mp */
*n*/
);
774 /* --- @mp_modsqrt@ --- *
776 * Arguments: @mp *d@ = destination integer
777 * @mp *a@ = source integer
778 * @mp *p@ = modulus (must be prime)
780 * Returns: If %$a$% is a quadratic residue, a square root of %$a$%; else
783 * Use: Returns an integer %$x$% such that %$x^2 \equiv a \pmod{p}$%,
784 * if one exists; else a null pointer. This function will not
785 * work if %$p$% is composite: you must factor the modulus, take
786 * a square root mod each factor, and recombine the results
787 * using the Chinese Remainder Theorem.
790 extern mp
*mp_modsqrt(mp */
*d*/
, mp */
*a*/
, mp */
*p*/
);
792 /*----- Test harness support ----------------------------------------------*/
794 #include <mLib/testrig.h>
796 #ifndef CATACOMB_MPTEXT_H
800 extern const test_type type_mp
;
802 /*----- That's all, folks -------------------------------------------------*/