projects / u / mdw / catacomb / blob
? search:


bc4b356ce58a74c4dc01db60dfeb20fecdcd0580
[u/mdw/catacomb] / mpx.h
1 /* -*-c-*-
2 *
3 * $Id: mpx.h,v 1.14 2002/10/09 00:36:03 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.14 2002/10/09 00:36:03 mdw
34 * Fix bounds on workspace for Karatsuba operations.
35 *
36 * Revision 1.13 2002/10/06 22:52:50 mdw
37 * Pile of changes for supporting two's complement properly.
38 *
39 * Revision 1.12 2001/04/03 19:36:05 mdw
40 * Add some simple bitwise operations so that Perl can use them.
41 *
42 * Revision 1.11 2000/10/08 15:48:35 mdw
43 * Rename Karatsuba constants now that we have @gfx_kmul@ too.
44 *
45 * Revision 1.10 2000/10/08 12:06:12 mdw
46 * Provide @mpx_ueq@ for rapidly testing equality of two integers.
47 *
48 * Revision 1.9 1999/12/22 15:49:07 mdw
49 * New function for division by a small integer.
50 *
51 * Revision 1.8 1999/12/11 10:57:43 mdw
52 * Karatsuba squaring algorithm.
53 *
54 * Revision 1.7 1999/12/11 01:51:28 mdw
55 * Change Karatsuba parameters slightly.
56 *
57 * Revision 1.6 1999/12/10 23:23:51 mdw
58 * Karatsuba-Ofman multiplication algorithm.
59 *
60 * Revision 1.5 1999/11/20 22:23:27 mdw
61 * Add function versions of some low-level macros with wider use.
62 *
63 * Revision 1.4 1999/11/17 18:04:43 mdw
64 * Add two's complement support. Fix a bug in MPX_UMLAN.
65 *
66 * Revision 1.3 1999/11/13 01:51:29 mdw
67 * Minor interface changes. Should be stable now.
68 *
69 * Revision 1.2 1999/11/11 17:47:55 mdw
70 * Minor changes for different `mptypes.h' format.
71 *
72 * Revision 1.1 1999/09/03 08:41:12 mdw
73 * Initial import.
74 *
75 */
76
77 #ifndef CATACOMB_MPX_H
78 #define CATACOMB_MPX_H
79
80 #ifdef __cplusplus
81 extern "C" {
82 #endif
83
84 /*----- The idea ----------------------------------------------------------*
85 *
86 * This file provides functions and macros which work on vectors of words as
87 * unsigned multiprecision integers. The interface works in terms of base
88 * and limit pointers (i.e., a pointer to the start of a vector, and a
89 * pointer just past its end) rather than base pointer and length, because
90 * that requires more arithmetic and state to work on.
91 *
92 * The interfaces are slightly bizarre in other ways. Try to use the
93 * higher-level functions where you can: they're rather better designed to
94 * actually be friendly and useful.
95 */
96
97 /*----- Header files ------------------------------------------------------*/
98
99 #include <string.h>
100
101 #ifndef CATACOMB_MPW_H
102 # include "mpw.h"
103 #endif
104
105 /*----- General manipulation ----------------------------------------------*/
106
107 /* --- @MPX_SHRINK@ --- *
108 *
109 * Arguments: @const mpw *v@ = pointer to vector of words
110 * @const mpw *vl@ = (updated) current limit of vector
111 *
112 * Use: Shrinks down the limit of a multiprecision integer vector.
113 */
114
115 #define MPX_SHRINK(v, vl) do { \
116 const mpw *_vv = (v), *_vvl = (vl); \
117 while (_vvl > _vv && !_vvl[-1]) \
118 _vvl--; \
119 (vl) = (mpw *)_vvl; \
120 } while (0)
121
122 /* --- @MPX_BITS@ --- *
123 *
124 * Arguments: @unsigned long b@ = result variable
125 * @const mpw *v@ = pointer to array of words
126 * @const mpw *vl@ = limit of vector (from @MPX_SHRINK@)
127 *
128 * Use: Calculates the number of bits in a multiprecision value.
129 */
130
131 #define MPX_BITS(b, v, vl) do { \
132 const mpw *_v = (v), *_vl = (vl); \
133 MPX_SHRINK(_v, _vl); \
134 if (_v == _vl) \
135 (b) = 0; \
136 else { \
137 unsigned long _b = MPW_BITS * (_vl - _v - 1) + 1; \
138 mpw _w = _vl[-1]; \
139 unsigned _k = MPW_BITS / 2; \
140 while (_k) { \
141 if (_w >> _k) { \
142 _w >>= _k; \
143 _b += _k; \
144 } \
145 _k >>= 1; \
146 } \
147 (b) = _b; \
148 } \
149 } while (0)
150
151 /* --- @MPX_OCTETS@ --- *
152 *
153 * Arguments: @size_t o@ = result variable
154 * @const mpw *v, *vl@ = pointer to array of words
155 *
156 * Use: Calculates the number of octets in a multiprecision value.
157 */
158
159 #define MPX_OCTETS(o, v, vl) do { \
160 unsigned long _bb; \
161 MPX_BITS(_bb, (v), (vl)); \
162 (o) = (_bb + 7) >> 3; \
163 } while (0)
164
165 /* --- @MPX_OCTETS2C@ --- *
166 *
167 * Arguments: @size_t o@ = result variable
168 * @const mpw *v, *vl@ = pointer to array of words
169 *
170 * Use: Calculates the number of octets in a multiprecision value, if
171 * you represent it as two's complement.
172 */
173
174 #define MPX_OCTETS2C(o, v, vl) do { \
175 unsigned long _bb; \
176 MPX_BITS(_bb, (v), (vl)); \
177 (o) = (_bb >> 3) + 1; \
178 } while (0)
179
180 /* --- @MPX_COPY@ --- *
181 *
182 * Arguments: @dv, dvl@ = destination vector base and limit
183 * @av, avl@ = source vector base and limit
184 *
185 * Use: Copies a multiprecision integer.
186 */
187
188 #define MPX_COPY(dv, dvl, av, avl) do { \
189 mpw *_dv = (dv), *_dvl = (dvl); \
190 size_t _dn = _dvl - _dv; \
191 const mpw *_av = (av), *_avl = (avl); \
192 size_t _an = _avl - _av; \
193 if (_av == _dv) { \
194 if (_dvl > _avl) \
195 memset(_dv, 0, MPWS(_dn - _an)); \
196 } else if (_an >= _dn) \
197 memmove(_dv, _av, MPWS(_dn)); \
198 else { \
199 memmove(_dv, _av, MPWS(_an)); \
200 memset(_dv + _an, 0, MPWS(_dn - _an)); \
201 } \
202 } while (0)
203
204 /* --- @MPX_ZERO@ --- *
205 *
206 * Arguments: @v, vl@ = base and limit of vector to clear
207 *
208 * Use: Zeroes the area between the two vector pointers.
209 */
210
211 #define MPX_ZERO(v, vl) do { \
212 mpw *_v = (v), *_vl = (vl); \
213 if (_v < _vl) \
214 memset(_v, 0, MPWS(_vl - _v)); \
215 } while (0)
216
217 /*----- Loading and storing -----------------------------------------------*/
218
219 /* --- @mpx_storel@ --- *
220 *
221 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
222 * @void *p@ = pointer to octet array
223 * @size_t sz@ = size of octet array
224 *
225 * Returns: ---
226 *
227 * Use: Stores an MP in an octet array, least significant octet
228 * first. High-end octets are silently discarded if there
229 * isn't enough space for them.
230 */
231
232 extern void mpx_storel(const mpw */*v*/, const mpw */*vl*/,
233 void */*p*/, size_t /*sz*/);
234
235 /* --- @mpx_loadl@ --- *
236 *
237 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
238 * @const void *p@ = pointer to octet array
239 * @size_t sz@ = size of octet array
240 *
241 * Returns: ---
242 *
243 * Use: Loads an MP in an octet array, least significant octet
244 * first. High-end octets are ignored if there isn't enough
245 * space for them.
246 */
247
248 extern void mpx_loadl(mpw */*v*/, mpw */*vl*/,
249 const void */*p*/, size_t /*sz*/);
250
251 /* --- @mpx_storeb@ --- *
252 *
253 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
254 * @void *p@ = pointer to octet array
255 * @size_t sz@ = size of octet array
256 *
257 * Returns: ---
258 *
259 * Use: Stores an MP in an octet array, most significant octet
260 * first. High-end octets are silently discarded if there
261 * isn't enough space for them.
262 */
263
264 extern void mpx_storeb(const mpw */*v*/, const mpw */*vl*/,
265 void */*p*/, size_t /*sz*/);
266
267 /* --- @mpx_loadb@ --- *
268 *
269 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
270 * @const void *p@ = pointer to octet array
271 * @size_t sz@ = size of octet array
272 *
273 * Returns: ---
274 *
275 * Use: Loads an MP in an octet array, most significant octet
276 * first. High-end octets are ignored if there isn't enough
277 * space for them.
278 */
279
280 extern void mpx_loadb(mpw */*v*/, mpw */*vl*/,
281 const void */*p*/, size_t /*sz*/);
282
283 /* --- @mpx_storel2cn@ --- *
284 *
285 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
286 * @void *pp@ = pointer to octet array
287 * @size_t sz@ = size of octet array
288 *
289 * Returns: ---
290 *
291 * Use: Stores a negative MP in an octet array, least significant
292 * octet first, as two's complement. High-end octets are
293 * silently discarded if there isn't enough space for them.
294 * This obviously makes the output bad.
295 */
296
297 extern void mpx_storel2cn(const mpw */*v*/, const mpw */*vl*/,
298 void */*p*/, size_t /*sz*/);
299
300 /* --- @mpx_loadl2cn@ --- *
301 *
302 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
303 * @const void *pp@ = pointer to octet array
304 * @size_t sz@ = size of octet array
305 *
306 * Returns: ---
307 *
308 * Use: Loads a negative MP in an octet array, least significant
309 * octet first, as two's complement. High-end octets are
310 * ignored if there isn't enough space for them. This probably
311 * means you made the wrong choice coming here.
312 */
313
314 extern void mpx_loadl2cn(mpw */*v*/, mpw */*vl*/,
315 const void */*p*/, size_t /*sz*/);
316
317 /* --- @mpx_storeb2cn@ --- *
318 *
319 * Arguments: @const mpw *v, *vl@ = base and limit of source vector
320 * @void *pp@ = pointer to octet array
321 * @size_t sz@ = size of octet array
322 *
323 * Returns: ---
324 *
325 * Use: Stores a negative MP in an octet array, most significant
326 * octet first, as two's complement. High-end octets are
327 * silently discarded if there isn't enough space for them,
328 * which probably isn't what you meant.
329 */
330
331 extern void mpx_storeb2cn(const mpw */*v*/, const mpw */*vl*/,
332 void */*p*/, size_t /*sz*/);
333
334 /* --- @mpx_loadb2cn@ --- *
335 *
336 * Arguments: @mpw *v, *vl@ = base and limit of destination vector
337 * @const void *pp@ = pointer to octet array
338 * @size_t sz@ = size of octet array
339 *
340 * Returns: ---
341 *
342 * Use: Loads a negative MP in an octet array, most significant octet
343 * first as two's complement. High-end octets are ignored if
344 * there isn't enough space for them. This probably means you
345 * chose this function wrongly.
346 */
347
348 extern void mpx_loadb2cn(mpw */*v*/, mpw */*vl*/,
349 const void */*p*/, size_t /*sz*/);
350
351
352 /*----- Logical shifting --------------------------------------------------*/
353
354 /* --- @mpx_lsl@ --- *
355 *
356 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
357 * @const mpw *av, *avl@ = source vector base and limit
358 * @size_t n@ = number of bit positions to shift by
359 *
360 * Returns: ---
361 *
362 * Use: Performs a logical shift left operation on an integer.
363 */
364
365 extern void mpx_lsl(mpw */*dv*/, mpw */*dvl*/,
366 const mpw */*av*/, const mpw */*avl*/,
367 size_t /*n*/);
368
369 /* --- @mpx_lsr@ --- *
370 *
371 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
372 * @const mpw *av, *avl@ = source vector base and limit
373 * @size_t n@ = number of bit positions to shift by
374 *
375 * Returns: ---
376 *
377 * Use: Performs a logical shift right operation on an integer.
378 */
379
380 extern void mpx_lsr(mpw */*dv*/, mpw */*dvl*/,
381 const mpw */*av*/, const mpw */*avl*/,
382 size_t /*n*/);
383
384 /*----- Bitwise operations ------------------------------------------------*/
385
386 /* --- How to implement them --- *
387 *
388 * x: 0011
389 * y: 0101
390 */
391
392 #define MPX_B0000(x, y) (0u)
393 #define MPX_B0001(x, y) ((x) & (y))
394 #define MPX_B0010(x, y) ((x) & ~(y))
395 #define MPX_B0011(x, y) (x)
396 #define MPX_B0100(x, y) (~(x) & ~(y))
397 #define MPX_B0101(x, y) (y)
398 #define MPX_B0110(x, y) ((x) ^ (y))
399 #define MPX_B0111(x, y) ((x) | (y))
400 #define MPX_B1000(x, y) (~((x) | (y)))
401 #define MPX_B1001(x, y) (~((x) ^ (y)))
402 #define MPX_B1010(x, y) (~(y))
403 #define MPX_B1011(x, y) ((x) | ~(y))
404 #define MPX_B1100(x, y) (~(x))
405 #define MPX_B1101(x, y) (~(x) | (y))
406 #define MPX_B1110(x, y) (~((x) & (y)))
407 #define MPX_B1111(x, y) (~0u)
408
409 /* --- @mpx_bitop@ --- *
410 *
411 * Arguments: @mpw *dv, *dvl@ = destination vector
412 * @const mpw *av, *avl@ = first source vector
413 * @const mpw *bv, *bvl@ = second source vector
414 *
415 * Returns: ---
416 *
417 * Use: Provide the dyadic boolean functions. The functions are
418 * named after the truth table they generate:
419 *
420 * a: 0011
421 * b: 0101
422 * @mpx_bitXXXX@
423 */
424
425 #define MPX_DOBIN(what) \
426 what(0000) what(0001) what(0010) what(0011) \
427 what(0100) what(0101) what(0110) what(0111) \
428 what(1000) what(1001) what(1010) what(1011) \
429 what(1100) what(1101) what(1110) what(1111)
430
431 #define MPX_BITDECL(string) \
432 extern void mpx_bit##string(mpw */*dv*/, mpw */*dvl*/, \
433 const mpw */*av*/, const mpw */*avl*/, \
434 const mpw */*bv*/, const mpw */*bvl*/);
435 MPX_DOBIN(MPX_BITDECL)
436
437 /* --- @mpx_[n]and@, @mpx_[n]or@, @mpx_xor@ --- *
438 *
439 * Synonyms for the commonly-used functions above.
440 */
441
442 #define mpx_and mpx_bit0001
443 #define mpx_or mpx_bit0111
444 #define mpx_nand mpx_bit1110
445 #define mpx_nor mpx_bit1000
446 #define mpx_xor mpx_bit0110
447
448 /* --- @mpx_not@ --- *
449 *
450 * Arguments: @mpw *dv, *dvl@ = destination vector
451 * @const mpw *av, *avl@ = first source vector
452 *
453 * Returns: ---
454 *
455 * Use; Bitwise NOT.
456 */
457
458 extern void mpx_not(mpw */*dv*/, mpw */*dvl*/,
459 const mpw */*av*/, const mpw */*avl*/);
460
461 /*----- Unsigned arithmetic -----------------------------------------------*/
462
463 /* --- @mpx_2c@ --- *
464 *
465 * Arguments: @mpw *dv, *dvl@ = destination vector
466 * @const mpw *v, *vl@ = source vector
467 *
468 * Returns: ---
469 *
470 * Use: Calculates the two's complement of @v@.
471 */
472
473 extern void mpx_2c(mpw */*dv*/, mpw */*dvl*/,
474 const mpw */*v*/, const mpw */*vl*/);
475
476 /* --- @mpx_ueq@ --- *
477 *
478 * Arguments: @const mpw *av, *avl@ = first argument vector base and limit
479 * @const mpw *bv, *bvl@ = second argument vector base and limit
480 *
481 * Returns: Nonzero if the two vectors are equal.
482 *
483 * Use: Performs an unsigned integer test for equality.
484 */
485
486 extern int mpx_ueq(const mpw */*av*/, const mpw */*avl*/,
487 const mpw */*bv*/, const mpw */*bvl*/);
488
489 /* --- @mpx_ucmp@ --- *
490 *
491 * Arguments: @const mpw *av, *avl@ = first argument vector base and limit
492 * @const mpw *bv, *bvl@ = second argument vector base and limit
493 *
494 * Returns: Less than, equal to, or greater than zero depending on
495 * whether @a@ is less than, equal to or greater than @b@,
496 * respectively.
497 *
498 * Use: Performs an unsigned integer comparison.
499 */
500
501 #define MPX_UCMP(av, avl, op, dv, dvl) \
502 (mpx_ucmp((av), (avl), (dv), (dvl)) op 0)
503
504 extern int mpx_ucmp(const mpw */*av*/, const mpw */*avl*/,
505 const mpw */*bv*/, const mpw */*bvl*/);
506
507 /* --- @mpx_uadd@ --- *
508 *
509 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
510 * @const mpw *av, *avl@ = first addend vector base and limit
511 * @const mpw *bv, *bvl@ = second addend vector base and limit
512 *
513 * Returns: ---
514 *
515 * Use: Performs unsigned integer addition. If the result overflows
516 * the destination vector, high-order bits are discarded. This
517 * means that two's complement addition happens more or less for
518 * free, although that's more a side-effect than anything else.
519 * The result vector may be equal to either or both source
520 * vectors, but may not otherwise overlap them.
521 */
522
523 extern void mpx_uadd(mpw */*dv*/, mpw */*dvl*/,
524 const mpw */*av*/, const mpw */*avl*/,
525 const mpw */*bv*/, const mpw */*bvl*/);
526
527 /* --- @mpx_uaddn@ --- *
528 *
529 * Arguments: @mpw *dv, *dvl@ = source and destination base and limit
530 * @mpw n@ = other addend
531 *
532 * Returns: ---
533 *
534 * Use: Adds a small integer to a multiprecision number.
535 */
536
537 #define MPX_UADDN(dv, dvl, n) do { \
538 mpw *_ddv = (dv), *_ddvl = (dvl); \
539 mpw _c = (n); \
540 \
541 while (_c && _ddv < _ddvl) { \
542 mpd _x = (mpd)*_ddv + (mpd)_c; \
543 *_ddv++ = MPW(_x); \
544 _c = _x >> MPW_BITS; \
545 } \
546 } while (0)
547
548 extern void mpx_uaddn(mpw */*dv*/, mpw */*dvl*/, mpw /*n*/);
549
550 /* --- @mpx_usub@ --- *
551 *
552 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
553 * @const mpw *av, *avl@ = first argument vector base and limit
554 * @const mpw *bv, *bvl@ = second argument vector base and limit
555 *
556 * Returns: ---
557 *
558 * Use: Performs unsigned integer subtraction. If the result
559 * overflows the destination vector, high-order bits are
560 * discarded. This means that two's complement subtraction
561 * happens more or less for free, although that's more a side-
562 * effect than anything else. The result vector may be equal to
563 * either or both source vectors, but may not otherwise overlap
564 * them.
565 */
566
567 extern void mpx_usub(mpw */*dv*/, mpw */*dvl*/,
568 const mpw */*av*/, const mpw */*avl*/,
569 const mpw */*bv*/, const mpw */*bvl*/);
570
571 /* --- @mpx_usubn@ --- *
572 *
573 * Arguments: @mpw *dv, *dvl@ = source and destination base and limit
574 * @n@ = subtrahend
575 *
576 * Returns: ---
577 *
578 * Use: Subtracts a small integer from a multiprecision number.
579 */
580
581 #define MPX_USUBN(dv, dvl, n) do { \
582 mpw *_ddv = (dv), *_ddvl = (dvl); \
583 mpw _c = (n); \
584 \
585 while (_ddv < _ddvl) { \
586 mpd _x = (mpd)*_ddv - (mpd)_c; \
587 *_ddv++ = MPW(_x); \
588 if (_x >> MPW_BITS) \
589 _c = 1; \
590 else \
591 break; \
592 } \
593 } while (0)
594
595 extern void mpx_usubn(mpw */*dv*/, mpw */*dvl*/, mpw /*n*/);
596
597 /* --- @mpx_umul@ --- *
598 *
599 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
600 * @const mpw *av, *avl@ = multiplicand vector base and limit
601 * @const mpw *bv, *bvl@ = multiplier vector base and limit
602 *
603 * Returns: ---
604 *
605 * Use: Performs unsigned integer multiplication. If the result
606 * overflows the desination vector, high-order bits are
607 * discarded. The result vector may not overlap the argument
608 * vectors in any way.
609 */
610
611 extern void mpx_umul(mpw */*dv*/, mpw */*dvl*/,
612 const mpw */*av*/, const mpw */*avl*/,
613 const mpw */*bv*/, const mpw */*bvl*/);
614
615 /* --- @mpx_umuln@ --- *
616 *
617 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
618 * @const mpw *av, *avl@ = multiplicand vector base and limit
619 * @mpw m@ = multiplier
620 *
621 * Returns: ---
622 *
623 * Use: Multiplies a multiprecision integer by a single-word value.
624 * The destination and source may be equal. The destination
625 * is completely cleared after use.
626 */
627
628 #define MPX_UMULN(dv, dvl, av, avl, m) do { \
629 mpw *_dv = (dv), *_dvl = (dvl); \
630 const mpw *_av = (av), *_avl = (avl); \
631 mpw _c = 0; \
632 mpd _m = (m); \
633 \
634 while (_av < _avl) { \
635 mpd _x; \
636 if (_dv >= _dvl) \
637 break; \
638 _x = (mpd)_m * (mpd)*_av++ + _c; \
639 *_dv++ = MPW(_x); \
640 _c = _x >> MPW_BITS; \
641 } \
642 if (_dv < _dvl) { \
643 *_dv++ = MPW(_c); \
644 MPX_ZERO(_dv, _dvl); \
645 } \
646 } while (0)
647
648 extern void mpx_umuln(mpw */*dv*/, mpw */*dvl*/,
649 const mpw */*av*/, const mpw */*avl*/, mpw m);
650
651 /* --- @mpx_umlan@ --- *
652 *
653 * Arguments: @mpw *dv, *dvl@ = destination/accumulator base and limit
654 * @const mpw *av, *avl@ = multiplicand vector base and limit
655 * @mpw m@ = multiplier
656 *
657 * Returns: ---
658 *
659 * Use: Multiplies a multiprecision integer by a single-word value
660 * and adds the result to an accumulator.
661 */
662
663 #define MPX_UMLAN(dv, dvl, av, avl, m) do { \
664 mpw *_dv = (dv), *_dvl = (dvl); \
665 const mpw *_av = (av), *_avl = (avl); \
666 mpw _cc = 0; \
667 mpd _m = (m); \
668 \
669 while (_dv < _dvl && _av < _avl) { \
670 mpd _x; \
671 _x = (mpd)*_dv + (mpd)_m * (mpd)*_av++ + _cc; \
672 *_dv++ = MPW(_x); \
673 _cc = _x >> MPW_BITS; \
674 } \
675 MPX_UADDN(_dv, _dvl, _cc); \
676 } while (0)
677
678 extern void mpx_umlan(mpw */*dv*/, mpw */*dvl*/,
679 const mpw */*av*/, const mpw */*avl*/, mpw m);
680
681 /* --- @mpx_usqr@ --- *
682 *
683 * Arguments: @mpw *dv, *dvl@ = destination vector base and limit
684 * @const mpw *av, *av@ = source vector base and limit
685 *
686 * Returns: ---
687 *
688 * Use: Performs unsigned integer squaring. The result vector must
689 * not overlap the source vector in any way.
690 */
691
692 extern void mpx_usqr(mpw */*dv*/, mpw */*dvl*/,
693 const mpw */*av*/, const mpw */*avl*/);
694
695 /* --- @mpx_udiv@ --- *
696 *
697 * Arguments: @mpw *qv, *qvl@ = quotient vector base and limit
698 * @mpw *rv, *rvl@ = dividend/remainder vector base and limit
699 * @const mpw *dv, *dvl@ = divisor vector base and limit
700 * @mpw *sv, *svl@ = scratch workspace
701 *
702 * Returns: ---
703 *
704 * Use: Performs unsigned integer division. If the result overflows
705 * the quotient vector, high-order bits are discarded. (Clearly
706 * the remainder vector can't overflow.) The various vectors
707 * may not overlap in any way. Yes, I know it's a bit odd
708 * requiring the dividend to be in the result position but it
709 * does make some sense really. The remainder must have
710 * headroom for at least two extra words. The scratch space
711 * must be at least one word larger than the divisor.
712 */
713
714 extern void mpx_udiv(mpw */*qv*/, mpw */*qvl*/, mpw */*rv*/, mpw */*rvl*/,
715 const mpw */*dv*/, const mpw */*dvl*/,
716 mpw */*sv*/, mpw */*svl*/);
717
718 /* --- @mpx_udivn@ --- *
719 *
720 * Arguments: @mpw *qv, *qvl@ = storage for the quotient (may overlap
721 * dividend)
722 * @const mpw *rv, *rvl@ = dividend
723 * @mpw d@ = single-precision divisor
724 *
725 * Returns: Remainder after divison.
726 *
727 * Use: Performs a single-precision division operation.
728 */
729
730 extern mpw mpx_udivn(mpw */*qv*/, mpw */*qvl*/,
731 const mpw */*rv*/, const mpw */*rvl*/, mpw /*d*/);
732
733 /*----- Karatsuba multiplication algorithms -------------------------------*/
734
735 /* --- @MPK_THRESH@ --- *
736 *
737 * This is the limiting length for using Karatsuba algorithms. It's best to
738 * use the simpler classical multiplication method on numbers smaller than
739 * this. It is unsafe to make this constant less than four (i.e., the
740 * algorithms will fail).
741 */
742
743 #define MPK_THRESH 16
744
745 /* --- @mpx_kmul@ --- *
746 *
747 * Arguments: @mpw *dv, *dvl@ = pointer to destination buffer
748 * @const mpw *av, *avl@ = pointer to first argument
749 * @const mpw *bv, *bvl@ = pointer to second argument
750 * @mpw *sv, *svl@ = pointer to scratch workspace
751 *
752 * Returns: ---
753 *
754 * Use: Multiplies two multiprecision integers using Karatsuba's
755 * algorithm. This is rather faster than traditional long
756 * multiplication (e.g., @mpx_umul@) on large numbers, although
757 * more expensive on small ones.
758 *
759 * The destination must be three times as large as the larger
760 * argument. The scratch space must be five times as large as
761 * the larger argument.
762 */
763
764 extern void mpx_kmul(mpw */*dv*/, mpw */*dvl*/,
765 const mpw */*av*/, const mpw */*avl*/,
766 const mpw */*bv*/, const mpw */*bvl*/,
767 mpw */*sv*/, mpw */*svl*/);
768
769 /* --- @mpx_ksqr@ --- *
770 *
771 * Arguments: @mpw *dv, *dvl@ = pointer to destination buffer
772 * @const mpw *av, *avl@ = pointer to first argument
773 * @mpw *sv, *svl@ = pointer to scratch workspace
774 *
775 * Returns: ---
776 *
777 * Use: Squares a multiprecision integers using something similar to
778 * Karatsuba's multiplication algorithm. This is rather faster
779 * than traditional long multiplication (e.g., @mpx_umul@) on
780 * large numbers, although more expensive on small ones, and
781 * rather simpler than full-blown Karatsuba multiplication.
782 *
783 * The destination must be three times as large as the larger
784 * argument. The scratch space must be five times as large as
785 * the larger argument.
786 */
787
788 extern void mpx_ksqr(mpw */*dv*/, mpw */*dvl*/,
789 const mpw */*av*/, const mpw */*avl*/,
790 mpw */*sv*/, mpw */*svl*/);
791
792 /*----- That's all, folks -------------------------------------------------*/
793
794 #ifdef __cplusplus
795 }
796 #endif
797
798 #endif
Catacomb cryptographic library (downstream of http://git.distorted.org.uk/~mdw/catacomb/)