|
|
1.1.1.2 root 1: /* C include file for MPI multiprecision integer math routines.
2:
3: Boulder Software Engineering
4: 3021 Eleventh Street
5: Boulder, CO 80304
6: (303) 541-0140
7:
8: (c) Copyright 1986-92 by Philip Zimmermann. All rights reserved.
9: The author assumes no liability for damages resulting from the use
10: of this software, even if the damage results from defects in this
11: software. No warranty is expressed or implied.
12:
13: These routines implement all of the multiprecision arithmetic necessary
14: for Rivest-Shamir-Adleman (RSA) public key cryptography, as well as
15: other number-theoretic algorithms such as ElGamal, Diffie-Hellman,
16: or Rabin.
17:
18: Although originally developed in Microsoft C for the IBM PC, this code
19: contains few machine dependencies. It assumes 2's complement
20: arithmetic. It can be adapted to 8-bit, 16-bit, or 32-bit machines,
21: lowbyte-highbyte order or highbyte-lowbyte order. This version
22: has been converted to ANSI C.
23:
24: Modified 8 Apr 92 - HAJK - Implement new VAX/VMS primitive support.
25: Modified 29 Nov 92 - Thad Smith
26: */
27:
28: #include <string.h>
29: #include "usuals.h" /* typedefs for byte, word16, boolean, etc. */
30: #include "platform.h" /* customization for different environments */
31:
32: /* Platform customization:
33: * A version which runs on almost any computer can be implemented by
34: * defining PORTABLE and MPORTABLE, preferably as a command line
35: * parameter. Faster versions can be generated by specifying specific
36: * parameters, such as size of unit and MULTUNIT, and by supplying some
37: * of the critical in assembly. See the file platform.h for more
38: * details on customization.
39: *
40: * The symbol HIGHFIRST, designating that integers and longs are stored
41: * with the most significant bit in the lowest address, should be defined
42: * on the command line for compiling all files, since it is used by files
43: * other than the mpilib routines.
44: */
45:
46: #ifndef ALIGN
47: #define ALIGN
48: #endif
49:
50: #ifndef PEASANT /* if not Russian peasant modulo multiply algorithm */
51: #ifndef MERRITT /* if not Merritt's modmult */
52: #ifndef UPTON /* if not Upton's modmult */
53: #ifndef SMITH
54: #define SMITH /* default: use Smith's modmult algorithm */
55: #endif
56: #endif
57: #endif
58: #endif
59:
60: #ifdef SMITH
61: #define UPTON_OR_SMITH /* enable common code */
62: #endif
63: #ifdef UPTON
64: #define UPTON_OR_SMITH /* enable common code */
65: #endif
66:
67: #ifndef UNIT32
68: #ifndef UNIT8
69: #define UNIT16 /* default--use 16-bit units */
70: #endif
71: #endif
72:
73: /*** CAUTION: If your machine has an unusual word size that is not a
74: power of 2 (8, 16, 32, or 64) bits wide, then the macros here that
75: use the symbol "LOG_UNITSIZE" must be changed.
76: ***/
77:
78: #ifdef UNIT8
79: typedef unsigned char unit;
80: typedef signed char signedunit;
81: #define UNITSIZE 8 /* number of bits in a unit */
82: #define LOG_UNITSIZE 3
83: #define uppermostbit ((unit) 0x80)
84: #define BYTES_PER_UNIT 1 /* number of bytes in a unit */
85: #define units2bits(n) ((n) << 3) /* fast multiply by UNITSIZE */
86: #define units2bytes(n) (n)
87: #define bits2units(n) (((n)+7) >> 3)
88: #define bytes2units(n) (n)
89: #endif
90:
91: #ifdef UNIT16
92: typedef word16 unit;
93: typedef short signedunit;
94: #define UNITSIZE 16 /* number of bits in a unit */
95: #define LOG_UNITSIZE 4
96: #define uppermostbit ((unit) 0x8000)
97: #define BYTES_PER_UNIT 2 /* number of bytes in a unit */
98: #define units2bits(n) ((n) << 4) /* fast multiply by UNITSIZE */
99: #define units2bytes(n) ((n) << 1)
100: #define bits2units(n) (((n)+15) >> 4)
101: #define bytes2units(n) (((n)+1) >> 1)
102: #endif
103:
104: #ifdef UNIT32
105: typedef word32 unit;
106: typedef long signedunit;
107: #define UNITSIZE 32 /* number of bits in a unit */
108: #define LOG_UNITSIZE 5
109: #define uppermostbit ((unit) 0x80000000L)
110: #define BYTES_PER_UNIT 4 /* number of bytes in a unit */
111: #define units2bits(n) ((n) << 5) /* fast multiply by UNITSIZE */
112: #define units2bytes(n) ((n) << 2)
113: #define bits2units(n) (((n)+31) >> 5)
114: #define bytes2units(n) (((n)+3) >> 2)
115: #endif
116:
117: #define power_of_2(b) ((unit) 1 << (b)) /* computes power-of-2 bit masks */
118: #define bits2bytes(n) (((n)+7) >> 3)
119: /* Some C compilers (like the ADSP2101) will not always collapse constant
120: expressions at compile time if the expressions contain shift operators. */
121: /* #define uppermostbit power_of_2(UNITSIZE-1) */
122: /* #define UNITSIZE units2bits(1) */ /* number of bits in a unit */
123: /* #define bytes2units(n) bits2units((n)<<3) */
124: /* #define BYTES_PER_UNIT (UNITSIZE >> 3) */
125: /* LOG_UNITSIZE is the log base 2 of UNITSIZE, ie: 4 for 16-bit units */
126: /* #define units2bits(n) ((n) << LOG_UNITSIZE) */ /* fast multiply by UNITSIZE */
127: /* #define units2bytes(n) ((n) << (LOG_UNITSIZE-3)) */
128: /* #define bits2units(n) (((n)+(UNITSIZE-1)) >> LOG_UNITSIZE) */
129: /* #define bytes2units(n) (((n)+(BYTES_PER_UNIT-1)) >> (LOG_UNITSIZE-3)) */
130:
131: typedef unit *unitptr;
132:
133:
134: /*--------------------- Byte ordering stuff -------------------*/
135: #ifdef HIGHFIRST
136:
137: /* these definitions assume MSB comes first */
138: #define tohigher(n) (-(n)) /* offset towards higher unit */
139: #define pre_higherunit(r) (--(r))
140: #define pre_lowerunit(r) (++(r))
141: #define post_higherunit(r) ((r)--)
142: #define post_lowerunit(r) ((r)++)
143: #define bit_index(n) (global_precision-bits2units((n)+1))
144: #define lsbptr(r,prec) ((r)+(prec)-1)
145: #define make_lsbptr(r,prec) (r) = lsbptr(r,prec)
146: #define unmake_lsbptr(r,prec) (r) = ((r)-(prec)+1)
147: #define msbptr(r,prec) (r)
148: #define make_msbptr(r,prec) /* (r) = msbptr(r,prec) */
149:
150: /* The macro rescale(r,current_precision,new_precision) rescales
151: a multiprecision integer by adjusting r and its precision to new values.
152: It can be used to reverse the effects of the normalize
153: routine given above. See the comments in normalize concerning
154: Intel vs. Motorola LSB/MSB conventions.
155: WARNING: You can only safely call rescale on registers that
156: you have previously normalized with the above normalize routine,
157: or are known to be big enough for the new precision. You may
158: specify a new precision that is smaller than the current precision.
159: You must be careful not to specify a new_precision value that is
160: too big, or which adjusts the r pointer out of range.
161: */
162: #define rescale(r,currentp,newp) r -= ((newp) - (currentp))
163:
164: /* The macro normalize(r,precision) "normalizes" a multiprecision integer
165: by adjusting r and precision to new values. For Motorola-style processors
166: (MSB-first), r is a pointer to the MSB of the register, and must
167: be adjusted to point to the first nonzero unit. For Intel/VAX-style
168: (LSB-first) processors, r is a pointer to the LSB of the register,
169: and must be left unchanged. The precision counter is always adjusted,
170: regardless of processor type. In the case of precision = 0,
171: r becomes undefined.
172: */
173: #define normalize(r,prec) \
174: { prec = significance(r); r += (global_precision-(prec)); }
175:
176: #else /* LOWFIRST byte order */
177:
178: /* these definitions assume LSB comes first */
179: #define tohigher(n) (n) /* offset towards higher unit */
180: #define pre_higherunit(r) (++(r))
181: #define pre_lowerunit(r) (--(r))
182: #define post_higherunit(r) ((r)++)
183: #define post_lowerunit(r) ((r)--)
184: #define bit_index(n) (bits2units((n)+1)-1)
185: #define lsbptr(r,prec) (r)
186: #define make_lsbptr(r,prec) /* (r) = lsbptr(r,prec) */
187: #define unmake_lsbptr(r,prec) /* (r) = (r) */
188: #define msbptr(r,prec) ((r)+(prec)-1)
189: #define make_msbptr(r,prec) (r) = msbptr(r,prec)
190:
191: #define rescale(r,currentp,newp) /* nil statement */
192: #define normalize(r,prec) prec = significance(r)
193:
194: #endif /* LOWFIRST byte order */
195: /*------------------ End byte ordering stuff -------------------*/
196:
197: /* Note that the address calculations require that lsbptr, msbptr,
198: make_lsbptr, make_msbptr, mp_tstbit, mp_setbit, mp_clrbit,
199: and bitptr all have unitptr arguments, not byte pointer arguments. */
200: #define bitptr(r,n) &((r)[bit_index(n)])
201: #define bitmsk(n) power_of_2((n) & (UNITSIZE-1))
202: /* bitmsk() assumes UNITSIZE is a power of 2 */
203: #define mp_tstbit(r,n) (*bitptr(r,n) & bitmsk(n))
204: #define mp_setbit(r,n) (*bitptr(r,n) |= bitmsk(n))
205: #define mp_clrbit(r,n) (*bitptr(r,n) &= ~bitmsk(n))
206: #define msunit(r) (*msbptr(r,global_precision))
207: #define lsunit(r) (*lsbptr(r,global_precision))
208: /* #define mp_tstminus(r) ((msunit(r) & uppermostbit)!=0) */
209: #define mp_tstminus(r) ((signedunit) msunit(r) < 0)
210:
211:
212: /* set working precision to specified number of bits. */
213: #ifdef mp_setp
214: void mp_setp(short nbits);
215: #define set_precision(prec) mp_setp(units2bits(global_precision=(prec)))
216: #else
217: #define set_precision(prec) (global_precision = (prec))
218: #endif
219:
220:
221: #ifdef PEASANT
222:
223: /* Define C names for Russian peasant modmult primitives. */
224: #define stage_modulus stage_peasant_modulus
225: #define mp_modmult peasant_modmult
226: #define modmult_burn peasant_burn
227: #define SLOP_BITS PEASANT_SLOP_BITS
228:
229: #else /* not PEASANT */
230: #ifdef MERRITT
231: /* Define C names for Merritt's modmult primitives. */
232: #define stage_modulus stage_merritt_modulus
233: #define mp_modmult merritt_modmult
234: #define modmult_burn merritt_burn
235: #define SLOP_BITS MERRITT_SLOP_BITS
236:
237: #else /* not PEASANT, MERRITT */
238: #ifdef UPTON
239: /* Define C names for Upton's modmult primitives. */
240: #define stage_modulus stage_upton_modulus
241: #define mp_modmult upton_modmult
242: #define modmult_burn upton_burn
243: #define SLOP_BITS UPTON_SLOP_BITS
244:
245: #else /* not PEASANT, MERRITT, UPTON */
246: #ifdef SMITH
247: /* Define C names for Smith's modmult primitives. */
248: #define stage_modulus stage_smith_modulus
249: #define mp_modmult smith_modmult
250: #define modmult_burn smith_burn
251: #define SLOP_BITS SMITH_SLOP_BITS
252:
253: #endif /* SMITH */
254: #endif /* UPTON */
255: #endif /* MERRITT */
256: #endif /* PEASANT */
257:
258:
259: #define mp_shift_left(r1) mp_rotate_left(r1,(boolean)0)
260: /* multiprecision shift left 1 bit */
261:
262: #define mp_add(r1,r2) mp_addc(r1,r2,(boolean)0)
263: /* multiprecision add with no carry */
264:
265: #define mp_sub(r1,r2) mp_subb(r1,r2,(boolean)0)
266: /* multiprecision subtract with no borrow */
267:
268: #define mp_abs(r) (mp_tstminus(r) ? (mp_neg(r),TRUE) : FALSE)
269:
270: #define msub(r,m) if (mp_compare(r,m) >= 0) mp_sub(r,m)
271: /* Prevents r from getting bigger than modulus m */
272:
273: #define testeq(r,i) \
274: ( (lsunit(r)==(i)) && (significance(r)<=1) )
275:
276: #define testne(r,i) \
277: ( (lsunit(r)!=(i)) || (significance(r)>1) )
278:
279: #define testge(r,i) \
280: ( (lsunit(r)>=(i)) || (significance(r)>1) )
281:
282: #define testle(r,i) \
283: ( (lsunit(r)<=(i)) && (significance(r)<=1) )
284:
285: #define mp_square(r1,r2) mp_mult(r1,r2,r2)
286: /* Square r2, returning product in r1 */
287:
288: #define mp_modsquare(r1,r2) mp_modmult(r1,r2,r2)
289: /* Square r2, returning modulo'ed product in r1 */
290:
291: #define countbytes(r) ((countbits(r)+7)>>3)
292:
293: /* SLOP_BITS is how many "carry bits" to allow for intermediate
294: calculation results to exceed the size of the modulus.
295: It is used by modexp to give some overflow elbow room for
296: modmult to use to perform modulo operations with the modulus.
297: The number of slop bits required is determined by the modmult
298: algorithm. The Russian peasant modmult algorithm only requires
299: 1 slop bit, for example. Note that if we use an external assembly
300: modmult routine, SLOP_BITS may be meaningless or may be defined in a
301: non-constant manner.
302: */
303: #define PEASANT_SLOP_BITS 1
304: #define MERRITT_SLOP_BITS UNITSIZE
305: #define UPTON_SLOP_BITS (UNITSIZE/2)
306: #ifdef mp_smul /* old version requires MS word = 0 */
307: #define SMITH_SLOP_BITS UNITSIZE
308: #else /* mp_smula or C version of mp_smul */
309: #define SMITH_SLOP_BITS 0
310: #endif /* mp_smul */
311:
1.1.1.3 root 312: /* MAX_BIT_PRECISION is upper limit that assembly primitives can handle.
313: It must be less than 32704 bits, or 4088 bytes. It should be an
314: integer multiple of UNITSIZE*2.
315: */
316: #if (SLOP_BITS > 0)
1.1.1.4 ! root 317: #define MAX_BIT_PRECISION (1280+(2*UNITSIZE))
1.1.1.3 root 318: #else
1.1.1.4 ! root 319: #define MAX_BIT_PRECISION 1280
1.1.1.3 root 320: #endif
321: #define MAX_BYTE_PRECISION (MAX_BIT_PRECISION/8)
322: #define MAX_UNIT_PRECISION (MAX_BIT_PRECISION/UNITSIZE)
323:
1.1.1.2 root 324:
325: /* global_precision is the unit precision last set by set_precision */
326: extern short global_precision;
327:
328:
329: /* The "bit sniffer" macros all begin sniffing at the MSB.
330: They are used internally by all the various multiply, divide,
331: modulo, exponentiation, and square root functions.
332: */
333: #define sniff_bit(bptr,bitmask) (*(bptr) & bitmask)
334:
335: #define init_bitsniffer(bptr,bitmask,prec,bits) \
336: { normalize(bptr,prec); \
337: if (!prec) \
338: return(0); \
339: bits = units2bits(prec); \
340: make_msbptr(bptr,prec); bitmask = uppermostbit; \
341: while (!sniff_bit(bptr,bitmask)) \
342: { bitmask >>= 1; bits--; \
343: } \
344: }
345:
346: #define bump_bitsniffer(bptr,bitmask) \
347: { if (!(bitmask >>= 1)) \
348: { bitmask = uppermostbit; \
349: post_lowerunit(bptr); \
350: } \
351: }
352:
353: /* bump_2bitsniffers is used internally by mp_udiv. */
354: #define bump_2bitsniffers(bptr,bptr2,bitmask) \
355: { if (!(bitmask >>= 1)) \
356: { bitmask = uppermostbit; \
357: post_lowerunit(bptr); \
358: post_lowerunit(bptr2); \
359: } \
360: }
361:
362: /* stuff_bit is used internally by mp_udiv and mp_sqrt. */
363: #define stuff_bit(bptr,bitmask) *(bptr) |= bitmask
364:
365:
366: boolean mp_addc
367: (register unitptr r1,register unitptr r2,register boolean carry);
368: /* multiprecision add with carry r2 to r1, result in r1 */
369:
370: boolean mp_subb
371: (register unitptr r1,register unitptr r2,register boolean borrow);
372: /* multiprecision subtract with borrow, r2 from r1, result in r1 */
373:
374: boolean mp_rotate_left(register unitptr r1,register boolean carry);
375: /* multiprecision rotate left 1 bit with carry, result in r1. */
376:
377: void mp_shift_right_bits(register unitptr r1,register short bits);
378: /* multiprecision shift right bits, result in r1. */
379:
380: short mp_compare(register unitptr r1,register unitptr r2);
381: /* Compares registers *r1, *r2, and returns -1, 0, or 1 */
382:
383: boolean mp_inc(register unitptr r);
384: /* Increment multiprecision integer r. */
385:
386: boolean mp_dec(register unitptr r);
387: /* Decrement multiprecision integer r. */
388:
389: void mp_neg(register unitptr r);
390: /* Compute 2's complement, the arithmetic negative, of r */
391:
392: #ifndef mp_move
393: #define mp_move(d,s) memcpy((void*)(d), (void*)(s), \
394: units2bytes(global_precision))
395: #endif
396: #ifndef unitfill0
397: #define unitfill0(r,ct) memset((void*)(r), 0, units2bytes(ct))
398: #endif
399:
400: #ifndef mp_burn
401: #define mp_burn(r) mp_init(r,0) /* for burning the evidence */
402: #define mp_init0(r) mp_init(r,0)
403: #endif
404:
405: #define empty_array(r) unitfill0(r, sizeof(r)/sizeof(r[0])/sizeof(unit))
406:
407: void mp_init(register unitptr r, word16 value);
408: /* Init multiprecision register r with short value. */
409:
410: short significance(register unitptr r);
411: /* Returns number of significant units in r */
412:
413: int mp_udiv(register unitptr remainder,register unitptr quotient,
414: register unitptr dividend,register unitptr divisor);
415: /* Unsigned divide, treats both operands as positive. */
416:
417: int mp_recip(register unitptr quotient,register unitptr divisor);
418: /* Compute reciprocal as 1/divisor. Used by faster modmult. */
419:
420: int mp_div(register unitptr remainder,register unitptr quotient,
421: register unitptr dividend,register unitptr divisor);
422: /* Signed divide, either or both operands may be negative. */
423:
424: word16 mp_shortdiv(register unitptr quotient,
425: register unitptr dividend,register word16 divisor);
426: /* Returns short remainder of unsigned divide. */
427:
428: int mp_mod(register unitptr remainder,
429: register unitptr dividend,register unitptr divisor);
430: /* Unsigned divide, treats both operands as positive. */
431:
432: word16 mp_shortmod(register unitptr dividend,register word16 divisor);
433: /* Just returns short remainder of unsigned divide. */
434:
435: int mp_mult(register unitptr prod,
436: register unitptr multiplicand,register unitptr multiplier);
437: /* Computes multiprecision prod = multiplicand * multiplier */
438:
439: int countbits(unitptr r);
440: /* Returns number of significant bits in r. */
441:
442: int stage_peasant_modulus(unitptr n);
443: int stage_merritt_modulus(unitptr n);
444: int stage_upton_modulus(unitptr n);
445: int stage_smith_modulus(unitptr n);
446: /* Must pass modulus to stage_modulus before calling modmult. */
447:
448: int peasant_modmult(register unitptr prod,
449: unitptr multiplicand,register unitptr multiplier);
450: int merritt_modmult(register unitptr prod,
451: unitptr multiplicand,register unitptr multiplier);
452: int upton_modmult(register unitptr prod,
453: unitptr multiplicand,register unitptr multiplier);
454: int smith_modmult(register unitptr prod,
455: unitptr multiplicand,register unitptr multiplier);
456: /* Performs combined multiply/modulo operation, with global modulus */
457:
458:
459:
460: int mp_modexp(register unitptr expout,register unitptr expin,
461: register unitptr exponent,register unitptr modulus);
462: /* Combined exponentiation/modulo algorithm. */
463:
1.1.1.4 ! root 464: int mp_modexp_crt(unitptr expout, unitptr expin,
! 465: unitptr p, unitptr q, unitptr ep, unitptr eq, unitptr u);
! 466: /* exponentiation and modulo using Chinese Remainder Theorem */
1.1.1.2 root 467:
468: /****************** end of MPI library ****************************/
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.