#include <math.h>
#include "..\h\config.h"
#include "..\h\rt.h"
#include "rproto.h"
#include <ctype.h>

#ifdef LargeInts

extern int over_flow;

/*
 *  Conventions:
 *
 *  Bignums entering this module and leaving it are too large to
 *  be represented with T_Integer.  So, externally, a given value
 *  is always T_Integer or always T_Bignum.
 *
 *  Routines outside this module operate on bignums by calling
 *  a routine like
 *
 *      bigadd(da, db, dx)
 *
 *  where da, db, and dx are pointers to tended descriptors.
 *  For the common case where one argument is a T_Integer, these
 *  call routines like
 *
 *      bigaddi(da, IntVal(*db), dx).
 *
 *  The bigxxxi routines can convert an integer to bignum form;
 *  they use itobig.
 *
 *  The routines that actually do the work take (length, address)
 *  pairs specifying unsigned base-B digit strings.  The sign handling
 *  is done in the bigxxx routines.
 */

/*
 * Type for doing arithmetic on (2 * NB)-bit nonnegative numbers.
 *  Normally unsigned but may be signed (with NB reduced appropriately)
 *  if unsigned arithmetic is slow.
 */

/* The bignum radix, B */

#define B            ((word)1 << NB)

/* Bignum digits in a word */

#define WORDLEN  (WordBits / NB + (WordBits % NB != 0))

/* size of a bignum block that will hold an integer */

#define INTBIGBLK  sizeof(struct b_bignum) + sizeof(DIGIT) * WORDLEN

/* lo(uword d) :            the low digit of a uword
   hi(uword d) :            the rest, d is unsigned
   signed_hi(uword d) :     the rest, d is signed
   dbl(DIGIT a, DIGIT b) : the two-digit uword [a,b] */

#define lo(d)        ((d) & (B - 1))
#define hi(d)        ((uword)(d) >> NB)
#define dbl(a,b)     (((uword)(a) << NB) + (b))

#if ((-1) >> 1) < 0
#define signed_hi(d) ((word)(d) >> NB)
#else
#define signbit      ((uword)1 << (WordBits - NB - 1))
#define signed_hi(d) ((word)((((uword)(d) >> NB) ^ signbit) - signbit))
#endif

/* BigNum(dptr dp) : the struct b_bignum pointed to by dp */

#define BigNum(dp)   ((struct b_bignum *)&BlkLoc(*dp)->bignumblk)

/* LEN(struct b_bignum *b) : number of significant digits */

#define LEN(b)       ((b)->lsd - (b)->msd + 1)

/* DIG(struct b_bignum *b, word i): pointer to ith most significant digit */

#define DIG(b,i)     (&(b)->digits[(b)->msd+(i)])

/* ceil, ln: ceil may be 1 too high in case ln is inaccurate */

#undef ceil
#define ceil(x)      ((word)((x) + 1.01))
#define ln(n)        (log((double)n))

/* determine the number of words needed for a bignum block with n digits */

#define BigNeed(n)   ( ((sizeof(struct b_bignum) + ((n) - 1) * sizeof(DIGIT)) \
		       + WordSize - 1) & -WordSize )

/* copied from rconv.c */

#if !EBCDIC
#define tonum(c)     (isdigit(c) ? (c)-'0' : 10+(((c)|(040))-'a'))
#else					/* !EBCDIC */
int tonum(c)
int c;
{
#ifdef SASC
   const static char *alphanum = "0123456789abcdefghijklmnopqrstuvwxyz" ;
   char *where;

   where = memchr(alphanum, tolower(c), 36);
   if (where == 0) return -1;
   return where - alphanum;
#else					/* SASC */
   if(isdigit(c)) return (c - '0');
   if( (c | ' ') >= 'A' & (c | ' ') <= 'I') return( c - 'A' );
   if( (c | ' ') >= 'J' & (c | ' ') <= 'R') return( (c - 'J') + 9 );
   if( (c | ' ') >= 'S' & (c | ' ') <= 'Z') return( (c - 'S') + 18 );
   return 0;
#endif					/* SASC */
   }
#endif					/* !EBCDIC */

/* copied from oref.c */

#define RandVal (RanScale*(k_random=(RandA*(long)k_random+RandC)&0x7fffffffL))


hidden int mkdesc	Params((struct b_bignum *x, dptr dx));
hidden novalue itobig	Params((word i, struct b_bignum *x, dptr dx));

hidden novalue decout	Params((FILE *f, DIGIT *n, word l));

hidden int bigaddi	Params((dptr da, word i, dptr dx));
hidden int bigsubi	Params((dptr da, word i, dptr dx));
hidden int bigmuli	Params((dptr da, word i, dptr dx));
hidden int bigdivi	Params((dptr da, word i, dptr dx));
hidden int bigmodi	Params((dptr da, word i, dptr dx));
hidden int bigpowi	Params((dptr da, word i, dptr dx));
hidden int bigpowii	Params((word a, word i, dptr dx));
hidden word bigcmpi	Params((dptr da, word i));

hidden DIGIT add1	Params((DIGIT *u, DIGIT *v, DIGIT *w, word n));
hidden word sub1	Params((DIGIT *u, DIGIT *v, DIGIT *w, word n));
hidden novalue mul1	Params((DIGIT *u, DIGIT *v, DIGIT *w, word n, word m));
hidden novalue div1	
		Params((DIGIT *a, DIGIT *b, DIGIT *q, DIGIT *r, word m, word n));
hidden novalue compl1	Params((DIGIT *u, DIGIT *w, word n));
hidden word cmp1	Params((DIGIT *u, DIGIT *v, word n));
hidden DIGIT addi1	Params((DIGIT *u, word k, DIGIT *w, word n));
hidden novalue subi1	Params((DIGIT *u, word k, DIGIT *w, word n));
hidden DIGIT muli1	Params((DIGIT *u, word k, int c, DIGIT *w, word n));
hidden DIGIT divi1	Params((DIGIT *u, word k, DIGIT *w, word n));
hidden DIGIT shifti1	Params((DIGIT *u, word k, DIGIT c, DIGIT *w, word n));
hidden word cmpi1	Params((DIGIT *u, word k, word n));

hidden novalue bdzero	Params((DIGIT *dest, word l));
hidden novalue bdcopy	Params((DIGIT *src, DIGIT *dest, word l));

/*
 * mkdesc -- put value into a descriptor
 */

static int mkdesc(x, dx)
struct b_bignum *x;
dptr dx;
{
   word xlen, cmp;
   static DIGIT maxword[WORDLEN] = { 1 << ((WordBits - 1) % NB) };

   /* suppress leading zero digits */

   while (x->msd != x->lsd && *DIG(x,0) == 0)
      x->msd++;

   /* put it into a word if it fits, otherwise return the bignum */

   xlen = LEN(x);

   if (xlen < WORDLEN ||
       (xlen == WORDLEN && ((cmp = cmp1(DIG(x,0), maxword, WORDLEN)) < 0 ||
        (cmp == (word)0 && x->sign)))) {
      word val = -(word)*DIG(x,0);
      word i;

      for (i = x->msd; ++i <= x->lsd; )
         val = (val << NB) - x->digits[i];
      if (!x->sign)
	 val = -val;
      dx->dword = D_Integer;
      IntVal(*dx) = val;
      }
   else {
      dx->dword = D_Bignum;
      BlkLoc(*dx) = (union block *)x;
      }
   return Success;
}

/*
 *  i -> big
 */

static novalue itobig(i, x, dx)
word i;
struct b_bignum *x;
dptr dx;
{
   x->lsd = WORDLEN - 1;
   x->msd = WORDLEN;
   x->sign = 0;

   if (i == 0) {
      x->msd--;
      *DIG(x,0) = 0;
      }
   else if (i < 0) {
      word d = lo(i);

      if (d != 0) {
         d = B - d;
         i += B;
         }
      i = - signed_hi(i);
      x->msd--;
      *DIG(x,0) = d;
      x->sign = 1;
      }
            
   while (i != 0) {
      x->msd--;
      *DIG(x,0) = lo(i);
      i = hi(i);
      }

   dx->dword = D_Bignum;
   BlkLoc(*dx) = (union block *)x;
}

/*
 *  string -> bignum 
 */

word bigradix(sign, r, s, dx)
int sign;                      /* '-' or not */
int r;                          /* radix 2 .. 36 */
char *s;                        /* input string */
dptr dx;                        /* output T_Integer or T_Bignum */
{
   struct b_bignum *b;
   DIGIT *bd;
   word len;
   int c;

   if (r == 0)
      return CvtFail;
   len = ceil(strlen(s) * ln(r) / ln(B));
   if (blkreq(BigNeed(len)) == Error)
      return CvtFail;
   b = alcbignum(len);
   bd = DIG(b,0);

   bdzero(bd, len);

   if (r < 2 || r > 36)
      return CvtFail;

   for (c = *s++; isalnum(c); c = *s++) {
      c = tonum(c);
      if (c >= r)
	 return CvtFail;
      muli1(bd, (word)r, c, bd, len);
      }

   while (isspace(c))
      c = *s++;

   if (c != '\0')
      return CvtFail;

   if (sign == '-')
      b->sign = 1;

   /* put value into dx and return the type */

   (void)mkdesc(b, dx);
   return Type(*dx);
}

/*
 *  bignum -> real
 */

double bigtoreal(da)
dptr da;
{
   word i;
   double r = 0;
   struct b_bignum *b = &BlkLoc(*da)->bignumblk;

   for (i = b->msd; i <= b->lsd; i++)
      r = r * B + b->digits[i];

   return (b->sign ? -r : r);
}
/*
 *  real -> bignum
 */

int realtobig(da, dx)
dptr da, dx;
{
   double x = BlkLoc(*da)->realblk.realval;
   struct b_bignum *b;
   word i, blen;
   word d;

   if (x > 0.9999 * MinLong && x < 0.9999 * MaxLong) {
      MakeInt((word)x, dx);
      return Success;		/* got lucky; a simple integer suffices */
      }

   x = x > 0 ? x : -x;
   blen = ln(x) / ln(B) + 0.99;
   for (i = 0; i < blen; i++)
      x /= B;
   if (x >= 1.0) {
      x /= B;
      blen += 1;
      }

   if (blkreq(BigNeed(blen)) == Error)
      return Error;
   b = alcbignum(blen);
   for (i = 0; i < blen; i++) {
      d = (x *= B);
      *DIG(b,i) = d;
      x -= d;
      }
     
   b->sign = x < 0;
   return mkdesc(b, dx);
}

/*
 *  bignum -> string
 */

int bigtos(da, dx)
dptr da, dx;
{
   struct b_bignum *a, *temp;
   word alen = LEN(BigNum(da));
   word slen = ceil(alen * ln(B) / ln(10));
   char *p, *q;

   if (strreq(slen) == Error || blkreq(BigNeed(alen)) == Error) 
      return CvtFail;
   a = BigNum(da);
   temp = alcbignum(alen);
   if (a->sign)
      slen++;
   q = alcstr("",slen);
   bdcopy(DIG(a,0), DIG(temp,0), alen);
   p = q += slen;
   while (cmpi1(DIG(temp,0), (word)0, alen))
      *--p = '0' + divi1(DIG(temp,0), (word)10, DIG(temp,0), alen);
   if (a->sign)
      *--p = '-';
   StrLen(*dx) = q - p;
   StrLoc(*dx) = p;
   return Cvt;
}

/*
 *  bignum -> file 
 */

novalue bigprint(f, da)
FILE *f;
dptr da;
{
   struct b_bignum *a, *temp;
   word alen = LEN(BigNum(da));
   word slen, dlen;

   slen = (BlkLoc(*da)->bignumblk.lsd - BlkLoc(*da)->bignumblk.msd + 1);
   dlen = slen * NB * 0.3010299956639812;	/* 1 / log2(10) */
   if (dlen > MaxDigits) {
      fprintf(f, "integer(~%ld)",dlen - 2);	/* center estimate */
      return;
      }

   if (blkreq(BigNeed(alen)) == Error) {
      fatalerr(0, NULL);		/* not worth passing this one back */
      }
   temp = alcbignum(alen);
   a = BigNum(da);
   bdcopy(DIG(a,0), DIG(temp,0), alen);
   if (a->sign)
      putc('-', f);
   decout(f, DIG(temp,0), alen);
}

static novalue decout(f, n, l)
FILE *f;
DIGIT *n;
word l;
{
   word i = divi1(n, (word)10, n, l);

   if (cmpi1(n, (word)0, l))
      decout(f, n, l);
   putc('0' + i, f);
}

/*
 *  da -> dx
 */

int cpbignum(da, dx)
dptr da, dx;
{
   struct b_bignum *a;
   word alen = LEN(BigNum(da));
   struct b_bignum *x;

   if (blkreq(BigNeed(alen)) == Error)
      return Error;
   x = alcbignum(alen);
   a = BigNum(da);
   bdcopy(DIG(a,0), DIG(x,0), alen);
   x->sign = a->sign;
   return mkdesc(x, dx);
}

/*
 *  da + db -> dx
 */

int bigadd(da, db, dx)
dptr da, db;
dptr dx;
{
   struct b_bignum *x, *a, *b;
   word alen, blen;
   word c;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
	 return Error;
      a = BigNum(da);
      b = BigNum(db);
      if (a->sign == b->sign) {
         if (alen > blen) {
            x = alcbignum(alen + 1);
            c = add1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen+1), blen);
            *DIG(x,0) = addi1(DIG(a,0), c, DIG(x,1), alen-blen);
            }
         else if (alen == blen) {
            x = alcbignum(alen + 1);
            *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
            }
         else {
            x = alcbignum(blen + 1);
            c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
            *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
            }
         x->sign = a->sign;
         }
      else {
         if (alen > blen) {
            x = alcbignum(alen);
            c = sub1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen), blen);
            subi1(DIG(a,0), -c, DIG(x,0), alen-blen);
            x->sign = a->sign;
            }
         else if (alen == blen) {
            x = alcbignum(alen);
            if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
               (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
               x->sign = a->sign;
               }
            else {
               (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
               x->sign = b->sign;
               }
            }
         else {
            x = alcbignum(blen);
            c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
            subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
            x->sign = b->sign;
            }
         }
      return mkdesc(x, dx);
      }
   else if (Type(*da) == T_Bignum)    /* bignum + integer */
      return bigaddi(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum)    /* integer + bignum */
      return bigaddi(db, IntVal(*da), dx);
   else {                             /* integer + integer */
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      return bigaddi(&td, IntVal(*db), dx);
      }
}

/*
 *  da - db -> dx
 */ 

int bigsub(da, db, dx)
dptr da, db, dx;
{
   struct descrip td;
   char tdigits[INTBIGBLK];
   struct b_bignum *a, *b, *x;
   word alen, blen;
   word c;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
	 return Error;
      a = BigNum(da);
      b = BigNum(db);
      if (a->sign != b->sign) {
         if (alen > blen) {
            x = alcbignum(alen + 1);
            c = add1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen+1), blen);
            *DIG(x,0) = addi1(DIG(a,0), c, DIG(x,1), alen-blen);
            }
         else if (alen == blen) {
            x = alcbignum(alen + 1);
            *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
            }
         else {
            x = alcbignum(blen + 1);
            c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
            *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
            }
         x->sign = a->sign;
         }
      else {
         if (alen > blen) {
            x = alcbignum(alen);
            c = sub1(DIG(a,alen-blen), DIG(b,0), DIG(x,alen-blen), blen);
            subi1(DIG(a,0), -c, DIG(x,0), alen-blen);
            x->sign = a->sign;
            }
         else if (alen == blen) {
            x = alcbignum(alen);
            if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
               (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
               x->sign = a->sign;
               }
            else {
               (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
               x->sign = 1 ^ b->sign;
               }
            }
         else {
            x = alcbignum(blen);
            c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
            subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
            x->sign = 1 ^ b->sign;
            }
         }
      return mkdesc(x, dx);
      }
   else if (Type(*da) == T_Bignum)     /* bignum - integer */
      return bigsubi(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum) {   /* integer - bignum */
      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(&td));
      blen = LEN(BigNum(db));
      if (blkreq(BigNeed(alen > blen ? alen + 1 : blen + 1)) == Error)
	 return Error;
      a = BigNum(&td);
      b = BigNum(db);
      if (a->sign != b->sign) {
         if (alen == blen) {
            x = alcbignum(alen + 1);
            *DIG(x,0) = add1(DIG(a,0), DIG(b,0), DIG(x,1), alen);
            }
         else {
            x = alcbignum(blen + 1);
            c = add1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen+1), alen);
            *DIG(x,0) = addi1(DIG(b,0), c, DIG(x,1), blen-alen);
            }
         x->sign = a->sign;
         }
      else {
         if (alen == blen) {
            x = alcbignum(alen);
            if (cmp1(DIG(a,0), DIG(b,0), alen) > 0) {
               (void)sub1(DIG(a,0), DIG(b,0), DIG(x,0), alen);
               x->sign = a->sign;
               }
            else {
               (void)sub1(DIG(b,0), DIG(a,0), DIG(x,0), alen);
               x->sign = 1 ^ b->sign;
               }
            }
         else {
            x = alcbignum(blen);
            c = sub1(DIG(b,blen-alen), DIG(a,0), DIG(x,blen-alen), alen);
            subi1(DIG(b,0), -c, DIG(x,0), blen-alen);
            x->sign = 1 ^ b->sign;
            }
         }
      return mkdesc(x, dx);
      }
   else {                              /* integer - integer */
      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      return bigsubi(&td, IntVal(*db), dx);
      }
      
}

/*
 *  da * db -> dx
 */

int bigmul(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x;
   word alen, blen;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      if (blkreq(BigNeed(alen + blen)) == Error)
	 return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(alen + blen);
      mul1(DIG(a,0), DIG(b,0), DIG(x,0), alen, blen);
      x->sign = a->sign ^ b->sign;
      return mkdesc(x, dx);
      }
   else if (Type(*da) == T_Bignum)    /* bignum * integer */
      return bigmuli(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum)    /* integer * bignum */
      return bigmuli(db, IntVal(*da), dx);
   else {                             /* integer * integer */
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      return bigmuli(&td, IntVal(*db), dx);
      }
}

/*
 *  da / db -> dx
 */
 
int bigdiv(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x;
   word alen, blen;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      if (alen < blen) {
         MakeInt(0, dx);
         return Success;
         }
      if (blkreq(BigNeed(alen-blen+1)+BigNeed(alen+1)+BigNeed(blen)) == Error)
	 return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(alen - blen + 1);
      if (blen == 1)
         divi1(DIG(a,0), (word)*DIG(b,0), DIG(x,0), alen);
      else
         div1(DIG(a,0), DIG(b,0), DIG(x,0), NULL, alen-blen, blen);
      x->sign = a->sign ^ b->sign;
      return mkdesc(x, dx);
      }
   else if (Type(*da) == T_Bignum)     /* bignum / integer */
      return bigdivi(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum) {   /* integer / bignum */
      MakeInt(0, dx);
      return Success;
      }
   /* not called for integer / integer */
}

/*
 *  da % db -> dx
 */

int bigmod(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x, *temp;
   word alen, blen;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      if (alen < blen) {
         cpbignum(da, dx);
         return Success;
         }
      if (blkreq(BigNeed(blen)+BigNeed(alen+1)+BigNeed(blen)) == Error)
	 return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(blen);
      if (blen == 1) {
	 temp = alcbignum(alen);
         *DIG(x,0) = divi1(DIG(a,0), (word)*DIG(b,0), DIG(temp,0), alen);
         }
      else
         div1(DIG(a,0), DIG(b,0), NULL, DIG(x,0), alen-blen, blen);
      x->sign = a->sign;
      return mkdesc(x, dx);
      }
   else if (Type(*da) == T_Bignum)     /* bignum % integer */
      return bigmodi(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum) {   /* integer % bignum */
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      cpbignum(&td, dx);
      return Success;
      }
   /* not called for integer % integer */
}

/*
 *  -i -> dx
 */

int bigneg(da, dx)
dptr da, dx;
{
   struct descrip td;
   char tdigits[INTBIGBLK];

   itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
   BigNum(&td)->sign ^= 1;
   return cpbignum(&td, dx);
}

/*
 *  da ^ db -> dx
 */

int bigpow(da, db, dx)
dptr da, db, dx;
{
   word n;

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      if (BigNum(db)->sign) {
         MakeInt(0, dx);
         }
      else {
         n = LEN(BigNum(db)) * NB;
         /* scan bits left to right.  skip leading 1. */
         while (--n >= 0)
            if ((*DIG(BigNum(db), n / NB) & (1 << (n % NB))))
               break;
         /* then, for each zero, square the partial result;
          *  for each one, square it and multiply it by a */
         *dx = *da;
         while (--n >= 0) {
            if (bigmul(dx, dx, dx) == Error)
	       return Error;
            if ((*DIG(BigNum(db), n / NB) & (1 << (n % NB))))
               if (bigmul(dx, da, dx) == Error)
		  return Error;
	    }
         }
      return Success;
      }
   else if (Type(*da) == T_Bignum)    /* bignum ^ integer */
      return bigpowi(da, IntVal(*db), dx);
   else if (Type(*db) == T_Bignum)    /* integer ^ bignum */
      return bigpowii(IntVal(*da), (word)bigtoreal(db), dx);
   else                               /* integer ^ integer */
      return bigpowii(IntVal(*da), IntVal(*db), dx);
}

/*
 *  iand(da, db) -> dx
 */

int bigand(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x;
   word alen, blen, xlen;
   word i;
   struct b_bignum *tad, *tbd;
   DIGIT *ad, *bd;
   struct descrip td;
   char tdigits[INTBIGBLK];

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(xlen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(xlen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] & bd[i];

      if (a->sign & b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*da) == T_Bignum) {   /* iand(bignum,integer) */
      itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(&td));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(alen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(&td);
      x = alcbignum(alen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] & bd[i];

      if (a->sign & b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*db) == T_Bignum) {   /* iand(integer,bignum) */
      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(&td));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(blen)) == Error)
         return Error;
      a = BigNum(&td);
      b = BigNum(db);
      x = alcbignum(blen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] & bd[i];

      if (a->sign & b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   /* not called for iand(integer,integer) */

   return mkdesc(x, dx);
}

/*
 *  ior(da, db) -> dx
 */

int bigor(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x;
   word alen, blen, xlen;
   word i;
   struct b_bignum *tad, *tbd;
   DIGIT *ad, *bd;
   struct descrip td;
   char tdigits[INTBIGBLK];

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(xlen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(xlen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] | bd[i];

      if (a->sign | b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*da) == T_Bignum) {   /* ior(bignum,integer) */
      itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(&td));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(alen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(&td);
      x = alcbignum(alen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] | bd[i];

      if (a->sign | b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*db) == T_Bignum) {   /* ior(integer,bignym) */
      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(&td));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(blen)) == Error)
         return Error;
      a = BigNum(&td);
      b = BigNum(db);
      x = alcbignum(blen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }
        
      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] | bd[i];

      if (a->sign | b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   /* not called for ior(integer,integer) */

   return mkdesc(x, dx);
}

/*
 *  xor(da, db) -> dx
 */

int bigxor(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *b, *x;
   word alen, blen, xlen;
   word i;
   struct b_bignum *tad, *tbd;
   DIGIT *ad, *bd;
   struct descrip td;
   char tdigits[INTBIGBLK];

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(xlen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(db);
      x = alcbignum(xlen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }

      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] ^ bd[i];

      if (a->sign ^ b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*da) == T_Bignum) {   /* ixor(bignum,integer) */
      itobig(IntVal(*db), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(da));
      blen = LEN(BigNum(&td));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(alen)) == Error)
         return Error;
      a = BigNum(da);
      b = BigNum(&td);
      x = alcbignum(alen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }

      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] ^ bd[i];

      if (a->sign ^ b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   else if (Type(*db) == T_Bignum) {   /* ixor(integer,bignum) */
      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(&td));
      blen = LEN(BigNum(db));
      xlen = alen > blen ? alen : blen;
      if (blkreq(3 * BigNeed(blen)) == Error)
         return Error;
      a = BigNum(&td);
      b = BigNum(db);
      x = alcbignum(blen);

      if (alen == xlen && !a->sign)
         ad = DIG(a,0);
      else {
         tad = alcbignum(xlen);
         ad = DIG(tad,0);
         bdzero(ad, xlen - alen);
         bdcopy(DIG(a,0), &ad[xlen-alen], alen);
         if (a->sign)
	    compl1(ad, ad, xlen);
         }

      if (blen == xlen && !b->sign)
         bd = DIG(b,0);
      else {
         tbd = alcbignum(xlen);
         bd = DIG(tbd,0);
         bdzero(bd, xlen - blen);
         bdcopy(DIG(b,0), &bd[xlen-blen], blen);
         if (b->sign)
	    compl1(bd, bd, xlen);
         }

      for (i = 0; i < xlen; i++)
         *DIG(x,i) = ad[i] ^ bd[i];

      if (a->sign ^ b->sign) {
         x->sign = 1;
         compl1(DIG(x,0), DIG(x,0), xlen);
         }
      }
   /* not called for ixor(integer,integer) */

   return mkdesc(x, dx);
}

/*
 *  bigshift(da, db) -> dx
 */

int bigshift(da, db, dx)
dptr da, db, dx;
{
   struct b_bignum *a, *x;
   word alen;
   word r = IntVal(*db) % NB;
   word q = (r >= 0 ? IntVal(*db) : (IntVal(*db) - (r += NB))) / NB;
   word xlen;
   struct b_bignum *tad;
   DIGIT *ad;

   if (Type(*da) == T_Bignum) {
      alen = LEN(BigNum(da));
      xlen = alen + q + 1;
      if (xlen <= 0) {
         MakeInt(0, dx);
         return Success;
         }
      else {
         if (blkreq(BigNeed(xlen) + BigNeed(alen)) == Error)
	    return Error;
         a = BigNum(da);
         x = alcbignum(xlen);
         ad = DIG(a,0);

         if (q >= 0) {
            *DIG(x,0) = shifti1(ad, r, (DIGIT)0, DIG(x,1), alen);
            bdzero(DIG(x,alen+1), q);
            x->sign = a->sign;
            }
         else  {
            *DIG(x,0) = shifti1(ad, r, ad[alen+q] >> (NB-r), DIG(x,1), alen+q);
            }

         return mkdesc(x, dx);
	 }
      }
   else {                              /* da is integer */
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(IntVal(*da), (struct b_bignum *)tdigits, &td);
      alen = LEN(BigNum(&td));
      xlen = alen + q + 1;
      if (xlen <= 0) {
         MakeInt(0, dx);
         return Success;
         }
      else {
         if (blkreq(BigNeed(xlen) + BigNeed(alen)) == Error)
	    return Error;
         a = BigNum(&td);
         x = alcbignum(xlen);
         ad = DIG(a,0);

         if (q >= 0) {
            *DIG(x,0) = shifti1(ad, r, (DIGIT)0, DIG(x,1), alen);
            bdzero(DIG(x,alen+1), q);
            x->sign = a->sign;
            }
         else 
            *DIG(x,0) = shifti1(ad, r, ad[alen+q] >> (NB-r), DIG(x,1), alen+q);

         return mkdesc(x, dx);
         }
      }
}

/*
 *  negative if da < db
 *  zero if da == db
 *  positive if da > db
 */

word bigcmp(da, db)
dptr da, db;
{
   struct b_bignum *a = BigNum(da);
   struct b_bignum *b = BigNum(db);
   word alen, blen; 

   if (Type(*da) == T_Bignum && Type(*db) == T_Bignum) {
      if (a->sign != b->sign)
         return (b->sign - a->sign);
      alen = LEN(a);
      blen = LEN(b);
      if (alen != blen)
         return (a->sign ? blen - alen : alen - blen);

      if (a->sign)
         return cmp1(DIG(b,0), DIG(a,0), alen);
      else
         return cmp1(DIG(a,0), DIG(b,0), alen);
      }
   else if (Type(*da) == T_Bignum)    /* cmp(bignum, integer) */
      return bigcmpi(da, IntVal(*db));
   else                               /* cmp(integer, bignum) */
      return -bigcmpi(db, IntVal(*da));
}

/*
 *  ?da -> dx
 */  

int bigrand(da, dx)
dptr da, dx;
{
   struct b_bignum *a;
   word alen = LEN(BigNum(da));
   struct b_bignum *x;
   struct b_bignum *td;
   DIGIT *d;
   word i;
   double rval;

   if (blkreq(4 * BigNeed(alen + 1) + 4) == Error)
      return Error;
   x = alcbignum(alen);
   td = alcbignum(alen + 1);
   d = DIG(td,0);
   a = BigNum(da);

   for (i = alen; i >= 0; i--) {
      rval = RandVal;
      d[i] = rval * B;
      }
    
   div1(d, DIG(a,0), NULL, DIG(x,0), (word)1, alen);
   addi1(DIG(x,0), (word)1, DIG(x,0), alen);
   return mkdesc(x, dx);
}

/*
 *  da + i -> dx
 */

static int bigaddi(da, i, dx)
dptr da, dx;
word i;
{
   struct b_bignum *a, *x; 
   word alen; 

   if (i < 0)
      return bigsubi(da, -i, dx);
   else if (i != (DIGIT)i) {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigadd(da, &td, dx);
      }
   else {
      alen = LEN(BigNum(da));
      if (blkreq(BigNeed(alen + 1)) == Error)
         return Error;
      a = BigNum(da);
      if (a->sign) {
	 x = alcbignum(alen);
         subi1(DIG(a,0), i, DIG(x,0), alen);
         }
      else {
         x = alcbignum(alen + 1);
         *DIG(x,0) = addi1(DIG(a,0), i, DIG(x,1), alen);
         }
      x->sign = a->sign;
      return mkdesc(x, dx);
      }
}

/*
 *  da - i -> dx
 */

static int bigsubi(da, i, dx)
dptr da, dx;
word i;
{
   struct b_bignum *a, *x; 
   word alen;

   if (i < 0)
      return bigaddi(da, -i, dx);
   else if (i != (DIGIT)i) {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigsub(da, &td, dx);
      }
   else {
      alen = LEN(BigNum(da));
      if (blkreq(BigNeed(alen + 1)) == Error)
	 return Error;
      a = BigNum(da);
      if (a->sign) {
         x = alcbignum(alen + 1);
         *DIG(x,0) = addi1(DIG(a,0), i, DIG(x,1), alen);
         }
      else {
         x = alcbignum(alen);
         subi1(DIG(a,0), i, DIG(x,0), alen);
         }
      x->sign = a->sign;
      return mkdesc(x, dx);
      }
}

/*
 *  da * i -> dx
 */

static int bigmuli(da, i, dx)
dptr da, dx;
word i;
{
   struct b_bignum *a, *x; 
   word alen;

   if (i <= -B || i >= B) {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigmul(da, &td, dx);
      }
   else {
      alen = LEN(BigNum(da));
      if (blkreq(BigNeed(alen + 1)) == Error)
	 return Error;
      a = BigNum(da);
      x = alcbignum(alen + 1);
      if (i >= 0)
         x->sign = a->sign;
      else {
         x->sign = 1 ^ a->sign;
         i = -i;
         }
      *DIG(x,0) = muli1(DIG(a,0), i, 0, DIG(x,1), alen);
      return mkdesc(x, dx);
      }
}

/*
 *  da / i -> dx
 */

static int bigdivi(da, i, dx)
dptr da, dx;
word i;
{
   struct b_bignum *a, *x; 
   word alen;

   if (i <= -B || i >= B) {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigdiv(da, &td, dx);
      }
   else {
      alen = LEN(BigNum(da));
      if (blkreq(BigNeed(alen)) == Error)
	 return Error;
      a = BigNum(da);
      x = alcbignum(alen);
      if (i >= 0)
         x->sign = a->sign;
      else {
         x->sign = 1 ^ a->sign;
         i = -i;
         }
      divi1(DIG(a,0), i, DIG(x,0), alen);
      return mkdesc(x, dx);
      }
}

/*
 *  da % i -> dx
 */

static int bigmodi(da, i, dx)
dptr da, dx;
word i;
{
   struct b_bignum *a, *temp;
   word alen;
   word x;

   if (i <= -B || i >= B) {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigmod(da, &td, dx);
      }
   else {
      alen = LEN(BigNum(da));
      if (blkreq(BigNeed(alen)) == Error)
	 return Error;
      a = BigNum(da);
      temp = alcbignum(alen);
      x = divi1(DIG(a,0), Abs(i), DIG(temp,0), alen);
      if (a->sign)
	 x = -x;
      MakeInt(x, dx);
      return Success;
      }
}

/*
 *  da ^ i -> dx
 */

static int bigpowi(da, i, dx)
dptr da, dx;
word i;
{
   int n = WordBits;

   if (i > 0) {
      /* scan bits left to right.  skip leading 1. */
      while (--n >= 0)
         if (i & ((word)1 << n))
	    break;
      /* then, for each zero, square the partial result;
         for each one, square it and multiply it by a */
      *dx = *da;
      while (--n >= 0) {
         if (bigmul(dx, dx, dx) == Error)
	    return Error;
         if (i & ((word)1 << n))
            if (bigmul(dx, da, dx) == Error)
	       return Error;
         }
      }
   else {
      MakeInt(0, dx);
      }
   return Success;
}

/*
 *  a ^ i -> dx
 */

static int bigpowii(a, i, dx)
word a, i;
dptr dx;
{
   word x, y;
   int n = WordBits;
   int isbig = 0;

   if (a == 0 || i <= 0) {              /* special cases */
      if (a == 0 && i <= 0)             /* 0 ^ negative -> error */
         RetError(-204, nulldesc);
      if (a == -1) {                    /* -1 ^ [odd,even] -> [-1,+1] */
         if (!(i & 1))
	    a = 1;
         }
      else if (a != 1) {                /* 1 ^ any -> 1 */
         a = 0;
         }                   /* others ^ negative -> 0 */
      MakeInt(a, dx);
      }
   else {
      struct descrip td;
      char tdigits[INTBIGBLK];

      /* scan bits left to right.  skip leading 1. */
      while (--n >= 0)
         if (i & ((word)1 << n))
	    break;
      /* then, for each zero, square the partial result;
         for each one, square it and multiply it by a */
      x = a;
      while (--n >= 0) {
         if (isbig) {
            if (bigmul(dx, dx, dx) == Error)
	       return Error;
	    }
         else {
            y = mul(x, x);
            if (!over_flow)
               x = y;
            else {
               itobig(x, (struct b_bignum *)tdigits, &td);
               if (bigmul(&td, &td, dx) == Error)
   	          return Error;
               isbig = (Type(*dx) == T_Bignum);
               } 
            }
         if (i & ((word)1 << n)) {
            if (isbig) {
               if (bigmuli(dx, a, dx) == Error)
		  return Error;
	       }
            else {
               y = mul(x, a);
               if (!over_flow)
                  x = y;
               else {
                  itobig(x, (struct b_bignum *)tdigits, &td);
                  if (bigmuli(&td, a, dx) == Error)
   		  return Error;
                  isbig = (Type(*dx) == T_Bignum);
                  }
               }
            }
         }
      if (!isbig) {
	 MakeInt(x, dx);
	 }
      }
   return Success;
}

/*
 *  negative if da < i
 *  zero if da == i
 *  positive if da > i
 */  
  
static word bigcmpi(da, i)
dptr da;
word i;
{
   struct b_bignum *a = BigNum(da);
   word alen = LEN(a);

   if (i > -B && i < B) {
      if (i >= 0)
         if (a->sign)
	    return -1;
         else
	    return cmpi1(DIG(a,0), i, alen);
      else
         if (a->sign)
	    return -cmpi1(DIG(a,0), -i, alen);
         else
	    return 1;
      }
   else {
      struct descrip td;
      char tdigits[INTBIGBLK];

      itobig(i, (struct b_bignum *)tdigits, &td);
      return bigcmp(da, &td);
      }
}


/* These are all straight out of Knuth vol. 2, Sec. 4.3.1. */

/*
 *  (u,n) + (v,n) -> (w,n)
 *
 *  returns carry, 0 or 1
 */

static DIGIT add1(u, v, w, n)
DIGIT *u, *v, *w;
word n;
{
   uword dig, carry; 
   word i;

   carry = 0;
   for (i = n; --i >= 0; ) {
      dig = (uword)u[i] + v[i] + carry;
      w[i] = lo(dig);
      carry = hi(dig);
      }
   return carry;
}

/*
 *  (u,n) - (v,n) -> (w,n)
 *
 *  returns carry, 0 or -1
 */

static word sub1(u, v, w, n)
DIGIT *u, *v, *w;
word n;
{
   uword dig, carry; 
   word i;

   carry = 0;
   for (i = n; --i >= 0; ) {
      dig = (uword)u[i] - v[i] + carry;
      w[i] = lo(dig);
      carry = signed_hi(dig);
      }
   return carry;
}

/*
 *  (u,n) * (v,m) -> (w,m+n)
 */

static novalue mul1(u, v, w, n, m)
DIGIT *u, *v, *w;
word n, m;
{
   word i, j;
   uword dig, carry;

   bdzero(&w[m], n);

   for (j = m; --j >= 0; ) {
      carry = 0;
      for (i = n; --i >= 0; ) {
         dig = (uword)u[i] * v[j] + w[i+j+1] + carry;
         w[i+j+1] = lo(dig);
         carry = hi(dig);
         }
      w[j] = carry;
      }
}

/*
 *  (a,m+n) / (b,n) -> (q,m+1) (r,n)
 *
 *  if q or r is NULL, the quotient or remainder is discarded
 */

static novalue div1(a, b, q, r, m, n)
DIGIT *a, *b, *q, *r;
word m, n;
{
   uword qhat, rhat;
   uword dig, carry;
   struct b_bignum *tu, *tv;
   DIGIT *u, *v;
   word d;
   word i, j;

   /* blkreq's done in calling routines */

   tu = alcbignum(m + n + 1);
   tv = alcbignum(n);
   u = DIG(tu,0);
   v = DIG(tv,0);

   /* D1 */
   for (d = 0; d < NB; d++)
      if (b[0] & (1 << (NB - 1 - d)))
         break;

   u[0] = shifti1(a, d, (DIGIT)0, &u[1], m+n);
   shifti1(b, d, (DIGIT)0, v, n);

   /* D2, D7 */
   for (j = 0; j <= m; j++) {
      /* D3 */
      if (u[j] == v[0]) {
         qhat = B - 1;
         rhat = (uword)v[0] + u[j+1];
         }
      else {
         uword numerator = dbl(u[j], u[j+1]);
         qhat = numerator / (uword)v[0];
         rhat = numerator % (uword)v[0];
         }

      while (rhat < B && qhat * v[1] > dbl(rhat, u[j+2])) {
         qhat -= 1;
         rhat += v[0];
         }
            
      /* D4 */
      carry = 0;
      for (i = n; i > 0; i--) {
         dig = u[i+j] - v[i-1] * qhat + carry;       /* -BSQ+B .. B-1 */
         u[i+j] = lo(dig);
         if ((uword)dig < B)
            carry = hi(dig);
         else carry = hi(dig) | -B;
         }
      carry = (word)(carry + u[j]) < 0;

      /* D5 */
      if (q)
	 q[j] = qhat;

      /* D6 */
      if (carry) {
         if (q)
	    q[j] -= 1;
         carry = 0;
         for (i = n; i > 0; i--) {
            dig = (uword)u[i+j] + v[i-1] + carry;
            u[i+j] = lo(dig);
            carry = hi(dig);
            }
         }
      }

   if (r) {
      if (d == 0)
         shifti1(&u[m+1], (word)d, (DIGIT)0, r, n);
      else
         r[0] = shifti1(&u[m+1], (word)(NB - d), u[m+n]>>d, &r[1], n - 1);
      }
}

/*
 *  - (u,n) -> (w,n)
 *
 */

static novalue compl1(u, w, n)
DIGIT *u, *w;
word n;
{
   uword dig, carry = 0;
   word i;

   for (i = n; --i >= 0; ) {
      dig = carry - u[i];
      w[i] = lo(dig);
      carry = signed_hi(dig);
      }
}

/*
 *  (u,n) : (v,n)
 */

static word cmp1(u, v, n)
DIGIT *u, *v;
word n;
{
   word i;

   for (i = 0; i < n; i++)
      if (u[i] != v[i])
         return u[i] > v[i] ? 1 : -1;
   return 0;
}

/*
 *  (u,n) + k -> (w,n)
 *
 *  k in 0 .. B-1
 *  returns carry, 0 or 1
 */

static DIGIT addi1(u, k, w, n)
DIGIT *u, *w;
word k;
word n;
{
   uword dig, carry;
   word i;
    
   carry = k;
   for (i = n; --i >= 0; ) {
      dig = (uword)u[i] + carry;
      w[i] = lo(dig);
      carry = hi(dig);
      }
   return carry;
}

/*
 *  (u,n) - k -> (w,n)
 *
 *  k in 0 .. B-1
 *  u must be greater than k
 */

static novalue subi1(u, k, w, n)
DIGIT *u, *w;
word k;
word n;
{
   uword dig, carry;
   word i;
    
   carry = -k;
   for (i = n; --i >= 0; ) {
      dig = (uword)u[i] + carry;
      w[i] = lo(dig);
      carry = signed_hi(dig);
      }
}

/*
 *  (u,n) * k + c -> (w,n)
 *
 *  k in 0 .. B-1
 *  returns carry, 0 .. B-1
 */

static DIGIT muli1(u, k, c, w, n)
DIGIT *u, *w;
word k;
int c;
word n;
{
   uword dig, carry;
   word i;

   carry = c;
   for (i = n; --i >= 0; ) {
      dig = (uword)k * u[i] + carry;
      w[i] = lo(dig);
      carry = hi(dig);
      }
   return carry;
}

/*
 *  (u,n) / k -> (w,n)
 *
 *  k in 0 .. B-1
 *  returns remainder, 0 .. B-1
 */

static DIGIT divi1(u, k, w, n)
DIGIT *u, *w;
word k;
word n;
{
   uword dig, remain;
   word i;

   remain = 0;
   for (i = 0; i < n; i++) {
      dig = dbl(remain, u[i]);
      w[i] = dig / k;
      remain = dig % k;
      }
   return remain;
}

/*
 *  ((u,n) << k) + c -> (w,n)
 *
 *  k in 0 .. NB-1
 *  c in 0 .. B-1 
 *  returns carry, 0 .. B-1
 */

static DIGIT shifti1(u, k, c, w, n)
DIGIT *u, c, *w;
word k;
word n;
{
   uword dig;
   word i;

   if (k == 0) {
      bdcopy(u, w, n);
      return 0;
      }
    
   for (i = n; --i >= 0; ) {
      dig = ((uword)u[i] << k) + c;
      w[i] = lo(dig);
      c = hi(dig);
      }
   return c;
}

/*
 *  (u,n) : k
 *
 *  k in 0 .. B-1
 */

static word cmpi1(u, k, n)
DIGIT *u;
word k;
word n;
{
   word i;

   for (i = 0; i < n-1; i++)
      if (u[i])
	 return 1;
   if (u[n - 1] == (DIGIT)k)
      return 0;
   return u[n - 1] > (DIGIT)k ? 1 : -1;
}

static novalue bdzero(dest, l)
DIGIT *dest;
word l;
{
   word i;

   for (i = 0; i < l; i++)
      dest[i] = 0;
}

static novalue bdcopy(src, dest, l)
DIGIT *src, *dest;
word l;
{
   word i;

   for (i = 0; i < l; i++)
      dest[i] = src[i];
}
#else					/* LargeInts */
static char x;			/* prevent empty module */
#endif					/* LargeInts */