Generic interface.
[u/mdw/catacomb] / mp.h
1 /* -*-c-*-
2 *
3 * $Id: mp.h,v 1.5 1999/11/22 20:50:37 mdw Exp $
4 *
5 * Simple 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: mp.h,v $
33 * Revision 1.5 1999/11/22 20:50:37 mdw
34 * Add support for computing Jacobi symbols.
35 *
36 * Revision 1.4 1999/11/21 22:13:02 mdw
37 * Add mp version of MPX_BITS.
38 *
39 * Revision 1.3 1999/11/19 13:19:14 mdw
40 * Fix const annotation.
41 *
42 * Revision 1.2 1999/11/17 18:02:16 mdw
43 * New multiprecision integer arithmetic suite.
44 *
45 */
46
47 #ifndef MP_H
48 #define MP_H
49
50 #ifdef __cplusplus
51 extern "C" {
52 #endif
53
54 /*----- Header files ------------------------------------------------------*/
55
56 #include <assert.h>
57 #include <string.h>
58
59 #include <mLib/sub.h>
60
61 #ifndef MPW_H
62 # include "mpw.h"
63 #endif
64
65 #ifndef MPX_H
66 # include "mpx.h"
67 #endif
68
69 /*----- Data structures ---------------------------------------------------*/
70
71 typedef struct mp {
72 mpw *v, *vl;
73 size_t sz;
74 unsigned f;
75 unsigned ref;
76 } mp;
77
78 #define MP_NEG 1u
79 #define MP_BURN 2u
80 #define MP_CONST 4u
81 #define MP_UNDEF 8u
82
83 /*----- Useful constants --------------------------------------------------*/
84
85 extern mp mp_const[];
86
87 #define MP_ZERO (&mp_const[0])
88 #define MP_ONE (&mp_const[1])
89 #define MP_TWO (&mp_const[2])
90 #define MP_THREE (&mp_const[3])
91 #define MP_FOUR (&mp_const[4])
92 #define MP_FIVE (&mp_const[5])
93 #define MP_TEN (&mp_const[6])
94 #define MP_MONE (&mp_const[7])
95
96 #define MP_NEW ((mp *)0)
97
98 /*----- Memory allocation hooks -------------------------------------------*/
99
100 #ifndef MPARENA_H
101 # include "mparena.h"
102 #endif
103
104 /* --- @MP_ARENA@ --- *
105 *
106 * This selects where memory is allocated from. Tweak to use more fancy
107 * things like custom arenas.
108 */
109
110 #ifndef MP_ARENA
111 # define MP_ARENA MPARENA_GLOBAL
112 #endif
113
114 /* --- @MP_ALLOC@ --- *
115 *
116 * Arguments: @size_t sz@ = size required
117 *
118 * Returns: Pointer to an allocated vector of the requested size.
119 *
120 * Use: Hook for vector allocation.
121 */
122
123 #ifndef MP_ALLOC
124 # define MP_ALLOC(sz) mpalloc(MP_ARENA, (sz))
125 #endif
126
127 /* --- @MP_FREE@ --- *
128 *
129 * Arguments: @mpw *v@ = pointer to vector
130 *
131 * Returns: ---
132 *
133 * Use: Hook for vector deallocation.
134 */
135
136 #ifndef MP_FREE
137 # define MP_FREE(v) mpfree(MP_ARENA, (v))
138 #endif
139
140 /*----- Paranoia management -----------------------------------------------*/
141
142 /* --- @mp_burn@ --- *
143 *
144 * Arguments: @mp *m@ = pointer to a multiprecision integer
145 *
146 * Returns: ---
147 *
148 * Use: Marks the integer as `burn-after-use'. When the integer's
149 * memory is deallocated, it is deleted so that traces can't
150 * remain in the swap file. In theory.
151 */
152
153 extern void mp_burn(mp */*m*/);
154
155 /*----- Trivial macros ----------------------------------------------------*/
156
157 /* --- @MP_LEN@ --- *
158 *
159 * Arguments: @mp *m@ = pointer to a multiprecision integer
160 *
161 * Returns: Length of the integer, in words.
162 */
163
164 #define MP_LEN(m) ((m)->vl - ((m)->v))
165
166 /*----- Memory management and reference counting --------------------------*/
167
168 /* --- @mp_create@ --- *
169 *
170 * Arguments: @size_t sz@ = size of vector required
171 *
172 * Returns: Pointer to pristine new MP structure with enough memory
173 * bolted onto it.
174 *
175 * Use: Creates a new multiprecision integer with indeterminate
176 * contents. The integer has a single reference.
177 */
178
179 extern mp *mp_create(size_t /*sz*/);
180
181 /* --- @mp_build@ --- *
182 *
183 * Arguments: @mp *m@ = pointer to an MP block to fill in
184 * @mpw *v@ = pointer to a word array
185 * @mpw *vl@ = pointer just past end of array
186 *
187 * Returns: ---
188 *
189 * Use: Creates a multiprecision integer representing some smallish
190 * number. You must provide storage for the number and dispose
191 * of it when you've finished with it. The number is marked as
192 * constant while it exists.
193 */
194
195 extern void mp_build(mp */*m*/, mpw */*v*/, mpw */*vl*/);
196
197 /* --- @mp_destroy@ --- *
198 *
199 * Arguments: @mp *m@ = pointer to a multiprecision integer
200 *
201 * Returns: ---
202 *
203 * Use: Destroys a multiprecision integer. The reference count isn't
204 * checked. Don't use this function if you don't know what
205 * you're doing: use @mp_drop@ instead.
206 */
207
208 extern void mp_destroy(mp */*m*/);
209
210 /* --- @mp_copy@ --- *
211 *
212 * Arguments: @mp *m@ = pointer to a multiprecision integer
213 *
214 * Returns: A copy of the given multiprecision integer.
215 *
216 * Use: Copies the given integer. In fact you just get another
217 * reference to the same old one again.
218 */
219
220 extern mp *mp_copy(mp */*m*/);
221
222 #define MP_COPY(m) ((m)->ref++, (m))
223
224 /* --- @mp_drop@ --- *
225 *
226 * Arguments: @mp *m@ = pointer to a multiprecision integer
227 *
228 * Returns: ---
229 *
230 * Use: Drops a reference to an integer which isn't wanted any more.
231 * If there are no more references, the integer is destroyed.
232 */
233
234 extern void mp_drop(mp */*m*/);
235
236 #define MP_DROP(m) do { \
237 mp *_mm = (m); \
238 if (_mm->ref > 1) \
239 _mm->ref--; \
240 else if (!(_mm->f & MP_CONST)) \
241 mp_destroy(_mm); \
242 } while (0)
243
244 /* --- @mp_split@ --- *
245 *
246 * Arguments: @mp *m@ = pointer to a multiprecision integer
247 *
248 * Returns: A reference to the same integer, possibly with a different
249 * address.
250 *
251 * Use: Splits off a modifiable version of the integer referred to.
252 */
253
254 extern mp *mp_split(mp */*m*/);
255
256 #define MP_SPLIT(m) do { \
257 mp *_mm = (m); \
258 if ((_mm->f & MP_CONST) || _mm->ref != 1) { \
259 mp *_dd = mp_create(_mm->sz); \
260 _dd->vl = _dd->v + MP_LEN(_mm); \
261 _dd->f = _mm->f & (MP_NEG | MP_BURN); \
262 memcpy(_dd->v, _mm->v, MPWS(MP_LEN(_mm))); \
263 _dd->ref = 1; \
264 _mm->ref--; \
265 (m) = _dd; \
266 } \
267 } while (0)
268
269 /* --- @mp_resize@ --- *
270 *
271 * Arguments: @mp *m@ = pointer to a multiprecision integer
272 * @size_t sz@ = new size
273 *
274 * Returns: ---
275 *
276 * Use: Resizes the vector containing the integer's digits. The new
277 * size must be at least as large as the current integer's
278 * length. The integer's length is increased and new digits are
279 * filled with zeroes. This isn't really intended for client
280 * use.
281 */
282
283 extern void mp_resize(mp */*m*/, size_t /*sz*/);
284
285 #define MP_RESIZE(m, ssz) do { \
286 mp *_m = (m); \
287 size_t _sz = (ssz); \
288 size_t _len = MP_LEN(_m); \
289 mpw *_v = MP_ALLOC(_sz); \
290 memcpy(_v, _m->v, MPWS(_len)); \
291 if (_m->f & MP_BURN) \
292 memset(_m->v, 0, MPWS(_m->sz)); \
293 MP_FREE(_m->v); \
294 _m->v = _v; \
295 _m->vl = _v + _len; \
296 _m->sz = _sz; \
297 } while (0)
298
299 /* --- @mp_ensure@ --- *
300 *
301 * Arguments: @mp *m@ = pointer to a multiprecision integer
302 * @size_t sz@ = required size
303 *
304 * Returns: ---
305 *
306 * Use: Ensures that the integer has enough space for @sz@ digits.
307 * The value is not changed.
308 */
309
310 extern void mp_ensure(mp */*m*/, size_t /*sz*/);
311
312 #define MP_ENSURE(m, ssz) do { \
313 mp *_mm = (m); \
314 size_t _ssz = (ssz); \
315 size_t _len = MP_LEN(_mm); \
316 if (_ssz > _mm->sz) \
317 MP_RESIZE(_mm, _ssz); \
318 if (!(_mm->f & MP_UNDEF) && _ssz > _len) { \
319 memset(_mm->vl, 0, MPWS(_ssz - _len)); \
320 _mm->vl = _mm->v + _ssz; \
321 } \
322 } while (0)
323
324 /* --- @mp_modify@ --- *
325 *
326 * Arguments: @mp *m@ = pointer to a multiprecision integer
327 * @size_t sz@ = size required
328 *
329 * Returns: Pointer to the integer (possibly different).
330 *
331 * Use: Prepares an integer to be overwritten. It's split off from
332 * other references to the same integer, and sufficient space is
333 * allocated.
334 */
335
336 extern mp *mp_modify(mp */*m*/, size_t /*sz*/);
337
338 #define MP_MODIFY(m, sz) do { \
339 size_t _rq = (sz); \
340 mp *_m = (m); \
341 if (_m == MP_NEW) \
342 _m = mp_create(_rq); \
343 else { \
344 MP_SPLIT(_m); \
345 MP_ENSURE(_m, _rq); \
346 } \
347 _m->vl = _m->v + _rq; \
348 (m) = _m; \
349 } while (0)
350
351 /*----- Size manipulation -------------------------------------------------*/
352
353 /* --- @mp_shrink@ --- *
354 *
355 * Arguments: @mp *m@ = pointer to a multiprecision integer
356 *
357 * Returns: ---
358 *
359 * Use: Reduces the recorded length of an integer. This doesn't
360 * reduce the amount of memory used, although it can improve
361 * performance a bit. To reduce memory, use @mp_minimize@
362 * instead. This can't change the value of an integer, and is
363 * therefore safe to use even when there are multiple
364 * references.
365 */
366
367 extern void mp_shrink(mp */*m*/);
368
369 #define MP_SHRINK(m) do { \
370 mp *_mm = (m); \
371 MPX_SHRINK(_mm->v, _mm->vl); \
372 if (!MP_LEN(_mm)) \
373 _mm->f &= ~MP_NEG; \
374 } while (0)
375
376 /* --- @mp_minimize@ --- *
377 *
378 * Arguments: @mp *m@ = pointer to a multiprecision integer
379 *
380 * Returns: ---
381 *
382 * Use: Reduces the amount of memory an integer uses. It's best to
383 * do this to numbers which aren't going to change in the
384 * future.
385 */
386
387 extern void mp_minimize(mp */*m*/);
388
389 /*----- Bit scanning ------------------------------------------------------*/
390
391 #ifndef MPSCAN_H
392 # include "mpscan.h"
393 #endif
394
395 /* --- @mp_scan@ --- *
396 *
397 * Arguments: @mpscan *sc@ = pointer to bitscanner block
398 * @const mp *m@ = pointer to a multiprecision integer
399 *
400 * Returns: ---
401 *
402 * Use: Initializes a bitscanner on a multiprecision integer.
403 */
404
405 extern void mp_scan(mpscan */*sc*/, const mp */*m*/);
406
407 #define MP_SCAN(sc, m) do { \
408 const mp *_mm = (m); \
409 mpscan *_sc = (sc); \
410 MPSCAN_INITX(_sc, _mm->v, _mm->vl); \
411 } while (0)
412
413 /* --- Other bitscanning aliases --- */
414
415 #define mp_step mpscan_step
416 #define mp_bit mpscan_bit
417
418 #define MP_STEP MPSCAN_STEP
419 #define MP_BIT MPSCAN_BIT
420
421 /*----- Loading and storing -----------------------------------------------*/
422
423 /* --- @mp_octets@ --- *
424 *
425 * Arguments: @const mp *m@ = a multiprecision integer
426 *
427 * Returns: The number of octets required to represent @m@.
428 *
429 * Use: Calculates the external storage required for a multiprecision
430 * integer.
431 */
432
433 extern size_t mp_octets(const mp */*m*/);
434
435 /* --- @mp_bits@ --- *
436 *
437 * Arguments: @const mp *m@ = a multiprecision integer
438 *
439 * Returns: The number of bits required to represent @m@.
440 *
441 * Use: Calculates the external storage required for a multiprecision
442 * integer.
443 */
444
445 extern unsigned long mp_bits(const mp */*m*/);
446
447 /* --- @mp_loadl@ --- *
448 *
449 * Arguments: @mp *d@ = destination
450 * @const void *pv@ = pointer to source data
451 * @size_t sz@ = size of the source data
452 *
453 * Returns: Resulting multiprecision number.
454 *
455 * Use: Loads a multiprecision number from an array of octets. The
456 * first byte in the array is the least significant. More
457 * formally, if the bytes are %$b_0, b_1, \ldots, b_{n-1}$%
458 * then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
459 */
460
461 extern mp *mp_loadl(mp */*d*/, const void */*pv*/, size_t /*sz*/);
462
463 /* --- @mp_storel@ --- *
464 *
465 * Arguments: @const mp *m@ = source
466 * @void *pv@ = pointer to output array
467 * @size_t sz@ = size of the output array
468 *
469 * Returns: ---
470 *
471 * Use: Stores a multiprecision number in an array of octets. The
472 * first byte in the array is the least significant. If the
473 * array is too small to represent the number, high-order bits
474 * are truncated; if the array is too large, high order bytes
475 * are filled with zeros. More formally, if the number is
476 * %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
477 * then the array is %$b_0, b_1, \ldots, b_{n-1}$%.
478 */
479
480 extern void mp_storel(const mp */*m*/, void */*pv*/, size_t /*sz*/);
481
482 /* --- @mp_loadb@ --- *
483 *
484 * Arguments: @mp *d@ = destination
485 * @const void *pv@ = pointer to source data
486 * @size_t sz@ = size of the source data
487 *
488 * Returns: Resulting multiprecision number.
489 *
490 * Use: Loads a multiprecision number from an array of octets. The
491 * last byte in the array is the least significant. More
492 * formally, if the bytes are %$b_{n-1}, b_{n-2}, \ldots, b_0$%
493 * then the result is %$N = \sum_{0 \le i < n} b_i 2^{8i}$%.
494 */
495
496 extern mp *mp_loadb(mp */*d*/, const void */*pv*/, size_t /*sz*/);
497
498 /* --- @mp_storeb@ --- *
499 *
500 * Arguments: @const mp *m@ = source
501 * @void *pv@ = pointer to output array
502 * @size_t sz@ = size of the output array
503 *
504 * Returns: ---
505 *
506 * Use: Stores a multiprecision number in an array of octets. The
507 * last byte in the array is the least significant. If the
508 * array is too small to represent the number, high-order bits
509 * are truncated; if the array is too large, high order bytes
510 * are filled with zeros. More formally, if the number is
511 * %$N = \sum{0 \le i} b_i 2^{8i}$% where %$0 \le b_i < 256$%,
512 * then the array is %$b_{n-1}, b_{n-2}, \ldots, b_0$%.
513 */
514
515 extern void mp_storeb(const mp */*m*/, void */*pv*/, size_t /*sz*/);
516
517 /*----- Simple arithmetic -------------------------------------------------*/
518
519 /* --- @mp_2c@ --- *
520 *
521 * Arguments: @mp *d@ = destination
522 * @mp *a@ = source
523 *
524 * Returns: Result, @a@ converted to two's complement notation.
525 */
526
527 extern mp *mp_2c(mp */*d*/, mp */*a*/);
528
529 /* --- @mp_sm@ --- *
530 *
531 * Arguments: @mp *d@ = destination
532 * @mp *a@ = source
533 *
534 * Returns: Result, @a@ converted to the native signed-magnitude
535 * notation.
536 */
537
538 extern mp *mp_sm(mp */*d*/, mp */*a*/);
539
540 /* --- @mp_lsl@ --- *
541 *
542 * Arguments: @mp *d@ = destination
543 * @const mp *a@ = source
544 * @size_t n@ = number of bits to move
545 *
546 * Returns: Result, @a@ shifted left by @n@.
547 */
548
549 extern mp *mp_lsl(mp */*d*/, const mp */*a*/, size_t /*n*/);
550
551 /* --- @mp_lsr@ --- *
552 *
553 * Arguments: @mp *d@ = destination
554 * @const mp *a@ = source
555 * @size_t n@ = number of bits to move
556 *
557 * Returns: Result, @a@ shifted left by @n@.
558 */
559
560 extern mp *mp_lsr(mp */*d*/, const mp */*a*/, size_t /*n*/);
561
562 /* --- @mp_cmp@ --- *
563 *
564 * Arguments: @const mp *a, *b@ = two numbers
565 *
566 * Returns: Less than, equal to or greater than zero, according to
567 * whether @a@ is less than, equal to or greater than @b@.
568 */
569
570 extern int mp_cmp(const mp */*a*/, const mp */*b*/);
571
572 #define MP_CMP(a, op, b) (mp_cmp((a), (b)) op 0)
573
574 /* --- @mp_add@ --- *
575 *
576 * Arguments: @mp *d@ = destination
577 * @const mp *a, *b@ = sources
578 *
579 * Returns: Result, @a@ added to @b@.
580 */
581
582 extern mp *mp_add(mp */*d*/, const mp */*a*/, const mp */*b*/);
583
584 /* --- @mp_sub@ --- *
585 *
586 * Arguments: @mp *d@ = destination
587 * @const mp *a, *b@ = sources
588 *
589 * Returns: Result, @b@ subtracted from @a@.
590 */
591
592 extern mp *mp_sub(mp */*d*/, const mp */*a*/, const mp */*b*/);
593
594 /* --- @mp_mul@ --- *
595 *
596 * Arguments: @mp *d@ = destination
597 * @const mp *a, *b@ = sources
598 *
599 * Returns: Result, @a@ multiplied by @b@.
600 */
601
602 extern mp *mp_mul(mp */*d*/, const mp */*a*/, const mp */*b*/);
603
604 /* --- @mp_sqr@ --- *
605 *
606 * Arguments: @mp *d@ = destination
607 * @const mp *a@ = source
608 *
609 * Returns: Result, @a@ squared.
610 */
611
612 extern mp *mp_sqr(mp */*d*/, const mp */*a*/);
613
614 /* --- @mp_div@ --- *
615 *
616 * Arguments: @mp **qq, **rr@ = destination, quotient and remainder
617 * @const mp *a, *b@ = sources
618 *
619 * Use: Calculates the quotient and remainder when @a@ is divided by
620 * @b@.
621 */
622
623 extern void mp_div(mp **/*qq*/, mp **/*rr*/,
624 const mp */*a*/, const mp */*b*/);
625
626 /*----- More advanced algorithms ------------------------------------------*/
627
628 /* --- @mp_gcd@ --- *
629 *
630 * Arguments: @mp **gcd, **xx, **yy@ = where to write the results
631 * @mp *a, *b@ = sources (must be nonzero)
632 *
633 * Returns: ---
634 *
635 * Use: Calculates @gcd(a, b)@, and two numbers @x@ and @y@ such that
636 * @ax + by = gcd(a, b)@. This is useful for computing modular
637 * inverses. Neither @a@ nor @b@ may be zero. Note that,
638 * unlike @mp_div@ for example, it is not possible to specify
639 * explicit destinations -- new MPs are always allocated.
640 */
641
642 extern void mp_gcd(mp **/*gcd*/, mp **/*xx*/, mp **/*yy*/,
643 mp */*a*/, mp */*b*/);
644
645 /* --- @mp_jacobi@ --- *
646 *
647 * Arguments: @mp *a@ = an integer less than @n@
648 * @mp *n@ = an odd integer
649 *
650 * Returns: @-1@, @0@ or @1@ -- the Jacobi symbol %$J(a, n)$%.
651 *
652 * Use: Computes the Jacobi symbol. If @n@ is prime, this is the
653 * Legendre symbol and is equal to 1 if and only if @a@ is a
654 * quadratic residue mod @n@. The result is zero if and only if
655 * @a@ and @n@ have a common factor greater than one.
656 */
657
658 int mp_jacobi(mp */*a*/, mp */*n*/);
659
660 /*----- Test harness support ----------------------------------------------*/
661
662 #include <mLib/testrig.h>
663
664 #ifndef MPTEXT_H
665 # include "mptext.h"
666 #endif
667
668 extern const test_type type_mp;
669
670 /*----- That's all, folks -------------------------------------------------*/
671
672 #ifdef __cplusplus
673 }
674 #endif
675
676 #endif