#include "mpi.h"
#include <stdlib.h>
#include <string.h>
#include <ctype.h>
#if MP_DEBUG
#include <stdio.h>
#define DIAG(T,V) {fprintf(stderr,T);mp_print(V,stderr);fputc('\n',stderr);}
#else
#define DIAG(T,V)
#endif
#if MP_LOGTAB
#include "logtab.h"
#define LOG_V_2(R) s_logv_2[(R)]
#else
#include <math.h>
#define LOG_V_2(R) (log(2.0)/log(R))
#endif
static unsigned int s_mp_defprec = MP_DEFPREC;
#define CARRYOUT(W) ((W)>>DIGIT_BIT)
#define ACCUM(W) ((W)&MP_DIGIT_MAX)
#define MP_LT -1
#define MP_EQ 0
#define MP_GT 1
static const char *mp_err_string[] = {
"unknown result code",
"boolean true",
"boolean false",
"out of memory",
"argument out of range",
"invalid input parameter",
"result is undefined"
};
static const char *s_dmap_1 =
"0123456789ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz+/";
#if 0
static const char *s_dmap_2 =
"ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/";
#endif
#if MP_MACRO == 0
void s_mp_setz(mp_digit *dp, mp_size count);
void s_mp_copy(mp_digit *sp, mp_digit *dp, mp_size count);
void *s_mp_alloc(size_t nb, size_t ni);
void s_mp_free(void *ptr);
#else
#if MP_MEMSET == 0
#define s_mp_setz(dp, count) \
{int ix;for(ix=0;ix<(count);ix++)(dp)[ix]=0;}
#else
#define s_mp_setz(dp, count) memset(dp, 0, (count) * sizeof(mp_digit))
#endif
#if MP_MEMCPY == 0
#define s_mp_copy(sp, dp, count) \
{int ix;for(ix=0;ix<(count);ix++)(dp)[ix]=(sp)[ix];}
#else
#define s_mp_copy(sp, dp, count) memcpy(dp, sp, (count) * sizeof(mp_digit))
#endif
#define s_mp_alloc(nb, ni) calloc(nb, ni)
#define s_mp_free(ptr) {if(ptr) free(ptr);}
#endif
mp_err s_mp_grow(mp_int *mp, mp_size min);
mp_err s_mp_pad(mp_int *mp, mp_size min);
void s_mp_clamp(mp_int *mp);
void s_mp_exch(mp_int *a, mp_int *b);
mp_err s_mp_lshd(mp_int *mp, mp_size p);
void s_mp_rshd(mp_int *mp, mp_size p);
void s_mp_div_2d(mp_int *mp, mp_digit d);
void s_mp_mod_2d(mp_int *mp, mp_digit d);
mp_err s_mp_mul_2d(mp_int *mp, mp_digit d);
void s_mp_div_2(mp_int *mp);
mp_err s_mp_mul_2(mp_int *mp);
mp_digit s_mp_norm(mp_int *a, mp_int *b);
mp_err s_mp_add_d(mp_int *mp, mp_digit d);
mp_err s_mp_sub_d(mp_int *mp, mp_digit d);
mp_err s_mp_mul_d(mp_int *mp, mp_digit d);
mp_err s_mp_div_d(mp_int *mp, mp_digit d, mp_digit *r);
mp_err s_mp_reduce(mp_int *x, mp_int *m, mp_int *mu);
mp_err s_mp_add(mp_int *a, mp_int *b);
mp_err s_mp_sub(mp_int *a, mp_int *b);
mp_err s_mp_mul(mp_int *a, mp_int *b);
#if 0
void s_mp_kmul(mp_digit *a, mp_digit *b, mp_digit *out, mp_size len);
#endif
#if MP_SQUARE
mp_err s_mp_sqr(mp_int *a);
#else
#define s_mp_sqr(a) s_mp_mul(a, a)
#endif
mp_err s_mp_div(mp_int *a, mp_int *b);
mp_err s_mp_2expt(mp_int *a, mp_digit k);
int s_mp_cmp(mp_int *a, mp_int *b);
int s_mp_cmp_d(mp_int *a, mp_digit d);
int s_mp_ispow2(mp_int *v);
int s_mp_ispow2d(mp_digit d);
int s_mp_tovalue(char ch, int r);
char s_mp_todigit(int val, int r, int low);
int s_mp_outlen(int bits, int r);
unsigned int mp_get_prec(void)
{
return s_mp_defprec;
}
void mp_set_prec(unsigned int prec)
{
if(prec == 0)
s_mp_defprec = MP_DEFPREC;
else
s_mp_defprec = prec;
}
mp_err mp_init(mp_int *mp)
{
return mp_init_size(mp, s_mp_defprec);
}
mp_err mp_init_array(mp_int mp[], int count)
{
mp_err res;
int pos;
ARGCHK(mp !=NULL && count > 0, MP_BADARG);
for(pos = 0; pos < count; ++pos) {
if((res = mp_init(&mp[pos])) != MP_OKAY)
goto CLEANUP;
}
return MP_OKAY;
CLEANUP:
while(--pos >= 0)
mp_clear(&mp[pos]);
return res;
}
mp_err mp_init_size(mp_int *mp, mp_size prec)
{
ARGCHK(mp != NULL && prec > 0, MP_BADARG);
if((DIGITS(mp) = s_mp_alloc(prec, sizeof(mp_digit))) == NULL)
return MP_MEM;
SIGN(mp) = MP_ZPOS;
USED(mp) = 1;
ALLOC(mp) = prec;
return MP_OKAY;
}
mp_err mp_init_copy(mp_int *mp, mp_int *from)
{
ARGCHK(mp != NULL && from != NULL, MP_BADARG);
if(mp == from)
return MP_OKAY;
if((DIGITS(mp) = s_mp_alloc(USED(from), sizeof(mp_digit))) == NULL)
return MP_MEM;
s_mp_copy(DIGITS(from), DIGITS(mp), USED(from));
USED(mp) = USED(from);
ALLOC(mp) = USED(from);
SIGN(mp) = SIGN(from);
return MP_OKAY;
}
mp_err mp_copy(mp_int *from, mp_int *to)
{
ARGCHK(from != NULL && to != NULL, MP_BADARG);
if(from == to)
return MP_OKAY;
{
mp_digit *tmp;
if(ALLOC(to) >= USED(from)) {
s_mp_setz(DIGITS(to) + USED(from), ALLOC(to) - USED(from));
s_mp_copy(DIGITS(from), DIGITS(to), USED(from));
} else {
if((tmp = s_mp_alloc(USED(from), sizeof(mp_digit))) == NULL)
return MP_MEM;
s_mp_copy(DIGITS(from), tmp, USED(from));
if(DIGITS(to) != NULL) {
#if MP_CRYPTO
s_mp_setz(DIGITS(to), ALLOC(to));
#endif
s_mp_free(DIGITS(to));
}
DIGITS(to) = tmp;
ALLOC(to) = USED(from);
}
USED(to) = USED(from);
SIGN(to) = SIGN(from);
}
return MP_OKAY;
}
void mp_exch(mp_int *mp1, mp_int *mp2)
{
#if MP_ARGCHK == 2
assert(mp1 != NULL && mp2 != NULL);
#else
if(mp1 == NULL || mp2 == NULL)
return;
#endif
s_mp_exch(mp1, mp2);
}
void mp_clear(mp_int *mp)
{
if(mp == NULL)
return;
if(DIGITS(mp) != NULL) {
#if MP_CRYPTO
s_mp_setz(DIGITS(mp), ALLOC(mp));
#endif
s_mp_free(DIGITS(mp));
DIGITS(mp) = NULL;
}
USED(mp) = 0;
ALLOC(mp) = 0;
}
void mp_clear_array(mp_int mp[], int count)
{
ARGCHK(mp != NULL && count > 0, MP_BADARG);
while(--count >= 0)
mp_clear(&mp[count]);
}
void mp_zero(mp_int *mp)
{
if(mp == NULL)
return;
s_mp_setz(DIGITS(mp), ALLOC(mp));
USED(mp) = 1;
SIGN(mp) = MP_ZPOS;
}
void mp_set(mp_int *mp, mp_digit d)
{
if(mp == NULL)
return;
mp_zero(mp);
DIGIT(mp, 0) = d;
}
mp_err mp_set_int(mp_int *mp, long z)
{
int ix;
unsigned long v = abs(z);
mp_err res;
ARGCHK(mp != NULL, MP_BADARG);
mp_zero(mp);
if(z == 0)
return MP_OKAY;
for(ix = sizeof(long) - 1; ix >= 0; ix--) {
if((res = s_mp_mul_2d(mp, CHAR_BIT)) != MP_OKAY)
return res;
res = s_mp_add_d(mp,
(mp_digit)((v >> (ix * CHAR_BIT)) & UCHAR_MAX));
if(res != MP_OKAY)
return res;
}
if(z < 0)
SIGN(mp) = MP_NEG;
return MP_OKAY;
}
mp_err mp_add_d(mp_int *a, mp_digit d, mp_int *b)
{
mp_err res = MP_OKAY;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
if(SIGN(b) == MP_ZPOS) {
res = s_mp_add_d(b, d);
} else if(s_mp_cmp_d(b, d) >= 0) {
res = s_mp_sub_d(b, d);
} else {
SIGN(b) = MP_ZPOS;
DIGIT(b, 0) = d - DIGIT(b, 0);
}
return res;
}
mp_err mp_sub_d(mp_int *a, mp_digit d, mp_int *b)
{
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
if(SIGN(b) == MP_NEG) {
if((res = s_mp_add_d(b, d)) != MP_OKAY)
return res;
} else if(s_mp_cmp_d(b, d) >= 0) {
if((res = s_mp_sub_d(b, d)) != MP_OKAY)
return res;
} else {
mp_neg(b, b);
DIGIT(b, 0) = d - DIGIT(b, 0);
SIGN(b) = MP_NEG;
}
if(s_mp_cmp_d(b, 0) == 0)
SIGN(b) = MP_ZPOS;
return MP_OKAY;
}
mp_err mp_mul_d(mp_int *a, mp_digit d, mp_int *b)
{
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if(d == 0) {
mp_zero(b);
return MP_OKAY;
}
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
res = s_mp_mul_d(b, d);
return res;
}
mp_err mp_mul_2(mp_int *a, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if((res = mp_copy(a, c)) != MP_OKAY)
return res;
return s_mp_mul_2(c);
}
mp_err mp_div_d(mp_int *a, mp_digit d, mp_int *q, mp_digit *r)
{
mp_err res;
mp_digit rem;
int pow;
ARGCHK(a != NULL, MP_BADARG);
if(d == 0)
return MP_RANGE;
if((pow = s_mp_ispow2d(d)) >= 0) {
mp_digit mask;
mask = (1 << pow) - 1;
rem = DIGIT(a, 0) & mask;
if(q) {
mp_copy(a, q);
s_mp_div_2d(q, pow);
}
if(r)
*r = rem;
return MP_OKAY;
}
if(q) {
if((res = mp_copy(a, q)) != MP_OKAY)
return res;
res = s_mp_div_d(q, d, &rem);
if(s_mp_cmp_d(q, 0) == MP_EQ)
SIGN(q) = MP_ZPOS;
} else {
mp_int qp;
if((res = mp_init_copy(&qp, a)) != MP_OKAY)
return res;
res = s_mp_div_d(&qp, d, &rem);
if(s_mp_cmp_d(&qp, 0) == 0)
SIGN(&qp) = MP_ZPOS;
mp_clear(&qp);
}
if(r)
*r = rem;
return res;
}
mp_err mp_div_2(mp_int *a, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if((res = mp_copy(a, c)) != MP_OKAY)
return res;
s_mp_div_2(c);
return MP_OKAY;
}
mp_err mp_expt_d(mp_int *a, mp_digit d, mp_int *c)
{
mp_int s, x;
mp_err res;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if((res = mp_init(&s)) != MP_OKAY)
return res;
if((res = mp_init_copy(&x, a)) != MP_OKAY)
goto X;
DIGIT(&s, 0) = 1;
while(d != 0) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY)
goto CLEANUP;
}
d >>= 1;
if((res = s_mp_sqr(&x)) != MP_OKAY)
goto CLEANUP;
}
s_mp_exch(&s, c);
CLEANUP:
mp_clear(&x);
X:
mp_clear(&s);
return res;
}
mp_err mp_abs(mp_int *a, mp_int *b)
{
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
SIGN(b) = MP_ZPOS;
return MP_OKAY;
}
mp_err mp_neg(mp_int *a, mp_int *b)
{
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
if(s_mp_cmp_d(b, 0) == MP_EQ)
SIGN(b) = MP_ZPOS;
else
SIGN(b) = (SIGN(b) == MP_NEG) ? MP_ZPOS : MP_NEG;
return MP_OKAY;
}
mp_err mp_add(mp_int *a, mp_int *b, mp_int *c)
{
mp_err res;
int cmp;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if(SIGN(a) == SIGN(b)) {
if(c == b) {
if((res = s_mp_add(c, a)) != MP_OKAY)
return res;
} else {
if(c != a && (res = mp_copy(a, c)) != MP_OKAY)
return res;
if((res = s_mp_add(c, b)) != MP_OKAY)
return res;
}
} else if((cmp = s_mp_cmp(a, b)) > 0) {
if(c == b) {
mp_int tmp;
if((res = mp_init_copy(&tmp, a)) != MP_OKAY)
return res;
if((res = s_mp_sub(&tmp, b)) != MP_OKAY) {
mp_clear(&tmp);
return res;
}
s_mp_exch(&tmp, c);
mp_clear(&tmp);
} else {
if(c != a && (res = mp_copy(a, c)) != MP_OKAY)
return res;
if((res = s_mp_sub(c, b)) != MP_OKAY)
return res;
}
} else if(cmp == 0) {
mp_zero(c);
return MP_OKAY;
} else {
if(c == a) {
mp_int tmp;
if((res = mp_init_copy(&tmp, b)) != MP_OKAY)
return res;
if((res = s_mp_sub(&tmp, a)) != MP_OKAY) {
mp_clear(&tmp);
return res;
}
s_mp_exch(&tmp, c);
mp_clear(&tmp);
} else {
if(c != b && (res = mp_copy(b, c)) != MP_OKAY)
return res;
if((res = s_mp_sub(c, a)) != MP_OKAY)
return res;
}
}
if(USED(c) == 1 && DIGIT(c, 0) == 0)
SIGN(c) = MP_ZPOS;
return MP_OKAY;
}
mp_err mp_sub(mp_int *a, mp_int *b, mp_int *c)
{
mp_err res;
int cmp;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if(SIGN(a) != SIGN(b)) {
if(c == a) {
if((res = s_mp_add(c, b)) != MP_OKAY)
return res;
} else {
if(c != b && ((res = mp_copy(b, c)) != MP_OKAY))
return res;
if((res = s_mp_add(c, a)) != MP_OKAY)
return res;
SIGN(c) = SIGN(a);
}
} else if((cmp = s_mp_cmp(a, b)) > 0) {
if(c == b) {
mp_int tmp;
if((res = mp_init_copy(&tmp, a)) != MP_OKAY)
return res;
if((res = s_mp_sub(&tmp, b)) != MP_OKAY) {
mp_clear(&tmp);
return res;
}
s_mp_exch(&tmp, c);
mp_clear(&tmp);
} else {
if(c != a && ((res = mp_copy(a, c)) != MP_OKAY))
return res;
if((res = s_mp_sub(c, b)) != MP_OKAY)
return res;
}
} else if(cmp == 0) {
mp_zero(c);
return MP_OKAY;
} else {
if(c == a) {
mp_int tmp;
if((res = mp_init_copy(&tmp, b)) != MP_OKAY)
return res;
if((res = s_mp_sub(&tmp, a)) != MP_OKAY) {
mp_clear(&tmp);
return res;
}
s_mp_exch(&tmp, c);
mp_clear(&tmp);
} else {
if(c != b && ((res = mp_copy(b, c)) != MP_OKAY))
return res;
if((res = s_mp_sub(c, a)) != MP_OKAY)
return res;
}
SIGN(c) = !SIGN(b);
}
if(USED(c) == 1 && DIGIT(c, 0) == 0)
SIGN(c) = MP_ZPOS;
return MP_OKAY;
}
mp_err mp_mul(mp_int *a, mp_int *b, mp_int *c)
{
mp_err res;
mp_sign sgn;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
sgn = (SIGN(a) == SIGN(b)) ? MP_ZPOS : MP_NEG;
if(c == b) {
if((res = s_mp_mul(c, a)) != MP_OKAY)
return res;
} else {
if((res = mp_copy(a, c)) != MP_OKAY)
return res;
if((res = s_mp_mul(c, b)) != MP_OKAY)
return res;
}
if(sgn == MP_ZPOS || s_mp_cmp_d(c, 0) == MP_EQ)
SIGN(c) = MP_ZPOS;
else
SIGN(c) = sgn;
return MP_OKAY;
}
mp_err mp_mul_2d(mp_int *a, mp_digit d, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if((res = mp_copy(a, c)) != MP_OKAY)
return res;
if(d == 0)
return MP_OKAY;
return s_mp_mul_2d(c, d);
}
#if MP_SQUARE
mp_err mp_sqr(mp_int *a, mp_int *b)
{
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if((res = mp_copy(a, b)) != MP_OKAY)
return res;
if((res = s_mp_sqr(b)) != MP_OKAY)
return res;
SIGN(b) = MP_ZPOS;
return MP_OKAY;
}
#endif
mp_err mp_div(mp_int *a, mp_int *b, mp_int *q, mp_int *r)
{
mp_err res;
mp_int qtmp, rtmp;
int cmp;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if(mp_cmp_z(b) == MP_EQ)
return MP_RANGE;
if((cmp = s_mp_cmp(a, b)) < 0) {
if(r) {
if((res = mp_copy(a, r)) != MP_OKAY)
return res;
}
if(q)
mp_zero(q);
return MP_OKAY;
} else if(cmp == 0) {
if(q) {
int qneg = (SIGN(a) != SIGN(b));
mp_set(q, 1);
if(qneg)
SIGN(q) = MP_NEG;
}
if(r)
mp_zero(r);
return MP_OKAY;
}
if((res = mp_init_copy(&qtmp, a)) != MP_OKAY)
return res;
if((res = mp_init_copy(&rtmp, b)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_div(&qtmp, &rtmp)) != MP_OKAY)
goto CLEANUP;
SIGN(&rtmp) = SIGN(a);
if(SIGN(a) == SIGN(b))
SIGN(&qtmp) = MP_ZPOS;
else
SIGN(&qtmp) = MP_NEG;
if(s_mp_cmp_d(&qtmp, 0) == MP_EQ)
SIGN(&qtmp) = MP_ZPOS;
if(s_mp_cmp_d(&rtmp, 0) == MP_EQ)
SIGN(&rtmp) = MP_ZPOS;
if(q)
s_mp_exch(&qtmp, q);
if(r)
s_mp_exch(&rtmp, r);
CLEANUP:
mp_clear(&rtmp);
mp_clear(&qtmp);
return res;
}
mp_err mp_div_2d(mp_int *a, mp_digit d, mp_int *q, mp_int *r)
{
mp_err res;
ARGCHK(a != NULL, MP_BADARG);
if(q) {
if((res = mp_copy(a, q)) != MP_OKAY)
return res;
s_mp_div_2d(q, d);
}
if(r) {
if((res = mp_copy(a, r)) != MP_OKAY)
return res;
s_mp_mod_2d(r, d);
}
return MP_OKAY;
}
mp_err mp_expt(mp_int *a, mp_int *b, mp_int *c)
{
mp_int s, x;
mp_err res;
mp_digit d;
int dig, bit;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if(mp_cmp_z(b) < 0)
return MP_RANGE;
if((res = mp_init(&s)) != MP_OKAY)
return res;
mp_set(&s, 1);
if((res = mp_init_copy(&x, a)) != MP_OKAY)
goto X;
for(dig = 0; dig < (USED(b) - 1); dig++) {
d = DIGIT(b, dig);
for(bit = 0; bit < DIGIT_BIT; bit++) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY)
goto CLEANUP;
}
d >>= 1;
if((res = s_mp_sqr(&x)) != MP_OKAY)
goto CLEANUP;
}
}
d = DIGIT(b, dig);
while(d) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY)
goto CLEANUP;
}
d >>= 1;
if((res = s_mp_sqr(&x)) != MP_OKAY)
goto CLEANUP;
}
if(mp_iseven(b))
SIGN(&s) = SIGN(a);
res = mp_copy(&s, c);
CLEANUP:
mp_clear(&x);
X:
mp_clear(&s);
return res;
}
mp_err mp_2expt(mp_int *a, mp_digit k)
{
ARGCHK(a != NULL, MP_BADARG);
return s_mp_2expt(a, k);
}
mp_err mp_mod(mp_int *a, mp_int *m, mp_int *c)
{
mp_err res;
int mag;
ARGCHK(a != NULL && m != NULL && c != NULL, MP_BADARG);
if(SIGN(m) == MP_NEG)
return MP_RANGE;
if((mag = s_mp_cmp(a, m)) > 0) {
if((res = mp_div(a, m, NULL, c)) != MP_OKAY)
return res;
if(SIGN(c) == MP_NEG) {
if((res = mp_add(c, m, c)) != MP_OKAY)
return res;
}
} else if(mag < 0) {
if((res = mp_copy(a, c)) != MP_OKAY)
return res;
if(mp_cmp_z(a) < 0) {
if((res = mp_add(c, m, c)) != MP_OKAY)
return res;
}
} else {
mp_zero(c);
}
return MP_OKAY;
}
mp_err mp_mod_d(mp_int *a, mp_digit d, mp_digit *c)
{
mp_err res;
mp_digit rem;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if(s_mp_cmp_d(a, d) > 0) {
if((res = mp_div_d(a, d, NULL, &rem)) != MP_OKAY)
return res;
} else {
if(SIGN(a) == MP_NEG)
rem = d - DIGIT(a, 0);
else
rem = DIGIT(a, 0);
}
if(c)
*c = rem;
return MP_OKAY;
}
mp_err mp_sqrt(mp_int *a, mp_int *b)
{
mp_int x, t;
mp_err res;
ARGCHK(a != NULL && b != NULL, MP_BADARG);
if(SIGN(a) == MP_NEG)
return MP_RANGE;
if(mp_cmp_d(a, 0) == MP_EQ || mp_cmp_d(a, 1) == MP_EQ)
return mp_copy(a, b);
if((res = mp_init_size(&t, USED(a))) != MP_OKAY)
return res;
if((res = mp_init_copy(&x, a)) != MP_OKAY)
goto X;
s_mp_rshd(&x, (USED(&x)/2)+1);
mp_add_d(&x, 1, &x);
for(;;) {
mp_copy(&x, &t);
if((res = mp_sqr(&t, &t)) != MP_OKAY ||
(res = mp_sub(&t, a, &t)) != MP_OKAY)
goto CLEANUP;
s_mp_mul_2(&x);
if((res = mp_div(&t, &x, &t, NULL)) != MP_OKAY)
goto CLEANUP;
s_mp_div_2(&x);
if(mp_cmp_z(&t) == MP_EQ)
break;
if((res = mp_sub(&x, &t, &x)) != MP_OKAY)
goto CLEANUP;
}
mp_sub_d(&x, 1, &x);
s_mp_exch(&x, b);
CLEANUP:
mp_clear(&x);
X:
mp_clear(&t);
return res;
}
#if MP_MODARITH
mp_err mp_addmod(mp_int *a, mp_int *b, mp_int *m, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && b != NULL && m != NULL && c != NULL, MP_BADARG);
if((res = mp_add(a, b, c)) != MP_OKAY)
return res;
if((res = mp_mod(c, m, c)) != MP_OKAY)
return res;
return MP_OKAY;
}
mp_err mp_submod(mp_int *a, mp_int *b, mp_int *m, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && b != NULL && m != NULL && c != NULL, MP_BADARG);
if((res = mp_sub(a, b, c)) != MP_OKAY)
return res;
if((res = mp_mod(c, m, c)) != MP_OKAY)
return res;
return MP_OKAY;
}
mp_err mp_mulmod(mp_int *a, mp_int *b, mp_int *m, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && b != NULL && m != NULL && c != NULL, MP_BADARG);
if((res = mp_mul(a, b, c)) != MP_OKAY)
return res;
if((res = mp_mod(c, m, c)) != MP_OKAY)
return res;
return MP_OKAY;
}
#if MP_SQUARE
mp_err mp_sqrmod(mp_int *a, mp_int *m, mp_int *c)
{
mp_err res;
ARGCHK(a != NULL && m != NULL && c != NULL, MP_BADARG);
if((res = mp_sqr(a, c)) != MP_OKAY)
return res;
if((res = mp_mod(c, m, c)) != MP_OKAY)
return res;
return MP_OKAY;
}
#endif
mp_err mp_exptmod(mp_int *a, mp_int *b, mp_int *m, mp_int *c)
{
mp_int s, x, mu;
mp_err res;
mp_digit d, *db = DIGITS(b);
mp_size ub = USED(b);
int dig, bit;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if(mp_cmp_z(b) < 0 || mp_cmp_z(m) <= 0)
return MP_RANGE;
if((res = mp_init(&s)) != MP_OKAY)
return res;
if((res = mp_init_copy(&x, a)) != MP_OKAY)
goto X;
if((res = mp_mod(&x, m, &x)) != MP_OKAY ||
(res = mp_init(&mu)) != MP_OKAY)
goto MU;
mp_set(&s, 1);
s_mp_add_d(&mu, 1);
s_mp_lshd(&mu, 2 * USED(m));
if((res = mp_div(&mu, m, &mu, NULL)) != MP_OKAY)
goto CLEANUP;
for(dig = 0; dig < (ub - 1); dig++) {
d = *db++;
for(bit = 0; bit < DIGIT_BIT; bit++) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_reduce(&s, m, &mu)) != MP_OKAY)
goto CLEANUP;
}
d >>= 1;
if((res = s_mp_sqr(&x)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_reduce(&x, m, &mu)) != MP_OKAY)
goto CLEANUP;
}
}
d = *db;
while(d) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_reduce(&s, m, &mu)) != MP_OKAY)
goto CLEANUP;
}
d >>= 1;
if((res = s_mp_sqr(&x)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_reduce(&x, m, &mu)) != MP_OKAY)
goto CLEANUP;
}
s_mp_exch(&s, c);
CLEANUP:
mp_clear(&mu);
MU:
mp_clear(&x);
X:
mp_clear(&s);
return res;
}
mp_err mp_exptmod_d(mp_int *a, mp_digit d, mp_int *m, mp_int *c)
{
mp_int s, x;
mp_err res;
ARGCHK(a != NULL && c != NULL, MP_BADARG);
if((res = mp_init(&s)) != MP_OKAY)
return res;
if((res = mp_init_copy(&x, a)) != MP_OKAY)
goto X;
mp_set(&s, 1);
while(d != 0) {
if(d & 1) {
if((res = s_mp_mul(&s, &x)) != MP_OKAY ||
(res = mp_mod(&s, m, &s)) != MP_OKAY)
goto CLEANUP;
}
d /= 2;
if((res = s_mp_sqr(&x)) != MP_OKAY ||
(res = mp_mod(&x, m, &x)) != MP_OKAY)
goto CLEANUP;
}
s_mp_exch(&s, c);
CLEANUP:
mp_clear(&x);
X:
mp_clear(&s);
return res;
}
#endif
int mp_cmp_z(mp_int *a)
{
if(SIGN(a) == MP_NEG)
return MP_LT;
else if(USED(a) == 1 && DIGIT(a, 0) == 0)
return MP_EQ;
else
return MP_GT;
}
int mp_cmp_d(mp_int *a, mp_digit d)
{
ARGCHK(a != NULL, MP_EQ);
if(SIGN(a) == MP_NEG)
return MP_LT;
return s_mp_cmp_d(a, d);
}
int mp_cmp(mp_int *a, mp_int *b)
{
ARGCHK(a != NULL && b != NULL, MP_EQ);
if(SIGN(a) == SIGN(b)) {
int mag;
if((mag = s_mp_cmp(a, b)) == MP_EQ)
return MP_EQ;
if(SIGN(a) == MP_ZPOS)
return mag;
else
return -mag;
} else if(SIGN(a) == MP_ZPOS) {
return MP_GT;
} else {
return MP_LT;
}
}
int mp_cmp_mag(mp_int *a, mp_int *b)
{
ARGCHK(a != NULL && b != NULL, MP_EQ);
return s_mp_cmp(a, b);
}
int mp_cmp_int(mp_int *a, long z)
{
mp_int tmp;
int out;
ARGCHK(a != NULL, MP_EQ);
mp_init(&tmp); mp_set_int(&tmp, z);
out = mp_cmp(a, &tmp);
mp_clear(&tmp);
return out;
}
int mp_isodd(mp_int *a)
{
ARGCHK(a != NULL, 0);
return (DIGIT(a, 0) & 1);
}
int mp_iseven(mp_int *a)
{
return !mp_isodd(a);
}
#if MP_NUMTH
mp_err mp_gcd(mp_int *a, mp_int *b, mp_int *c)
{
mp_err res;
mp_int u, v, t;
mp_size k = 0;
ARGCHK(a != NULL && b != NULL && c != NULL, MP_BADARG);
if(mp_cmp_z(a) == MP_EQ && mp_cmp_z(b) == MP_EQ)
return MP_RANGE;
if(mp_cmp_z(a) == MP_EQ) {
return mp_copy(b, c);
} else if(mp_cmp_z(b) == MP_EQ) {
return mp_copy(a, c);
}
if((res = mp_init(&t)) != MP_OKAY)
return res;
if((res = mp_init_copy(&u, a)) != MP_OKAY)
goto U;
if((res = mp_init_copy(&v, b)) != MP_OKAY)
goto V;
SIGN(&u) = MP_ZPOS;
SIGN(&v) = MP_ZPOS;
while(mp_iseven(&u) && mp_iseven(&v)) {
s_mp_div_2(&u);
s_mp_div_2(&v);
++k;
}
if(mp_isodd(&u)) {
if((res = mp_copy(&v, &t)) != MP_OKAY)
goto CLEANUP;
if(SIGN(&v) == MP_ZPOS)
SIGN(&t) = MP_NEG;
else
SIGN(&t) = MP_ZPOS;
} 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;
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);
res = mp_mul(&u, &v, c);
CLEANUP:
mp_clear(&v);
V:
mp_clear(&u);
U:
mp_clear(&t);
return res;
}
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);
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;
}
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;
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);
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);
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(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;
}
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;
}
#endif
#if MP_IOFUNC
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));
}
}
#endif
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) {
if(str[0])
SIGN(mp) = MP_NEG;
else
SIGN(mp) = MP_ZPOS;
}
return res;
}
int mp_signed_bin_size(mp_int *mp)
{
ARGCHK(mp != NULL, 0);
return mp_unsigned_bin_size(mp) + 1;
}
mp_err mp_to_signed_bin(mp_int *mp, unsigned char *str)
{
ARGCHK(mp != NULL && str != NULL, MP_BADARG);
str[0] = (char)SIGN(mp);
return mp_to_unsigned_bin(mp, str + 1);
}
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;
}
int mp_unsigned_bin_size(mp_int *mp)
{
mp_digit topdig;
int count;
ARGCHK(mp != NULL, 0);
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;
}
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;
if(dp == end && *dp == 0) {
*str = '\0';
return MP_OKAY;
}
while(dp < end) {
int ix;
d = *dp;
for(ix = 0; ix < sizeof(mp_digit); ++ix) {
*spos = d & UCHAR_MAX;
d >>= CHAR_BIT;
++spos;
}
++dp;
}
d = *end;
while(d != 0) {
*spos = d & UCHAR_MAX;
d >>= CHAR_BIT;
++spos;
}
while(--spos > str) {
unsigned char t = *str;
*str = *spos;
*spos = t;
++str;
}
return MP_OKAY;
}
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;
}
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);
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;
++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;
}
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;
if(mp_cmp_z(mp) < 0)
++len;
return len;
}
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);
}
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;
sgn = SIGN(&tmp); SIGN(&tmp) = MP_ZPOS;
while(mp_cmp_z(&tmp) != 0) {
if((res = s_mp_div_d(&tmp, rdx, &rem)) != MP_OKAY) {
mp_clear(&tmp);
return res;
}
ch = s_mp_todigit(rem, radix, 0);
str[pos++] = ch;
}
if(sgn == MP_NEG)
str[pos++] = '-';
str[pos--] = '\0';
ix = 0;
while(ix < pos) {
char tmp = str[ix];
str[ix] = str[pos];
str[pos] = tmp;
++ix;
--pos;
}
mp_clear(&tmp);
}
return MP_OKAY;
}
int mp_char2value(char ch, int r)
{
return s_mp_tovalue(ch, r);
}
const char *mp_strerror(mp_err ec)
{
int aec = (ec < 0) ? -ec : ec;
if(ec < MP_LAST_CODE || ec > MP_OKAY) {
return mp_err_string[0];
} else {
return mp_err_string[aec + 1];
}
}
mp_err s_mp_grow(mp_int *mp, mp_size min)
{
if(min > ALLOC(mp)) {
mp_digit *tmp;
min = ((min + (s_mp_defprec - 1)) / s_mp_defprec) * s_mp_defprec;
if((tmp = s_mp_alloc(min, sizeof(mp_digit))) == NULL)
return MP_MEM;
s_mp_copy(DIGITS(mp), tmp, USED(mp));
#if MP_CRYPTO
s_mp_setz(DIGITS(mp), ALLOC(mp));
#endif
s_mp_free(DIGITS(mp));
DIGITS(mp) = tmp;
ALLOC(mp) = min;
}
return MP_OKAY;
}
mp_err s_mp_pad(mp_int *mp, mp_size min)
{
if(min > USED(mp)) {
mp_err res;
if(min > ALLOC(mp) && (res = s_mp_grow(mp, min)) != MP_OKAY)
return res;
USED(mp) = min;
}
return MP_OKAY;
}
#if MP_MACRO == 0
void s_mp_setz(mp_digit *dp, mp_size count)
{
#if MP_MEMSET == 0
int ix;
for(ix = 0; ix < count; ix++)
dp[ix] = 0;
#else
memset(dp, 0, count * sizeof(mp_digit));
#endif
}
#endif
#if MP_MACRO == 0
void s_mp_copy(mp_digit *sp, mp_digit *dp, mp_size count)
{
#if MP_MEMCPY == 0
int ix;
for(ix = 0; ix < count; ix++)
dp[ix] = sp[ix];
#else
memcpy(dp, sp, count * sizeof(mp_digit));
#endif
}
#endif
#if MP_MACRO == 0
void *s_mp_alloc(size_t nb, size_t ni)
{
return calloc(nb, ni);
}
#endif
#if MP_MACRO == 0
void s_mp_free(void *ptr)
{
if(ptr)
free(ptr);
}
#endif
void s_mp_clamp(mp_int *mp)
{
mp_size du = USED(mp);
mp_digit *zp = DIGITS(mp) + du - 1;
while(du > 1 && !*zp--)
--du;
USED(mp) = du;
}
void s_mp_exch(mp_int *a, mp_int *b)
{
mp_int tmp;
tmp = *a;
*a = *b;
*b = tmp;
}
mp_err s_mp_lshd(mp_int *mp, mp_size p)
{
mp_err res;
mp_size pos;
mp_digit *dp;
int ix;
if(p == 0)
return MP_OKAY;
if((res = s_mp_pad(mp, USED(mp) + p)) != MP_OKAY)
return res;
pos = USED(mp) - 1;
dp = DIGITS(mp);
for(ix = pos - p; ix >= 0; ix--)
dp[ix + p] = dp[ix];
for(ix = 0; ix < p; ix++)
dp[ix] = 0;
return MP_OKAY;
}
void s_mp_rshd(mp_int *mp, mp_size p)
{
mp_size ix;
mp_digit *dp;
if(p == 0)
return;
if(p >= USED(mp)) {
s_mp_setz(DIGITS(mp), ALLOC(mp));
USED(mp) = 1;
SIGN(mp) = MP_ZPOS;
return;
}
dp = DIGITS(mp);
for(ix = p; ix < USED(mp); ix++)
dp[ix - p] = dp[ix];
ix -= p;
while(ix < USED(mp))
dp[ix++] = 0;
s_mp_clamp(mp);
}
void s_mp_div_2(mp_int *mp)
{
s_mp_div_2d(mp, 1);
}
mp_err s_mp_mul_2(mp_int *mp)
{
int ix;
mp_digit kin = 0, kout, *dp = DIGITS(mp);
mp_err res;
for(ix = 0; ix < USED(mp); ix++) {
kout = (dp[ix] >> (DIGIT_BIT - 1)) & 1;
dp[ix] = (dp[ix] << 1) | kin;
kin = kout;
}
if(kin) {
if(ix >= ALLOC(mp)) {
if((res = s_mp_grow(mp, ALLOC(mp) + 1)) != MP_OKAY)
return res;
dp = DIGITS(mp);
}
dp[ix] = kin;
USED(mp) += 1;
}
return MP_OKAY;
}
void s_mp_mod_2d(mp_int *mp, mp_digit d)
{
unsigned int ndig = (d / DIGIT_BIT), nbit = (d % DIGIT_BIT);
unsigned int ix;
mp_digit dmask, *dp = DIGITS(mp);
if(ndig >= USED(mp))
return;
dmask = (1 << nbit) - 1;
dp[ndig] &= dmask;
for(ix = ndig + 1; ix < USED(mp); ix++)
dp[ix] = 0;
s_mp_clamp(mp);
}
mp_err s_mp_mul_2d(mp_int *mp, mp_digit d)
{
mp_err res;
mp_digit save, next, mask, *dp;
mp_size used;
int ix;
if((res = s_mp_lshd(mp, d / DIGIT_BIT)) != MP_OKAY)
return res;
dp = DIGITS(mp); used = USED(mp);
d %= DIGIT_BIT;
mask = (1 << d) - 1;
if((dp[used - 1] >> (DIGIT_BIT - d)) & mask) {
if((res = s_mp_grow(mp, used + 1)) != MP_OKAY)
return res;
dp = DIGITS(mp);
}
save = 0;
for(ix = 0; ix < used; ix++) {
next = (dp[ix] >> (DIGIT_BIT - d)) & mask;
dp[ix] = (dp[ix] << d) | save;
save = next;
}
if(save) {
dp[used] = save;
USED(mp) += 1;
}
s_mp_clamp(mp);
return MP_OKAY;
}
void s_mp_div_2d(mp_int *mp, mp_digit d)
{
int ix;
mp_digit save, next, mask, *dp = DIGITS(mp);
s_mp_rshd(mp, d / DIGIT_BIT);
d %= DIGIT_BIT;
mask = (1 << d) - 1;
save = 0;
for(ix = USED(mp) - 1; ix >= 0; ix--) {
next = dp[ix] & mask;
dp[ix] = (dp[ix] >> d) | (save << (DIGIT_BIT - d));
save = next;
}
s_mp_clamp(mp);
}
mp_digit s_mp_norm(mp_int *a, mp_int *b)
{
mp_digit t, d = 0;
t = DIGIT(b, USED(b) - 1);
while(t < (RADIX / 2)) {
t <<= 1;
++d;
}
if(d != 0) {
s_mp_mul_2d(a, d);
s_mp_mul_2d(b, d);
}
return d;
}
mp_err s_mp_add_d(mp_int *mp, mp_digit d)
{
mp_word w, k = 0;
mp_size ix = 1, used = USED(mp);
mp_digit *dp = DIGITS(mp);
w = dp[0] + d;
dp[0] = ACCUM(w);
k = CARRYOUT(w);
while(ix < used && k) {
w = dp[ix] + k;
dp[ix] = ACCUM(w);
k = CARRYOUT(w);
++ix;
}
if(k != 0) {
mp_err res;
if((res = s_mp_pad(mp, USED(mp) + 1)) != MP_OKAY)
return res;
DIGIT(mp, ix) = k;
}
return MP_OKAY;
}
mp_err s_mp_sub_d(mp_int *mp, mp_digit d)
{
mp_word w, b = 0;
mp_size ix = 1, used = USED(mp);
mp_digit *dp = DIGITS(mp);
w = (RADIX + dp[0]) - d;
b = CARRYOUT(w) ? 0 : 1;
dp[0] = ACCUM(w);
while(b && ix < used) {
w = (RADIX + dp[ix]) - b;
b = CARRYOUT(w) ? 0 : 1;
dp[ix] = ACCUM(w);
++ix;
}
s_mp_clamp(mp);
if(b)
return MP_RANGE;
else
return MP_OKAY;
}
mp_err s_mp_mul_d(mp_int *a, mp_digit d)
{
mp_word w, k = 0;
mp_size ix, max;
mp_err res;
mp_digit *dp = DIGITS(a);
max = USED(a);
w = dp[max - 1] * d;
if(CARRYOUT(w) != 0) {
if((res = s_mp_pad(a, max + 1)) != MP_OKAY)
return res;
dp = DIGITS(a);
}
for(ix = 0; ix < max; ix++) {
w = (dp[ix] * d) + k;
dp[ix] = ACCUM(w);
k = CARRYOUT(w);
}
if(k) {
dp[max] = k;
USED(a) = max + 1;
}
s_mp_clamp(a);
return MP_OKAY;
}
mp_err s_mp_div_d(mp_int *mp, mp_digit d, mp_digit *r)
{
mp_word w = 0, t;
mp_int quot;
mp_err res;
mp_digit *dp = DIGITS(mp), *qp;
int ix;
if(d == 0)
return MP_RANGE;
if((res = mp_init_size(", USED(mp))) != MP_OKAY)
return res;
USED(") = USED(mp);
qp = DIGITS(");
for(ix = USED(mp) - 1; ix >= 0; ix--) {
w = (w << DIGIT_BIT) | dp[ix];
if(w >= d) {
t = w / d;
w = w % d;
} else {
t = 0;
}
qp[ix] = t;
}
if(r)
*r = w;
s_mp_clamp(");
mp_exch(", mp);
mp_clear(");
return MP_OKAY;
}
mp_err s_mp_add(mp_int *a, mp_int *b)
{
mp_word w = 0;
mp_digit *pa, *pb;
mp_size ix, used = USED(b);
mp_err res;
if((used > USED(a)) && (res = s_mp_pad(a, used)) != MP_OKAY)
return res;
pa = DIGITS(a);
pb = DIGITS(b);
for(ix = 0; ix < used; ++ix) {
w += *pa + *pb++;
*pa++ = ACCUM(w);
w = CARRYOUT(w);
}
used = USED(a);
while(w && ix < used) {
w += *pa;
*pa++ = ACCUM(w);
w = CARRYOUT(w);
++ix;
}
if(w) {
if((res = s_mp_pad(a, used + 1)) != MP_OKAY)
return res;
DIGIT(a, ix) = w;
}
return MP_OKAY;
}
mp_err s_mp_sub(mp_int *a, mp_int *b)
{
mp_word w = 0;
mp_digit *pa, *pb;
mp_size ix, used = USED(b);
pa = DIGITS(a);
pb = DIGITS(b);
for(ix = 0; ix < used; ++ix) {
w = (RADIX + *pa) - w - *pb++;
*pa++ = ACCUM(w);
w = CARRYOUT(w) ? 0 : 1;
}
used = USED(a);
while(ix < used) {
w = RADIX + *pa - w;
*pa++ = ACCUM(w);
w = CARRYOUT(w) ? 0 : 1;
++ix;
}
s_mp_clamp(a);
if(w)
return MP_RANGE;
else
return MP_OKAY;
}
mp_err s_mp_reduce(mp_int *x, mp_int *m, mp_int *mu)
{
mp_int q;
mp_err res;
mp_size um = USED(m);
if((res = mp_init_copy(&q, x)) != MP_OKAY)
return res;
s_mp_rshd(&q, um - 1);
s_mp_mul(&q, mu);
s_mp_rshd(&q, um + 1);
s_mp_mod_2d(x, (mp_digit)(DIGIT_BIT * (um + 1)));
#ifndef SHRT_MUL
s_mp_mul(&q, m);
s_mp_mod_2d(&q, (mp_digit)(DIGIT_BIT * (um + 1)));
#else
s_mp_mul_dig(&q, m, um + 1);
#endif
if((res = mp_sub(x, &q, x)) != MP_OKAY)
goto CLEANUP;
if(mp_cmp_z(x) < 0) {
mp_set(&q, 1);
if((res = s_mp_lshd(&q, um + 1)) != MP_OKAY)
goto CLEANUP;
if((res = mp_add(x, &q, x)) != MP_OKAY)
goto CLEANUP;
}
while(mp_cmp(x, m) >= 0) {
if((res = s_mp_sub(x, m)) != MP_OKAY)
break;
}
CLEANUP:
mp_clear(&q);
return res;
}
mp_err s_mp_mul(mp_int *a, mp_int *b)
{
mp_word w, k = 0;
mp_int tmp;
mp_err res;
mp_size ix, jx, ua = USED(a), ub = USED(b);
mp_digit *pa, *pb, *pt, *pbt;
if((res = mp_init_size(&tmp, ua + ub)) != MP_OKAY)
return res;
USED(&tmp) = ua + ub;
pbt = DIGITS(&tmp);
pb = DIGITS(b);
for(ix = 0; ix < ub; ++ix, ++pb) {
if(*pb == 0)
continue;
pa = DIGITS(a);
for(jx = 0; jx < ua; ++jx, ++pa) {
pt = pbt + ix + jx;
w = *pb * *pa + k + *pt;
*pt = ACCUM(w);
k = CARRYOUT(w);
}
pbt[ix + jx] = k;
k = 0;
}
s_mp_clamp(&tmp);
s_mp_exch(&tmp, a);
mp_clear(&tmp);
return MP_OKAY;
}
#if 0
void s_mp_kmul(mp_digit *a, mp_digit *b, mp_digit *out, mp_size len)
{
mp_word w, k = 0;
mp_size ix, jx;
mp_digit *pa, *pt;
for(ix = 0; ix < len; ++ix, ++b) {
if(*b == 0)
continue;
pa = a;
for(jx = 0; jx < len; ++jx, ++pa) {
pt = out + ix + jx;
w = *b * *pa + k + *pt;
*pt = ACCUM(w);
k = CARRYOUT(w);
}
out[ix + jx] = k;
k = 0;
}
}
#endif
#if MP_SQUARE
mp_err s_mp_sqr(mp_int *a)
{
mp_word w, k = 0;
mp_int tmp;
mp_err res;
mp_size ix, jx, kx, used = USED(a);
mp_digit *pa1, *pa2, *pt, *pbt;
if((res = mp_init_size(&tmp, 2 * used)) != MP_OKAY)
return res;
USED(&tmp) = 2 * used;
pbt = DIGITS(&tmp);
pa1 = DIGITS(a);
for(ix = 0; ix < used; ++ix, ++pa1) {
if(*pa1 == 0)
continue;
w = DIGIT(&tmp, ix + ix) + (*pa1 * *pa1);
pbt[ix + ix] = ACCUM(w);
k = CARRYOUT(w);
for(jx = ix + 1, pa2 = DIGITS(a) + jx; jx < used; ++jx, ++pa2) {
mp_word u = 0, v;
pt = pbt + ix + jx;
w = *pa1 * *pa2;
u = (w >> (MP_WORD_BIT - 1)) & 1;
w *= 2;
v = *pt + k;
u |= ((MP_WORD_MAX - v) < w);
w += v;
*pt = ACCUM(w);
k = CARRYOUT(w) | (u << DIGIT_BIT);
}
k = DIGIT(&tmp, ix + jx) + k;
pbt[ix + jx] = ACCUM(k);
k = CARRYOUT(k);
kx = 1;
while(k) {
k = pbt[ix + jx + kx] + 1;
pbt[ix + jx + kx] = ACCUM(k);
k = CARRYOUT(k);
++kx;
}
}
s_mp_clamp(&tmp);
s_mp_exch(&tmp, a);
mp_clear(&tmp);
return MP_OKAY;
}
#endif
mp_err s_mp_div(mp_int *a, mp_int *b)
{
mp_int quot, rem, t;
mp_word q;
mp_err res;
mp_digit d;
int ix;
if(mp_cmp_z(b) == 0)
return MP_RANGE;
if((ix = s_mp_ispow2(b)) >= 0) {
mp_copy(a, b);
s_mp_div_2d(a, (mp_digit)ix);
s_mp_mod_2d(b, (mp_digit)ix);
return MP_OKAY;
}
if((res = mp_init_size(", USED(a))) != MP_OKAY)
return res;
if((res = mp_init_size(&t, USED(a))) != MP_OKAY)
goto T;
if((res = mp_init_size(&rem, USED(a))) != MP_OKAY)
goto REM;
d = s_mp_norm(a, b);
ix = USED(a) - 1;
while(ix >= 0) {
while(s_mp_cmp(&rem, b) < 0 && ix >= 0) {
if((res = s_mp_lshd(&rem, 1)) != MP_OKAY)
goto CLEANUP;
if((res = s_mp_lshd(", 1)) != MP_OKAY)
goto CLEANUP;
DIGIT(&rem, 0) = DIGIT(a, ix);
s_mp_clamp(&rem);
--ix;
}
if(s_mp_cmp(&rem, b) < 0)
break;
q = DIGIT(&rem, USED(&rem) - 1);
if(q <= DIGIT(b, USED(b) - 1) && USED(&rem) > 1)
q = (q << DIGIT_BIT) | DIGIT(&rem, USED(&rem) - 2);
q /= DIGIT(b, USED(b) - 1);
if(q >= RADIX)
q = RADIX - 1;
mp_copy(b, &t);
if((res = s_mp_mul_d(&t, q)) != MP_OKAY)
goto CLEANUP;
while(s_mp_cmp(&t, &rem) > 0) {
--q;
s_mp_sub(&t, b);
}
if((res = s_mp_sub(&rem, &t)) != MP_OKAY)
goto CLEANUP;
DIGIT(", 0) = q;
}
if(d != 0)
s_mp_div_2d(&rem, d);
s_mp_clamp(");
s_mp_clamp(&rem);
s_mp_exch(", a);
s_mp_exch(&rem, b);
CLEANUP:
mp_clear(&rem);
REM:
mp_clear(&t);
T:
mp_clear(");
return res;
}
mp_err s_mp_2expt(mp_int *a, mp_digit k)
{
mp_err res;
mp_size dig, bit;
dig = k / DIGIT_BIT;
bit = k % DIGIT_BIT;
mp_zero(a);
if((res = s_mp_pad(a, dig + 1)) != MP_OKAY)
return res;
DIGIT(a, dig) |= (1 << bit);
return MP_OKAY;
}
int s_mp_cmp(mp_int *a, mp_int *b)
{
mp_size ua = USED(a), ub = USED(b);
if(ua > ub)
return MP_GT;
else if(ua < ub)
return MP_LT;
else {
int ix = ua - 1;
mp_digit *ap = DIGITS(a) + ix, *bp = DIGITS(b) + ix;
while(ix >= 0) {
if(*ap > *bp)
return MP_GT;
else if(*ap < *bp)
return MP_LT;
--ap; --bp; --ix;
}
return MP_EQ;
}
}
int s_mp_cmp_d(mp_int *a, mp_digit d)
{
mp_size ua = USED(a);
mp_digit *ap = DIGITS(a);
if(ua > 1)
return MP_GT;
if(*ap < d)
return MP_LT;
else if(*ap > d)
return MP_GT;
else
return MP_EQ;
}
int s_mp_ispow2(mp_int *v)
{
mp_digit d, *dp;
mp_size uv = USED(v);
int extra = 0, ix;
d = DIGIT(v, uv - 1);
while(d && ((d & 1) == 0)) {
d >>= 1;
++extra;
}
if(d == 1) {
ix = uv - 2;
dp = DIGITS(v) + ix;
while(ix >= 0) {
if(*dp)
return -1;
--dp; --ix;
}
return ((uv - 1) * DIGIT_BIT) + extra;
}
return -1;
}
int s_mp_ispow2d(mp_digit d)
{
int pow = 0;
while((d & 1) == 0) {
++pow; d >>= 1;
}
if(d == 1)
return pow;
return -1;
}
int s_mp_tovalue(char ch, int r)
{
int val, xch;
if(r > 36)
xch = ch;
else
xch = toupper(ch);
if(isdigit(xch))
val = xch - '0';
else if(isupper(xch))
val = xch - 'A' + 10;
else if(islower(xch))
val = xch - 'a' + 36;
else if(xch == '+')
val = 62;
else if(xch == '/')
val = 63;
else
return -1;
if(val < 0 || val >= r)
return -1;
return val;
}
char s_mp_todigit(int val, int r, int low)
{
char ch;
if(val < 0 || val >= r)
return 0;
ch = s_dmap_1[val];
if(r <= 36 && low)
ch = tolower(ch);
return ch;
}
int s_mp_outlen(int bits, int r)
{
return (int)((double)bits * LOG_V_2(r));
}