mpilib.cpp

RSA算法的VC源码

/* C source code for multiprecision arithmetic library routines.


   These routines implement all of the multiprecision arithmetic
   necessary for number-theoretic cryptographic algorithms such as
   ElGamal, Diffie-Hellman, Rabin, or factoring studies for large
   composite numbers, as well as Rivest-Shamir-Adleman (RSA) public
   key cryptography.
 */

/* #define COUNTMULTS *//* count modmults for performance studies */
#include "stdafx.h"

#ifndef poll_for_break
#define poll_for_break()	/* stub */
#endif

#include "mpilib.h"

#ifdef mp_smula
#ifdef mp_smul
Error:Both mp_smula and mp_smul cannot be defined.
#else
#define mp_smul	mp_smula
#endif
#endif

/* set macros for MULTUNIT */
#ifdef MUNIT8
#define MULTUNITSIZE   8
typedef unsigned char MULTUNIT;
#ifdef UNIT8
#define MULTUNIT_SIZE_SAME
#endif
#else				/* not MUNIT8 */
#ifdef MUNIT32
#define MULTUNITSIZE   32
typedef unsigned long MULTUNIT;
#ifdef UNIT32
#define MULTUNIT_SIZE_SAME
#else
/* #error is not portable, this has the same effect */
#include "UNITSIZE cannot be smaller than MULTUNITSIZE"
#endif
#else				/* assume MUNIT16 */
#define MULTUNITSIZE   16
typedef unsigned short MULTUNIT;
#ifdef UNIT16
#define MULTUNIT_SIZE_SAME
#endif				/* UNIT16 */
#ifdef UNIT8
#include "UNITSIZE cannot be smaller than MULTUNITSIZE"
#endif				/* UNIT8 */
#endif				/* MUNIT32 */
#endif				/* MUNIT8 */

#define MULTUNIT_msb    ((MULTUNIT)1 << (MULTUNITSIZE-1))	/* msb of
								   MULTUNIT */
#define DMULTUNIT_msb   (1L << (2*MULTUNITSIZE-1))
#define MULTUNIT_mask   ((MULTUNIT)((1L << MULTUNITSIZE)-1))
#define MULTUNITs_perunit   (UNITSIZE/MULTUNITSIZE)


void mp_smul(MULTUNIT * prod, MULTUNIT * multiplicand, MULTUNIT multiplier);
void mp_dmul(unitptr prod, unitptr multiplicand, unitptr multiplier);

short global_precision = 0;	/* units of precision for all routines */
/*      global_precision is the unit precision last set by set_precision.
   Initially, set_precision() should be called to define global_precision
   before using any of these other multiprecision library routines.
   i.e.:   set_precision(MAX_UNIT_PRECISION);
 */

/*************** multiprecision library primitives ****************/
/*      The following portable C primitives should be recoded in assembly.
   The entry point name should be defined, in "mpilib.h" to the external
   entry point name.  If undefined, the C version will be used.
 */

typedef unsigned long int ulint;

UINT random(UINT nMax)//0-->nMax-1
{
	int nRand=rand();
	float fMap=(float)(nMax-1)/RAND_MAX;
	float fRetVal=(float)nRand*fMap+0.5F;
	return (UINT)fRetVal;
}

#ifndef mp_addc
/* 
   multiprecision add with carry r2 to r1, result in r1
   carry is incoming carry flag-- value should be 0 or 1
*/
boolean mp_addc
 (register unitptr r1, register unitptr r2, register boolean carry)
{
    register unit x;
    short precision;		/* number of units to add */
    precision = global_precision;
    make_lsbptr(r1, precision);
    make_lsbptr(r2, precision);
    while (precision--) {
	if (carry) {
	    x = *r1 + *r2 + 1;
	    carry = (*r2 >= (unit) (~*r1));
	} else {
	    x = *r1 + *r2;
	    carry = (x < *r1);
	}
	post_higherunit(r2);
	*post_higherunit(r1) = x;
    }
    return carry;		/* return the final carry flag bit */
}				/* mp_addc */
#endif				/* mp_addc */

#ifndef mp_subb

/* 
   multiprecision subtract with borrow, r2 from r1, result in r1
   borrow is incoming borrow flag-- value should be 0 or 1
 */
boolean mp_subb
 (register unitptr r1, register unitptr r2, register boolean borrow)
{
    register unit x;
    short precision;		/* number of units to subtract */
    precision = global_precision;
    make_lsbptr(r1, precision);
    make_lsbptr(r2, precision);
    while (precision--) {
	if (borrow) {
	    x = *r1 - *r2 - 1;
	    borrow = (*r1 <= *r2);
	} else {
	    x = *r1 - *r2;
	    borrow = (*r1 < *r2);
	}
	post_higherunit(r2);
	*post_higherunit(r1) = x;
    }
    return borrow;		/* return the final carry/borrow flag bit */
}				/* mp_subb */
#endif				/* mp_subb */

#ifndef mp_rotate_left

/*
   multiprecision rotate left 1 bit with carry, result in r1.
   carry is incoming carry flag-- value should be 0 or 1
*/
boolean mp_rotate_left(register unitptr r1, register boolean carry)
{
    register int precision;	/* number of units to rotate */
    unsigned int mcarry = carry, nextcarry;	/* int is supposed to be
						 * the efficient size for ops*/
    precision = global_precision;
    make_lsbptr(r1, precision);
    while (precision--) {
	nextcarry = (((signedunit) * r1) < 0);
	*r1 = (*r1 << 1) | mcarry;
	mcarry = nextcarry;
	pre_higherunit(r1);
    }
    return nextcarry;		/* return the final carry flag bit */
}				/* mp_rotate_left */
#endif				/* mp_rotate_left */

/************** end of primitives that should be in assembly *************/

/* The mp_shift_right_bits function is not called in any time-critical
   situations in public-key cryptographic functions, so it doesn't
   need to be coded in assembly language.
 */

/*
   multiprecision shift right bits, result in r1.
   bits is how many bits to shift, must be <= UNITSIZE.
*/
void mp_shift_right_bits(register unitptr r1, register short bits)
{
    register short precision;	/* number of units to shift */
    register unit carry, nextcarry, bitmask;
    register short unbits;
    if (bits == 0)
	return;			/* shift zero bits is a no-op */
    carry = 0;
    bitmask = power_of_2(bits) - 1;
    unbits = UNITSIZE - bits;	/* shift bits must be <= UNITSIZE */
    precision = global_precision;
    make_msbptr(r1, precision);
    if (bits == UNITSIZE) {
	while (precision--) {
	    nextcarry = *r1;
	    *r1 = carry;
	    carry = nextcarry;
	    pre_lowerunit(r1);
	}
    } else {
	while (precision--) {
	    nextcarry = *r1 & bitmask;
	    *r1 >>= bits;
	    *r1 |= carry << unbits;
	    carry = nextcarry;
	    pre_lowerunit(r1);
	}
    }
}				/* mp_shift_right_bits */

#ifndef mp_compare

/*
 * Compares multiprecision integers *r1, *r2, and returns:
 * -1 iff *r1 < *r2
 *  0 iff *r1 == *r2
 * +1 iff *r1 > *r2
 */
short mp_compare(register unitptr r1, register unitptr r2)
{
    register short precision;	/* number of units to compare */

    precision = global_precision;
    make_msbptr(r1, precision);
    make_msbptr(r2, precision);
    do {
	if (*r1 < *r2)
	    return -1;
	if (*post_lowerunit(r1) > *post_lowerunit(r2))
	    return 1;
    } while (--precision);
    return 0;			/*  *r1 == *r2  */
}				/* mp_compare */
#endif				/* mp_compare */

/* Increment multiprecision integer r. */
boolean mp_inc(register unitptr r)
{
    register short precision;
    precision = global_precision;
    make_lsbptr(r, precision);
    do {
	if (++(*r))
	    return 0;		/* no carry */
	post_higherunit(r);
    } while (--precision);
    return 1;			/* return carry set */
}				/* mp_inc */

/* Decrement multiprecision integer r. */
boolean mp_dec(register unitptr r)
{
    register short precision;
    precision = global_precision;
    make_lsbptr(r, precision);
    do {
	if ((signedunit) (--(*r)) != (signedunit) - 1)
	    return 0;		/* no borrow */
	post_higherunit(r);
    } while (--precision);
    return 1;			/* return borrow set */
}				/* mp_dec */

/* Compute 2's complement, the arithmetic negative, of r */
void mp_neg(register unitptr r)
{
    register short precision;	/* number of units to negate */
    precision = global_precision;
    mp_dec(r);			/* 2's complement is 1's complement plus 1 */
    do {			/* now do 1's complement */
	*r = ~(*r);
	r++;
    } while (--precision);
}				/* mp_neg */

#ifndef mp_move

void mp_move(register unitptr dst, register unitptr src)
{
    register short precision;	/* number of units to move */
    precision = global_precision;
    do {
	*dst++ = *src++;
    } while (--precision);
}				/* mp_move */
#endif				/* mp_move */

/* Init multiprecision register r with short value. */
void mp_init(register unitptr r, word16 value)
{	/* Note that mp_init doesn't extend sign bit for >32767 */

    unitfill0(r, global_precision);
    make_lsbptr(r, global_precision);
    *post_higherunit(r) = value;
#ifdef UNIT8
    *post_higherunit(r) = value >> UNITSIZE;
#endif
}				/* mp_init */

/* Returns number of significant units in r */
short significance(register unitptr r)
{
    register short precision;
    precision = global_precision;
    make_msbptr(r, precision);
    do {
	if (*post_lowerunit(r))
	    return precision;
    } while (--precision);
    return precision;
}				/* significance */


#ifndef unitfill0

/* Zero-fill the unit buffer r. */
void unitfill0(unitptr r, word16 unitcount)
{
    while (unitcount--)
	*r++ = 0;
}				/* unitfill0 */
#endif				/* unitfill0 */

/* Unsigned divide, treats both operands as positive. */
int mp_udiv(register unitptr remainder, register unitptr quotient,
	    register unitptr dividend, register unitptr divisor)
{
    int bits;
    short dprec;
    register unit bitmask;

    if (testeq(divisor, 0))
	return -1;		/* zero divisor means divide error */
    mp_init0(remainder);
    mp_init0(quotient);
    /* normalize and compute number of bits in dividend first */
    init_bitsniffer(dividend, bitmask, dprec, bits);
    /* rescale quotient to same precision (dprec) as dividend */
    rescale(quotient, global_precision, dprec);
    make_msbptr(quotient, dprec);

    while (bits--) {
	mp_rotate_left(remainder,
		       (boolean) (sniff_bit(dividend, bitmask) != 0));
	if (mp_compare(remainder, divisor) >= 0) {
	    mp_sub(remainder, divisor);
	    stuff_bit(quotient, bitmask);
	}
	bump_2bitsniffers(dividend, quotient, bitmask);
    }
    return 0;
}				/* mp_udiv */

#ifdef UPTON_OR_SMITH

#define RECIPMARGIN 0		/* extra margin bits used by mp_recip() */

/* Compute reciprocal (quotient) as 1/divisor.  Used by faster modmult. */
int mp_recip(register unitptr quotient, register unitptr divisor)
{
    int bits;
    short qprec;
    register unit bitmask;
    unit remainder[MAX_UNIT_PRECISION];
    if (testeq(divisor, 0))
	return -1;		/* zero divisor means divide error */
    mp_init0(remainder);
    mp_init0(quotient);

    /* normalize and compute number of bits in quotient first */
    bits = countbits(divisor) + RECIPMARGIN;
    bitmask = bitmsk(bits);	/* bitmask within a single unit */
    qprec = bits2units(bits + 1);
    mp_setbit(remainder, (bits - RECIPMARGIN) - 1);
    /* rescale quotient to precision of divisor + RECIPMARGIN bits */
    rescale(quotient, global_precision, qprec);
    make_msbptr(quotient, qprec);

    while (bits--) {
	mp_shift_left(remainder);
	if (mp_compare(remainder, divisor) >= 0) {
	    mp_sub(remainder, divisor);