source: trunk/src/bignum.c@ 592

/* Implementation of bignums by Henrik.Johansson@Nexus.Comm.SE in 1991
 * version 1.2
 */

/* Modifications in this file (Nov 1999, Rainer Wichmann): 
 * - a few compiler warnings fixed:
 *   gcc insists on 'if ( ( b_eq_qr = ((b == q) || (b == r)) ) )'
 *   instead of     'if (   b_eq_qr = ((b == q) || (b == r))   )'
 * - code that is not used in samhain is #defined out by:
 *    #if (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER))
 *   or
 *    #ifdef UNUSED_CODE
 * - the error message (error_string[]) is enclosed in N_(),_() macros
 * - the following four #includes have been added
 */
#include "config_xor.h"


#include <string.h>

#ifdef HAVE_MEMORY_H
#include <memory.h>
#endif

#include "samhain.h"
#include "bignum.h"


#if (!defined(HAVE_LIBGMP) || !defined(HAVE_GMP_H)) && (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER))

#define DIGIT BIGNUM_DIGIT
#define DIGIT2 BIGNUM_TWO_DIGITS

#ifndef TRUE
#define TRUE (1 == 1)
#endif
#ifndef FALSE
#define FALSE (1 != 1)
#endif


#define MIN_ALLOC ((sizeof(long) / sizeof(DIGIT)) << 1)

/* Don't expect the "sign" field to have the right value (BIG_SIGN_0) at
 * all times, so compare real content of bignum with zerop.
 * Maybe it would be better to be sure at all times that the sign is right?
 */
#define zerop(BIG) ((*(BIG)->dp == 0) && ((BIG)->dgs_used == 1))
#define uonep(BIG) ((*(BIG)->dp == 1) && ((BIG)->dgs_used == 1))
#define NEGATE_SIGN(SGN) (-(SGN))

#define DIGIT_BITS BIGNUM_DIGIT_BITS
#define DIGIT2_BITS BIGNUM_TWO_DIGITS_BITS

#define DIGIT_PART(N) ((DIGIT)((N) & ((((DIGIT2)1) << DIGIT_BITS) - 1)))
#define RESULT_MINUSP(RES) (((RES) & ((DIGIT2)1 << (DIGIT2_BITS - 1))) != 0)
#define SHIFT_DIGIT_DOWN(N) (DIGIT_PART((N) >> DIGIT_BITS))
#define SHIFT_DIGIT_UP(N) ((DIGIT2)((N) << DIGIT_BITS))

#define DIGIT_BASE (((DIGIT2)1) << DIGIT_BITS)
#define MAX_DIGIT ((((DIGIT2)1) << DIGIT_BITS) - 1)

#define uint unsigned int
#define ulong unsigned long

char *last_string    = NULL;
char *big_end_string = NULL;

ulong length_last_string;
static char big_base_digits1[] =
N_("00112233445566778899AaBbCcDdEeFfGgHhIiJjKkLlMmNnOoPpQqRrSsTtUuVvWwXxYyZz");
static char error_string[] = N_("(big_errno != 0, something has gone wrong)");
static char big_base_digits[73] = "\0";

bignum big_one;
bignum tmp_string, tmp_add, tmp_mul, tmp_q, tmp_r, tmp_rand, tmp_round;
DIGIT *tmp_digits;
ulong  tmp_ulong;

/* #define COUNT_ALLOC */
#ifdef COUNT_ALLOC
ulong dgs_alloc = 0, digits_freed = 0;
#endif

int big_errno = BIG_OK;

double log2tbl[] =
{
    -1.0, -1.0,
    1.000000000000000, 1.584961891174317, 2.000000000000002, 2.321928024291994,
    2.584960937500003, 2.807353973388675, 3.000000000000004, 3.169923782348637,
    3.321928024291997, 3.459430694580083, 3.584960937500005, 3.700439453125006,
    3.807353973388678, 3.906888961791999, 4.000000000000014, 4.087459564208999,
    4.169921875000016, 4.247924804687517, 4.321926116943377, 4.392314910888691,
    4.459430694580098, 4.523559570312520, 4.584960937500021, 4.643856048584007,
    4.700439453125023, 4.754886627197290, 4.807353973388697, 4.857978820800807,
    4.906887054443386, 4.954193115234403, 5.000000000000028, 5.044391632080107,
    5.087459564209015, 5.129280090332062, 5.169921875000032
};

struct digit_blocks
{
    int digCnt;
    DIGIT dig_base;
} big_block[37];                /* The code uses index 2 to 36 */


/* ----------------------------------------------------------------------
 * All static utility functions are placed first in this file.
 */

#ifndef HAVE_STRCHR
static char * strchr (char * str_ptr, char ch)      
{
    do
    {
        if (*str_ptr == ch)
        {
            return str_ptr;
        }
    } while (*str_ptr++ != '\0');
    return NULL;
}
#endif

#if 0

extern int printf();

static void
print_digits(prefix, a)
char *prefix;
bignum *a;
{
    unsigned long i;
    ulong tmp;

    printf("%s", prefix);
    if (a->sign == BIG_SIGN_MINUS)
    {
        printf("- ");
    }
    else
    {
        if (a->sign == BIG_SIGN_PLUS)
        {
            printf("+ ");
        }
        else
        {
            printf("+/- ");
        }
    }
    for (i = a->dgs_used - 1; i > 0; i--)
    {
        tmp = a->dp[i];
        printf("%lu, ", tmp);
    }
    tmp = *a->dp;
    printf("%lu\n", tmp);
}

#else

#define print_digits(prefix, big) /* */

#endif

static void
init_digit_blocks(void)
{
    uint digcnt, base;
    DIGIT maxdigit, curdigit, tmp;

    for (base = 2; base <= 36; base++)
    {
        tmp = ((1U << (DIGIT_BITS - 1)) / base);
        maxdigit = tmp * 2;
        curdigit = 1;
        digcnt = 0;
        while (curdigit < maxdigit)
        {
            digcnt++;
            curdigit *= base;
        }
        big_block[base].digCnt = digcnt;
        big_block[base].dig_base = curdigit;
    }
}

#ifdef COUNT_ALLOC
static void
free_digits(DIGIT * dp, ulong count)
{
  /* Mar 4 2000 R. Wichmann: check for dp == NULL */
  if (!dp) 
    return;
  digits_freed += count;
  free((char *)dp);
  dp = NULL;
}
#else
/* June 5 2000 R. Wichmann: check for dp == NULL */
/* #define free_digits(DP,CNT) free((char *)DP) */
static void
free_digits(DIGIT * dp, ulong  count)
{
  /* Mar 4 2000 R. Wichmann: check for dp == NULL */
  if (!dp) 
    return;
  (void) count;
  free((char *)dp);
  dp = NULL;
}
#endif

static DIGIT *
allocate_digits(ulong alloclen)
{
    DIGIT *digit_ptr;

    if (big_errno != BIG_OK)
    {
        return NULL;
    }
#ifdef COUNT_ALLOC
    dgs_alloc += alloclen;
#endif
    digit_ptr = calloc(alloclen, sizeof(DIGIT));
    if (digit_ptr == NULL)
    {
        big_errno = BIG_MEMERR;
        return NULL;
    }
    return digit_ptr;
}

static bigerr_t
newsize(DIGIT **dp_p, ulong *cursize_p, ulong minsz, ulong newsz)
{
    if (*cursize_p >= minsz)
    {
        return big_errno;
    }
    free_digits(*dp_p, *cursize_p);
    if (newsz < MIN_ALLOC)
    {
        newsz = MIN_ALLOC;
    }
    if ((*dp_p = allocate_digits(newsz)) == NULL)
    {
        return big_errno;
    }
    *cursize_p = newsz;
    return big_errno;
}

/* `memcpy' uses an `int' as counter.  If an int can't hold a big enough
 * counter, then `digits_cpy' becomes a function (and a little slower,
 * probably.  Unfortunately `memcpy' takes a signed counter, but _that_
 * size of numbers are almost too big!  Or isn't it?
 */
#ifdef MEMCPY_LONG_COUNTER
/* extern char *memcpy(); */
#define digits_cpy(dst, src, count) memcpy((char *)(dst), (char *)(src), \
                                           (count) * sizeof(DIGIT))
#else
static void
digits_cpy(DIGIT *dst, DIGIT *src, ulong count)
{
    while (count-- > 0)
    {
        *dst++ = *src++;
    }
}
#endif

static DIGIT *
copy_digits(DIGIT *src, ulong cp_count, ulong alloclen)
{
    DIGIT *dst;

    if (big_errno != BIG_OK)
    {
        return NULL;
    }
    if ((dst = allocate_digits(alloclen)) == NULL)
    {
        return NULL;
    }
    digits_cpy(dst, src, cp_count);
    return dst;
}

static void
add_digit(bignum * a, DIGIT d)
{
    DIGIT2 res = d;
    DIGIT *last_digit, *digit_ptr = a->dp, *digvect;
    
    last_digit = digit_ptr + a->dgs_used;
    while (digit_ptr < last_digit)
    {
        res = *digit_ptr + res;
        *digit_ptr = DIGIT_PART(res);
        res = SHIFT_DIGIT_DOWN(res);
        digit_ptr++;
        if (res == 0)
        {
            break;
        }
    }
    if (res != 0)
    {
        if (a->dgs_used < a->dgs_alloc)
        {
            *digit_ptr = DIGIT_PART(res);
        }
        else
        {
            if ((digvect = copy_digits(a->dp, a->dgs_used,
                                       a->dgs_used + 4)) == NULL)
            {
                return;         /* big_errno will be set to last error */
            }
            digvect[a->dgs_used] = DIGIT_PART(res);
            free_digits(a->dp, a->dgs_alloc);
            a->dgs_alloc = a->dgs_used + 4;
            a->dp = digvect;
        }
        a->dgs_used++;
    }
}

static DIGIT
vect_div_digit(DIGIT *digvect, ulong *len, DIGIT d)
{
    DIGIT *digit_ptr = digvect + *len - 1, newval;
    DIGIT2 rest = 0;

    if (d == 0)
    {
        big_errno = BIG_DIV_ZERO;
        return -1;
    }
    if (d == 1)
    {
        return 0;
    }
    while (digit_ptr >= digvect)
    {
        rest = (DIGIT2)SHIFT_DIGIT_UP(rest);
        newval = DIGIT_PART((*digit_ptr + rest) / d);
        rest = (*digit_ptr + rest) % d;
        *digit_ptr = newval;
        digit_ptr--;
    }
    if (*len > 1)
    {
        if (digvect[*len - 1] == 0)
        {
            *len -= 1;
        }
    }
    return DIGIT_PART(rest);
}

static DIGIT
udiv_digit(bignum *a, DIGIT d)
{
    DIGIT res;

    res = vect_div_digit(a->dp, &a->dgs_used, d);
    if (zerop(a))
    {
        a->sign = BIG_SIGN_0;
    }
    return res;
}

static DIGIT
vect_mul_digit(DIGIT *digvect, ulong len, DIGIT x)
{
    DIGIT *digit_ptr = digvect + len;
    DIGIT2 res = 0;

    while (digvect < digit_ptr)
    {
        res += *digvect * (DIGIT2)x;
        *digvect = DIGIT_PART(res);
        res = SHIFT_DIGIT_DOWN(res);
        digvect++;
    }
    return DIGIT_PART(res);
}

static void
umul_digit(bignum *a, DIGIT x)
{
    DIGIT ovfl, *digvect;

    ovfl = vect_mul_digit(a->dp, a->dgs_used, x);
    if (ovfl != 0)
    {
        if (a->dgs_used < a->dgs_alloc)
        {
            a->dp[a->dgs_used] = ovfl;
        }
        else
        {
            digvect = copy_digits(a->dp,
                                  a->dgs_used,
                                  a->dgs_used + 4);
            digvect[a->dgs_used] = ovfl;
            free_digits(a->dp, a->dgs_alloc);
            a->dgs_alloc = a->dgs_used + 4;
            a->dp = digvect;
        }
        a->dgs_used++;
    }
}

static int
ucompare_digits(bignum *a, bignum *b)
{
    DIGIT *a_ptr, *b_ptr;
    int retval = 0;

    if (a->dgs_used  == b->dgs_used)
    {
        a_ptr = a->dp + a->dgs_used - 1;
        b_ptr = b->dp + b->dgs_used - 1;
        while ((*a_ptr == *b_ptr) && (a_ptr >= a->dp))
        {
            a_ptr--;
            b_ptr--;
        }
        if (a_ptr < a->dp)
        {
            return retval;
        }
        else
        {
            if (retval == 0)
                retval = (*a_ptr > *b_ptr) ? 1 : -1;
        }
        return retval;
    }
    return (a->dgs_used > b->dgs_used) ? 1 : -1;
}

static void
uadd_digits(bignum *a, bignum *b, bignum *c)
{
    DIGIT *dp_x, *dp_y, *dp_z, *res_dp, *end_x, *end_y;
    ulong len_x, len_y;
    DIGIT2 res = 0;

    if (a->dgs_used > b->dgs_used)
    {
        dp_x = a->dp;
        len_x = a->dgs_used;
        dp_y = b->dp;
        len_y = b->dgs_used;
    }
    else
    {
        dp_x = b->dp;
        len_x = b->dgs_used;
        dp_y = a->dp;
        len_y = a->dgs_used;
    }
    end_x = dp_x + len_x;
    end_y = dp_y + len_y;
    if (c->dgs_alloc >= len_x)
    {
        dp_z = c->dp;
    }
    else
    {
        if (newsize(&tmp_add.dp, &tmp_add.dgs_alloc,
                    len_x, len_x + 4) != BIG_OK)
        {
            return;
        }
        dp_z = tmp_add.dp;
    }
    res_dp = dp_z;
    while (dp_y < end_y)
    {   
        res += ((DIGIT2)*dp_x++) + *dp_y++;
        *res_dp++ = DIGIT_PART(res);
        res = SHIFT_DIGIT_DOWN(res);
    }
    while (dp_x < end_x)
    {
        res += *dp_x++;
        *res_dp++ = DIGIT_PART(res);
        res = SHIFT_DIGIT_DOWN(res);
    }
    if (res != 0)
    {
        *res_dp++ = DIGIT_PART(res);
    }

    if (dp_z != c->dp)
    {
        tmp_digits = c->dp;
        c->dp = tmp_add.dp;
        tmp_add.dp = tmp_digits;

        tmp_ulong = c->dgs_alloc;
        c->dgs_alloc = tmp_add.dgs_alloc;
        tmp_add.dgs_alloc = tmp_ulong;
    }
    c->dgs_used = res_dp - c->dp;
}

static void
usub_digits(bignum *a, bignum *b, bignum *c)
{
    DIGIT *dp_x, *dp_y, *dp_z, *res_dp, *end_x, *end_y;
    ulong len_x, len_y;
    DIGIT2 res = 0;

    dp_x = a->dp;
    len_x = a->dgs_used;
    dp_y = b->dp;
    len_y = b->dgs_used;

    end_x = dp_x + len_x;
    end_y = dp_y + len_y;
    if (c->dgs_alloc >= len_x)
    {
        dp_z = c->dp;
    }
    else
    {
        if (newsize(&tmp_add.dp, &tmp_add.dgs_alloc,
                    len_x, len_x + 2) != BIG_OK)
        {
            return;
        }
        dp_z = tmp_add.dp;
    }
    res_dp = dp_z;
    while (dp_y < end_y)
    {   
        res = ((DIGIT2)*dp_x++) - *dp_y++ - RESULT_MINUSP(res);
        *res_dp++ = DIGIT_PART(res);
    }
    while (dp_x < end_x)
    {
        res = *dp_x++ - RESULT_MINUSP(res);
        *res_dp++ = DIGIT_PART(res);
    }
#ifdef BIG_CHECK_LIMITS
    if (RESULT_MINUSP(res) != 0)
    {
        big_errno = BIG_ALGERR;
        return;
    }
#endif
    while ((*--res_dp == 0) && (res_dp > dp_z))
    {
        /* Move pointer down until we reach a non-zero */
    }
    if (dp_z != c->dp)
    {
        tmp_digits = c->dp;
        c->dp = tmp_add.dp;
        tmp_add.dp = tmp_digits;

        tmp_ulong = tmp_add.dgs_alloc;
        tmp_add.dgs_alloc = c->dgs_alloc;
        c->dgs_alloc = tmp_ulong;
    }
    c->dgs_used = res_dp - dp_z + 1;
}

/*
 *   This (pseudo) random number generator is not very good.  It has a long
 * period [ 2 ^ (DIGIT_BITS * 2) (?)] before it starts with the same sequence
 * again, but the lower bits are "less random" than they should be.  I've
 * solved it by using word of double length, and returning the upper bits.
 * Please mail me if you know how to make it into a "more random" generator.
 *   One important thing though: it will have to reasonably fast.
 *   The good thing with this one is that it is very portable, but that doesn't
 * help much when you want _good_ random numbers.
 *                                              Henrik.Johansson@Nexus.Comm.SE
 */
#ifdef UNUSED_CODE
static DIGIT
rand(void)
{
    static DIGIT2 x1 = 33, x2 = 45; /* Just give them two starting values */

    if (x2 = BIG_RAND_A2 * x2 + BIG_RAND_C2, x2 == 45)
    {
        x1 = BIG_RAND_A1 * x1 + BIG_RAND_C1;
    }
    return SHIFT_DIGIT_DOWN(x1 + x2); /* Skip least significant bits */
}
/* UNUSED_CODE */
#endif

/* ----------------------------------------------------------------------
 * All external functions comes here.
 */

bigerr_t
big_init_pkg()
{
    strcpy (big_base_digits, _(big_base_digits1));      /* known to fit  */
    init_digit_blocks();
    big_create(&tmp_string);
    big_create(&tmp_add);
    big_create(&tmp_mul);
    big_create(&tmp_q);
    big_create(&tmp_r);
    big_create(&tmp_rand);
    big_create(&big_one);
    big_set_long((long)1, &big_one);
    length_last_string = 10;
    if ((last_string = calloc(1,length_last_string)) == NULL)
    {
        big_errno = BIG_MEMERR;
    }
    return big_errno;
}

void
big_release_pkg()
{
    big_destroy(&tmp_string);
    big_destroy(&tmp_add);
    big_destroy(&tmp_mul);
    big_destroy(&tmp_q);
    big_destroy(&tmp_r);
    big_destroy(&tmp_round);
    big_destroy(&big_one);
    if (last_string != NULL)
      free(last_string);
#ifdef COUNT_ALLOC
    printf("Allocated digits: %lu\n", dgs_alloc);
    printf("Freed digits:     %lu\n", digits_freed);
#endif
}

bigerr_t
big_create(bignum * a)
{
  /* Make sure that big_create alway returns a NULL bignum
   * that can safely be destroyed.
   */
    if (a != NULL)
      {
        a->dp = NULL;
      }
    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    if (a == NULL)
    {
      big_errno = BIG_ARGERR;
      return big_errno;
    }

    a->sign = BIG_SIGN_0;
    a->dgs_used = 1;
    if ((a->dp = allocate_digits((long)sizeof(long))) == NULL)
    {
        return big_errno;
    }
    a->dgs_alloc = sizeof(long);
    *a->dp = 0;
    return big_errno;
}

void
big_destroy(bignum *a)
{
  /* Mar 4, 2000 R. Wichmann: check for a == NULL */
  /* free_digits() checks for a->dp == NULL       */
  if (a != NULL)
    free_digits(a->dp, a->dgs_alloc);
}

ulong
big_bitcount(bignum *a)
{
    int bits = 0;
    DIGIT high_digit;

    if (big_errno != BIG_OK)
    {
        return 0;
    }
    high_digit = a->dp[a->dgs_used - 1];
    while (high_digit != 0)
    {
        bits++;
        high_digit >>= 1;
    }
    return DIGIT_BITS * (a->dgs_used - 1) + bits;
}

bigerr_t
big_set_big(bignum *a, bignum *b)
{
    if ((big_errno != BIG_OK) || (a == b))
    {
        return big_errno;
    }
    if (newsize(&b->dp, &b->dgs_alloc, a->dgs_used, a->dgs_used) != BIG_OK)
    {
        return big_errno;
    }
    b->dgs_used = a->dgs_used;
    b->sign = a->sign;
    digits_cpy(b->dp, a->dp, a->dgs_used);
    return big_errno;
}


void
big_set_ulong(ulong n, bignum *a)
{
    unsigned int i;

    if (big_errno != BIG_OK)
    {
        return;
    }
    if (n == 0)
    {
        *a->dp = 0;
        a->dgs_used = 1;
        a->sign = BIG_SIGN_0;
    }
    else
    {
        a->dgs_used = 0;
        for (i = 0; i < sizeof(long) / sizeof(DIGIT); i++)
        {
            if (n == 0)
            {
                break;
            }
            a->dgs_used++;
            a->dp[i] = DIGIT_PART(n);
            n = SHIFT_DIGIT_DOWN(n);
        }
        a->sign = BIG_SIGN_PLUS;
    }
    print_digits("set_(u)long: a = ", a);
}

void
big_set_long(long n, bignum *a)
{
    ulong m;

    m = (n < 0) ? - n : n;
    big_set_ulong(m, a);
    if (!(n >= 0))
    {
        a->sign = BIG_SIGN_MINUS;
    }
}


bigerr_t
big_set_string(char * numstr, int base, bignum *a)
{
    char *src, *maxdigit, *chrptr;
    DIGIT dig_base, dig_sum, last_base;
    int cnt, maxcnt;

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }

    big_end_string = numstr;
    if ((base < 2) || (base > 36) || (numstr == NULL) || (a == NULL))
    {
        big_errno = BIG_ARGERR;
        return big_errno;
    }

    src            = numstr;

    maxdigit = big_base_digits + (base << 1);

    maxcnt = big_block[base].digCnt;
    dig_base = big_block[base].dig_base;
    a->dgs_used = 1;
    *a->dp = 0;
    a->sign = BIG_SIGN_PLUS;    /* Assume it will be positive */
    while (strchr(" \t\n\r", *src) != NULL) /* Skip whitespaces */
    {
        src++;
    }
    if ((*src == '-') || (*src == '+')) /* First non-white is a sign? */
    {
        a->sign = (*src == '-') ? BIG_SIGN_MINUS : BIG_SIGN_PLUS;
        src++;
    }
    chrptr = strchr(big_base_digits, *src);
    if ((chrptr == NULL) || (chrptr >= maxdigit)) /* Next chr not a digit? */
    {
        big_end_string = src;
        big_errno = BIG_ARGERR;
        return big_errno;
    }
    while (*src == '0')         /* Next chr a '0'? */
    {
        src++;                  /* Read past all '0'es */
    }
    chrptr = strchr(big_base_digits, *src);
    if ((chrptr == NULL) || (chrptr >= maxdigit)) /* Next char not a digit */
    {
        big_end_string = src;
        a->sign = BIG_SIGN_0;   /* It was just a plain 0 */
        return big_errno;
    }
    dig_sum = 0;
    cnt = 0;
    while ((chrptr = strchr(big_base_digits, *src)),
           (chrptr != NULL) && (chrptr < maxdigit))
    {
        dig_sum = dig_sum * base + ((chrptr - big_base_digits) >> 1);
        if (++cnt >= maxcnt)
        {
            umul_digit(a, dig_base);
            add_digit(a, dig_sum);
            dig_sum = 0;
            cnt = 0;
        }
        src++;
    }
    if (cnt > 0)
    {
        last_base = base;
        while (--cnt > 0)
        {
            last_base *= base;
        }
        umul_digit(a, last_base);
        add_digit(a, dig_sum);
    }
    big_end_string = src;
    return big_errno;
}


int
big_ulong(bignum *a, ulong *n)
{
    ulong old_n;
    DIGIT *dp;

    if (big_errno != BIG_OK)
    {
        return FALSE;
    }
    if (a->dgs_used > sizeof(ulong) / sizeof(DIGIT))
    {
        return FALSE;
    }
    dp = a->dp + a->dgs_used - 1;
    old_n = *n = *dp--;
    while ((dp >= a->dp) && (old_n < *n))
    {
        old_n = *n;
        *n = SHIFT_DIGIT_UP(*n) + *dp--;
    } 
    if (old_n >= *n)
    {
        return FALSE;
    }
    return FALSE;
}

int
big_long(bignum *a, long *n)
{
    long old_n;
    DIGIT *dp;

    if (a->dgs_used > sizeof(ulong) / sizeof(DIGIT))
    {
        return FALSE;
    }
    dp = a->dp + a->dgs_used - 1;
    old_n = *n = *dp--;
    while ((dp >= a->dp) && (old_n <= *n))
    {
        old_n = *n;
        *n = SHIFT_DIGIT_UP(*n) + *dp--;
    } 
    if (old_n >= *n)
    {
        return FALSE;
    }
    if (a->sign == BIG_SIGN_MINUS)
    {
      *n = -*n;
    }
    return FALSE;
}


char *
big_string(bignum *a, int base)
{
    ulong bit_length, str_length;
    char *dst;
    DIGIT dig_base, rem;
    int cnt, maxcnt;

    if (big_errno != BIG_OK)
    {
        return _(error_string);
    }
    big_set_big(a, &tmp_string);
                                /* Need more room than minimum here... */
    bit_length = tmp_string.dgs_used * DIGIT_BITS;
    /*    bit_length = big_bitcount(&tmp_string); */
    str_length = (ulong)(bit_length / log2tbl[base] + 4);
    if (str_length > length_last_string)
    {
        if (last_string != NULL)
          free(last_string);
        if ((last_string = calloc(1,str_length)) == NULL)
        {
            big_errno = BIG_MEMERR;
            return NULL;
        }
        length_last_string = str_length;
    }
            
    dst = last_string + length_last_string - 1;
    *dst-- = '\0';
    maxcnt = big_block[base].digCnt;
    dig_base = big_block[base].dig_base;
    while (tmp_string.dgs_used > 1)
    {
        rem = udiv_digit(&tmp_string, dig_base);
        for (cnt = 0; cnt < maxcnt; cnt++)
        {
            *dst-- = big_base_digits[(rem % base) << 1];
            rem /= base;
        }
    }
    rem = *tmp_string.dp;
    do
    {
        *dst-- = big_base_digits[(rem % base) << 1];
        rem /= base;
    } while (rem != 0);

    if (a->sign == BIG_SIGN_MINUS)
    {
        *dst = '-';
    }
    else
    {
        dst++;
    }
    return dst;
}

bigerr_t
big_negate(bignum *a, bignum *b)
{
    big_set_big(a, b);
    b->sign = NEGATE_SIGN(a->sign);
    return big_errno;
}

int
big_sign(bignum *a)
{
    return a->sign;
}

bigerr_t
big_abs(bignum *a, bignum *b)
{
    big_set_big(a, b);
    if (a->sign == BIG_SIGN_MINUS)
    {
        b->sign = NEGATE_SIGN(a->sign);
    }
    return big_errno;
}

int
big_compare(bignum *a, bignum *b)
{
    if (a->sign == b->sign)
    {
        if (a->sign == 0)
        {
            return 0;
        }
        return
            (a->sign == BIG_SIGN_MINUS)
            ? -ucompare_digits(a, b)
            : ucompare_digits(a, b);
    }
    return b->sign - a->sign;
}

int
big_lessp(bignum *a, bignum *b)
{
    return big_compare(a, b) < 0;
}

int
big_leqp(bignum *a, bignum *b)
{
    return !(big_compare(a, b) > 0);
}

int
big_equalp(bignum *a, bignum *b)
{
    return big_compare(a, b) == 0;
}

int
big_geqp(bignum *a, bignum *b)
{
    return !(big_compare(a, b) < 0);
}

int
big_greaterp(bignum *a, bignum *b)
{
    return big_compare(a, b) > 0;
}

int
big_zerop(bignum *a)
{
    return a->sign == BIG_SIGN_0;
}

int
big_evenp(bignum *a)
{
    return ((*a->dp & 0x01) == 0);
}

int
big_oddp(bignum *a)
{
    return ((*a->dp & 0x01) == 1);
}

bigerr_t
big_add(bignum *a, bignum *b, bignum *c)
{
    int cmp;

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    print_digits("add:\ta = ", a);
    print_digits("\tb = ", b);
    if (a->sign == b->sign)
    {
        uadd_digits(a, b, c);
        c->sign = a->sign;
    }
    else
    {
        cmp = ucompare_digits(a, b);
        if (cmp < 0)
        {
            usub_digits(b, a, c);
            if (zerop(c))
            {
                c->sign = BIG_SIGN_0;
            }
            else
            {
                c->sign = b->sign;
            }
        }
        else if (cmp > 0)
        {
            usub_digits(a, b, c);
            if (zerop(c))
            {
                c->sign = BIG_SIGN_0;
            }
            else
            {
                c->sign = a->sign;
            }
        }
        else
        {
            c->dgs_used = 1;
            *c->dp = 0;
            c->sign = BIG_SIGN_0;
        }
    }
    print_digits("\tc = ", c);
    return big_errno;
}

bigerr_t
big_sub(bignum *a, bignum *b, bignum *c)
{
    int cmp;

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    print_digits("sub:\ta = ", a);
    print_digits("\tb = ", b);
    if (a->sign == BIG_SIGN_0)
    {
        big_set_big(b, c);
        big_negate(c, c);
        print_digits("\tc = ", c);
        return big_errno;
    }
    if (b->sign == BIG_SIGN_0)
    {
        big_set_big(a, c);
        print_digits("\tc = ", c);
        return big_errno;
    }

    cmp = ucompare_digits(a, b);
    if (cmp <= 0)
    {
        if (a->sign != b->sign)
        {
            uadd_digits(a, b, c);
            c->sign = a->sign;
        }
        else
        {
            usub_digits(b, a, c);
            c->sign = (zerop(c) ? BIG_SIGN_0 : NEGATE_SIGN(a->sign));
        }
    }
    else if (cmp > 0)
    {
        if (a->sign != b->sign)
        {
            uadd_digits(a, b, c);
            c->sign = a->sign;
        }
        else
        {
            usub_digits(a, b, c);
            c->sign = (zerop(c) ? BIG_SIGN_0 : b->sign);
        }
    }
    else
    {
        c->dgs_used = 1;
        *c->dp = 0;
        c->sign = BIG_SIGN_0;
    }
    print_digits("\tc = ", c);
    return big_errno;
}

bigerr_t
big_mul(bignum *a, bignum *b, bignum *c)
{
    DIGIT *dp_x, *dp_xstart, *dp_xend;
    DIGIT *dp_y, *dp_ystart, *dp_yend;
    DIGIT *dp_z, *dp_zstart, *dp_zend, *dp_zsumstart;
    ulong len_x, len_y /* , len_z */;
    DIGIT2 res;
    DIGIT tmp_res;           /* Without use of this, gcc (v 1.39) generates */
                             /* erroneous code on a sun386 machine */
                             /* Should be removed with #ifdef's, when */
                             /* not on a sun386, but... */

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    print_digits("mul:\ta = ", a);
    print_digits("\tb = ", b);

    if (zerop(a) || zerop(b))
    {
        c->sign = BIG_SIGN_0;
        c->dgs_used = 1;
        *c->dp = 0;
        print_digits("(a=0 || b=0)c = ", c);
        return big_errno;
    }
    if (uonep(a))
    {
        big_set_big(b, c);
        c->sign = (a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS;
        print_digits("(abs(a)=1)c = ", c);
        return big_errno;
    }
    if (uonep(b))
    {
        big_set_big(a, c);
        c->sign = (a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS;
        print_digits("(abs(b)=1)c = ", c);
        return big_errno;
    }

    if (a->dgs_used < b->dgs_used)
    {
        dp_xstart = a->dp;
        len_x = a->dgs_used;
        dp_ystart = b->dp;
        len_y = b->dgs_used;
    }
    else
    {
        dp_xstart = b->dp;
        len_x = b->dgs_used;
        dp_ystart = a->dp;
        len_y = a->dgs_used;
    }
    if ((c == a) || (c == b))
    {
        if (newsize(&tmp_mul.dp, &tmp_mul.dgs_alloc,
                    len_x + len_y, len_x + len_y + 2) != BIG_OK)
        {
            return big_errno;
        }
        dp_zsumstart = tmp_mul.dp;
        /* len_z = tmp_mul.dgs_alloc; */
    }   
    else
    {
        if (newsize(&c->dp, &c->dgs_alloc,
                    len_x + len_y, len_x + len_y + 2) != BIG_OK)
        {
            return big_errno;
        }
        dp_zsumstart = c->dp;
        /* len_z = c->dgs_alloc; */
    }

    dp_xend = dp_xstart + len_x;
    dp_yend = dp_ystart + len_y;
    dp_zend = dp_zsumstart + len_y;

    for (dp_z = dp_zsumstart; dp_z < dp_zend; dp_z++)
    {
        *dp_z = 0;              /* Zero out rightmost digits */
    }

    dp_zstart = dp_zsumstart;

    for (dp_x = dp_xstart; dp_x < dp_xend; dp_x++)
    {
        dp_z = dp_zstart;
        tmp_res = 0;
        for (dp_y = dp_ystart; dp_y < dp_yend; dp_y++)
        {
            res = (DIGIT2)(*dp_x) * (*dp_y) + (*dp_z) + tmp_res;
            *dp_z++ = DIGIT_PART(res);
            tmp_res = SHIFT_DIGIT_DOWN(res);
        }
        *dp_z = tmp_res;
        dp_zstart++;
    }
    if (dp_zsumstart != c->dp)
    {
        tmp_digits = c->dp;
        c->dp = tmp_mul.dp;
        tmp_mul.dp = tmp_digits;

        tmp_ulong = c->dgs_alloc;
        c->dgs_alloc = tmp_mul.dgs_alloc;
        tmp_mul.dgs_alloc = tmp_ulong;
    }
    if (*dp_z == 0)
    {
        dp_z--;
    }
    c->dgs_used = dp_z - dp_zsumstart + 1;
    c->sign = a->sign * b->sign;
    print_digits("\tc = ", c);
    return big_errno;
}

bigerr_t
big_trunc(bignum *a, bignum *b, bignum *q, bignum *r)
{
    DIGIT *v_end, *q_end, *r_end, *src, *dst;
    DIGIT norm, qhat, t1, t2, t3, v1, v2, u1, u2, u3;
    DIGIT2 temp, res, carry;
    long a_l, b_l, q_l, r_l;
#if 0
    ulong i; 
#endif
    ulong j, m, n;
    int cmp, q_eq_a, q_eq_b, r_eq_a, r_eq_b;

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    print_digits("div:\ta = ", a);
    print_digits("\tb = ", b);
    if (zerop(b))
    {
        big_errno = BIG_DIV_ZERO;
        return big_errno;
    }
    
    if (q == r)
    {
        big_errno = BIG_ARGERR;
        return big_errno;
    }

    if (b->dgs_used == 1)
    {
        big_set_big(a, q);
        q->sign = ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
        *r->dp = udiv_digit(q, *b->dp);
        r->dgs_used = 1;
        r->sign = (zerop(r) ? BIG_SIGN_0 : a->sign);
        print_digits("\t3:q = ", q);
        print_digits("\t  r = ", r);
        return big_errno;
    }
    
    if (big_long(a, &a_l))      /* Pretend it is a signed value */
    {
        big_long(b, &b_l);      /* |a| < |b| so this will succeed */
        q_l = a_l / b_l;        /* Compute with unsigned operators */
        r_l = a_l % b_l;
        big_set_long((long)q_l, q);
        big_set_long((long)r_l, r);
        print_digits("\t4:q = ", q);
        print_digits("\t  r = ", r);
        return big_errno;
    }

    cmp = ucompare_digits(a, b); /* Unsigned compare, that is... */
    if (cmp < 0)
    {
        big_set_big(a, r);      /* r = a (don't care about big_errno here) */
        q->sign = BIG_SIGN_0;   /* q = 0 */
        *q->dp = 0;
        q->dgs_used = 1;
        print_digits("\t1:q = ", q);
        print_digits("\t  r = ", r);
        return big_errno;
    }
    else
    if (cmp == 0)
    {
        q->sign = ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
        *q->dp = 1;     /* q = 1 */
        q->dgs_used = 1;
        r->sign = BIG_SIGN_0;   /* r = 0 */
        *r->dp = 0;
        r->dgs_used = 1;
        print_digits("\t2:q = ", q);
        print_digits("\t  r = ", r);
        return big_errno;
    }

    q_eq_a = (q == a);
    q_eq_b = (q == b);
    if (q_eq_a || q_eq_b)
    {
        q = &tmp_q;
    }

    r_eq_a = (r == a);
    r_eq_b = (r == b);
    if (r_eq_a || r_eq_b)
    {
        r = &tmp_r;
    }

    if (newsize(&r->dp, &r->dgs_alloc, /* At least one more dig in r */
                a->dgs_used + 1, a->dgs_used + 2) != BIG_OK)
    {
        return big_errno;
    }
    big_set_big(a, r);
    r->dp[a->dgs_used] = 0; /* In case no overflow in mult. */
    
    n = b->dgs_used;
    v_end = b->dp + n - 1;
    norm = DIGIT_PART((DIGIT_BASE / ((DIGIT2)*v_end + 1)));
    if (norm != 1)
    {
        umul_digit(r, norm);
        umul_digit(b, norm);
        print_digits("r = ", r);
        print_digits("b = ", b);
    }
    m = a->dgs_used + 1 - b->dgs_used;
    r_end = r->dp + a->dgs_used;
    
    if (newsize(&q->dp, &q->dgs_alloc, m, m + 2) != BIG_OK)
    {
        return big_errno;
    }
    q_end = q->dp + m - 1;
    
    v1 = *v_end;
    v2 = *(v_end - 1);
    
    for (j = 0; j < m; j++)                   /* m steps through division */
    {
        u1 = *r_end;                            /*  routine */
        u2 = *(r_end - 1);
        u3 = *(r_end - 2);
        qhat = ((u1 == v1) ?
                MAX_DIGIT :
                DIGIT_PART(((DIGIT2)u1 * DIGIT_BASE + u2) / v1));
        while (1)
        {
            t3 = DIGIT_PART(temp = (DIGIT2)qhat * v2);
            t2 = DIGIT_PART(temp = SHIFT_DIGIT_DOWN(temp) + v1 * (DIGIT2)qhat);
            t1 = DIGIT_PART(SHIFT_DIGIT_DOWN(temp));
#if 0
            printf("t1 = %lu, ", (ulong)t1);
            printf("t2 = %lu, ", (ulong)t2);
            printf("t3 = %lu\n", (ulong)t3);
#endif
            if (t1 < u1) break;
            if (t1 > u1) {--qhat; continue; }
            if (t2 < u2) break;
            if (t2 > u2) { --qhat; continue; }
            if (t3 <= u3) break;
            qhat--;
        }
        
        /* This is a tricky one - multiply and subtract simultaneously */
        carry = 1;
        res = 0;
        src = b->dp;
        dst = r->dp + m - j - 1;
        while (src <= v_end)
        {
            res = (DIGIT2)qhat * *(src++) + SHIFT_DIGIT_DOWN(res);
            carry += (DIGIT2)(*dst) + MAX_DIGIT - DIGIT_PART(res);
            *(dst++) = DIGIT_PART(carry);
            carry = DIGIT_PART(SHIFT_DIGIT_DOWN(carry));
        }
        carry += (DIGIT2)(*dst) + MAX_DIGIT - SHIFT_DIGIT_DOWN(res);
        *dst = DIGIT_PART(carry);
        carry = DIGIT_PART(SHIFT_DIGIT_DOWN(carry));
        
        if (carry == 0)
        {
            qhat--;
            src = b->dp;
            dst = r->dp + m - j - 1;
            while (dst <= r_end)
            {
                carry = (DIGIT2)(*dst) + *src++ + SHIFT_DIGIT_DOWN(carry);
                *dst++ = DIGIT_PART(carry);
            }
            *dst = 0;
        }
        *(q_end - j) = DIGIT_PART(qhat);
        r_end--;
    }
    r->sign = a->sign;
#if 0
    i = r->dgs_used;
#endif
    while ((*r_end == 0) && (r_end > r->dp))
    {
        r_end--;
    }
    if (r_end == r->dp)
    {
        r->dgs_used = 1;
        r->sign = BIG_SIGN_0;
    }
    else
    {
        r->dgs_used = r_end - r->dp + 1;
        r->sign = a->sign;
    }
    if (norm != 1)
    {
        udiv_digit(b, norm);
        udiv_digit(r, norm);
    }
    while ((*q_end == 0) && (q_end > q->dp)) /* This is not needed!(?) */
    {
        q_end--;
    }
    /* was ifdef'd out already in original */
#if 0
    i = m - 1;
    while ((i > 0) && (q->dp[i--] == 0))
    {
        /* Loop through all zeroes */
    }
#endif
    q->dgs_used = q_end - q->dp + 1;
    q->sign =  ((a->sign == b->sign) ? BIG_SIGN_PLUS : BIG_SIGN_MINUS);
    
    if (q_eq_a)
    {
        big_set_big(q, a);
    }
    else
    if (q_eq_b)
    {
        big_set_big(q, b);
    }

    if (r_eq_b)
    {
        big_set_big(r, b);
    }
    else
    if (r_eq_a)
    {
        big_set_big(r, a);
    }

    print_digits("\t5:q = ", q);
    print_digits("\t  r = ", r);
    return big_errno;
}

bigerr_t
big_floor(bignum *a, bignum *b, bignum *q, bignum *r)
{
    int b_eq_qr, sign_eq;

    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
    {
        big_set_big(b, &tmp_mul);
    }
    sign_eq = a->sign == b->sign;
    big_trunc(a, b, q, r);
    if (sign_eq)
    {
        return big_errno;
    }
    if (!zerop(r))
    {
        if (b_eq_qr)
        {
            big_add(r, &tmp_mul, r);
        }
        else
        {
            big_add(r, b, r);
        }
        print_digits("big_one = ", &big_one);
        big_sub(q, &big_one, q);
    }
    return big_errno;
}

bigerr_t
big_ceil(bignum *a, bignum *b, bignum *q, bignum *r)
{
    int b_eq_qr, sign_diff;

    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
    {
        big_set_big(b, &tmp_mul);
    }
    sign_diff = a->sign != b->sign;
    big_trunc(a, b, q, r);
    if (sign_diff)
    {
        return big_errno;
    }
    if (!zerop(r))
    {
        if (b_eq_qr)
        {
            big_sub(r, &tmp_mul, r);
        }
        else
        {
            big_sub(r, b, r);
        }
        big_add(q, &big_one, q);
    }
    return big_errno;
}

#ifdef UNUSED_CODE
/* This one doesn't work to 100%.  I was a little braindamaged when I wrote
 * this, but I'll eventually fix it.   Zzzzzzz.
 */
bigerr_t
big_round(bignum *a, bignum *b, bignum *q, bignum *r)
{
    int b_eq_qr, b_neg_p, a_sgn_neq_b_sgn;

    if ( ( b_eq_qr = ((b == q) || (b == r)) ) )
    {
        big_set_big(b, &tmp_round);
    }
    b_neg_p = b->sign == BIG_SIGN_MINUS;
    a_sgn_neq_b_sgn = a->sign != b->sign;
    big_trunc(a, b, q, r);

    big_set_big(r, &tmp_add);
    umul_digit(&tmp_add, 2);
    if (ucompare_digits(&tmp_add, b) > 0) /* |2 * r| > |b| */
    {
        if (q->sign == BIG_SIGN_0)
        {
            if (a_sgn_neq_b_sgn)
            {
                big_sub(q, &big_one, q);
            }
            else
            {
                big_add(q, &big_one, q);
            }
        }
        else
        {
            if (q->sign == BIG_SIGN_MINUS)
            {
                big_sub(q, &big_one, q);
            }
            else
            {
                big_add(q, &big_one, q);
            }
        }
        if (b_eq_qr)
        {
            if (q->sign == BIG_SIGN_PLUS)
            {
                big_sub(r, &tmp_round, r);
            }
            else
            {
                big_add(r, &tmp_round, r);
            }
        }
        else
        {
            if (q->sign == BIG_SIGN_PLUS)
            {
                big_sub(r, b, r);
            }
            else
            {
                big_add(r, b, r);
            }
        }
    }
    return big_errno;
}

bigerr_t
big_random(bignum *a, bignum *b)
{
    unsigned long i;
    int a_sgn = a->sign;

    if (big_errno != BIG_OK)
    {
        return big_errno;
    }
    if (zerop(a))               /* a = 0 -> big_random => 0 (special case) */
    {
        *b->dp = 0;
        b->dgs_used = 1;
        b->sign = a_sgn;
        return big_errno;
    }
    if (newsize(&tmp_rand.dp, &tmp_rand.dgs_alloc,
                a->dgs_used + 1, a->dgs_used + 1) != BIG_OK)
    {
        return big_errno;
    }
    for (i = 0; i <= a->dgs_used; i++)
    {
        tmp_rand.dp[i] = rand();
    }
    while (tmp_rand.dp[a->dgs_used] == 0) /* Make sure high digit is non-0 */
    {
        tmp_rand.dp[a->dgs_used] = rand();
    }
    tmp_rand.dgs_used = a->dgs_used + 1;
    tmp_rand.sign = BIG_SIGN_PLUS;
    a->sign = BIG_SIGN_PLUS;
    big_trunc(&tmp_rand, a, &tmp_q, b); /* Dangerous to use tmp_q here... */
    a->sign = a_sgn;
    b->sign = zerop(&tmp_rand) ? BIG_SIGN_0 : a_sgn;
    return big_errno;
}

/* UNUSED_CODE */
#endif

/* ----------------------------------------------------------------------
 * External functions that do not need to know anything about the internal
 * representation of a bignum.
 */

#ifdef UNUSED_CODE

int
big_expt(bignum *a, unsigned long z, bignum *x)
{
    bignum b;

    big_create(&b);
    big_set_big(a, &b);
    big_set_long((long)1, x);
    
    while ((z != 0) && (big_errno == BIG_OK))
    {
        while ((z & 0x01) == 0)
        {
            z >>= 1;
            big_mul(&b, &b, &b);
        }
        z -= 1;
        big_mul(x, &b, x);
    }

    big_destroy(&b);
    return big_errno;
}

/* UNUSED_CODE */
#endif

int
big_exptmod(bignum *a_in, bignum *z_in, bignum *n, bignum *x)
{
    bignum a, z, b0, b1, b2, dummy;

    big_create(&a);
    big_create(&z);
    big_create(&b0);
    big_create(&b1);
    big_create(&b2);
    big_create(&dummy);

    big_set_big(a_in, &a);
    big_set_big(z_in, &z);
    big_set_long((long)1, x);
    big_set_long((long)0, &b0);
    big_set_long((long)1, &b1);
    big_set_long((long)2, &b2);

    /* No foolproof testing on big_errno - it really ought to be done */
    while ((big_compare(&z, &b0) != 0) && (big_errno == BIG_OK))
    {
        while (big_evenp(&z) && (big_errno == BIG_OK))
        {
            big_trunc(&z, &b2, &z, &dummy);
            big_mul(&a, &a, &a);
            big_trunc(&a, n, &dummy, &a);
        }
        big_sub(&z, &b1, &z);
        big_mul(x, &a, x);
        big_trunc(x, n, &dummy, x);
    }

    big_destroy(&dummy);
    big_destroy(&b2);
    big_destroy(&b1);
    big_destroy(&b0);
    big_destroy(&z);
    big_destroy(&a);

    return big_errno;
}

#ifdef UNUSED_CODE

int
big_gcd(bignum *a, bignum *b, bignum *g)
{
    bignum a1, b1, tmp;

    if (big_zerop(b))
    {
        big_abs(a, g);
        return big_errno;
    }

    big_create(&a1);
    big_create(&b1);
    big_create(&tmp);

    big_abs(a, &a1);
    big_abs(b, &b1);

    while (!big_zerop(&b1) && (big_errno == BIG_OK))
    {
        big_floor(&a1, &b1, &tmp, &a1);
        if (big_zerop(&a1))
        {
            break;
        }
        big_floor(&b1, &a1, &tmp, &b1);
    }

    if (big_zerop(&a1))
    {
        big_set_big(&b1, g);
    }
    else
    {
        big_set_big(&a1, g);
    }

    big_destroy(&tmp);
    big_destroy(&b1);
    big_destroy(&a1);

    return big_errno;
}
/* UNUSED_CODE */
#endif

/* #if (defined (SH_WITH_CLIENT) || defined (SH_WITH_SERVER)) */
#endif