mpi.c

PeerSec Networks MatrixSSL?is an embedded SSL implementation designed for small footprint applicatio

/*	
 *	mpi.c
 *	Release $Name: MATRIXSSL_1_7_3_OPEN $
 *
 *	multiple-precision integer library
 */
/*
 *	Copyright (c) PeerSec Networks, 2002-2005. All Rights Reserved.
 *	The latest version of this code is available at http://www.matrixssl.org
 *
 *	This software is open source; you can redistribute it and/or modify
 *	it under the terms of the GNU General Public License as published by
 *	the Free Software Foundation; either version 2 of the License, or
 *	(at your option) any later version.
 *
 *	This General Public License does NOT permit incorporating this software 
 *	into proprietary programs.  If you are unable to comply with the GPL, a 
 *	commercial license for this software may be purchased from PeerSec Networks
 *	at http://www.peersec.com
 *	
 *	This program is distributed in WITHOUT ANY WARRANTY; without even the 
 *	implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. 
 *	See the GNU General Public License for more details.
 *	
 *	You should have received a copy of the GNU General Public License
 *	along with this program; if not, write to the Free Software
 *	Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 *	http://www.gnu.org/copyleft/gpl.html
 */
/******************************************************************************/

#include "../cryptoLayer.h"
#include <stdarg.h>

static int32 mp_exptmod_fast (psPool_t *pool, mp_int * G, mp_int * X,
				mp_int * P, mp_int * Y, int32 redmode);

/******************************************************************************/
/*
	FUTURE
	1. Convert the mp_init and mp_clear functions to not use malloc + free,
	but to use static storage within the bignum variable instead - but
	how to handle grow()?  Maybe use a simple memory allocator
	2. verify stack usage of all functions and use of MP_LOW_MEM:
		fast_mp_montgomery_reduce
		fast_s_mp_mul_digs
		fast_s_mp_sqr
		fast_s_mp_mul_high_digs
	3. HAC stands for Handbook of Applied Cryptography
		http://www.cacr.math.uwaterloo.ca/hac/
*/
/******************************************************************************/
/*
	Utility functions
*/
void psZeromem(void *dst, size_t len)
{
	unsigned char *mem = (unsigned char *)dst;
	
	if (dst == NULL) {
		return;
	}
	while (len-- > 0) {
		*mem++ = 0;
	}
}

void psBurnStack(unsigned long len)
{
	unsigned char buf[32];
	
	psZeromem(buf, sizeof(buf));
	if (len > (unsigned long)sizeof(buf)) {
		psBurnStack(len - sizeof(buf));
	}
}

/******************************************************************************/
/*
	Multiple precision integer functions
	Note: we don't use va_args here to prevent portability issues.
*/
int32 _mp_init_multi(psPool_t *pool, mp_int *mp0, mp_int *mp1, mp_int *mp2, 
					mp_int *mp3, mp_int *mp4, mp_int *mp5, 
					mp_int *mp6, mp_int *mp7)
{
	mp_err	res		= MP_OKAY;		/* Assume ok until proven otherwise */
	int32		n		= 0;			/* Number of ok inits */
	mp_int	*tempArray[9] = {mp0, mp1, mp2, mp3, mp4, mp5, mp6, mp7, NULL};

	while (tempArray[n] != NULL) {
		if (mp_init(pool, tempArray[n]) != MP_OKAY) {
			res = MP_MEM;
			break;
		}
		n++;
	}

	if (res == MP_MEM) {
		n = 0;
		while (tempArray[n] != NULL) {
			mp_clear(tempArray[n]);
			n++;
		}
	}
	return res;		/* Assumed ok, if error flagged above. */
}
/******************************************************************************/
/*
	Reads a unsigned char array, assumes the msb is stored first [big endian]
 */
int32 mp_read_unsigned_bin (mp_int * a, unsigned char *b, int32 c)
{
	int32		res;

/*
	Make sure there are at least two digits.
 */
	if (a->alloc < 2) {
		if ((res = mp_grow(a, 2)) != MP_OKAY) {
			return res;
		}
	}

/*
	Zero the int32.
 */
	mp_zero (a);

/*
	read the bytes in
 */
	while (c-- > 0) {
		if ((res = mp_mul_2d (a, 8, a)) != MP_OKAY) {
			return res;
		}

#ifndef MP_8BIT
		a->dp[0] |= *b++;
		a->used += 1;
#else
		a->dp[0] = (*b & MP_MASK);
		a->dp[1] |= ((*b++ >> 7U) & 1);
		a->used += 2;
#endif /* MP_8BIT */
	}
	mp_clamp (a);
	return MP_OKAY;
}

/******************************************************************************/
/* 
	Compare two ints (signed)
 */
int32 mp_cmp (mp_int * a, mp_int * b)
{
/*
	compare based on sign
 */
	if (a->sign != b->sign) {
		if (a->sign == MP_NEG) {
			return MP_LT;
		} else {
			return MP_GT;
		}
	}

/*
	compare digits
 */
	if (a->sign == MP_NEG) {
		/* if negative compare opposite direction */
		return mp_cmp_mag(b, a);
	} else {
		return mp_cmp_mag(a, b);
	}
}

/******************************************************************************/
/*
	Store in unsigned [big endian] format.
*/
int32 mp_to_unsigned_bin(psPool_t *pool, mp_int * a, unsigned char *b)
{
	int32			x, res;
	mp_int		t;

	if ((res = mp_init_copy(pool, &t, a)) != MP_OKAY) {
		return res;
	}

	x = 0;
	while (mp_iszero (&t) == 0) {
#ifndef MP_8BIT
		b[x++] = (unsigned char) (t.dp[0] & 255);
#else
		b[x++] = (unsigned char) (t.dp[0] | ((t.dp[1] & 0x01) << 7));
#endif /* MP_8BIT */
		if ((res = mp_div_2d (pool, &t, 8, &t, NULL)) != MP_OKAY) {
			mp_clear (&t);
			return res;
		}
	}
	bn_reverse (b, x);
	mp_clear (&t);
	return MP_OKAY;
}

void _mp_clear_multi(mp_int *mp0, mp_int *mp1, mp_int *mp2, mp_int *mp3,
				  mp_int *mp4, mp_int *mp5, mp_int *mp6, mp_int *mp7)
{
	int32		n		= 0;		/* Number of ok inits */

	mp_int	*tempArray[9] = {mp0, mp1, mp2, mp3, mp4, mp5, mp6, mp7, NULL};

	for (n = 0; tempArray[n] != NULL; n++) {
		mp_clear(tempArray[n]);
	}
}

/******************************************************************************/
/*
	Init a new mp_int.
*/
int32 mp_init (psPool_t *pool, mp_int * a)
{
	int32		i;
/*
	allocate memory required and clear it
 */
	a->dp = OPT_CAST(mp_digit) psMalloc(pool, sizeof (mp_digit) * MP_PREC);
	if (a->dp == NULL) {
		return MP_MEM;
	}

/*
	set the digits to zero
 */
	for (i = 0; i < MP_PREC; i++) {
		a->dp[i] = 0;
	}
/*
	set the used to zero, allocated digits to the default precision and sign
	to positive
 */
	a->used  = 0;
	a->alloc = MP_PREC;
	a->sign  = MP_ZPOS;

	return MP_OKAY;
}

/******************************************************************************/
/*
	clear one (frees).
 */
void mp_clear (mp_int * a)
{
	int32		i;
/*
	only do anything if a hasn't been freed previously
 */
	if (a->dp != NULL) {
/*
		first zero the digits
 */
		for (i = 0; i < a->used; i++) {
			a->dp[i] = 0;
		}

		/* free ram */
		psFree (a->dp);

/*
		reset members to make debugging easier
 */
		a->dp		= NULL;
		a->alloc	= a->used = 0;
		a->sign		= MP_ZPOS;
	}
}

/******************************************************************************/
/*
	Get the size for an unsigned equivalent.
 */
int32 mp_unsigned_bin_size (mp_int * a)
{
	int32		size = mp_count_bits (a);

	return	(size / 8 + ((size & 7) != 0 ? 1 : 0));
}

/******************************************************************************/
/*
	Trim unused digits 

	This is used to ensure that leading zero digits are trimed and the 
	leading "used" digit will be non-zero. Typically very fast.  Also fixes 
	the sign if there are no more leading digits
*/
void mp_clamp (mp_int * a)
{
/*
	decrease used while the most significant digit is zero.
 */
	while (a->used > 0 && a->dp[a->used - 1] == 0) {
		--(a->used);
	}

/*
	reset the sign flag if used == 0
 */
	if (a->used == 0) {
		a->sign = MP_ZPOS;
	}
}

/******************************************************************************/
/*
	Shift left by a certain bit count.
 */
int32 mp_mul_2d (mp_int * a, int32 b, mp_int * c)
{
	mp_digit	d;
	int32			res;

/*
	Copy
 */
	if (a != c) {
		if ((res = mp_copy (a, c)) != MP_OKAY) {
			return res;
		}
	}

	if (c->alloc < (int32)(c->used + b/DIGIT_BIT + 1)) {
		if ((res = mp_grow (c, c->used + b / DIGIT_BIT + 1)) != MP_OKAY) {
			return res;
		}
	}

/*
	Shift by as many digits in the bit count
 */
	if (b >= (int32)DIGIT_BIT) {
		if ((res = mp_lshd (c, b / DIGIT_BIT)) != MP_OKAY) {
			return res;
		}
	}

/*
	shift any bit count < DIGIT_BIT
 */
	d = (mp_digit) (b % DIGIT_BIT);
	if (d != 0) {
		register mp_digit *tmpc, shift, mask, r, rr;
		register int32 x;

/*
		bitmask for carries
 */
		mask = (((mp_digit)1) << d) - 1;

/*
		shift for msbs
 */
		shift = DIGIT_BIT - d;

		/* alias */
		tmpc = c->dp;

		/* carry */
		r = 0;
		for (x = 0; x < c->used; x++) {
/*
			get the higher bits of the current word
 */
			rr = (*tmpc >> shift) & mask;

/*
			shift the current word and OR in the carry
 */
			*tmpc = ((*tmpc << d) | r) & MP_MASK;
			++tmpc;

/*
			set the carry to the carry bits of the current word
 */
			r = rr;
		}

/*
		set final carry
 */
		if (r != 0) {
			c->dp[(c->used)++] = r;
		}
	}
	mp_clamp (c);
	return MP_OKAY;
}

/******************************************************************************/
/* 
	Set to zero.
 */
void mp_zero (mp_int * a)
{
	int			n;
	mp_digit	*tmp;

	a->sign = MP_ZPOS;
	a->used = 0;

	tmp = a->dp;
	for (n = 0; n < a->alloc; n++) {
		*tmp++ = 0;
	}
}

#ifdef MP_LOW_MEM
#define TAB_SIZE 32
#else
#define TAB_SIZE 256
#endif /* MP_LOW_MEM */

/******************************************************************************/
/*
	Compare maginitude of two ints (unsigned).
 */
int32 mp_cmp_mag (mp_int * a, mp_int * b)
{
	int32			n;
	mp_digit	*tmpa, *tmpb;

/*
	compare based on # of non-zero digits
 */
	if (a->used > b->used) {
		return MP_GT;
	}

	if (a->used < b->used) {
		return MP_LT;
	}

	/* alias for a */
	tmpa = a->dp + (a->used - 1);

	/* alias for b */
	tmpb = b->dp + (a->used - 1);

/*
	compare based on digits
 */
	for (n = 0; n < a->used; ++n, --tmpa, --tmpb) {
		if (*tmpa > *tmpb) {
			return MP_GT;
		}

		if (*tmpa < *tmpb) {
			return MP_LT;
		}
	}
	return MP_EQ;
}

/******************************************************************************/
/*
	computes Y == G**X mod P, HAC pp.616, Algorithm 14.85

	Uses a left-to-right k-ary sliding window to compute the modular 
	exponentiation. The value of k changes based on the size of the exponent.

	Uses Montgomery or Diminished Radix reduction [whichever appropriate]
*/
int32 mp_exptmod(psPool_t *pool, mp_int * G, mp_int * X, mp_int * P, mp_int * Y)
{

/*
	modulus P must be positive
 */
	if (P->sign == MP_NEG) {
		return MP_VAL;
	}

/*
	if exponent X is negative we have to recurse
 */
	if (X->sign == MP_NEG) {
		mp_int tmpG, tmpX;
		int32 err;

/*
		first compute 1/G mod P
 */
		if ((err = mp_init(pool, &tmpG)) != MP_OKAY) {
			return err;
		}
		if ((err = mp_invmod(pool, G, P, &tmpG)) != MP_OKAY) {
			mp_clear(&tmpG);
			return err;
		}

/*
		now get |X|
 */
		if ((err = mp_init(pool, &tmpX)) != MP_OKAY) {
			mp_clear(&tmpG);
			return err;
		}
		if ((err = mp_abs(X, &tmpX)) != MP_OKAY) {
			mp_clear(&tmpG);
			mp_clear(&tmpX);
			return err;
		}

/*
		and now compute (1/G)**|X| instead of G**X [X < 0]
 */
		err = mp_exptmod(pool, &tmpG, &tmpX, P, Y);
		mp_clear(&tmpG);
		mp_clear(&tmpX);
		return err;
	}

/*
	if the modulus is odd or dr != 0 use the fast method
 */
	if (mp_isodd (P) == 1) {
		return mp_exptmod_fast (pool, G, X, P, Y, 0);
	} else {
/*
		no exptmod for evens
 */
		return MP_VAL;
	}
}

/******************************************************************************/
/*
	Call only from mp_exptmod to make sure this fast version qualifies
*/
static int32 mp_exptmod_fast(psPool_t *pool, mp_int * G, mp_int * X,
				mp_int * P, mp_int * Y, int32 redmode)
{
	mp_int		M[TAB_SIZE], res;
	mp_digit	buf, mp;
	int32		err, bitbuf, bitcpy, bitcnt, mode, digidx, x, y, winsize;


/*
	use a pointer to the reduction algorithm.  This allows us to use
	one of many reduction algorithms without modding the guts of
	the code with if statements everywhere.
 */
	int32	(*redux)(mp_int*,mp_int*,mp_digit);

/*
	find window size
 */
	x = mp_count_bits (X);
	if (x <= 7) {
		winsize = 2;
	} else if (x <= 36) {
		winsize = 3;
	} else if (x <= 140) {
		winsize = 4;
	} else if (x <= 450) {
		winsize = 5;
	} else if (x <= 1303) {
		winsize = 6;
	} else if (x <= 3529) {
		winsize = 7;
	} else {
		winsize = 8;
	}

#ifdef MP_LOW_MEM
	if (winsize > 5) {
		winsize = 5;
	}
#endif

/*
	init M array
	init first cell
 */
	if ((err = mp_init(pool, &M[1])) != MP_OKAY) {
		return err;
	}

/*
	now init the second half of the array
 */
	for (x = 1<<(winsize-1); x < (1 << winsize); x++) {
		if ((err = mp_init(pool, &M[x])) != MP_OKAY) {
			for (y = 1<<(winsize-1); y < x; y++) {
				mp_clear(&M[y]);
			}
			mp_clear(&M[1]);
			return err;
		}
	}


/*
	now setup montgomery
 */
	if ((err = mp_montgomery_setup(P, &mp)) != MP_OKAY) {
		goto LBL_M;
	}

/*
	automatically pick the comba one if available
 */
	if (((P->used * 2 + 1) < MP_WARRAY) &&
			P->used < (1 << ((CHAR_BIT * sizeof (mp_word)) - (2 * DIGIT_BIT)))) {
		redux = fast_mp_montgomery_reduce;
	} else {
/*
		use slower baseline Montgomery method
 */
		redux = mp_montgomery_reduce;
	}

/*
	setup result
 */
	if ((err = mp_init(pool, &res)) != MP_OKAY) {
		goto LBL_M;
	}

/*
	create M table. The first half of the table is not computed
	though accept for M[0] and M[1]
*/

/*
	now we need R mod m
 */
	if ((err = mp_montgomery_calc_normalization(&res, P)) != MP_OKAY) {
		goto LBL_RES;
	}

/*
	now set M[1] to G * R mod m
 */
	if ((err = mp_mulmod(pool, G, &res, P, &M[1])) != MP_OKAY) {
		goto LBL_RES;
	}

/*
	compute the value at M[1<<(winsize-1)] by squaring
	M[1] (winsize-1) times
*/
	if ((err = mp_copy(&M[1], &M[1 << (winsize - 1)])) != MP_OKAY) {
		goto LBL_RES;
	}

	for (x = 0; x < (winsize - 1); x++) {
		if ((err = mp_sqr(pool, &M[1 << (winsize - 1)],
				&M[1 << (winsize - 1)])) != MP_OKAY) {
			goto LBL_RES;
		}
		if ((err = redux(&M[1 << (winsize - 1)], P, mp)) != MP_OKAY) {
			goto LBL_RES;
		}
	}

/*
	create upper table
 */
	for (x = (1 << (winsize - 1)) + 1; x < (1 << winsize); x++) {
		if ((err = mp_mul(pool, &M[x - 1], &M[1], &M[x])) != MP_OKAY) {
			goto LBL_RES;
		}
		if ((err = redux(&M[x], P, mp)) != MP_OKAY) {
			goto LBL_RES;
		}
	}

/*
	set initial mode and bit cnt
 */
	mode   = 0;
	bitcnt = 1;