/*
 *
 *	KerN describe in C.
 *	[Bepaul] 18 Juin 1988
 *
 */

/*
 *	Description of types and constants.
 *
 * Several conventions are used in the commentary:  
 *    A "BigNum" is the name for an infinite-precision number.  
 *    Capital letters (e.g., "N") are used to refer to the value of BigNums.
 *    The notation "N[d]" is used to refer to a particular digit of N.
 *    The notation "Size(N)" refers to the number of digits in N,
 *	 which is typically passed to the subroutine as "nl".
 *    The notation "Length(N)" refers to the number of digits in N,
 *       not including any leading zeros.
 *    The word "Base" is used for the number 2 ** BN←DIGIT←SIZE, where
 *       BN←DIGIT←SIZE is the number of bits in a single BigNum digit.
 *    The expression "BBase(N)" is used for Base ** NumDigits(N).
 *    The term "leading zeros" refers to any zeros before the most 
 *       significant digit of a number.
 */
#include "bignum.h"
#include "bignumintern.h"


/*
 *	Creation and access to type and length fields.
 */
extern char *malloc();
/* Allocates a BigNum structure and returns a pointer to it */
BigNum BnAlloc(size) int size; {
	register BigNum n;

	n = (BigNum) (malloc(sizeof(struct BigNumHeader) +
				size * sizeof(BigNumDigit))
			+ sizeof(struct BigNumHeader));
	BN←LENGTH(n) = size;
	return(n);
}

/* Allocates a BigNum, inserts its Type, and returns a pointer to it */
BigNum BnCreate(type, size) BigNumType type; int size; {
	register BigNum n;

	n = BnAlloc(size);
	BN←TYPE(n) = type;
	BnSetToZero(n, 0, size);
	return(n);
}

/* Frees a BigNum structure */
BnFree(n) BigNum n; {
	return(free(((struct BigNumHeader *) n) - 1));
}

/* Returns the BigNum's Type */
BigNumType BnGetType(n) BigNum n; {
        return(BN←TYPE(n));
}

/* Sets the BigNum's Type */
BnSetType(n, type) BigNum n; BigNumType type; {
        BN←TYPE(n) = type;
}

/* Returns the number of digits allocated for the BigNum */
BnGetSize(n) BigNum n; {
	return(BN←LENGTH(n));
}

/*
 *	non arithmetic access to digits
 */
/* Sets all the specified digits of the BigNum to 0 */
BnSetToZero(n, nd, nl) BigNum n; int nd, nl; {
	register BigNum rn = &n[nd];
	register int rnl = nl;

	while(rnl) {
		rnl--;
		*(rn++) = 0;
}	}

/* Copies N => M */
BnAssign(m, md, n, nd, nl) register BigNum m, n; register int md, nd, nl; {

	if((n != m) || (md <= nd) || (nd + nl <= md)) {
		n = &n[nd];
		m = &m[md];
		while(nl != 0) {
			nl--;
			*m++ = *n++;
		}
	} else {
		n = &n[nd + nl];
		m = &m[md + nl];
		while(nl != 0) {
			nl--;
			*--m = *--n;
}	}	}

/* Sets a single digit of N to the passed value */
BnSetDigit(n, nd, digit) BigNum n; int nd, digit; {
	n[nd] = digit;
}

/* Returns a single digit of N, if the digit can be represented in just 
 * BN←WORD←SIZE bits.
 */
BnGetDigit(n, nd) register BigNum n; int nd; {
	return(n[nd]);
}

/* Returns the number of digits of N, not counting leading zeros */
BnNumDigits(n, nd, nl) register BigNum n; register int nd, nl; {

	n = &n[nd + nl];
	while((nl != 0) && (*(--n) == 0)) nl--;
        if(nl == 0) return(1);
	return(nl);
}

/* Returns the number of leading zero bits in a digit of N */
BnNumLeadingZeroBitsInDigit(n, nd) BigNum n; int nd; {
	int register p = 0;
	register BigNumDigit c = n[nd], mask = 1 << (BN←DIGIT←SIZE - 1);

	if (c == 0) return(BN←DIGIT←SIZE);
	while((c & mask) == 0) {
		p++;
		mask >>= 1;
	}
	return(p);
}
	
/*
 *	Predicates on one digit.
 */
 
/* Returns TRUE iff the digit can be represented in just BN←WORD←SIZE bits */
Boolean BnDoesDigitFitInWord (n, nd) BigNum n; int nd; {
	/* The C compiler must evaluate the predicate at compile time */
	if(BN←DIGIT←SIZE > BN←WORD←SIZE) {
		if(n[nd] >= 1 << BN←WORD←SIZE)
			return(0);
		else	return(1);
	} else {
		return(1);
}	}

/* Returns TRUE iff N[d] = 0 */
Boolean BnIsDigitZero(n, nd) BigNum n; int nd; {
	return((n[nd] == 0));
}

/* Returns TRUE iff Base/2 <= N[d] < Base; i.e., if N[d]'s leading bit is 1 */
Boolean BnIsDigitNormalized(n, nd) BigNum n; int nd; {
	if(n[nd] & (1 << (BN←DIGIT←SIZE - 1)))
		return(1);
	else	return(0);
}

/* Returns TRUE iff N[d] is odd */
Boolean BnIsDigitOdd(n, nd) BigNum n; int nd; {
	return(n[nd] & 1);
}

/* Returns 1 if M[d] > N[d], -1 if M[d] < N[d], 0 otherwise. */
Sign BnCompareDigits(m, md, n, nd) BigNum m, n; int md, nd; {
	register BigNumDigit cm = m[md];
	register BigNumDigit cn = n[nd];

	if(cm < cn) return(-1);
	if(cm == cn) return(0);
	return(1);
}

/*
 *	Logical operations.
 */
/* Performs the computation BBase(N) - N - 1 => N */
BnComplement(n, nd, nl) register BigNum n; register int nd, nl; {
	n = &n[nd];
	while(nl != 0) {
		nl--;
		*(n++) ↑= -1;
}	}

/* Performs the logical computation M And N => M */
BnAndDigits(m, md, n, nd) BigNum m, n; int md, nd; {
	m[md] &= n[nd];
}

/* Performs the logical computation M Or N => M */
BnOrDigits(m, md, n, nd) BigNum m, n; int md, nd; {
	m[md] |= n[nd];
}

/* Performs the logical computation M Xor N => M */
BnXorDigits(m, md, n, nd) BigNum m, n; int md, nd; {
	m[md] ↑= n[nd];
}

/*
 *	Shift operations
 */
/* Shifts M left by "nbits", filling with 0s.  
 * The leftmost "nbits" of M replace the digit N[d].
 * Assumes 0 <= nbits < BN←DIGIT←SIZE. 
 */
BnShiftLeft(m, md, ml, n, nd, nbits)
		register BigNum m; register int ml, nbits;
			 BigNum n;	    int md, nd; {

	if(nbits == 0) {
		n[nd] = 0;
	} else {
		register int rnbits = BN←DIGIT←SIZE - nbits;
		register BigNumDigit res = 0, save;

		m = &m[md];
		while(ml != 0) {
			ml--;
			save = *m;
			*m++ = (save << nbits) | res;
			res = save >> rnbits;
		}
		n[nd] = res;
}	}

/* Shifts M right by "nbits", filling with 0s.  
 * The rightmost "nbits" of M replace the digit N[d].
 * Assumes 0 <= nbits < BN←DIGIT←SIZE. 
 */
BnShiftRight(m, md, ml, n, nd, nbits)
		register BigNum m; register int ml, nbits;
			 BigNum n;	    int md, nd; {

	if(nbits == 0) {
		n[nd] = 0;
	} else {
		register int lnbits = BN←DIGIT←SIZE - nbits;
		register BigNumDigit res = 0, save;

		m = &m[md + ml];
		while(ml != 0) {
			ml--;
			save = *(--m);
			*m = (save >> nbits) | res;
			res = save << lnbits;
		}
		n[nd] = res;
}	}

/*
 *	Additions.
 */
/* Performs the sum N + CarryIn => N.  Returns the CarryOut. */
BigNumCarry BnAddCarry(n, nd, nl, carryin)
		register BigNum n;
		register int nd, nl;
		BigNumCarry carryin; {

	if(carryin == 0) return(0);
	if(nl == 0) return(1);
	n = &n[nd];
	while((nl != 0) && !(++(*n++))) nl--;
	if(nl != 0) return(0);
	return(1);
}

/* Performs the sum M + N + CarryIn => M.  Returns the CarryOut. 
 * Assumes Size(M) >= Size(N).
 */
BigNumCarry BnAdd(m, md, ml, n, nd, nl, carryin)
		register BigNum m, n;
		register int ml, nl;
		int md, nd;
		BigNumCarry carryin; {
	register BigNumProduct d = carryin;

	m = &m[md]; n = &n[nd];
	ml -= nl;
	while(nl != 0) {
		nl--;
		if(sizeof(BigNumProduct) > sizeof(BigNumDigit)) {
			d += *m + *(n++);
			*(m++) = d;
			d >>= BN←DIGIT←SIZE;
		} else {
				register BigNumProduct save;
			save = *m;
			d += save;
			if(d < save)  {
				*(m++) = *(n++);
				d = 1;
			} else {
				save = *(n++);
				d += save;
				*(m++) = d;
				d = (d < save) ? 1 : 0;
	}	}	}
	if(d == 0) return(0);
	if(ml == 0) return(1);
	while((ml != 0) && !(++(*m++))) ml--;
	if(ml != 0) return(0);
	return(1);
}

/*
 *	Subtraction.
 */
/* Performs the difference N + CarryIn - 1 => N.  Returns the CarryOut. */
BigNumCarry BnSubtractBorrow(n, nd, nl, carryin)
		register BigNum n;
		register int nd, nl;
		BigNumCarry carryin; {

	if(carryin == 1) return(1);
	if(nl == 0) return(0);
	n = &n[nd];
	while((nl != 0) && !((*n++)--)) nl--;
	if(nl != 0) return(1);
	return(0);
}

/* Performs the difference M - N + CarryIn - 1 => M.  Returns the CarryOut.
 * Assumes Size(M) >= Size(N).
 */
BigNumCarry BnSubtract(m, md, ml, n, nd, nl, carryin)
		register BigNum m, n;
		register int ml, nl;
		int md, nd;
		BigNumCarry carryin; {
	register BigNumProduct d = carryin;

	ml -= nl;
	m = &m[md]; n = &n[nd];
	while(nl != 0) {
		nl--;
		if(sizeof(BigNumProduct) > sizeof(BigNumDigit)) {
				register BigNumDigit invn;
			invn = *(n++) ↑ -1;
			d += *m + invn;
			*(m++) = d;
			d >>= BN←DIGIT←SIZE;
		} else {
				register BigNumProduct save;
				register BigNumDigit invn;
			save = *m;
			invn = *(n++) ↑ -1;
			d += save;
			if(d < save)  {
				*(m++) = invn;
				d = 1;
			} else {
				d += invn;
				*(m++) = d;
				d = (d < invn) ? 1 : 0;
	}	}	}
	if(d == 1) return(1);
	if(ml == 0) return(0);
	while((ml != 0) && !((*m++)--)) ml--;
	if(ml != 0) return(1);
	return(0);
}

/*
 *	Multiplication.
 */
/* Performs the product:
 *    Q = P + M * N[d]
 *    BB = BBase(P)
 *    Q mod BB => P
 *    Q div BB => CarryOut
 * Returns the CarryOut.  Assumes Size(P) >= Size(M) + 1.
 */
#define LOW(x) (x & ((1 << (BN←DIGIT←SIZE / 2)) - 1))
#define HIGHT(x) (x >> (BN←DIGIT←SIZE / 2))
#define L2H(x) (x << (BN←DIGIT←SIZE / 2))
BigNumCarry BnMultiplyDigit(p, pd, pl, m, md, ml, n, nd)
		register BigNum p, m;
		register int pl, ml;
		BigNum n;
		int pd, md, nd; {
	register BigNumDigit c = n[nd];
	register BigNumProduct d = 0;


	if(c == 0) return(0);
	if(c == 1) return(BnAdd(p, pd, pl, m, md, ml, 0));
	p = &p[pd];
	m = &m[md];
	pl -= ml;
	while(ml != 0) {
		ml--;
		if(sizeof(BigNumProduct) > sizeof(BigNumDigit)) {
			d += *p + (c * (*(m++)));
			*(p++) = d;
			d >>= BN←DIGIT←SIZE;
		} else {
			BigNumDigit save, Lm, Hm, Lc, Hc, X0, X1, X2, X3;

			save = *(m++);
			Lc = LOW(c); Hc = HIGHT(c);
			Lm = LOW(save); Hm = HIGHT(save);
			X0 = Lc * Lm; X1 = Lc * Hm;
			X2 = Hc * Lm; X3 = Hc * Hm;
			d += X0;	if(d < X0) X3++;
			d += L2H(X1);	if(d < L2H(X1)) X3++;
			d += L2H(X2);	if(d < L2H(X2)) X3++;
			d += *p;	if(d < *p) X3++;
			*(p++) = d;
			d = X3 + HIGHT(X1) + HIGHT(X2);
	}	}
	if(sizeof(BigNumProduct) > sizeof(BigNumDigit)) {
		while(pl != 0) {
			pl--;
			d += *p;
			*(p++) = d;
			d >>= BN←DIGIT←SIZE;
		}
		return(d);
	} else {
		BigNumDigit save;

		save = *p;
		d += save;
		*(p++) = d;
		if(d >= save) 	return(0);
		pl--;
		while((pl != 0) && !(++(*p++))) pl--;
		if(pl != 0) return(0);
		return(1);
}	}

/*
 *	Division.
 */
/* Performs the quotient:
 *    N div D[d] => Q
 *    N mod D[d] => R
 * Assumes leading digit of N < D[d], and D[d] > 0.
 */
extern double pow();
BnDivideDigit(q, qd, r, rd, n, nd, nl, d, dd)
		register BigNum q, n;
		BigNum r, d;
		register int nl;
		int qd, rd, nd, dd; {
	register BigNumDigit c = d[dd];
	register BigNumProduct quad;

	n = &n[nd + nl];
	nl--;
	q = &q[qd + nl];
	quad = *(--n);
	while(nl != 0) {
		nl--;
		if(sizeof(BigNumProduct) > sizeof(BigNumDigit)) {
			quad = (quad << BN←DIGIT←SIZE) | *(--n);
			*(--q) = quad / c;
			quad = quad % c;
		} else {
			printf("Not yet implemented\n");
			exit(0);
	}	}
	r[rd] = quad;
}