mpi.c

Dropbear is an SSH 2 server, designed to be usable in small memory environments. It supports:

  } else {
    if((res = mp_copy(&u, &t)) != MP_OKAY)
      goto CLEANUP;

  }

  for(;;) {
    while(mp_iseven(&t)) {
      s_mp_div_2(&t);
    }

    if(mp_cmp_z(&t) == MP_GT) {
      if((res = mp_copy(&t, &u)) != MP_OKAY)
	goto CLEANUP;

    } else {
      if((res = mp_copy(&t, &v)) != MP_OKAY)
	goto CLEANUP;

      /* v = -t */
      if(SIGN(&t) == MP_ZPOS)
	SIGN(&v) = MP_NEG;
      else
	SIGN(&v) = MP_ZPOS;
    }

    if((res = mp_sub(&u, &v, &t)) != MP_OKAY)
      goto CLEANUP;

    if(s_mp_cmp_d(&t, 0) == MP_EQ)
      break;
  }

  s_mp_2expt(&v, k);       /* v = 2^k   */
  res = mp_mul(&u, &v, c); /* c = u * v */

 CLEANUP:
  mp_clear(&v);
 V:
  mp_clear(&u);
 U:
  mp_clear(&t);

  return res;

} /* end mp_bgcd() */

/* }}} */

/* {{{ mp_lcm(a, b, c) */

/* We compute the least common multiple using the rule:

   ab = [a, b](a, b)

   ... by computing the product, and dividing out the gcd.
 */

mp_err mp_lcm(mp_int *a, mp_int *b, mp_int *c)
{
  mp_int  gcd, prod;
  mp_err  res;

  ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);

  /* Set up temporaries */
  if((res = mp_init(&gcd)) != MP_OKAY)
    return res;
  if((res = mp_init(&prod)) != MP_OKAY)
    goto GCD;

  if((res = mp_mul(a, b, &prod)) != MP_OKAY)
    goto CLEANUP;
  if((res = mp_gcd(a, b, &gcd)) != MP_OKAY)
    goto CLEANUP;

  res = mp_div(&prod, &gcd, c, NULL);

 CLEANUP:
  mp_clear(&prod);
 GCD:
  mp_clear(&gcd);

  return res;

} /* end mp_lcm() */

/* }}} */

/* {{{ mp_xgcd(a, b, g, x, y) */

/*
  mp_xgcd(a, b, g, x, y)

  Compute g = (a, b) and values x and y satisfying Bezout's identity
  (that is, ax + by = g).  This uses the extended binary GCD algorithm
  based on the Stein algorithm used for mp_gcd()
 */

mp_err mp_xgcd(mp_int *a, mp_int *b, mp_int *g, mp_int *x, mp_int *y)
{
  mp_int   gx, xc, yc, u, v, A, B, C, D;
  mp_int  *clean[9];
  mp_err   res;
  int      last = -1;

  if(mp_cmp_z(b) == 0)
    return MP_RANGE;

  /* Initialize all these variables we need */
  if((res = mp_init(&u)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &u;
  if((res = mp_init(&v)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &v;
  if((res = mp_init(&gx)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &gx;
  if((res = mp_init(&A)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &A;
  if((res = mp_init(&B)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &B;
  if((res = mp_init(&C)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &C;
  if((res = mp_init(&D)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &D;
  if((res = mp_init_copy(&xc, a)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &xc;
  mp_abs(&xc, &xc);
  if((res = mp_init_copy(&yc, b)) != MP_OKAY) goto CLEANUP;
  clean[++last] = &yc;
  mp_abs(&yc, &yc);

  mp_set(&gx, 1);

  /* Divide by two until at least one of them is even */
  while(mp_iseven(&xc) && mp_iseven(&yc)) {
    s_mp_div_2(&xc);
    s_mp_div_2(&yc);
    if((res = s_mp_mul_2(&gx)) != MP_OKAY)
      goto CLEANUP;
  }

  mp_copy(&xc, &u);
  mp_copy(&yc, &v);
  mp_set(&A, 1); mp_set(&D, 1);

  /* Loop through binary GCD algorithm */
  for(;;) {
    while(mp_iseven(&u)) {
      s_mp_div_2(&u);

      if(mp_iseven(&A) && mp_iseven(&B)) {
	s_mp_div_2(&A); s_mp_div_2(&B);
      } else {
	if((res = mp_add(&A, &yc, &A)) != MP_OKAY) goto CLEANUP;
	s_mp_div_2(&A);
	if((res = mp_sub(&B, &xc, &B)) != MP_OKAY) goto CLEANUP;
	s_mp_div_2(&B);
      }
    }

    while(mp_iseven(&v)) {
      s_mp_div_2(&v);

      if(mp_iseven(&C) && mp_iseven(&D)) {
	s_mp_div_2(&C); s_mp_div_2(&D);
      } else {
	if((res = mp_add(&C, &yc, &C)) != MP_OKAY) goto CLEANUP;
	s_mp_div_2(&C);
	if((res = mp_sub(&D, &xc, &D)) != MP_OKAY) goto CLEANUP;
	s_mp_div_2(&D);
      }
    }

    if(mp_cmp(&u, &v) >= 0) {
      if((res = mp_sub(&u, &v, &u)) != MP_OKAY) goto CLEANUP;
      if((res = mp_sub(&A, &C, &A)) != MP_OKAY) goto CLEANUP;
      if((res = mp_sub(&B, &D, &B)) != MP_OKAY) goto CLEANUP;

    } else {
      if((res = mp_sub(&v, &u, &v)) != MP_OKAY) goto CLEANUP;
      if((res = mp_sub(&C, &A, &C)) != MP_OKAY) goto CLEANUP;
      if((res = mp_sub(&D, &B, &D)) != MP_OKAY) goto CLEANUP;

    }

    /* If we're done, copy results to output */
    if(mp_cmp_z(&u) == 0) {
      if(x)
	if((res = mp_copy(&C, x)) != MP_OKAY) goto CLEANUP;

      if(y)
	if((res = mp_copy(&D, y)) != MP_OKAY) goto CLEANUP;
      
      if(g)
	if((res = mp_mul(&gx, &v, g)) != MP_OKAY) goto CLEANUP;

      break;
    }
  }

 CLEANUP:
  while(last >= 0)
    mp_clear(clean[last--]);

  return res;

} /* end mp_xgcd() */

/* }}} */

/* {{{ mp_invmod(a, m, c) */

/*
  mp_invmod(a, m, c)

  Compute c = a^-1 (mod m), if there is an inverse for a (mod m).
  This is equivalent to the question of whether (a, m) = 1.  If not,
  MP_UNDEF is returned, and there is no inverse.
 */

mp_err mp_invmod(mp_int *a, mp_int *m, mp_int *c)
{
  mp_int  g, x;
  mp_err  res;

  ARGCHK(a && m && c, MP_BADARG);

  if(mp_cmp_z(a) == 0 || mp_cmp_z(m) == 0)
    return MP_RANGE;

  if((res = mp_init(&g)) != MP_OKAY)
    return res;
  if((res = mp_init(&x)) != MP_OKAY)
    goto X;

  if((res = mp_xgcd(a, m, &g, &x, NULL)) != MP_OKAY)
    goto CLEANUP;

  if(mp_cmp_d(&g, 1) != MP_EQ) {
    res = MP_UNDEF;
    goto CLEANUP;
  }

  res = mp_mod(&x, m, c);
  SIGN(c) = SIGN(a);

CLEANUP:
  mp_clear(&x);
X:
  mp_clear(&g);

  return res;

} /* end mp_invmod() */

/* }}} */
#endif /* if MP_NUMTH */

/* }}} */

/*------------------------------------------------------------------------*/
/* {{{ mp_print(mp, ofp) */

#if MP_IOFUNC
/*
  mp_print(mp, ofp)

  Print a textual representation of the given mp_int on the output
  stream 'ofp'.  Output is generated using the internal radix.
 */

void   mp_print(mp_int *mp, FILE *ofp)
{
  int   ix;

  if(mp == NULL || ofp == NULL)
    return;

  fputc((SIGN(mp) == MP_NEG) ? '-' : '+', ofp);

  for(ix = USED(mp) - 1; ix >= 0; ix--) {
    fprintf(ofp, DIGIT_FMT, DIGIT(mp, ix));
  }

} /* end mp_print() */

#endif /* if MP_IOFUNC */

/* }}} */

/*------------------------------------------------------------------------*/
/* {{{ More I/O Functions */

/* {{{ mp_read_signed_bin(mp, str, len) */

/* 
   mp_read_signed_bin(mp, str, len)

   Read in a raw value (base 256) into the given mp_int
 */

mp_err  mp_read_signed_bin(mp_int *mp, unsigned char *str, int len)
{
  mp_err         res;

  ARGCHK(mp != NULL && str != NULL && len > 0, MP_BADARG);

  if((res = mp_read_unsigned_bin(mp, str + 1, len - 1)) == MP_OKAY) {
    /* Get sign from first byte */
    if(str[0])
      SIGN(mp) = MP_NEG;
    else
      SIGN(mp) = MP_ZPOS;
  }

  return res;

} /* end mp_read_signed_bin() */

/* }}} */

/* {{{ mp_signed_bin_size(mp) */

int    mp_signed_bin_size(mp_int *mp)
{
  ARGCHK(mp != NULL, 0);

  return mp_unsigned_bin_size(mp) + 1;

} /* end mp_signed_bin_size() */

/* }}} */

/* {{{ mp_to_signed_bin(mp, str) */

mp_err mp_to_signed_bin(mp_int *mp, unsigned char *str)
{
  ARGCHK(mp != NULL && str != NULL, MP_BADARG);

  /* Caller responsible for allocating enough memory (use mp_raw_size(mp)) */
  str[0] = (char)SIGN(mp);

  return mp_to_unsigned_bin(mp, str + 1);

} /* end mp_to_signed_bin() */

/* }}} */

/* {{{ mp_read_unsigned_bin(mp, str, len) */

/*
  mp_read_unsigned_bin(mp, str, len)

  Read in an unsigned value (base 256) into the given mp_int
 */

mp_err  mp_read_unsigned_bin(mp_int *mp, unsigned char *str, int len)
{
  int     ix;
  mp_err  res;

  ARGCHK(mp != NULL && str != NULL && len > 0, MP_BADARG);

  mp_zero(mp);

  for(ix = 0; ix < len; ix++) {
    if((res = s_mp_mul_2d(mp, CHAR_BIT)) != MP_OKAY)
      return res;

    if((res = mp_add_d(mp, str[ix], mp)) != MP_OKAY)
      return res;
  }
  
  return MP_OKAY;
  
} /* end mp_read_unsigned_bin() */

/* }}} */

/* {{{ mp_unsigned_bin_size(mp) */

int     mp_unsigned_bin_size(mp_int *mp) 
{
  mp_digit   topdig;
  int        count;

  ARGCHK(mp != NULL, 0);

  /* Special case for the value zero */
  if(USED(mp) == 1 && DIGIT(mp, 0) == 0)
    return 1;

  count = (USED(mp) - 1) * sizeof(mp_digit);
  topdig = DIGIT(mp, USED(mp) - 1);

  while(topdig != 0) {
    ++count;
    topdig >>= CHAR_BIT;
  }

  return count;

} /* end mp_unsigned_bin_size() */

/* }}} */

/* {{{ mp_to_unsigned_bin(mp, str) */

mp_err mp_to_unsigned_bin(mp_int *mp, unsigned char *str)
{
  mp_digit      *dp, *end, d;
  unsigned char *spos;

  ARGCHK(mp != NULL && str != NULL, MP_BADARG);

  dp = DIGITS(mp);
  end = dp + USED(mp) - 1;
  spos = str;

  /* Special case for zero, quick test */
  if(dp == end && *dp == 0) {
    *str = '\0';
    return MP_OKAY;
  }

  /* Generate digits in reverse order */
  while(dp < end) {
    int      ix;

    d = *dp;
    for(ix = 0; ix < sizeof(mp_digit); ++ix) {
      *spos = d & UCHAR_MAX;
      d >>= CHAR_BIT;
      ++spos;
    }

    ++dp;
  }

  /* Now handle last digit specially, high order zeroes are not written */
  d = *end;
  while(d != 0) {
    *spos = d & UCHAR_MAX;
    d >>= CHAR_BIT;
    ++spos;
  }

  /* Reverse everything to get digits in the correct order */
  while(--spos > str) {
    unsigned char t = *str;
    *str = *spos;
    *spos = t;

    ++str;
  }

  return MP_OKAY;

} /* end mp_to_unsigned_bin() */

/* }}} */

/* {{{ mp_count_bits(mp) */

int    mp_count_bits(mp_int *mp)
{
  int      len;
  mp_digit d;

  ARGCHK(mp != NULL, MP_BADARG);

  len = DIGIT_BIT * (USED(mp) - 1);
  d = DIGIT(mp, USED(mp) - 1);

  while(d != 0) {
    ++len;
    d >>= 1;
  }

  return len;
  
} /* end mp_count_bits() */

/* }}} */

/* {{{ mp_read_radix(mp, str, radix) */

/*
  mp_read_radix(mp, str, radix)

  Read an integer from the given string, and set mp to the resulting
  value.  The input is presumed to be in base 10.  Leading non-digit
  characters are ignored, and the function reads until a non-digit
  character or the end of the string.
 */

mp_err  mp_read_radix(mp_int *mp, unsigned char *str, int radix)
{
  int     ix = 0, val = 0;
  mp_err  res;
  mp_sign sig = MP_ZPOS;

  ARGCHK(mp != NULL && str != NULL && radix >= 2 && radix <= MAX_RADIX, 
	 MP_BADARG);

  mp_zero(mp);

  /* Skip leading non-digit characters until a digit or '-' or '+' */
  while(str[ix] && 
	(s_mp_tovalue(str[ix], radix) < 0) && 
	str[ix] != '-' &&
	str[ix] != '+') {
    ++ix;
  }

  if(str[ix] == '-') {
    sig = MP_NEG;
    ++ix;
  } else if(str[ix] == '+') {
    sig = MP_ZPOS; /* this is the default anyway... */
    ++ix;
  }

  while((val = s_mp_tovalue(str[ix], radix)) >= 0) {
    if((res = s_mp_mul_d(mp, radix)) != MP_OKAY)
      return res;
    if((res = s_mp_add_d(mp, val)) != MP_OKAY)
      return res;
    ++ix;
  }

  if(s_mp_cmp_d(mp, 0) == MP_EQ)
    SIGN(mp) = MP_ZPOS;
  else
    SIGN(mp) = sig;

  return MP_OKAY;

} /* end mp_read_radix() */

/* }}} */

/* {{{ mp_radix_size(mp, radix) */

int    mp_radix_size(mp_int *mp, int radix)
{
  int  len;
  ARGCHK(mp != NULL, 0);

  len = s_mp_outlen(mp_count_bits(mp), radix) + 1; /* for NUL terminator */

  if(mp_cmp_z(mp) < 0)
    ++len; /* for sign */

  return len;

} /* end mp_radix_size() */

/* }}} */

/* {{{ mp_value_radix_size(num, qty, radix) */

/* num = number of digits
   qty = number of bits per digit
   radix = target base
   
   Return the number of digits in the specified radix that would be
   needed to express 'num' digits of 'qty' bits each.
 */
int    mp_value_radix_size(int num, int qty, int radix)
{
  ARGCHK(num >= 0 && qty > 0 && radix >= 2 && radix <= MAX_RADIX, 0);

  return s_mp_outlen(num * qty, radix);

} /* end mp_value_radix_size() */

/* }}} */

/* {{{ mp_toradix(mp, str, radix) */

mp_err mp_toradix(mp_int *mp, unsigned char *str, int radix)
{
  int  ix, pos = 0;

  ARGCHK(mp != NULL && str != NULL, MP_BADARG);
  ARGCHK(radix > 1 && radix <= MAX_RADIX, MP_RANGE);

  if(mp_cmp_z(mp) == MP_EQ) {
    str[0] = '0';
    str[1] = '\0';
  } else {
    mp_err   res;
    mp_int   tmp;
    mp_sign  sgn;
    mp_digit rem, rdx = (mp_digit)radix;
    char     ch;

    if((res = mp_init_copy(&tmp, mp)) != MP_OKAY)
      return res;

    /* Save sign for later, and take absolute value */
    sgn = SIGN(&tmp); SIGN(&tmp) = MP_ZPOS;

    /* Generate output digits in reverse order      */
    while(mp_cmp_z(&tmp) != 0) {
      if((res = s_mp_div_d(&tmp, rdx, &rem)) != MP_OKAY) {
	mp_clear(&tmp);
	return res;
      }

      /* Generate digits, use capital letters */
      ch = s_mp_todigit(rem, radix, 0);

      str[pos++] = ch;
    }

    /* Add - sign if original value was negative */
    if(sgn == MP_NEG)
      str[pos++] = '-';

    /* Add trailing NUL to end the string        */
    str[pos--] = '\0';

    /* Reverse the digits and sign indicator     */
    ix = 0;
    while(ix < pos) {
      char tmp = str[ix];

      str[ix] = str[pos];
      str[pos] = tmp;
      ++ix;
      --pos;
    }
    
    mp_clear(&tmp);
  }

  return MP_OKAY;

} /* end mp_toradix() */

/* }}} */

/* {{{ mp_char2value(ch, r) */

int    mp_char2value(char ch, int r)
{
  return s_mp_tovalue(ch, r);

} /* end mp_tovalue() */

/* }}} */