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