2006-07-13 06:15:25 +02:00
|
|
|
/* imath version 1.3 */
|
|
|
|
/*
|
|
|
|
Name: imath.c
|
|
|
|
Purpose: Arbitrary precision integer arithmetic routines.
|
2010-04-02 17:21:20 +02:00
|
|
|
Author: M. J. Fromberger <http://spinning-yarns.org/michael/sw/>
|
2006-07-19 19:05:50 +02:00
|
|
|
Info: Id: imath.c 21 2006-04-02 18:58:36Z sting
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
Copyright (C) 2002 Michael J. Fromberger, All Rights Reserved.
|
|
|
|
|
|
|
|
Permission is hereby granted, free of charge, to any person
|
|
|
|
obtaining a copy of this software and associated documentation files
|
|
|
|
(the "Software"), to deal in the Software without restriction,
|
|
|
|
including without limitation the rights to use, copy, modify, merge,
|
|
|
|
publish, distribute, sublicense, and/or sell copies of the Software,
|
|
|
|
and to permit persons to whom the Software is furnished to do so,
|
|
|
|
subject to the following conditions:
|
|
|
|
|
|
|
|
The above copyright notice and this permission notice shall be
|
|
|
|
included in all copies or substantial portions of the Software.
|
|
|
|
|
|
|
|
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
|
|
|
|
EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
|
|
|
|
MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
|
|
|
|
NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS
|
|
|
|
BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN
|
|
|
|
ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
|
|
|
|
CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
|
|
|
|
SOFTWARE.
|
|
|
|
*/
|
2010-09-20 22:08:53 +02:00
|
|
|
/* contrib/pgcrypto/imath.c */
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
#include "postgres.h"
|
|
|
|
#include "px.h"
|
|
|
|
#include "imath.h"
|
|
|
|
|
|
|
|
#undef assert
|
2006-07-15 17:27:14 +02:00
|
|
|
#define assert(TEST) Assert(TEST)
|
2006-07-13 06:15:25 +02:00
|
|
|
#define TRACEABLE_CLAMP 0
|
|
|
|
#define TRACEABLE_FREE 0
|
|
|
|
|
|
|
|
/* {{{ Constants */
|
|
|
|
|
|
|
|
const mp_result MP_OK = 0; /* no error, all is well */
|
|
|
|
const mp_result MP_FALSE = 0; /* boolean false */
|
|
|
|
const mp_result MP_TRUE = -1; /* boolean true */
|
|
|
|
const mp_result MP_MEMORY = -2; /* out of memory */
|
|
|
|
const mp_result MP_RANGE = -3; /* argument out of range */
|
|
|
|
const mp_result MP_UNDEF = -4; /* result undefined */
|
|
|
|
const mp_result MP_TRUNC = -5; /* output truncated */
|
|
|
|
const mp_result MP_BADARG = -6; /* invalid null argument */
|
|
|
|
|
|
|
|
const mp_sign MP_NEG = 1; /* value is strictly negative */
|
|
|
|
const mp_sign MP_ZPOS = 0; /* value is non-negative */
|
|
|
|
|
|
|
|
static const char *s_unknown_err = "unknown result code";
|
|
|
|
static const char *s_error_msg[] = {
|
|
|
|
"error code 0",
|
|
|
|
"boolean true",
|
|
|
|
"out of memory",
|
|
|
|
"argument out of range",
|
|
|
|
"result undefined",
|
|
|
|
"output truncated",
|
|
|
|
"invalid null argument",
|
|
|
|
NULL
|
|
|
|
};
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* Optional library flags */
|
|
|
|
#define MP_CAP_DIGITS 1 /* flag bit to capitalize letter digits */
|
|
|
|
|
|
|
|
/* Argument checking macros
|
|
|
|
Use CHECK() where a return value is required; NRCHECK() elsewhere */
|
|
|
|
#define CHECK(TEST) assert(TEST)
|
|
|
|
#define NRCHECK(TEST) assert(TEST)
|
|
|
|
|
|
|
|
/* {{{ Logarithm table for computing output sizes */
|
|
|
|
|
|
|
|
/* The ith entry of this table gives the value of log_i(2).
|
|
|
|
|
|
|
|
An integer value n requires ceil(log_i(n)) digits to be represented
|
|
|
|
in base i. Since it is easy to compute lg(n), by counting bits, we
|
|
|
|
can compute log_i(n) = lg(n) * log_i(2).
|
|
|
|
*/
|
|
|
|
static const double s_log2[] = {
|
|
|
|
0.000000000, 0.000000000, 1.000000000, 0.630929754, /* 0 1 2 3 */
|
|
|
|
0.500000000, 0.430676558, 0.386852807, 0.356207187, /* 4 5 6 7 */
|
|
|
|
0.333333333, 0.315464877, 0.301029996, 0.289064826, /* 8 9 10 11 */
|
|
|
|
0.278942946, 0.270238154, 0.262649535, 0.255958025, /* 12 13 14 15 */
|
|
|
|
0.250000000, 0.244650542, 0.239812467, 0.235408913, /* 16 17 18 19 */
|
|
|
|
0.231378213, 0.227670249, 0.224243824, 0.221064729, /* 20 21 22 23 */
|
|
|
|
0.218104292, 0.215338279, 0.212746054, 0.210309918, /* 24 25 26 27 */
|
|
|
|
0.208014598, 0.205846832, 0.203795047, 0.201849087, /* 28 29 30 31 */
|
|
|
|
0.200000000, 0.198239863, 0.196561632, 0.194959022, /* 32 33 34 35 */
|
|
|
|
0.193426404, 0.191958720, 0.190551412, 0.189200360, /* 36 37 38 39 */
|
|
|
|
0.187901825, 0.186652411, 0.185449023, 0.184288833, /* 40 41 42 43 */
|
|
|
|
0.183169251, 0.182087900, 0.181042597, 0.180031327, /* 44 45 46 47 */
|
|
|
|
0.179052232, 0.178103594, 0.177183820, 0.176291434, /* 48 49 50 51 */
|
|
|
|
0.175425064, 0.174583430, 0.173765343, 0.172969690, /* 52 53 54 55 */
|
|
|
|
0.172195434, 0.171441601, 0.170707280, 0.169991616, /* 56 57 58 59 */
|
|
|
|
0.169293808, 0.168613099, 0.167948779, 0.167300179, /* 60 61 62 63 */
|
|
|
|
0.166666667
|
|
|
|
};
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
/* {{{ Various macros */
|
|
|
|
|
|
|
|
/* Return the number of digits needed to represent a static value */
|
|
|
|
#define MP_VALUE_DIGITS(V) \
|
|
|
|
((sizeof(V)+(sizeof(mp_digit)-1))/sizeof(mp_digit))
|
|
|
|
|
|
|
|
/* Round precision P to nearest word boundary */
|
|
|
|
#define ROUND_PREC(P) ((mp_size)(2*(((P)+1)/2)))
|
|
|
|
|
|
|
|
/* Set array P of S digits to zero */
|
|
|
|
#define ZERO(P, S) \
|
|
|
|
do{mp_size i__=(S)*sizeof(mp_digit);mp_digit *p__=(P);memset(p__,0,i__);}while(0)
|
|
|
|
|
|
|
|
/* Copy S digits from array P to array Q */
|
|
|
|
#define COPY(P, Q, S) \
|
|
|
|
do{mp_size i__=(S)*sizeof(mp_digit);mp_digit *p__=(P),*q__=(Q);\
|
|
|
|
memcpy(q__,p__,i__);}while(0)
|
|
|
|
|
|
|
|
/* Reverse N elements of type T in array A */
|
|
|
|
#define REV(T, A, N) \
|
|
|
|
do{T *u_=(A),*v_=u_+(N)-1;while(u_<v_){T xch=*u_;*u_++=*v_;*v_--=xch;}}while(0)
|
|
|
|
|
|
|
|
#if TRACEABLE_CLAMP
|
|
|
|
#define CLAMP(Z) s_clamp(Z)
|
|
|
|
#else
|
|
|
|
#define CLAMP(Z) \
|
2021-12-13 01:12:00 +01:00
|
|
|
do{ \
|
|
|
|
mp_int z_=(Z); \
|
|
|
|
mp_size uz_=MP_USED(z_); \
|
|
|
|
mp_digit *dz_=MP_DIGITS(z_)+uz_-1; \
|
|
|
|
while(uz_ > 1 && (*dz_-- == 0)) --uz_; \
|
|
|
|
MP_USED(z_)=uz_; \
|
|
|
|
}while(0)
|
2006-07-13 06:15:25 +02:00
|
|
|
#endif
|
|
|
|
|
|
|
|
#undef MIN
|
|
|
|
#undef MAX
|
|
|
|
#define MIN(A, B) ((B)<(A)?(B):(A))
|
|
|
|
#define MAX(A, B) ((B)>(A)?(B):(A))
|
|
|
|
#define SWAP(T, A, B) do{T t_=(A);A=(B);B=t_;}while(0)
|
|
|
|
|
|
|
|
#define TEMP(K) (temp + (K))
|
|
|
|
#define SETUP(E, C) \
|
|
|
|
do{if((res = (E)) != MP_OK) goto CLEANUP; ++(C);}while(0)
|
|
|
|
|
|
|
|
#define CMPZ(Z) \
|
|
|
|
(((Z)->used==1&&(Z)->digits[0]==0)?0:((Z)->sign==MP_NEG)?-1:1)
|
|
|
|
|
|
|
|
#define UMUL(X, Y, Z) \
|
|
|
|
do{mp_size ua_=MP_USED(X),ub_=MP_USED(Y);mp_size o_=ua_+ub_;\
|
|
|
|
ZERO(MP_DIGITS(Z),o_);\
|
|
|
|
(void) s_kmul(MP_DIGITS(X),MP_DIGITS(Y),MP_DIGITS(Z),ua_,ub_);\
|
|
|
|
MP_USED(Z)=o_;CLAMP(Z);}while(0)
|
|
|
|
|
|
|
|
#define USQR(X, Z) \
|
|
|
|
do{mp_size ua_=MP_USED(X),o_=ua_+ua_;ZERO(MP_DIGITS(Z),o_);\
|
|
|
|
(void) s_ksqr(MP_DIGITS(X),MP_DIGITS(Z),ua_);MP_USED(Z)=o_;CLAMP(Z);}while(0)
|
|
|
|
|
|
|
|
#define UPPER_HALF(W) ((mp_word)((W) >> MP_DIGIT_BIT))
|
|
|
|
#define LOWER_HALF(W) ((mp_digit)(W))
|
|
|
|
#define HIGH_BIT_SET(W) ((W) >> (MP_WORD_BIT - 1))
|
|
|
|
#define ADD_WILL_OVERFLOW(W, V) ((MP_WORD_MAX - (V)) < (W))
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* Default number of digits allocated to a new mp_int */
|
|
|
|
static mp_size default_precision = 64;
|
|
|
|
|
|
|
|
/* Minimum number of digits to invoke recursive multiply */
|
|
|
|
static mp_size multiply_threshold = 32;
|
|
|
|
|
|
|
|
/* Default library configuration flags */
|
|
|
|
static mp_word mp_flags = MP_CAP_DIGITS;
|
|
|
|
|
|
|
|
/* Allocate a buffer of (at least) num digits, or return
|
|
|
|
NULL if that couldn't be done. */
|
|
|
|
static mp_digit *s_alloc(mp_size num);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
#if TRACEABLE_FREE
|
|
|
|
static void s_free(void *ptr);
|
|
|
|
#else
|
|
|
|
#define s_free(P) px_free(P)
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* Insure that z has at least min digits allocated, resizing if
|
|
|
|
necessary. Returns true if successful, false if out of memory. */
|
|
|
|
static int s_pad(mp_int z, mp_size min);
|
|
|
|
|
|
|
|
/* Normalize by removing leading zeroes (except when z = 0) */
|
|
|
|
#if TRACEABLE_CLAMP
|
|
|
|
static void s_clamp(mp_int z);
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* Fill in a "fake" mp_int on the stack with a given value */
|
|
|
|
static void s_fake(mp_int z, int value, mp_digit vbuf[]);
|
|
|
|
|
|
|
|
/* Compare two runs of digits of given length, returns <0, 0, >0 */
|
|
|
|
static int s_cdig(mp_digit *da, mp_digit *db, mp_size len);
|
|
|
|
|
|
|
|
/* Pack the unsigned digits of v into array t */
|
|
|
|
static int s_vpack(int v, mp_digit t[]);
|
|
|
|
|
|
|
|
/* Compare magnitudes of a and b, returns <0, 0, >0 */
|
|
|
|
static int s_ucmp(mp_int a, mp_int b);
|
|
|
|
|
|
|
|
/* Compare magnitudes of a and v, returns <0, 0, >0 */
|
|
|
|
static int s_vcmp(mp_int a, int v);
|
|
|
|
|
|
|
|
/* Unsigned magnitude addition; assumes dc is big enough.
|
|
|
|
Carry out is returned (no memory allocated). */
|
|
|
|
static mp_digit s_uadd(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b);
|
|
|
|
|
|
|
|
/* Unsigned magnitude subtraction. Assumes dc is big enough. */
|
|
|
|
static void s_usub(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b);
|
|
|
|
|
|
|
|
/* Unsigned recursive multiplication. Assumes dc is big enough. */
|
|
|
|
static int s_kmul(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b);
|
|
|
|
|
|
|
|
/* Unsigned magnitude multiplication. Assumes dc is big enough. */
|
|
|
|
static void s_umul(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b);
|
|
|
|
|
|
|
|
/* Unsigned recursive squaring. Assumes dc is big enough. */
|
|
|
|
static int s_ksqr(mp_digit *da, mp_digit *dc, mp_size size_a);
|
|
|
|
|
|
|
|
/* Unsigned magnitude squaring. Assumes dc is big enough. */
|
|
|
|
static void s_usqr(mp_digit *da, mp_digit *dc, mp_size size_a);
|
|
|
|
|
|
|
|
/* Single digit addition. Assumes a is big enough. */
|
|
|
|
static void s_dadd(mp_int a, mp_digit b);
|
|
|
|
|
|
|
|
/* Single digit multiplication. Assumes a is big enough. */
|
|
|
|
static void s_dmul(mp_int a, mp_digit b);
|
|
|
|
|
|
|
|
/* Single digit multiplication on buffers; assumes dc is big enough. */
|
|
|
|
static void s_dbmul(mp_digit *da, mp_digit b, mp_digit *dc,
|
|
|
|
mp_size size_a);
|
|
|
|
|
|
|
|
/* Single digit division. Replaces a with the quotient,
|
|
|
|
returns the remainder. */
|
|
|
|
static mp_digit s_ddiv(mp_int a, mp_digit b);
|
|
|
|
|
|
|
|
/* Quick division by a power of 2, replaces z (no allocation) */
|
|
|
|
static void s_qdiv(mp_int z, mp_size p2);
|
|
|
|
|
|
|
|
/* Quick remainder by a power of 2, replaces z (no allocation) */
|
|
|
|
static void s_qmod(mp_int z, mp_size p2);
|
|
|
|
|
|
|
|
/* Quick multiplication by a power of 2, replaces z.
|
|
|
|
Allocates if necessary; returns false in case this fails. */
|
|
|
|
static int s_qmul(mp_int z, mp_size p2);
|
|
|
|
|
|
|
|
/* Quick subtraction from a power of 2, replaces z.
|
|
|
|
Allocates if necessary; returns false in case this fails. */
|
|
|
|
static int s_qsub(mp_int z, mp_size p2);
|
|
|
|
|
|
|
|
/* Return maximum k such that 2^k divides z. */
|
|
|
|
static int s_dp2k(mp_int z);
|
|
|
|
|
|
|
|
/* Return k >= 0 such that z = 2^k, or -1 if there is no such k. */
|
|
|
|
static int s_isp2(mp_int z);
|
|
|
|
|
|
|
|
/* Set z to 2^k. May allocate; returns false in case this fails. */
|
|
|
|
static int s_2expt(mp_int z, int k);
|
|
|
|
|
|
|
|
/* Normalize a and b for division, returns normalization constant */
|
|
|
|
static int s_norm(mp_int a, mp_int b);
|
|
|
|
|
|
|
|
/* Compute constant mu for Barrett reduction, given modulus m, result
|
|
|
|
replaces z, m is untouched. */
|
|
|
|
static mp_result s_brmu(mp_int z, mp_int m);
|
|
|
|
|
|
|
|
/* Reduce a modulo m, using Barrett's algorithm. */
|
|
|
|
static int s_reduce(mp_int x, mp_int m, mp_int mu, mp_int q1, mp_int q2);
|
|
|
|
|
|
|
|
/* Modular exponentiation, using Barrett reduction */
|
|
|
|
static mp_result s_embar(mp_int a, mp_int b, mp_int m, mp_int mu, mp_int c);
|
|
|
|
|
|
|
|
/* Unsigned magnitude division. Assumes |a| > |b|. Allocates
|
|
|
|
temporaries; overwrites a with quotient, b with remainder. */
|
|
|
|
static mp_result s_udiv(mp_int a, mp_int b);
|
|
|
|
|
|
|
|
/* Compute the number of digits in radix r required to represent the
|
|
|
|
given value. Does not account for sign flags, terminators, etc. */
|
|
|
|
static int s_outlen(mp_int z, mp_size r);
|
|
|
|
|
|
|
|
/* Guess how many digits of precision will be needed to represent a
|
|
|
|
radix r value of the specified number of digits. Returns a value
|
|
|
|
guaranteed to be no smaller than the actual number required. */
|
|
|
|
static mp_size s_inlen(int len, mp_size r);
|
|
|
|
|
|
|
|
/* Convert a character to a digit value in radix r, or
|
|
|
|
-1 if out of range */
|
|
|
|
static int s_ch2val(char c, int r);
|
|
|
|
|
|
|
|
/* Convert a digit value to a character */
|
|
|
|
static char s_val2ch(int v, int caps);
|
|
|
|
|
|
|
|
/* Take 2's complement of a buffer in place */
|
|
|
|
static void s_2comp(unsigned char *buf, int len);
|
|
|
|
|
|
|
|
/* Convert a value to binary, ignoring sign. On input, *limpos is the
|
|
|
|
bound on how many bytes should be written to buf; on output, *limpos
|
|
|
|
is set to the number of bytes actually written. */
|
|
|
|
static mp_result s_tobin(mp_int z, unsigned char *buf, int *limpos, int pad);
|
|
|
|
|
|
|
|
#if 0
|
|
|
|
/* Dump a representation of the mp_int to standard output */
|
|
|
|
void s_print(char *tag, mp_int z);
|
|
|
|
void s_print_buf(char *tag, mp_digit *buf, mp_size num);
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* {{{ get_default_precision() */
|
|
|
|
|
|
|
|
mp_size
|
|
|
|
mp_get_default_precision(void)
|
|
|
|
{
|
|
|
|
return default_precision;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_set_default_precision(s) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_set_default_precision(mp_size s)
|
|
|
|
{
|
|
|
|
NRCHECK(s > 0);
|
|
|
|
|
|
|
|
default_precision = (mp_size) ROUND_PREC(s);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_get_multiply_threshold() */
|
|
|
|
|
|
|
|
mp_size
|
|
|
|
mp_get_multiply_threshold(void)
|
|
|
|
{
|
|
|
|
return multiply_threshold;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_set_multiply_threshold(s) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_set_multiply_threshold(mp_size s)
|
|
|
|
{
|
|
|
|
multiply_threshold = s;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_init(z) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_init(mp_int z)
|
|
|
|
{
|
|
|
|
return mp_int_init_size(z, default_precision);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_alloc() */
|
|
|
|
|
|
|
|
mp_int
|
|
|
|
mp_int_alloc(void)
|
|
|
|
{
|
|
|
|
mp_int out = px_alloc(sizeof(mpz_t));
|
|
|
|
|
|
|
|
assert(out != NULL);
|
|
|
|
out->digits = NULL;
|
|
|
|
out->used = 0;
|
|
|
|
out->alloc = 0;
|
|
|
|
out->sign = 0;
|
|
|
|
|
|
|
|
return out;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_init_size(z, prec) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_init_size(mp_int z, mp_size prec)
|
|
|
|
{
|
|
|
|
CHECK(z != NULL);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
prec = (mp_size) ROUND_PREC(prec);
|
|
|
|
prec = MAX(prec, default_precision);
|
|
|
|
|
|
|
|
if ((MP_DIGITS(z) = s_alloc(prec)) == NULL)
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
z->digits[0] = 0;
|
|
|
|
MP_USED(z) = 1;
|
|
|
|
MP_ALLOC(z) = prec;
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_init_copy(z, old) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_init_copy(mp_int z, mp_int old)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
mp_size uold,
|
|
|
|
target;
|
|
|
|
|
|
|
|
CHECK(z != NULL && old != NULL);
|
|
|
|
|
|
|
|
uold = MP_USED(old);
|
|
|
|
target = MAX(uold, default_precision);
|
|
|
|
|
|
|
|
if ((res = mp_int_init_size(z, target)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
MP_USED(z) = uold;
|
|
|
|
MP_SIGN(z) = MP_SIGN(old);
|
|
|
|
COPY(MP_DIGITS(old), MP_DIGITS(z), uold);
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_init_value(z, value) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_init_value(mp_int z, int value)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
if ((res = mp_int_init(z)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
return mp_int_set_value(z, value);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_set_value(z, value) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_set_value(mp_int z, int value)
|
|
|
|
{
|
|
|
|
mp_size ndig;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
/* How many digits to copy */
|
|
|
|
ndig = (mp_size) MP_VALUE_DIGITS(value);
|
|
|
|
|
|
|
|
if (!s_pad(z, ndig))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
MP_USED(z) = (mp_size) s_vpack(value, MP_DIGITS(z));
|
|
|
|
MP_SIGN(z) = (value < 0) ? MP_NEG : MP_ZPOS;
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_clear(z) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_int_clear(mp_int z)
|
|
|
|
{
|
|
|
|
if (z == NULL)
|
|
|
|
return;
|
|
|
|
|
|
|
|
if (MP_DIGITS(z) != NULL)
|
|
|
|
{
|
|
|
|
s_free(MP_DIGITS(z));
|
|
|
|
MP_DIGITS(z) = NULL;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_free(z) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_int_free(mp_int z)
|
|
|
|
{
|
|
|
|
NRCHECK(z != NULL);
|
|
|
|
|
|
|
|
if (z->digits != NULL)
|
|
|
|
mp_int_clear(z);
|
|
|
|
|
|
|
|
px_free(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_copy(a, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_copy(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
CHECK(a != NULL && c != NULL);
|
|
|
|
|
|
|
|
if (a != c)
|
|
|
|
{
|
|
|
|
mp_size ua = MP_USED(a);
|
|
|
|
mp_digit *da,
|
|
|
|
*dc;
|
|
|
|
|
|
|
|
if (!s_pad(c, ua))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
da = MP_DIGITS(a);
|
|
|
|
dc = MP_DIGITS(c);
|
|
|
|
COPY(da, dc, ua);
|
|
|
|
|
|
|
|
MP_USED(c) = ua;
|
|
|
|
MP_SIGN(c) = MP_SIGN(a);
|
|
|
|
}
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_swap(a, c) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_int_swap(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
if (a != c)
|
|
|
|
{
|
|
|
|
mpz_t tmp = *a;
|
|
|
|
|
|
|
|
*a = *c;
|
|
|
|
*c = tmp;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_zero(z) */
|
|
|
|
|
|
|
|
void
|
|
|
|
mp_int_zero(mp_int z)
|
|
|
|
{
|
|
|
|
NRCHECK(z != NULL);
|
|
|
|
|
|
|
|
z->digits[0] = 0;
|
|
|
|
MP_USED(z) = 1;
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_abs(a, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_abs(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
|
|
|
|
CHECK(a != NULL && c != NULL);
|
|
|
|
|
|
|
|
if ((res = mp_int_copy(a, c)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
MP_SIGN(c) = MP_ZPOS;
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_neg(a, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_neg(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
|
|
|
|
CHECK(a != NULL && c != NULL);
|
|
|
|
|
|
|
|
if ((res = mp_int_copy(a, c)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
if (CMPZ(c) != 0)
|
|
|
|
MP_SIGN(c) = 1 - MP_SIGN(a);
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_add(a, b, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_add(mp_int a, mp_int b, mp_int c)
|
|
|
|
{
|
|
|
|
mp_size ua,
|
|
|
|
ub,
|
|
|
|
uc,
|
|
|
|
max;
|
|
|
|
|
|
|
|
CHECK(a != NULL && b != NULL && c != NULL);
|
|
|
|
|
|
|
|
ua = MP_USED(a);
|
|
|
|
ub = MP_USED(b);
|
|
|
|
uc = MP_USED(c);
|
|
|
|
max = MAX(ua, ub);
|
|
|
|
|
|
|
|
if (MP_SIGN(a) == MP_SIGN(b))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Same sign -- add magnitudes, preserve sign of addends */
|
|
|
|
mp_digit carry;
|
|
|
|
|
|
|
|
if (!s_pad(c, max))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
carry = s_uadd(MP_DIGITS(a), MP_DIGITS(b), MP_DIGITS(c), ua, ub);
|
|
|
|
uc = max;
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
if (carry)
|
2006-07-13 06:15:25 +02:00
|
|
|
{
|
|
|
|
if (!s_pad(c, max + 1))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
c->digits[max] = carry;
|
|
|
|
++uc;
|
|
|
|
}
|
|
|
|
|
|
|
|
MP_USED(c) = uc;
|
|
|
|
MP_SIGN(c) = MP_SIGN(a);
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Different signs -- subtract magnitudes, preserve sign of greater */
|
|
|
|
mp_int x,
|
2006-10-04 02:30:14 +02:00
|
|
|
y;
|
2012-04-24 04:43:09 +02:00
|
|
|
int cmp = s_ucmp(a, b); /* magnitude comparison, sign ignored */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Set x to max(a, b), y to min(a, b) to simplify later code */
|
|
|
|
if (cmp >= 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
|
|
|
x = a;
|
|
|
|
y = b;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
x = b;
|
|
|
|
y = a;
|
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (!s_pad(c, MP_USED(x)))
|
|
|
|
return MP_MEMORY;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Subtract smaller from larger */
|
|
|
|
s_usub(MP_DIGITS(x), MP_DIGITS(y), MP_DIGITS(c), MP_USED(x), MP_USED(y));
|
|
|
|
MP_USED(c) = MP_USED(x);
|
|
|
|
CLAMP(c);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Give result the sign of the larger */
|
|
|
|
MP_SIGN(c) = MP_SIGN(x);
|
|
|
|
}
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_add_value(a, value, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_add_value(mp_int a, int value, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t vtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
return mp_int_add(a, &vtmp, c);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_sub(a, b, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_sub(mp_int a, mp_int b, mp_int c)
|
|
|
|
{
|
|
|
|
mp_size ua,
|
|
|
|
ub,
|
|
|
|
uc,
|
|
|
|
max;
|
|
|
|
|
|
|
|
CHECK(a != NULL && b != NULL && c != NULL);
|
|
|
|
|
|
|
|
ua = MP_USED(a);
|
|
|
|
ub = MP_USED(b);
|
|
|
|
uc = MP_USED(c);
|
|
|
|
max = MAX(ua, ub);
|
|
|
|
|
|
|
|
if (MP_SIGN(a) != MP_SIGN(b))
|
|
|
|
{
|
|
|
|
/* Different signs -- add magnitudes and keep sign of a */
|
|
|
|
mp_digit carry;
|
|
|
|
|
|
|
|
if (!s_pad(c, max))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
carry = s_uadd(MP_DIGITS(a), MP_DIGITS(b), MP_DIGITS(c), ua, ub);
|
|
|
|
uc = max;
|
|
|
|
|
|
|
|
if (carry)
|
|
|
|
{
|
|
|
|
if (!s_pad(c, max + 1))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
c->digits[max] = carry;
|
|
|
|
++uc;
|
|
|
|
}
|
|
|
|
|
|
|
|
MP_USED(c) = uc;
|
|
|
|
MP_SIGN(c) = MP_SIGN(a);
|
|
|
|
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
/* Same signs -- subtract magnitudes */
|
|
|
|
mp_int x,
|
|
|
|
y;
|
|
|
|
mp_sign osign;
|
|
|
|
int cmp = s_ucmp(a, b);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (!s_pad(c, max))
|
|
|
|
return MP_MEMORY;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (cmp >= 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
|
|
|
x = a;
|
|
|
|
y = b;
|
2006-07-13 06:15:25 +02:00
|
|
|
osign = MP_ZPOS;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
x = b;
|
|
|
|
y = a;
|
2006-07-13 06:15:25 +02:00
|
|
|
osign = MP_NEG;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (MP_SIGN(a) == MP_NEG && cmp != 0)
|
|
|
|
osign = 1 - osign;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_usub(MP_DIGITS(x), MP_DIGITS(y), MP_DIGITS(c), MP_USED(x), MP_USED(y));
|
|
|
|
MP_USED(c) = MP_USED(x);
|
|
|
|
CLAMP(c);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
MP_SIGN(c) = osign;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
return MP_OK;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_sub_value(a, value, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_sub_value(mp_int a, int value, mp_int c)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
mpz_t vtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
return mp_int_sub(a, &vtmp, c);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_mul(a, b, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_mul(mp_int a, mp_int b, mp_int c)
|
|
|
|
{
|
|
|
|
mp_digit *out;
|
|
|
|
mp_size osize,
|
2006-10-04 02:30:14 +02:00
|
|
|
ua,
|
|
|
|
ub,
|
|
|
|
p = 0;
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_sign osign;
|
|
|
|
|
|
|
|
CHECK(a != NULL && b != NULL && c != NULL);
|
|
|
|
|
|
|
|
/* If either input is zero, we can shortcut multiplication */
|
|
|
|
if (mp_int_compare_zero(a) == 0 || mp_int_compare_zero(b) == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_zero(c);
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Output is positive if inputs have same sign, otherwise negative */
|
|
|
|
osign = (MP_SIGN(a) == MP_SIGN(b)) ? MP_ZPOS : MP_NEG;
|
|
|
|
|
|
|
|
/*
|
|
|
|
* If the output is not equal to any of the inputs, we'll write the
|
|
|
|
* results there directly; otherwise, allocate a temporary space.
|
|
|
|
*/
|
|
|
|
ua = MP_USED(a);
|
|
|
|
ub = MP_USED(b);
|
2015-02-02 16:00:45 +01:00
|
|
|
osize = MAX(ua, ub);
|
|
|
|
osize = 4 * ((osize + 1) / 2);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (c == a || c == b)
|
|
|
|
{
|
|
|
|
p = ROUND_PREC(osize);
|
|
|
|
p = MAX(p, default_precision);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((out = s_alloc(p)) == NULL)
|
|
|
|
return MP_MEMORY;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
if (!s_pad(c, osize))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
out = MP_DIGITS(c);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
ZERO(out, osize);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (!s_kmul(MP_DIGITS(a), MP_DIGITS(b), out, ua, ub))
|
|
|
|
return MP_MEMORY;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/*
|
|
|
|
* If we allocated a new buffer, get rid of whatever memory c was already
|
|
|
|
* using, and fix up its fields to reflect that.
|
|
|
|
*/
|
|
|
|
if (out != MP_DIGITS(c))
|
|
|
|
{
|
|
|
|
s_free(MP_DIGITS(c));
|
|
|
|
MP_DIGITS(c) = out;
|
|
|
|
MP_ALLOC(c) = p;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
MP_USED(c) = osize; /* might not be true, but we'll fix it ... */
|
|
|
|
CLAMP(c); /* ... right here */
|
|
|
|
MP_SIGN(c) = osign;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_mul_value(a, value, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_mul_value(mp_int a, int value, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t vtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
return mp_int_mul(a, &vtmp, c);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_mul_pow2(a, p2, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_mul_pow2(mp_int a, int p2, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
|
|
|
|
CHECK(a != NULL && c != NULL && p2 >= 0);
|
|
|
|
|
|
|
|
if ((res = mp_int_copy(a, c)) != MP_OK)
|
|
|
|
return res;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (s_qmul(c, (mp_size) p2))
|
|
|
|
return MP_OK;
|
|
|
|
else
|
|
|
|
return MP_MEMORY;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_sqr(a, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_sqr(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
mp_digit *out;
|
|
|
|
mp_size osize,
|
|
|
|
p = 0;
|
|
|
|
|
|
|
|
CHECK(a != NULL && c != NULL);
|
|
|
|
|
|
|
|
/* Get a temporary buffer big enough to hold the result */
|
2015-02-02 16:00:45 +01:00
|
|
|
osize = (mp_size) 4 * ((MP_USED(a) + 1) / 2);
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
if (a == c)
|
|
|
|
{
|
|
|
|
p = ROUND_PREC(osize);
|
|
|
|
p = MAX(p, default_precision);
|
|
|
|
|
|
|
|
if ((out = s_alloc(p)) == NULL)
|
|
|
|
return MP_MEMORY;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
else
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
if (!s_pad(c, osize))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
out = MP_DIGITS(c);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
ZERO(out, osize);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_ksqr(MP_DIGITS(a), out, MP_USED(a));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* Get rid of whatever memory c was already using, and fix up its fields
|
|
|
|
* to reflect the new digit array it's using
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
if (out != MP_DIGITS(c))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
s_free(MP_DIGITS(c));
|
|
|
|
MP_DIGITS(c) = out;
|
|
|
|
MP_ALLOC(c) = p;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
MP_USED(c) = osize; /* might not be true, but we'll fix it ... */
|
|
|
|
CLAMP(c); /* ... right here */
|
|
|
|
MP_SIGN(c) = MP_ZPOS;
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_div(a, b, q, r) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_div(mp_int a, mp_int b, mp_int q, mp_int r)
|
|
|
|
{
|
|
|
|
int cmp,
|
|
|
|
last = 0,
|
|
|
|
lg;
|
|
|
|
mp_result res = MP_OK;
|
|
|
|
mpz_t temp[2];
|
|
|
|
mp_int qout,
|
|
|
|
rout;
|
|
|
|
mp_sign sa = MP_SIGN(a),
|
|
|
|
sb = MP_SIGN(b);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(a != NULL && b != NULL && q != r);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (CMPZ(b) == 0)
|
|
|
|
return MP_UNDEF;
|
|
|
|
else if ((cmp = s_ucmp(a, b)) < 0)
|
|
|
|
{
|
|
|
|
/*
|
|
|
|
* If |a| < |b|, no division is required: q = 0, r = a
|
|
|
|
*/
|
|
|
|
if (r && (res = mp_int_copy(a, r)) != MP_OK)
|
|
|
|
return res;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (q)
|
|
|
|
mp_int_zero(q);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
else if (cmp == 0)
|
|
|
|
{
|
|
|
|
/*
|
|
|
|
* If |a| = |b|, no division is required: q = 1 or -1, r = 0
|
|
|
|
*/
|
|
|
|
if (r)
|
|
|
|
mp_int_zero(r);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (q)
|
|
|
|
{
|
|
|
|
mp_int_zero(q);
|
|
|
|
q->digits[0] = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (sa != sb)
|
|
|
|
MP_SIGN(q) = MP_NEG;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return MP_OK;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/*
|
|
|
|
* When |a| > |b|, real division is required. We need someplace to store
|
|
|
|
* quotient and remainder, but q and r are allowed to be NULL or to
|
|
|
|
* overlap with the inputs.
|
|
|
|
*/
|
|
|
|
if ((lg = s_isp2(b)) < 0)
|
|
|
|
{
|
|
|
|
if (q && b != q && (res = mp_int_copy(a, q)) == MP_OK)
|
|
|
|
{
|
|
|
|
qout = q;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
qout = TEMP(last);
|
|
|
|
SETUP(mp_int_init_copy(TEMP(last), a), last);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (r && a != r && (res = mp_int_copy(b, r)) == MP_OK)
|
|
|
|
{
|
|
|
|
rout = r;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
rout = TEMP(last);
|
|
|
|
SETUP(mp_int_init_copy(TEMP(last), b), last);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = s_udiv(qout, rout)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if (q && (res = mp_int_copy(a, q)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (r && (res = mp_int_copy(a, r)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (q)
|
|
|
|
s_qdiv(q, (mp_size) lg);
|
|
|
|
qout = q;
|
|
|
|
if (r)
|
|
|
|
s_qmod(r, (mp_size) lg);
|
|
|
|
rout = r;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Recompute signs for output */
|
|
|
|
if (rout)
|
|
|
|
{
|
|
|
|
MP_SIGN(rout) = sa;
|
|
|
|
if (CMPZ(rout) == 0)
|
|
|
|
MP_SIGN(rout) = MP_ZPOS;
|
|
|
|
}
|
|
|
|
if (qout)
|
|
|
|
{
|
|
|
|
MP_SIGN(qout) = (sa == sb) ? MP_ZPOS : MP_NEG;
|
|
|
|
if (CMPZ(qout) == 0)
|
|
|
|
MP_SIGN(qout) = MP_ZPOS;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (q && (res = mp_int_copy(qout, q)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (r && (res = mp_int_copy(rout, r)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_mod(a, m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_mod(mp_int a, mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
mpz_t tmp;
|
|
|
|
mp_int out;
|
|
|
|
|
|
|
|
if (m == c)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_init(&tmp)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
out = &tmp;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
out = c;
|
|
|
|
}
|
|
|
|
|
|
|
|
if ((res = mp_int_div(a, m, NULL, out)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if (CMPZ(out) < 0)
|
|
|
|
res = mp_int_add(out, m, c);
|
|
|
|
else
|
|
|
|
res = mp_int_copy(out, c);
|
|
|
|
|
|
|
|
CLEANUP:
|
|
|
|
if (out != c)
|
|
|
|
mp_int_clear(&tmp);
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
|
|
|
|
/* {{{ mp_int_div_value(a, value, q, r) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_div_value(mp_int a, int value, mp_int q, int *r)
|
|
|
|
{
|
|
|
|
mpz_t vtmp,
|
|
|
|
rtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
mp_result res;
|
|
|
|
|
|
|
|
if ((res = mp_int_init(&rtmp)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
if ((res = mp_int_div(a, &vtmp, q, &rtmp)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if (r)
|
|
|
|
(void) mp_int_to_int(&rtmp, r); /* can't fail */
|
|
|
|
|
|
|
|
CLEANUP:
|
|
|
|
mp_int_clear(&rtmp);
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_div_pow2(a, p2, q, r) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_div_pow2(mp_int a, int p2, mp_int q, mp_int r)
|
|
|
|
{
|
|
|
|
mp_result res = MP_OK;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(a != NULL && p2 >= 0 && q != r);
|
|
|
|
|
|
|
|
if (q != NULL && (res = mp_int_copy(a, q)) == MP_OK)
|
|
|
|
s_qdiv(q, (mp_size) p2);
|
|
|
|
|
|
|
|
if (res == MP_OK && r != NULL && (res = mp_int_copy(a, r)) == MP_OK)
|
|
|
|
s_qmod(r, (mp_size) p2);
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_expt(a, b, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_expt(mp_int a, int b, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t t;
|
|
|
|
mp_result res;
|
|
|
|
unsigned int v = abs(b);
|
|
|
|
|
|
|
|
CHECK(b >= 0 && c != NULL);
|
|
|
|
|
|
|
|
if ((res = mp_int_init_copy(&t, a)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
(void) mp_int_set_value(c, 1);
|
|
|
|
while (v != 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
|
|
|
if (v & 1)
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_mul(c, &t, c)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
v >>= 1;
|
2006-07-13 06:15:25 +02:00
|
|
|
if (v == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
break;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_sqr(&t, &t)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
mp_int_clear(&t);
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_expt_value(a, b, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_expt_value(int a, int b, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t t;
|
|
|
|
mp_result res;
|
|
|
|
unsigned int v = abs(b);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(b >= 0 && c != NULL);
|
|
|
|
|
|
|
|
if ((res = mp_int_init_value(&t, a)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
(void) mp_int_set_value(c, 1);
|
|
|
|
while (v != 0)
|
|
|
|
{
|
|
|
|
if (v & 1)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_mul(c, &t, c)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
|
|
|
|
v >>= 1;
|
|
|
|
if (v == 0)
|
|
|
|
break;
|
|
|
|
|
|
|
|
if ((res = mp_int_sqr(&t, &t)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
mp_int_clear(&t);
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_compare(a, b) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_compare(mp_int a, mp_int b)
|
|
|
|
{
|
|
|
|
mp_sign sa;
|
|
|
|
|
|
|
|
CHECK(a != NULL && b != NULL);
|
|
|
|
|
|
|
|
sa = MP_SIGN(a);
|
|
|
|
if (sa == MP_SIGN(b))
|
|
|
|
{
|
|
|
|
int cmp = s_ucmp(a, b);
|
|
|
|
|
|
|
|
/*
|
|
|
|
* If they're both zero or positive, the normal comparison applies; if
|
|
|
|
* both negative, the sense is reversed.
|
|
|
|
*/
|
Allow btree comparison functions to return INT_MIN.
Historically we forbade datatype-specific comparison functions from
returning INT_MIN, so that it would be safe to invert the sort order
just by negating the comparison result. However, this was never
really safe for comparison functions that directly return the result
of memcmp(), strcmp(), etc, as POSIX doesn't place any such restriction
on those library functions. Buildfarm results show that at least on
recent Linux on s390x, memcmp() actually does return INT_MIN sometimes,
causing sort failures.
The agreed-on answer is to remove this restriction and fix relevant
call sites to not make such an assumption; code such as "res = -res"
should be replaced by "INVERT_COMPARE_RESULT(res)". The same is needed
in a few places that just directly negated the result of memcmp or
strcmp.
To help find places having this problem, I've also added a compile option
to nbtcompare.c that causes some of the commonly used comparators to
return INT_MIN/INT_MAX instead of their usual -1/+1. It'd likely be
a good idea to have at least one buildfarm member running with
"-DSTRESS_SORT_INT_MIN". That's far from a complete test of course,
but it should help to prevent fresh introductions of such bugs.
This is a longstanding portability hazard, so back-patch to all supported
branches.
Discussion: https://postgr.es/m/20180928185215.ffoq2xrq5d3pafna@alap3.anarazel.de
2018-10-05 22:01:29 +02:00
|
|
|
if (sa != MP_ZPOS)
|
|
|
|
INVERT_COMPARE_RESULT(cmp);
|
|
|
|
return cmp;
|
2006-07-13 06:15:25 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if (sa == MP_ZPOS)
|
|
|
|
return 1;
|
|
|
|
else
|
|
|
|
return -1;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_compare_unsigned(a, b) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_compare_unsigned(mp_int a, mp_int b)
|
|
|
|
{
|
|
|
|
NRCHECK(a != NULL && b != NULL);
|
|
|
|
|
|
|
|
return s_ucmp(a, b);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_compare_zero(z) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_compare_zero(mp_int z)
|
|
|
|
{
|
|
|
|
NRCHECK(z != NULL);
|
|
|
|
|
|
|
|
if (MP_USED(z) == 1 && z->digits[0] == 0)
|
|
|
|
return 0;
|
|
|
|
else if (MP_SIGN(z) == MP_ZPOS)
|
|
|
|
return 1;
|
|
|
|
else
|
|
|
|
return -1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_compare_value(z, value) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_compare_value(mp_int z, int value)
|
|
|
|
{
|
|
|
|
mp_sign vsign = (value < 0) ? MP_NEG : MP_ZPOS;
|
|
|
|
int cmp;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
if (vsign == MP_SIGN(z))
|
|
|
|
{
|
|
|
|
cmp = s_vcmp(z, value);
|
|
|
|
|
Allow btree comparison functions to return INT_MIN.
Historically we forbade datatype-specific comparison functions from
returning INT_MIN, so that it would be safe to invert the sort order
just by negating the comparison result. However, this was never
really safe for comparison functions that directly return the result
of memcmp(), strcmp(), etc, as POSIX doesn't place any such restriction
on those library functions. Buildfarm results show that at least on
recent Linux on s390x, memcmp() actually does return INT_MIN sometimes,
causing sort failures.
The agreed-on answer is to remove this restriction and fix relevant
call sites to not make such an assumption; code such as "res = -res"
should be replaced by "INVERT_COMPARE_RESULT(res)". The same is needed
in a few places that just directly negated the result of memcmp or
strcmp.
To help find places having this problem, I've also added a compile option
to nbtcompare.c that causes some of the commonly used comparators to
return INT_MIN/INT_MAX instead of their usual -1/+1. It'd likely be
a good idea to have at least one buildfarm member running with
"-DSTRESS_SORT_INT_MIN". That's far from a complete test of course,
but it should help to prevent fresh introductions of such bugs.
This is a longstanding portability hazard, so back-patch to all supported
branches.
Discussion: https://postgr.es/m/20180928185215.ffoq2xrq5d3pafna@alap3.anarazel.de
2018-10-05 22:01:29 +02:00
|
|
|
if (vsign != MP_ZPOS)
|
|
|
|
INVERT_COMPARE_RESULT(cmp);
|
|
|
|
return cmp;
|
2006-07-13 06:15:25 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if (value < 0)
|
|
|
|
return 1;
|
|
|
|
else
|
|
|
|
return -1;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_exptmod(a, b, m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_exptmod(mp_int a, mp_int b, mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
mp_size um;
|
|
|
|
mpz_t temp[3];
|
|
|
|
mp_int s;
|
|
|
|
int last = 0;
|
|
|
|
|
|
|
|
CHECK(a != NULL && b != NULL && c != NULL && m != NULL);
|
|
|
|
|
|
|
|
/* Zero moduli and negative exponents are not considered. */
|
|
|
|
if (CMPZ(m) == 0)
|
|
|
|
return MP_UNDEF;
|
|
|
|
if (CMPZ(b) < 0)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
um = MP_USED(m);
|
|
|
|
SETUP(mp_int_init_size(TEMP(0), 2 * um), last);
|
|
|
|
SETUP(mp_int_init_size(TEMP(1), 2 * um), last);
|
|
|
|
|
|
|
|
if (c == b || c == m)
|
|
|
|
{
|
|
|
|
SETUP(mp_int_init_size(TEMP(2), 2 * um), last);
|
|
|
|
s = TEMP(2);
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
s = c;
|
|
|
|
}
|
|
|
|
|
|
|
|
if ((res = mp_int_mod(a, m, TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if ((res = s_brmu(TEMP(1), m)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if ((res = s_embar(TEMP(0), b, m, TEMP(1), s)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
res = mp_int_copy(s, c);
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
CLEANUP:
|
2006-07-13 06:15:25 +02:00
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_exptmod_evalue(a, value, m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_exptmod_evalue(mp_int a, int value, mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t vtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
return mp_int_exptmod(a, &vtmp, m, c);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_exptmod_bvalue(v, b, m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_exptmod_bvalue(int value, mp_int b,
|
|
|
|
mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
mpz_t vtmp;
|
|
|
|
mp_digit vbuf[MP_VALUE_DIGITS(value)];
|
|
|
|
|
|
|
|
s_fake(&vtmp, value, vbuf);
|
|
|
|
|
|
|
|
return mp_int_exptmod(&vtmp, b, m, c);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_exptmod_known(a, b, m, mu, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_exptmod_known(mp_int a, mp_int b, mp_int m, mp_int mu, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
mp_size um;
|
|
|
|
mpz_t temp[2];
|
|
|
|
mp_int s;
|
|
|
|
int last = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(a && b && m && c);
|
|
|
|
|
|
|
|
/* Zero moduli and negative exponents are not considered. */
|
|
|
|
if (CMPZ(m) == 0)
|
|
|
|
return MP_UNDEF;
|
|
|
|
if (CMPZ(b) < 0)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
um = MP_USED(m);
|
|
|
|
SETUP(mp_int_init_size(TEMP(0), 2 * um), last);
|
|
|
|
|
|
|
|
if (c == b || c == m)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
SETUP(mp_int_init_size(TEMP(1), 2 * um), last);
|
|
|
|
s = TEMP(1);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
s = c;
|
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
if ((res = mp_int_mod(a, m, TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if ((res = s_embar(TEMP(0), b, m, mu, s)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
res = mp_int_copy(s, c);
|
|
|
|
|
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_redux_const(m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_redux_const(mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
CHECK(m != NULL && c != NULL && m != c);
|
|
|
|
|
|
|
|
return s_brmu(c, m);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_invmod(a, m, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_invmod(mp_int a, mp_int m, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
mp_sign sa;
|
|
|
|
int last = 0;
|
|
|
|
mpz_t temp[2];
|
|
|
|
|
|
|
|
CHECK(a != NULL && m != NULL && c != NULL);
|
|
|
|
|
|
|
|
if (CMPZ(a) == 0 || CMPZ(m) <= 0)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
sa = MP_SIGN(a); /* need this for the result later */
|
|
|
|
|
|
|
|
for (last = 0; last < 2; ++last)
|
|
|
|
if ((res = mp_int_init(TEMP(last))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if ((res = mp_int_egcd(a, m, TEMP(0), TEMP(1), NULL)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if (mp_int_compare_value(TEMP(0), 1) != 0)
|
|
|
|
{
|
|
|
|
res = MP_UNDEF;
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* It is first necessary to constrain the value to the proper range */
|
|
|
|
if ((res = mp_int_mod(TEMP(1), m, TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
/*
|
|
|
|
* Now, if 'a' was originally negative, the value we have is actually the
|
|
|
|
* magnitude of the negative representative; to get the positive value we
|
|
|
|
* have to subtract from the modulus. Otherwise, the value is okay as it
|
|
|
|
* stands.
|
|
|
|
*/
|
|
|
|
if (sa == MP_NEG)
|
|
|
|
res = mp_int_sub(m, TEMP(1), c);
|
|
|
|
else
|
|
|
|
res = mp_int_copy(TEMP(1), c);
|
|
|
|
|
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_gcd(a, b, c) */
|
|
|
|
|
|
|
|
/* Binary GCD algorithm due to Josef Stein, 1961 */
|
|
|
|
mp_result
|
|
|
|
mp_int_gcd(mp_int a, mp_int b, mp_int c)
|
|
|
|
{
|
|
|
|
int ca,
|
|
|
|
cb,
|
|
|
|
k = 0;
|
|
|
|
mpz_t u,
|
|
|
|
v,
|
|
|
|
t;
|
|
|
|
mp_result res;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(a != NULL && b != NULL && c != NULL);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
ca = CMPZ(a);
|
|
|
|
cb = CMPZ(b);
|
|
|
|
if (ca == 0 && cb == 0)
|
|
|
|
return MP_UNDEF;
|
|
|
|
else if (ca == 0)
|
|
|
|
return mp_int_abs(b, c);
|
|
|
|
else if (cb == 0)
|
|
|
|
return mp_int_abs(a, c);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_init(&t)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
if ((res = mp_int_init_copy(&u, a)) != MP_OK)
|
|
|
|
goto U;
|
|
|
|
if ((res = mp_int_init_copy(&v, b)) != MP_OK)
|
|
|
|
goto V;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
MP_SIGN(&u) = MP_ZPOS;
|
|
|
|
MP_SIGN(&v) = MP_ZPOS;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
{ /* Divide out common factors of 2 from u and v */
|
|
|
|
int div2_u = s_dp2k(&u),
|
|
|
|
div2_v = s_dp2k(&v);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
k = MIN(div2_u, div2_v);
|
|
|
|
s_qdiv(&u, (mp_size) k);
|
|
|
|
s_qdiv(&v, (mp_size) k);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (mp_int_is_odd(&u))
|
|
|
|
{
|
|
|
|
if ((res = mp_int_neg(&v, &t)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if ((res = mp_int_copy(&u, &t)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (;;)
|
|
|
|
{
|
|
|
|
s_qdiv(&t, s_dp2k(&t));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (CMPZ(&t) > 0)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_copy(&t, &u)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if ((res = mp_int_neg(&t, &v)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_sub(&u, &v, &t)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (CMPZ(&t) == 0)
|
|
|
|
break;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_abs(&u, c)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (!s_qmul(c, (mp_size) k))
|
|
|
|
res = MP_MEMORY;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
mp_int_clear(&v);
|
|
|
|
V: mp_int_clear(&u);
|
|
|
|
U: mp_int_clear(&t);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_egcd(a, b, c, x, y) */
|
|
|
|
|
|
|
|
/* This is the binary GCD algorithm again, but this time we keep track
|
|
|
|
of the elementary matrix operations as we go, so we can get values
|
|
|
|
x and y satisfying c = ax + by.
|
|
|
|
*/
|
|
|
|
mp_result
|
|
|
|
mp_int_egcd(mp_int a, mp_int b, mp_int c,
|
|
|
|
mp_int x, mp_int y)
|
|
|
|
{
|
|
|
|
int k,
|
|
|
|
last = 0,
|
|
|
|
ca,
|
|
|
|
cb;
|
|
|
|
mpz_t temp[8];
|
|
|
|
mp_result res;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(a != NULL && b != NULL && c != NULL &&
|
|
|
|
(x != NULL || y != NULL));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
ca = CMPZ(a);
|
|
|
|
cb = CMPZ(b);
|
|
|
|
if (ca == 0 && cb == 0)
|
|
|
|
return MP_UNDEF;
|
|
|
|
else if (ca == 0)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_abs(b, c)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
mp_int_zero(x);
|
|
|
|
(void) mp_int_set_value(y, 1);
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
else if (cb == 0)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_abs(a, c)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
(void) mp_int_set_value(x, 1);
|
|
|
|
mp_int_zero(y);
|
|
|
|
return MP_OK;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/*
|
|
|
|
* Initialize temporaries: A:0, B:1, C:2, D:3, u:4, v:5, ou:6, ov:7
|
|
|
|
*/
|
|
|
|
for (last = 0; last < 4; ++last)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_init(TEMP(last))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
TEMP(0)->digits[0] = 1;
|
|
|
|
TEMP(3)->digits[0] = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
SETUP(mp_int_init_copy(TEMP(4), a), last);
|
|
|
|
SETUP(mp_int_init_copy(TEMP(5), b), last);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* We will work with absolute values here */
|
|
|
|
MP_SIGN(TEMP(4)) = MP_ZPOS;
|
|
|
|
MP_SIGN(TEMP(5)) = MP_ZPOS;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
{ /* Divide out common factors of 2 from u and v */
|
|
|
|
int div2_u = s_dp2k(TEMP(4)),
|
|
|
|
div2_v = s_dp2k(TEMP(5));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
k = MIN(div2_u, div2_v);
|
|
|
|
s_qdiv(TEMP(4), k);
|
|
|
|
s_qdiv(TEMP(5), k);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
SETUP(mp_int_init_copy(TEMP(6), TEMP(4)), last);
|
|
|
|
SETUP(mp_int_init_copy(TEMP(7), TEMP(5)), last);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (;;)
|
|
|
|
{
|
|
|
|
while (mp_int_is_even(TEMP(4)))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
s_qdiv(TEMP(4), 1);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (mp_int_is_odd(TEMP(0)) || mp_int_is_odd(TEMP(1)))
|
|
|
|
{
|
|
|
|
if ((res = mp_int_add(TEMP(0), TEMP(7), TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(1), TEMP(6), TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_qdiv(TEMP(0), 1);
|
|
|
|
s_qdiv(TEMP(1), 1);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
while (mp_int_is_even(TEMP(5)))
|
|
|
|
{
|
|
|
|
s_qdiv(TEMP(5), 1);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (mp_int_is_odd(TEMP(2)) || mp_int_is_odd(TEMP(3)))
|
|
|
|
{
|
|
|
|
if ((res = mp_int_add(TEMP(2), TEMP(7), TEMP(2))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(3), TEMP(6), TEMP(3))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_qdiv(TEMP(2), 1);
|
|
|
|
s_qdiv(TEMP(3), 1);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (mp_int_compare(TEMP(4), TEMP(5)) >= 0)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_sub(TEMP(4), TEMP(5), TEMP(4))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(0), TEMP(2), TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(1), TEMP(3), TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
if ((res = mp_int_sub(TEMP(5), TEMP(4), TEMP(5))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(2), TEMP(0), TEMP(2))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_sub(TEMP(3), TEMP(1), TEMP(3))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (CMPZ(TEMP(4)) == 0)
|
|
|
|
{
|
|
|
|
if (x && (res = mp_int_copy(TEMP(2), x)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (y && (res = mp_int_copy(TEMP(3), y)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (c)
|
|
|
|
{
|
|
|
|
if (!s_qmul(TEMP(5), k))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
res = MP_MEMORY;
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
res = mp_int_copy(TEMP(5), c);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
break;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_divisible_value(a, v) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_divisible_value(mp_int a, int v)
|
|
|
|
{
|
|
|
|
int rem = 0;
|
|
|
|
|
|
|
|
if (mp_int_div_value(a, v, NULL, &rem) != MP_OK)
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
return rem == 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_is_pow2(z) */
|
|
|
|
|
|
|
|
int
|
|
|
|
mp_int_is_pow2(mp_int z)
|
|
|
|
{
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
return s_isp2(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_sqrt(a, c) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_sqrt(mp_int a, mp_int c)
|
|
|
|
{
|
|
|
|
mp_result res = MP_OK;
|
|
|
|
mpz_t temp[2];
|
|
|
|
int last = 0;
|
|
|
|
|
|
|
|
CHECK(a != NULL && c != NULL);
|
|
|
|
|
|
|
|
/* The square root of a negative value does not exist in the integers. */
|
|
|
|
if (MP_SIGN(a) == MP_NEG)
|
|
|
|
return MP_UNDEF;
|
|
|
|
|
|
|
|
SETUP(mp_int_init_copy(TEMP(last), a), last);
|
|
|
|
SETUP(mp_int_init(TEMP(last)), last);
|
|
|
|
|
|
|
|
for (;;)
|
|
|
|
{
|
|
|
|
if ((res = mp_int_sqr(TEMP(0), TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if (mp_int_compare_unsigned(a, TEMP(1)) == 0)
|
|
|
|
break;
|
|
|
|
|
|
|
|
if ((res = mp_int_copy(a, TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_div(TEMP(1), TEMP(0), TEMP(1), NULL)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_add(TEMP(0), TEMP(1), TEMP(1))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if ((res = mp_int_div_pow2(TEMP(1), 1, TEMP(1), NULL)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
|
|
|
|
if (mp_int_compare_unsigned(TEMP(0), TEMP(1)) == 0)
|
|
|
|
break;
|
|
|
|
if ((res = mp_int_sub_value(TEMP(0), 1, TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
if (mp_int_compare_unsigned(TEMP(0), TEMP(1)) == 0)
|
|
|
|
break;
|
|
|
|
|
|
|
|
if ((res = mp_int_copy(TEMP(1), TEMP(0))) != MP_OK)
|
|
|
|
goto CLEANUP;
|
|
|
|
}
|
|
|
|
|
|
|
|
res = mp_int_copy(TEMP(0), c);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_to_int(z, out) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_to_int(mp_int z, int *out)
|
|
|
|
{
|
|
|
|
unsigned int uv = 0;
|
|
|
|
mp_size uz;
|
|
|
|
mp_digit *dz;
|
|
|
|
mp_sign sz;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
/* Make sure the value is representable as an int */
|
|
|
|
sz = MP_SIGN(z);
|
|
|
|
if ((sz == MP_ZPOS && mp_int_compare_value(z, INT_MAX) > 0) ||
|
|
|
|
mp_int_compare_value(z, INT_MIN) < 0)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
uz = MP_USED(z);
|
|
|
|
dz = MP_DIGITS(z) + uz - 1;
|
|
|
|
|
|
|
|
while (uz > 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
uv <<= MP_DIGIT_BIT / 2;
|
|
|
|
uv = (uv << (MP_DIGIT_BIT / 2)) | *dz--;
|
2006-10-04 02:30:14 +02:00
|
|
|
--uz;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (out)
|
2006-07-13 06:15:25 +02:00
|
|
|
*out = (sz == MP_NEG) ? -(int) uv : (int) uv;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_to_string(z, radix, str, limit) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_to_string(mp_int z, mp_size radix,
|
|
|
|
char *str, int limit)
|
|
|
|
{
|
|
|
|
mp_result res;
|
|
|
|
int cmp = 0;
|
|
|
|
|
|
|
|
CHECK(z != NULL && str != NULL && limit >= 2);
|
|
|
|
|
|
|
|
if (radix < MP_MIN_RADIX || radix > MP_MAX_RADIX)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
if (CMPZ(z) == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
*str++ = s_val2ch(0, mp_flags & MP_CAP_DIGITS);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
mpz_t tmp;
|
2006-10-04 02:30:14 +02:00
|
|
|
char *h,
|
|
|
|
*t;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((res = mp_int_init_copy(&tmp, z)) != MP_OK)
|
|
|
|
return res;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (MP_SIGN(z) == MP_NEG)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
*str++ = '-';
|
2006-10-04 02:30:14 +02:00
|
|
|
--limit;
|
|
|
|
}
|
|
|
|
h = str;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Generate digits in reverse order until finished or limit reached */
|
|
|
|
for ( /* */ ; limit > 0; --limit)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_digit d;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((cmp = CMPZ(&tmp)) == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
break;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
d = s_ddiv(&tmp, (mp_digit) radix);
|
|
|
|
*str++ = s_val2ch(d, mp_flags & MP_CAP_DIGITS);
|
|
|
|
}
|
|
|
|
t = str - 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Put digits back in correct output order */
|
|
|
|
while (h < t)
|
|
|
|
{
|
|
|
|
char tc = *h;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
*h++ = *t;
|
|
|
|
*t-- = tc;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_clear(&tmp);
|
|
|
|
}
|
|
|
|
|
|
|
|
*str = '\0';
|
|
|
|
if (cmp == 0)
|
|
|
|
return MP_OK;
|
|
|
|
else
|
|
|
|
return MP_TRUNC;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_string_len(z, radix) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_string_len(mp_int z, mp_size radix)
|
|
|
|
{
|
|
|
|
int len;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
if (radix < MP_MIN_RADIX || radix > MP_MAX_RADIX)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
len = s_outlen(z, radix) + 1; /* for terminator */
|
|
|
|
|
|
|
|
/* Allow for sign marker on negatives */
|
|
|
|
if (MP_SIGN(z) == MP_NEG)
|
|
|
|
len += 1;
|
|
|
|
|
|
|
|
return len;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_read_string(z, radix, *str) */
|
|
|
|
|
|
|
|
/* Read zero-terminated string into z */
|
|
|
|
mp_result
|
|
|
|
mp_int_read_string(mp_int z, mp_size radix, const char *str)
|
|
|
|
{
|
|
|
|
return mp_int_read_cstring(z, radix, str, NULL);
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_read_cstring(z, radix, *str, **end) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_read_cstring(mp_int z, mp_size radix, const char *str, char **end)
|
|
|
|
{
|
|
|
|
int ch;
|
|
|
|
|
|
|
|
CHECK(z != NULL && str != NULL);
|
|
|
|
|
|
|
|
if (radix < MP_MIN_RADIX || radix > MP_MAX_RADIX)
|
|
|
|
return MP_RANGE;
|
|
|
|
|
|
|
|
/* Skip leading whitespace */
|
2006-09-22 23:39:58 +02:00
|
|
|
while (isspace((unsigned char) *str))
|
2006-07-13 06:15:25 +02:00
|
|
|
++str;
|
|
|
|
|
|
|
|
/* Handle leading sign tag (+/-, positive default) */
|
|
|
|
switch (*str)
|
|
|
|
{
|
|
|
|
case '-':
|
|
|
|
MP_SIGN(z) = MP_NEG;
|
|
|
|
++str;
|
|
|
|
break;
|
|
|
|
case '+':
|
|
|
|
++str; /* fallthrough */
|
|
|
|
default:
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
break;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Skip leading zeroes */
|
|
|
|
while ((ch = s_ch2val(*str, radix)) == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
++str;
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
/* Make sure there is enough space for the value */
|
|
|
|
if (!s_pad(z, s_inlen(strlen(str), radix)))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
MP_USED(z) = 1;
|
|
|
|
z->digits[0] = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
while (*str != '\0' && ((ch = s_ch2val(*str, radix)) >= 0))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
s_dmul(z, (mp_digit) radix);
|
|
|
|
s_dadd(z, (mp_digit) ch);
|
|
|
|
++str;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
CLAMP(z);
|
|
|
|
|
|
|
|
/* Override sign for zero, even if negative specified. */
|
|
|
|
if (CMPZ(z) == 0)
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
|
|
|
|
if (end != NULL)
|
|
|
|
*end = (char *) str;
|
|
|
|
|
|
|
|
/*
|
|
|
|
* Return a truncation error if the string has unprocessed characters
|
|
|
|
* remaining, so the caller can tell if the whole string was done
|
|
|
|
*/
|
|
|
|
if (*str != '\0')
|
|
|
|
return MP_TRUNC;
|
|
|
|
else
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_count_bits(z) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_count_bits(mp_int z)
|
|
|
|
{
|
|
|
|
mp_size nbits = 0,
|
|
|
|
uz;
|
|
|
|
mp_digit d;
|
|
|
|
|
|
|
|
CHECK(z != NULL);
|
|
|
|
|
|
|
|
uz = MP_USED(z);
|
|
|
|
if (uz == 1 && z->digits[0] == 0)
|
|
|
|
return 1;
|
|
|
|
|
|
|
|
--uz;
|
|
|
|
nbits = uz * MP_DIGIT_BIT;
|
|
|
|
d = z->digits[uz];
|
|
|
|
|
|
|
|
while (d != 0)
|
|
|
|
{
|
|
|
|
d >>= 1;
|
|
|
|
++nbits;
|
|
|
|
}
|
|
|
|
|
|
|
|
return nbits;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_to_binary(z, buf, limit) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_to_binary(mp_int z, unsigned char *buf, int limit)
|
|
|
|
{
|
|
|
|
static const int PAD_FOR_2C = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_result res;
|
|
|
|
int limpos = limit;
|
|
|
|
|
|
|
|
CHECK(z != NULL && buf != NULL);
|
|
|
|
|
|
|
|
res = s_tobin(z, buf, &limpos, PAD_FOR_2C);
|
|
|
|
|
|
|
|
if (MP_SIGN(z) == MP_NEG)
|
|
|
|
s_2comp(buf, limpos);
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_read_binary(z, buf, len) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_read_binary(mp_int z, unsigned char *buf, int len)
|
|
|
|
{
|
|
|
|
mp_size need,
|
|
|
|
i;
|
|
|
|
unsigned char *tmp;
|
|
|
|
mp_digit *dz;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CHECK(z != NULL && buf != NULL && len > 0);
|
|
|
|
|
|
|
|
/* Figure out how many digits are needed to represent this value */
|
|
|
|
need = ((len * CHAR_BIT) + (MP_DIGIT_BIT - 1)) / MP_DIGIT_BIT;
|
|
|
|
if (!s_pad(z, need))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
mp_int_zero(z);
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* If the high-order bit is set, take the 2's complement before reading
|
|
|
|
* the value (it will be restored afterward)
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
if (buf[0] >> (CHAR_BIT - 1))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
MP_SIGN(z) = MP_NEG;
|
|
|
|
s_2comp(buf, len);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
dz = MP_DIGITS(z);
|
|
|
|
for (tmp = buf, i = len; i > 0; --i, ++tmp)
|
|
|
|
{
|
|
|
|
s_qmul(z, (mp_size) CHAR_BIT);
|
|
|
|
*dz |= *tmp;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Restore 2's complement if we took it before */
|
|
|
|
if (MP_SIGN(z) == MP_NEG)
|
|
|
|
s_2comp(buf, len);
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_binary_len(z) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_binary_len(mp_int z)
|
|
|
|
{
|
|
|
|
mp_result res = mp_int_count_bits(z);
|
|
|
|
int bytes = mp_int_unsigned_len(z);
|
|
|
|
|
|
|
|
if (res <= 0)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
bytes = (res + (CHAR_BIT - 1)) / CHAR_BIT;
|
|
|
|
|
|
|
|
/*
|
|
|
|
* If the highest-order bit falls exactly on a byte boundary, we need to
|
|
|
|
* pad with an extra byte so that the sign will be read correctly when
|
|
|
|
* reading it back in.
|
|
|
|
*/
|
|
|
|
if (bytes * CHAR_BIT == res)
|
|
|
|
++bytes;
|
|
|
|
|
|
|
|
return bytes;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_to_unsigned(z, buf, limit) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_to_unsigned(mp_int z, unsigned char *buf, int limit)
|
|
|
|
{
|
|
|
|
static const int NO_PADDING = 0;
|
|
|
|
|
|
|
|
CHECK(z != NULL && buf != NULL);
|
|
|
|
|
|
|
|
return s_tobin(z, buf, &limit, NO_PADDING);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_read_unsigned(z, buf, len) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_read_unsigned(mp_int z, unsigned char *buf, int len)
|
|
|
|
{
|
|
|
|
mp_size need,
|
|
|
|
i;
|
|
|
|
unsigned char *tmp;
|
|
|
|
mp_digit *dz;
|
|
|
|
|
|
|
|
CHECK(z != NULL && buf != NULL && len > 0);
|
|
|
|
|
|
|
|
/* Figure out how many digits are needed to represent this value */
|
|
|
|
need = ((len * CHAR_BIT) + (MP_DIGIT_BIT - 1)) / MP_DIGIT_BIT;
|
|
|
|
if (!s_pad(z, need))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
mp_int_zero(z);
|
|
|
|
|
|
|
|
dz = MP_DIGITS(z);
|
|
|
|
for (tmp = buf, i = len; i > 0; --i, ++tmp)
|
|
|
|
{
|
|
|
|
(void) s_qmul(z, CHAR_BIT);
|
|
|
|
*dz |= *tmp;
|
|
|
|
}
|
|
|
|
|
|
|
|
return MP_OK;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_int_unsigned_len(z) */
|
|
|
|
|
|
|
|
mp_result
|
|
|
|
mp_int_unsigned_len(mp_int z)
|
|
|
|
{
|
|
|
|
mp_result res = mp_int_count_bits(z);
|
|
|
|
int bytes;
|
|
|
|
|
|
|
|
if (res <= 0)
|
|
|
|
return res;
|
|
|
|
|
|
|
|
bytes = (res + (CHAR_BIT - 1)) / CHAR_BIT;
|
|
|
|
|
|
|
|
return bytes;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ mp_error_string(res) */
|
|
|
|
|
|
|
|
const char *
|
|
|
|
mp_error_string(mp_result res)
|
|
|
|
{
|
|
|
|
int ix;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (res > 0)
|
|
|
|
return s_unknown_err;
|
|
|
|
|
|
|
|
res = -res;
|
|
|
|
for (ix = 0; ix < res && s_error_msg[ix] != NULL; ++ix)
|
|
|
|
;
|
|
|
|
|
|
|
|
if (s_error_msg[ix] != NULL)
|
|
|
|
return s_error_msg[ix];
|
|
|
|
else
|
|
|
|
return s_unknown_err;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/*------------------------------------------------------------------------*/
|
|
|
|
/* Private functions for internal use. These make assumptions. */
|
|
|
|
|
|
|
|
/* {{{ s_alloc(num) */
|
|
|
|
|
|
|
|
static mp_digit *
|
|
|
|
s_alloc(mp_size num)
|
|
|
|
{
|
|
|
|
mp_digit *out = px_alloc(num * sizeof(mp_digit));
|
|
|
|
|
|
|
|
assert(out != NULL); /* for debugging */
|
|
|
|
|
|
|
|
return out;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_realloc(old, num) */
|
|
|
|
|
|
|
|
static mp_digit *
|
|
|
|
s_realloc(mp_digit *old, mp_size num)
|
|
|
|
{
|
|
|
|
mp_digit *new = px_realloc(old, num * sizeof(mp_digit));
|
|
|
|
|
|
|
|
assert(new != NULL); /* for debugging */
|
|
|
|
|
|
|
|
return new;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_free(ptr) */
|
|
|
|
|
|
|
|
#if TRACEABLE_FREE
|
|
|
|
static void
|
|
|
|
s_free(void *ptr)
|
|
|
|
{
|
|
|
|
px_free(ptr);
|
|
|
|
}
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_pad(z, min) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_pad(mp_int z, mp_size min)
|
|
|
|
{
|
|
|
|
if (MP_ALLOC(z) < min)
|
|
|
|
{
|
|
|
|
mp_size nsize = ROUND_PREC(min);
|
|
|
|
mp_digit *tmp = s_realloc(MP_DIGITS(z), nsize);
|
|
|
|
|
|
|
|
if (tmp == NULL)
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
MP_DIGITS(z) = tmp;
|
|
|
|
MP_ALLOC(z) = nsize;
|
|
|
|
}
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_clamp(z) */
|
|
|
|
|
|
|
|
#if TRACEABLE_CLAMP
|
|
|
|
static void
|
|
|
|
s_clamp(mp_int z)
|
|
|
|
{
|
|
|
|
mp_size uz = MP_USED(z);
|
|
|
|
mp_digit *zd = MP_DIGITS(z) + uz - 1;
|
|
|
|
|
|
|
|
while (uz > 1 && (*zd-- == 0))
|
|
|
|
--uz;
|
|
|
|
|
|
|
|
MP_USED(z) = uz;
|
|
|
|
}
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_fake(z, value, vbuf) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_fake(mp_int z, int value, mp_digit vbuf[])
|
|
|
|
{
|
|
|
|
mp_size uv = (mp_size) s_vpack(value, vbuf);
|
|
|
|
|
|
|
|
z->used = uv;
|
|
|
|
z->alloc = MP_VALUE_DIGITS(value);
|
|
|
|
z->sign = (value < 0) ? MP_NEG : MP_ZPOS;
|
|
|
|
z->digits = vbuf;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_cdig(da, db, len) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_cdig(mp_digit *da, mp_digit *db, mp_size len)
|
|
|
|
{
|
|
|
|
mp_digit *dat = da + len - 1,
|
|
|
|
*dbt = db + len - 1;
|
|
|
|
|
|
|
|
for ( /* */ ; len != 0; --len, --dat, --dbt)
|
|
|
|
{
|
|
|
|
if (*dat > *dbt)
|
|
|
|
return 1;
|
|
|
|
else if (*dat < *dbt)
|
|
|
|
return -1;
|
|
|
|
}
|
|
|
|
|
|
|
|
return 0;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_vpack(v, t[]) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_vpack(int v, mp_digit t[])
|
|
|
|
{
|
|
|
|
unsigned int uv = (unsigned int) ((v < 0) ? -v : v);
|
|
|
|
int ndig = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (uv == 0)
|
|
|
|
t[ndig++] = 0;
|
|
|
|
else
|
|
|
|
{
|
|
|
|
while (uv != 0)
|
|
|
|
{
|
|
|
|
t[ndig++] = (mp_digit) uv;
|
|
|
|
uv >>= MP_DIGIT_BIT / 2;
|
|
|
|
uv >>= MP_DIGIT_BIT / 2;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
return ndig;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_ucmp(a, b) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_ucmp(mp_int a, mp_int b)
|
|
|
|
{
|
|
|
|
mp_size ua = MP_USED(a),
|
|
|
|
ub = MP_USED(b);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (ua > ub)
|
|
|
|
return 1;
|
|
|
|
else if (ub > ua)
|
|
|
|
return -1;
|
|
|
|
else
|
|
|
|
return s_cdig(MP_DIGITS(a), MP_DIGITS(b), ua);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_vcmp(a, v) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_vcmp(mp_int a, int v)
|
|
|
|
{
|
|
|
|
mp_digit vdig[MP_VALUE_DIGITS(v)];
|
|
|
|
int ndig = 0;
|
|
|
|
mp_size ua = MP_USED(a);
|
|
|
|
|
|
|
|
ndig = s_vpack(v, vdig);
|
|
|
|
|
|
|
|
if (ua > ndig)
|
|
|
|
return 1;
|
|
|
|
else if (ua < ndig)
|
|
|
|
return -1;
|
|
|
|
else
|
|
|
|
return s_cdig(MP_DIGITS(a), vdig, ndig);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_uadd(da, db, dc, size_a, size_b) */
|
|
|
|
|
|
|
|
static mp_digit
|
|
|
|
s_uadd(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b)
|
|
|
|
{
|
|
|
|
mp_size pos;
|
|
|
|
mp_word w = 0;
|
|
|
|
|
|
|
|
/* Insure that da is the longer of the two to simplify later code */
|
|
|
|
if (size_b > size_a)
|
|
|
|
{
|
|
|
|
SWAP(mp_digit *, da, db);
|
|
|
|
SWAP(mp_size, size_a, size_b);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Add corresponding digits until the shorter number runs out */
|
|
|
|
for (pos = 0; pos < size_b; ++pos, ++da, ++db, ++dc)
|
|
|
|
{
|
|
|
|
w = w + (mp_word) *da + (mp_word) *db;
|
|
|
|
*dc = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Propagate carries as far as necessary */
|
|
|
|
for ( /* */ ; pos < size_a; ++pos, ++da, ++dc)
|
|
|
|
{
|
|
|
|
w = w + *da;
|
|
|
|
|
|
|
|
*dc = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Return carry out */
|
|
|
|
return (mp_digit) w;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_usub(da, db, dc, size_a, size_b) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_usub(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b)
|
|
|
|
{
|
|
|
|
mp_size pos;
|
|
|
|
mp_word w = 0;
|
|
|
|
|
|
|
|
/* We assume that |a| >= |b| so this should definitely hold */
|
|
|
|
assert(size_a >= size_b);
|
|
|
|
|
|
|
|
/* Subtract corresponding digits and propagate borrow */
|
|
|
|
for (pos = 0; pos < size_b; ++pos, ++da, ++db, ++dc)
|
|
|
|
{
|
|
|
|
w = ((mp_word) MP_DIGIT_MAX + 1 + /* MP_RADIX */
|
|
|
|
(mp_word) *da) - w - (mp_word) *db;
|
|
|
|
|
|
|
|
*dc = LOWER_HALF(w);
|
|
|
|
w = (UPPER_HALF(w) == 0);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Finish the subtraction for remaining upper digits of da */
|
|
|
|
for ( /* */ ; pos < size_a; ++pos, ++da, ++dc)
|
|
|
|
{
|
|
|
|
w = ((mp_word) MP_DIGIT_MAX + 1 + /* MP_RADIX */
|
|
|
|
(mp_word) *da) - w;
|
|
|
|
|
|
|
|
*dc = LOWER_HALF(w);
|
|
|
|
w = (UPPER_HALF(w) == 0);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* If there is a borrow out at the end, it violates the precondition */
|
|
|
|
assert(w == 0);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_kmul(da, db, dc, size_a, size_b) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_kmul(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b)
|
|
|
|
{
|
|
|
|
mp_size bot_size;
|
|
|
|
|
|
|
|
/* Make sure b is the smaller of the two input values */
|
|
|
|
if (size_b > size_a)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
SWAP(mp_digit *, da, db);
|
|
|
|
SWAP(mp_size, size_a, size_b);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* Insure that the bottom is the larger half in an odd-length split; the
|
|
|
|
* code below relies on this being true.
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
bot_size = (size_a + 1) / 2;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* If the values are big enough to bother with recursion, use the
|
|
|
|
* Karatsuba algorithm to compute the product; otherwise use the normal
|
|
|
|
* multiplication algorithm
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
if (multiply_threshold &&
|
|
|
|
size_a >= multiply_threshold &&
|
|
|
|
size_b > bot_size)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_digit *t1,
|
2006-10-04 02:30:14 +02:00
|
|
|
*t2,
|
|
|
|
*t3,
|
|
|
|
carry;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_digit *a_top = da + bot_size;
|
|
|
|
mp_digit *b_top = db + bot_size;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_size at_size = size_a - bot_size;
|
|
|
|
mp_size bt_size = size_b - bot_size;
|
|
|
|
mp_size buf_size = 2 * bot_size;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* Do a single allocation for all three temporary buffers needed; each
|
|
|
|
* buffer must be big enough to hold the product of two bottom halves,
|
|
|
|
* and one buffer needs space for the completed product; twice the
|
|
|
|
* space is plenty.
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((t1 = s_alloc(4 * buf_size)) == NULL)
|
|
|
|
return 0;
|
|
|
|
t2 = t1 + buf_size;
|
|
|
|
t3 = t2 + buf_size;
|
|
|
|
ZERO(t1, 4 * buf_size);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* t1 and t2 are initially used as temporaries to compute the inner
|
|
|
|
* product (a1 + a0)(b1 + b0) = a1b1 + a1b0 + a0b1 + a0b0
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
carry = s_uadd(da, a_top, t1, bot_size, at_size); /* t1 = a1 + a0 */
|
|
|
|
t1[bot_size] = carry;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
carry = s_uadd(db, b_top, t2, bot_size, bt_size); /* t2 = b1 + b0 */
|
|
|
|
t2[bot_size] = carry;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) s_kmul(t1, t2, t3, bot_size + 1, bot_size + 1); /* t3 = t1 * t2 */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* Now we'll get t1 = a0b0 and t2 = a1b1, and subtract them out so
|
|
|
|
* that we're left with only the pieces we want: t3 = a1b0 + a0b1
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2015-02-02 16:00:45 +01:00
|
|
|
ZERO(t1, buf_size);
|
|
|
|
ZERO(t2, buf_size);
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) s_kmul(da, db, t1, bot_size, bot_size); /* t1 = a0 * b0 */
|
|
|
|
(void) s_kmul(a_top, b_top, t2, at_size, bt_size); /* t2 = a1 * b1 */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Subtract out t1 and t2 to get the inner product */
|
|
|
|
s_usub(t3, t1, t3, buf_size + 2, buf_size);
|
|
|
|
s_usub(t3, t2, t3, buf_size + 2, buf_size);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Assemble the output value */
|
|
|
|
COPY(t1, dc, buf_size);
|
2015-02-02 16:00:45 +01:00
|
|
|
carry = s_uadd(t3, dc + bot_size, dc + bot_size,
|
|
|
|
buf_size + 1, buf_size);
|
|
|
|
assert(carry == 0);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2015-02-02 16:00:45 +01:00
|
|
|
carry = s_uadd(t2, dc + 2 * bot_size, dc + 2 * bot_size,
|
|
|
|
buf_size, buf_size);
|
|
|
|
assert(carry == 0);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_free(t1); /* note t2 and t3 are just internal pointers
|
|
|
|
* to t1 */
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
s_umul(da, db, dc, size_a, size_b);
|
|
|
|
}
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_umul(da, db, dc, size_a, size_b) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_umul(mp_digit *da, mp_digit *db, mp_digit *dc,
|
|
|
|
mp_size size_a, mp_size size_b)
|
|
|
|
{
|
|
|
|
mp_size a,
|
|
|
|
b;
|
|
|
|
mp_word w;
|
|
|
|
|
|
|
|
for (a = 0; a < size_a; ++a, ++dc, ++da)
|
|
|
|
{
|
|
|
|
mp_digit *dct = dc;
|
|
|
|
mp_digit *dbt = db;
|
|
|
|
|
|
|
|
if (*da == 0)
|
|
|
|
continue;
|
|
|
|
|
|
|
|
w = 0;
|
|
|
|
for (b = 0; b < size_b; ++b, ++dbt, ++dct)
|
|
|
|
{
|
|
|
|
w = (mp_word) *da * (mp_word) *dbt + w + (mp_word) *dct;
|
|
|
|
|
|
|
|
*dct = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
}
|
|
|
|
|
|
|
|
*dct = (mp_digit) w;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_ksqr(da, dc, size_a) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_ksqr(mp_digit *da, mp_digit *dc, mp_size size_a)
|
|
|
|
{
|
|
|
|
if (multiply_threshold && size_a > multiply_threshold)
|
|
|
|
{
|
|
|
|
mp_size bot_size = (size_a + 1) / 2;
|
|
|
|
mp_digit *a_top = da + bot_size;
|
|
|
|
mp_digit *t1,
|
|
|
|
*t2,
|
|
|
|
*t3;
|
|
|
|
mp_size at_size = size_a - bot_size;
|
|
|
|
mp_size buf_size = 2 * bot_size;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((t1 = s_alloc(4 * buf_size)) == NULL)
|
|
|
|
return 0;
|
|
|
|
t2 = t1 + buf_size;
|
|
|
|
t3 = t2 + buf_size;
|
|
|
|
ZERO(t1, 4 * buf_size);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) s_ksqr(da, t1, bot_size); /* t1 = a0 ^ 2 */
|
|
|
|
(void) s_ksqr(a_top, t2, at_size); /* t2 = a1 ^ 2 */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) s_kmul(da, a_top, t3, bot_size, at_size); /* t3 = a0 * a1 */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Quick multiply t3 by 2, shifting left (can't overflow) */
|
|
|
|
{
|
|
|
|
int i,
|
|
|
|
top = bot_size + at_size;
|
|
|
|
mp_word w,
|
|
|
|
save = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (i = 0; i < top; ++i)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
w = t3[i];
|
|
|
|
w = (w << 1) | save;
|
|
|
|
t3[i] = LOWER_HALF(w);
|
|
|
|
save = UPPER_HALF(w);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
t3[i] = LOWER_HALF(save);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Assemble the output value */
|
|
|
|
COPY(t1, dc, 2 * bot_size);
|
|
|
|
(void) s_uadd(t3, dc + bot_size, dc + bot_size,
|
|
|
|
buf_size + 1, buf_size + 1);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) s_uadd(t2, dc + 2 * bot_size, dc + 2 * bot_size,
|
|
|
|
buf_size, buf_size);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
px_free(t1); /* note that t2 and t2 are internal pointers
|
2006-10-04 02:30:14 +02:00
|
|
|
* only */
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
s_usqr(da, dc, size_a);
|
|
|
|
}
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_usqr(da, dc, size_a) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_usqr(mp_digit *da, mp_digit *dc, mp_size size_a)
|
|
|
|
{
|
|
|
|
mp_size i,
|
|
|
|
j;
|
|
|
|
mp_word w;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (i = 0; i < size_a; ++i, dc += 2, ++da)
|
|
|
|
{
|
|
|
|
mp_digit *dct = dc,
|
|
|
|
*dat = da;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (*da == 0)
|
|
|
|
continue;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Take care of the first digit, no rollover */
|
|
|
|
w = (mp_word) *dat * (mp_word) *dat + (mp_word) *dct;
|
|
|
|
*dct = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
++dat;
|
|
|
|
++dct;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (j = i + 1; j < size_a; ++j, ++dat, ++dct)
|
|
|
|
{
|
|
|
|
mp_word t = (mp_word) *da * (mp_word) *dat;
|
|
|
|
mp_word u = w + (mp_word) *dct,
|
|
|
|
ov = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Check if doubling t will overflow a word */
|
|
|
|
if (HIGH_BIT_SET(t))
|
|
|
|
ov = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
w = t + t;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Check if adding u to w will overflow a word */
|
|
|
|
if (ADD_WILL_OVERFLOW(w, u))
|
|
|
|
ov = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
w += u;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
*dct = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
if (ov)
|
|
|
|
{
|
|
|
|
w += MP_DIGIT_MAX; /* MP_RADIX */
|
|
|
|
++w;
|
|
|
|
}
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
w = w + *dct;
|
|
|
|
*dct = (mp_digit) w;
|
|
|
|
while ((w = UPPER_HALF(w)) != 0)
|
|
|
|
{
|
|
|
|
++dct;
|
|
|
|
w = w + *dct;
|
|
|
|
*dct = LOWER_HALF(w);
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
assert(w == 0);
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_dadd(a, b) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_dadd(mp_int a, mp_digit b)
|
|
|
|
{
|
|
|
|
mp_word w = 0;
|
|
|
|
mp_digit *da = MP_DIGITS(a);
|
|
|
|
mp_size ua = MP_USED(a);
|
|
|
|
|
|
|
|
w = (mp_word) *da + b;
|
|
|
|
*da++ = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
|
|
|
|
for (ua -= 1; ua > 0; --ua, ++da)
|
|
|
|
{
|
|
|
|
w = (mp_word) *da + w;
|
|
|
|
|
|
|
|
*da = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
}
|
|
|
|
|
|
|
|
if (w)
|
|
|
|
{
|
|
|
|
*da = (mp_digit) w;
|
|
|
|
MP_USED(a) += 1;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_dmul(a, b) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_dmul(mp_int a, mp_digit b)
|
|
|
|
{
|
|
|
|
mp_word w = 0;
|
|
|
|
mp_digit *da = MP_DIGITS(a);
|
|
|
|
mp_size ua = MP_USED(a);
|
|
|
|
|
|
|
|
while (ua > 0)
|
|
|
|
{
|
|
|
|
w = (mp_word) *da * b + w;
|
|
|
|
*da++ = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
--ua;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (w)
|
|
|
|
{
|
|
|
|
*da = (mp_digit) w;
|
|
|
|
MP_USED(a) += 1;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_dbmul(da, b, dc, size_a) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_dbmul(mp_digit *da, mp_digit b, mp_digit *dc, mp_size size_a)
|
|
|
|
{
|
|
|
|
mp_word w = 0;
|
|
|
|
|
|
|
|
while (size_a > 0)
|
|
|
|
{
|
|
|
|
w = (mp_word) *da++ * (mp_word) b + w;
|
|
|
|
|
|
|
|
*dc++ = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w);
|
|
|
|
--size_a;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (w)
|
|
|
|
*dc = LOWER_HALF(w);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_ddiv(da, d, dc, size_a) */
|
|
|
|
|
|
|
|
static mp_digit
|
|
|
|
s_ddiv(mp_int a, mp_digit b)
|
|
|
|
{
|
|
|
|
mp_word w = 0,
|
|
|
|
qdigit;
|
|
|
|
mp_size ua = MP_USED(a);
|
|
|
|
mp_digit *da = MP_DIGITS(a) + ua - 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for ( /* */ ; ua > 0; --ua, --da)
|
|
|
|
{
|
|
|
|
w = (w << MP_DIGIT_BIT) | *da;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (w >= b)
|
|
|
|
{
|
|
|
|
qdigit = w / b;
|
|
|
|
w = w % b;
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
qdigit = 0;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
*da = (mp_digit) qdigit;
|
|
|
|
}
|
|
|
|
|
|
|
|
CLAMP(a);
|
|
|
|
return (mp_digit) w;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_qdiv(z, p2) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_qdiv(mp_int z, mp_size p2)
|
|
|
|
{
|
|
|
|
mp_size ndig = p2 / MP_DIGIT_BIT,
|
|
|
|
nbits = p2 % MP_DIGIT_BIT;
|
|
|
|
mp_size uz = MP_USED(z);
|
|
|
|
|
|
|
|
if (ndig)
|
|
|
|
{
|
|
|
|
mp_size mark;
|
|
|
|
mp_digit *to,
|
|
|
|
*from;
|
|
|
|
|
|
|
|
if (ndig >= uz)
|
|
|
|
{
|
|
|
|
mp_int_zero(z);
|
|
|
|
return;
|
|
|
|
}
|
|
|
|
|
|
|
|
to = MP_DIGITS(z);
|
|
|
|
from = to + ndig;
|
|
|
|
|
|
|
|
for (mark = ndig; mark < uz; ++mark)
|
|
|
|
*to++ = *from++;
|
|
|
|
|
|
|
|
MP_USED(z) = uz - ndig;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (nbits)
|
|
|
|
{
|
|
|
|
mp_digit d = 0,
|
|
|
|
*dz,
|
|
|
|
save;
|
|
|
|
mp_size up = MP_DIGIT_BIT - nbits;
|
|
|
|
|
|
|
|
uz = MP_USED(z);
|
|
|
|
dz = MP_DIGITS(z) + uz - 1;
|
|
|
|
|
|
|
|
for ( /* */ ; uz > 0; --uz, --dz)
|
|
|
|
{
|
|
|
|
save = *dz;
|
|
|
|
|
|
|
|
*dz = (*dz >> nbits) | (d << up);
|
|
|
|
d = save;
|
|
|
|
}
|
|
|
|
|
|
|
|
CLAMP(z);
|
|
|
|
}
|
|
|
|
|
|
|
|
if (MP_USED(z) == 1 && z->digits[0] == 0)
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_qmod(z, p2) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_qmod(mp_int z, mp_size p2)
|
|
|
|
{
|
|
|
|
mp_size start = p2 / MP_DIGIT_BIT + 1,
|
|
|
|
rest = p2 % MP_DIGIT_BIT;
|
|
|
|
mp_size uz = MP_USED(z);
|
|
|
|
mp_digit mask = (1 << rest) - 1;
|
|
|
|
|
|
|
|
if (start <= uz)
|
|
|
|
{
|
|
|
|
MP_USED(z) = start;
|
|
|
|
z->digits[start - 1] &= mask;
|
|
|
|
CLAMP(z);
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_qmul(z, p2) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_qmul(mp_int z, mp_size p2)
|
|
|
|
{
|
|
|
|
mp_size uz,
|
|
|
|
need,
|
|
|
|
rest,
|
|
|
|
extra,
|
|
|
|
i;
|
|
|
|
mp_digit *from,
|
|
|
|
*to,
|
|
|
|
d;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (p2 == 0)
|
|
|
|
return 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
uz = MP_USED(z);
|
|
|
|
need = p2 / MP_DIGIT_BIT;
|
|
|
|
rest = p2 % MP_DIGIT_BIT;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/*
|
|
|
|
* Figure out if we need an extra digit at the top end; this occurs if the
|
|
|
|
* topmost `rest' bits of the high-order digit of z are not zero, meaning
|
|
|
|
* they will be shifted off the end if not preserved
|
|
|
|
*/
|
|
|
|
extra = 0;
|
|
|
|
if (rest != 0)
|
|
|
|
{
|
|
|
|
mp_digit *dz = MP_DIGITS(z) + uz - 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if ((*dz >> (MP_DIGIT_BIT - rest)) != 0)
|
|
|
|
extra = 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (!s_pad(z, uz + need + extra))
|
|
|
|
return 0;
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
/*
|
2006-07-13 06:15:25 +02:00
|
|
|
* If we need to shift by whole digits, do that in one pass, then to back
|
|
|
|
* and shift by partial digits.
|
2006-10-04 02:30:14 +02:00
|
|
|
*/
|
2006-07-13 06:15:25 +02:00
|
|
|
if (need > 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
from = MP_DIGITS(z) + uz - 1;
|
|
|
|
to = from + need;
|
|
|
|
|
|
|
|
for (i = 0; i < uz; ++i)
|
|
|
|
*to-- = *from--;
|
|
|
|
|
|
|
|
ZERO(MP_DIGITS(z), need);
|
|
|
|
uz += need;
|
|
|
|
}
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
if (rest)
|
|
|
|
{
|
|
|
|
d = 0;
|
2006-07-13 06:15:25 +02:00
|
|
|
for (i = need, from = MP_DIGITS(z) + need; i < uz; ++i, ++from)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_digit save = *from;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
*from = (*from << rest) | (d >> (MP_DIGIT_BIT - rest));
|
|
|
|
d = save;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
d >>= (MP_DIGIT_BIT - rest);
|
|
|
|
if (d != 0)
|
|
|
|
{
|
|
|
|
*from = d;
|
|
|
|
uz += extra;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
MP_USED(z) = uz;
|
|
|
|
CLAMP(z);
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_qsub(z, p2) */
|
|
|
|
|
|
|
|
/* Subtract |z| from 2^p2, assuming 2^p2 > |z|, and set z to be positive */
|
|
|
|
static int
|
|
|
|
s_qsub(mp_int z, mp_size p2)
|
|
|
|
{
|
|
|
|
mp_digit hi = (1 << (p2 % MP_DIGIT_BIT)),
|
|
|
|
*zp;
|
|
|
|
mp_size tdig = (p2 / MP_DIGIT_BIT),
|
|
|
|
pos;
|
|
|
|
mp_word w = 0;
|
|
|
|
|
|
|
|
if (!s_pad(z, tdig + 1))
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
for (pos = 0, zp = MP_DIGITS(z); pos < tdig; ++pos, ++zp)
|
|
|
|
{
|
|
|
|
w = ((mp_word) MP_DIGIT_MAX + 1) - w - (mp_word) *zp;
|
|
|
|
|
|
|
|
*zp = LOWER_HALF(w);
|
|
|
|
w = UPPER_HALF(w) ? 0 : 1;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
w = ((mp_word) MP_DIGIT_MAX + 1 + hi) - w - (mp_word) *zp;
|
|
|
|
*zp = LOWER_HALF(w);
|
|
|
|
|
|
|
|
assert(UPPER_HALF(w) != 0); /* no borrow out should be possible */
|
|
|
|
|
|
|
|
MP_SIGN(z) = MP_ZPOS;
|
|
|
|
CLAMP(z);
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_dp2k(z) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_dp2k(mp_int z)
|
|
|
|
{
|
|
|
|
int k = 0;
|
|
|
|
mp_digit *dp = MP_DIGITS(z),
|
|
|
|
d;
|
|
|
|
|
|
|
|
if (MP_USED(z) == 1 && *dp == 0)
|
|
|
|
return 1;
|
|
|
|
|
|
|
|
while (*dp == 0)
|
|
|
|
{
|
|
|
|
k += MP_DIGIT_BIT;
|
|
|
|
++dp;
|
|
|
|
}
|
|
|
|
|
2006-10-04 02:30:14 +02:00
|
|
|
d = *dp;
|
2006-07-13 06:15:25 +02:00
|
|
|
while ((d & 1) == 0)
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
|
|
|
d >>= 1;
|
|
|
|
++k;
|
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
return k;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_isp2(z) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_isp2(mp_int z)
|
|
|
|
{
|
|
|
|
mp_size uz = MP_USED(z),
|
|
|
|
k = 0;
|
|
|
|
mp_digit *dz = MP_DIGITS(z),
|
|
|
|
d;
|
|
|
|
|
|
|
|
while (uz > 1)
|
|
|
|
{
|
|
|
|
if (*dz++ != 0)
|
|
|
|
return -1;
|
|
|
|
k += MP_DIGIT_BIT;
|
|
|
|
--uz;
|
|
|
|
}
|
|
|
|
|
|
|
|
d = *dz;
|
|
|
|
while (d > 1)
|
|
|
|
{
|
|
|
|
if (d & 1)
|
|
|
|
return -1;
|
|
|
|
++k;
|
|
|
|
d >>= 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
return (int) k;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_2expt(z, k) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_2expt(mp_int z, int k)
|
|
|
|
{
|
|
|
|
mp_size ndig,
|
|
|
|
rest;
|
|
|
|
mp_digit *dz;
|
|
|
|
|
|
|
|
ndig = (k + MP_DIGIT_BIT) / MP_DIGIT_BIT;
|
|
|
|
rest = k % MP_DIGIT_BIT;
|
|
|
|
|
|
|
|
if (!s_pad(z, ndig))
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
dz = MP_DIGITS(z);
|
|
|
|
ZERO(dz, ndig);
|
|
|
|
*(dz + ndig - 1) = (1 << rest);
|
|
|
|
MP_USED(z) = ndig;
|
|
|
|
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_norm(a, b) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_norm(mp_int a, mp_int b)
|
|
|
|
{
|
|
|
|
mp_digit d = b->digits[MP_USED(b) - 1];
|
|
|
|
int k = 0;
|
|
|
|
|
2007-07-16 00:43:40 +02:00
|
|
|
while (d < (mp_digit) ((mp_digit) 1 << (MP_DIGIT_BIT - 1)))
|
2006-07-13 06:15:25 +02:00
|
|
|
{ /* d < (MP_RADIX / 2) */
|
|
|
|
d <<= 1;
|
|
|
|
++k;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* These multiplications can't fail */
|
|
|
|
if (k != 0)
|
|
|
|
{
|
|
|
|
(void) s_qmul(a, (mp_size) k);
|
|
|
|
(void) s_qmul(b, (mp_size) k);
|
|
|
|
}
|
|
|
|
|
|
|
|
return k;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_brmu(z, m) */
|
|
|
|
|
|
|
|
static mp_result
|
|
|
|
s_brmu(mp_int z, mp_int m)
|
|
|
|
{
|
|
|
|
mp_size um = MP_USED(m) * 2;
|
|
|
|
|
|
|
|
if (!s_pad(z, um))
|
|
|
|
return MP_MEMORY;
|
|
|
|
|
|
|
|
s_2expt(z, MP_DIGIT_BIT * um);
|
|
|
|
return mp_int_div(z, m, z, NULL);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_reduce(x, m, mu, q1, q2) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_reduce(mp_int x, mp_int m, mp_int mu, mp_int q1, mp_int q2)
|
|
|
|
{
|
|
|
|
mp_size um = MP_USED(m),
|
|
|
|
umb_p1,
|
|
|
|
umb_m1;
|
|
|
|
|
|
|
|
umb_p1 = (um + 1) * MP_DIGIT_BIT;
|
|
|
|
umb_m1 = (um - 1) * MP_DIGIT_BIT;
|
|
|
|
|
|
|
|
if (mp_int_copy(x, q1) != MP_OK)
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
/* Compute q2 = floor((floor(x / b^(k-1)) * mu) / b^(k+1)) */
|
|
|
|
s_qdiv(q1, umb_m1);
|
|
|
|
UMUL(q1, mu, q2);
|
|
|
|
s_qdiv(q2, umb_p1);
|
|
|
|
|
|
|
|
/* Set x = x mod b^(k+1) */
|
|
|
|
s_qmod(x, umb_p1);
|
|
|
|
|
|
|
|
/*
|
|
|
|
* Now, q is a guess for the quotient a / m. Compute x - q * m mod
|
|
|
|
* b^(k+1), replacing x. This may be off by a factor of 2m, but no more
|
|
|
|
* than that.
|
|
|
|
*/
|
|
|
|
UMUL(q2, m, q1);
|
|
|
|
s_qmod(q1, umb_p1);
|
|
|
|
(void) mp_int_sub(x, q1, x); /* can't fail */
|
|
|
|
|
|
|
|
/*
|
|
|
|
* The result may be < 0; if it is, add b^(k+1) to pin it in the proper
|
|
|
|
* range.
|
|
|
|
*/
|
|
|
|
if ((CMPZ(x) < 0) && !s_qsub(x, umb_p1))
|
|
|
|
return 0;
|
|
|
|
|
|
|
|
/*
|
|
|
|
* If x > m, we need to back it off until it is in range. This will be
|
|
|
|
* required at most twice.
|
|
|
|
*/
|
|
|
|
if (mp_int_compare(x, m) >= 0)
|
|
|
|
(void) mp_int_sub(x, m, x);
|
|
|
|
if (mp_int_compare(x, m) >= 0)
|
|
|
|
(void) mp_int_sub(x, m, x);
|
|
|
|
|
|
|
|
/* At this point, x has been properly reduced. */
|
|
|
|
return 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_embar(a, b, m, mu, c) */
|
|
|
|
|
|
|
|
/* Perform modular exponentiation using Barrett's method, where mu is
|
|
|
|
the reduction constant for m. Assumes a < m, b > 0. */
|
|
|
|
static mp_result
|
|
|
|
s_embar(mp_int a, mp_int b, mp_int m, mp_int mu, mp_int c)
|
|
|
|
{
|
|
|
|
mp_digit *db,
|
|
|
|
*dbt,
|
|
|
|
umu,
|
|
|
|
d;
|
|
|
|
mpz_t temp[3];
|
|
|
|
mp_result res;
|
|
|
|
int last = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
umu = MP_USED(mu);
|
|
|
|
db = MP_DIGITS(b);
|
|
|
|
dbt = db + MP_USED(b) - 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
while (last < 3)
|
2015-02-02 16:00:45 +01:00
|
|
|
{
|
|
|
|
SETUP(mp_int_init_size(TEMP(last), 4 * umu), last);
|
|
|
|
ZERO(MP_DIGITS(TEMP(last - 1)), MP_ALLOC(TEMP(last - 1)));
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) mp_int_set_value(c, 1);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Take care of low-order digits */
|
|
|
|
while (db < dbt)
|
|
|
|
{
|
|
|
|
int i;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
for (d = *db, i = MP_DIGIT_BIT; i > 0; --i, d >>= 1)
|
|
|
|
{
|
|
|
|
if (d & 1)
|
|
|
|
{
|
|
|
|
/* The use of a second temporary avoids allocation */
|
|
|
|
UMUL(c, a, TEMP(0));
|
|
|
|
if (!s_reduce(TEMP(0), m, mu, TEMP(1), TEMP(2)))
|
|
|
|
{
|
|
|
|
res = MP_MEMORY;
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_copy(TEMP(0), c);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
USQR(a, TEMP(0));
|
|
|
|
assert(MP_SIGN(TEMP(0)) == MP_ZPOS);
|
|
|
|
if (!s_reduce(TEMP(0), m, mu, TEMP(1), TEMP(2)))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
res = MP_MEMORY;
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
assert(MP_SIGN(TEMP(0)) == MP_ZPOS);
|
|
|
|
mp_int_copy(TEMP(0), a);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
++db;
|
2006-07-13 06:15:25 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
/* Take care of highest-order digit */
|
|
|
|
d = *dbt;
|
2006-10-04 02:30:14 +02:00
|
|
|
for (;;)
|
|
|
|
{
|
|
|
|
if (d & 1)
|
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
UMUL(c, a, TEMP(0));
|
|
|
|
if (!s_reduce(TEMP(0), m, mu, TEMP(1), TEMP(2)))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
res = MP_MEMORY;
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_copy(TEMP(0), c);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
|
|
|
d >>= 1;
|
|
|
|
if (!d)
|
|
|
|
break;
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
USQR(a, TEMP(0));
|
|
|
|
if (!s_reduce(TEMP(0), m, mu, TEMP(1), TEMP(2)))
|
2006-10-04 02:30:14 +02:00
|
|
|
{
|
2006-07-13 06:15:25 +02:00
|
|
|
res = MP_MEMORY;
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
(void) mp_int_copy(TEMP(0), a);
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
2006-07-13 06:15:25 +02:00
|
|
|
|
|
|
|
CLEANUP:
|
|
|
|
while (--last >= 0)
|
|
|
|
mp_int_clear(TEMP(last));
|
|
|
|
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_udiv(a, b) */
|
|
|
|
|
|
|
|
/* Precondition: a >= b and b > 0
|
|
|
|
Postcondition: a' = a / b, b' = a % b
|
|
|
|
*/
|
|
|
|
static mp_result
|
|
|
|
s_udiv(mp_int a, mp_int b)
|
|
|
|
{
|
|
|
|
mpz_t q,
|
|
|
|
r,
|
|
|
|
t;
|
|
|
|
mp_size ua,
|
|
|
|
ub,
|
|
|
|
qpos = 0;
|
|
|
|
mp_digit *da,
|
|
|
|
btop;
|
|
|
|
mp_result res = MP_OK;
|
|
|
|
int k,
|
|
|
|
skip = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Force signs to positive */
|
|
|
|
MP_SIGN(a) = MP_ZPOS;
|
|
|
|
MP_SIGN(b) = MP_ZPOS;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Normalize, per Knuth */
|
|
|
|
k = s_norm(a, b);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
ua = MP_USED(a);
|
|
|
|
ub = MP_USED(b);
|
|
|
|
btop = b->digits[ub - 1];
|
|
|
|
if ((res = mp_int_init_size(&q, ua)) != MP_OK)
|
|
|
|
return res;
|
|
|
|
if ((res = mp_int_init_size(&t, ua + 1)) != MP_OK)
|
|
|
|
goto CLEANUP;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
da = MP_DIGITS(a);
|
|
|
|
r.digits = da + ua - 1; /* The contents of r are shared with a */
|
|
|
|
r.used = 1;
|
|
|
|
r.sign = MP_ZPOS;
|
|
|
|
r.alloc = MP_ALLOC(a);
|
|
|
|
ZERO(t.digits, t.alloc);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Solve for quotient digits, store in q.digits in reverse order */
|
|
|
|
while (r.digits >= da)
|
|
|
|
{
|
|
|
|
assert(qpos <= q.alloc);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (s_ucmp(b, &r) > 0)
|
|
|
|
{
|
|
|
|
r.digits -= 1;
|
|
|
|
r.used += 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (++skip > 1)
|
|
|
|
q.digits[qpos++] = 0;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
CLAMP(&r);
|
|
|
|
}
|
|
|
|
else
|
|
|
|
{
|
|
|
|
mp_word pfx = r.digits[r.used - 1];
|
|
|
|
mp_word qdigit;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
if (r.used > 1 && (pfx < btop || r.digits[r.used - 2] == 0))
|
|
|
|
{
|
|
|
|
pfx <<= MP_DIGIT_BIT / 2;
|
|
|
|
pfx <<= MP_DIGIT_BIT / 2;
|
|
|
|
pfx |= r.digits[r.used - 2];
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
qdigit = pfx / btop;
|
|
|
|
if (qdigit > MP_DIGIT_MAX)
|
|
|
|
qdigit = 1;
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_dbmul(MP_DIGITS(b), (mp_digit) qdigit, t.digits, ub);
|
|
|
|
t.used = ub + 1;
|
|
|
|
CLAMP(&t);
|
|
|
|
while (s_ucmp(&t, &r) > 0)
|
|
|
|
{
|
|
|
|
--qdigit;
|
|
|
|
(void) mp_int_sub(&t, b, &t); /* cannot fail */
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
s_usub(r.digits, t.digits, r.digits, r.used, t.used);
|
|
|
|
CLAMP(&r);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
q.digits[qpos++] = (mp_digit) qdigit;
|
|
|
|
ZERO(t.digits, t.used);
|
|
|
|
skip = 0;
|
|
|
|
}
|
2006-10-04 02:30:14 +02:00
|
|
|
}
|
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Put quotient digits in the correct order, and discard extra zeroes */
|
|
|
|
q.used = qpos;
|
|
|
|
REV(mp_digit, q.digits, qpos);
|
|
|
|
CLAMP(&q);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
/* Denormalize the remainder */
|
|
|
|
CLAMP(a);
|
|
|
|
if (k != 0)
|
|
|
|
s_qdiv(a, k);
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_copy(a, b); /* ok: 0 <= r < b */
|
|
|
|
mp_int_copy(&q, a); /* ok: q <= a */
|
2006-10-04 02:30:14 +02:00
|
|
|
|
2006-07-13 06:15:25 +02:00
|
|
|
mp_int_clear(&t);
|
|
|
|
CLEANUP:
|
|
|
|
mp_int_clear(&q);
|
|
|
|
return res;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_outlen(z, r) */
|
|
|
|
|
|
|
|
/* Precondition: 2 <= r < 64 */
|
|
|
|
static int
|
|
|
|
s_outlen(mp_int z, mp_size r)
|
|
|
|
{
|
|
|
|
mp_result bits;
|
|
|
|
double raw;
|
|
|
|
|
|
|
|
bits = mp_int_count_bits(z);
|
|
|
|
raw = (double) bits * s_log2[r];
|
|
|
|
|
|
|
|
return (int) (raw + 0.999999);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_inlen(len, r) */
|
|
|
|
|
|
|
|
static mp_size
|
|
|
|
s_inlen(int len, mp_size r)
|
|
|
|
{
|
|
|
|
double raw = (double) len / s_log2[r];
|
|
|
|
mp_size bits = (mp_size) (raw + 0.5);
|
|
|
|
|
|
|
|
return (mp_size) ((bits + (MP_DIGIT_BIT - 1)) / MP_DIGIT_BIT);
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_ch2val(c, r) */
|
|
|
|
|
|
|
|
static int
|
|
|
|
s_ch2val(char c, int r)
|
|
|
|
{
|
|
|
|
int out;
|
|
|
|
|
2006-09-22 23:39:58 +02:00
|
|
|
if (isdigit((unsigned char) c))
|
2006-07-13 06:15:25 +02:00
|
|
|
out = c - '0';
|
2006-09-22 23:39:58 +02:00
|
|
|
else if (r > 10 && isalpha((unsigned char) c))
|
|
|
|
out = toupper((unsigned char) c) - 'A' + 10;
|
2006-07-13 06:15:25 +02:00
|
|
|
else
|
|
|
|
return -1;
|
|
|
|
|
|
|
|
return (out >= r) ? -1 : out;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_val2ch(v, caps) */
|
|
|
|
|
|
|
|
static char
|
|
|
|
s_val2ch(int v, int caps)
|
|
|
|
{
|
|
|
|
assert(v >= 0);
|
|
|
|
|
|
|
|
if (v < 10)
|
|
|
|
return v + '0';
|
|
|
|
else
|
|
|
|
{
|
|
|
|
char out = (v - 10) + 'a';
|
|
|
|
|
|
|
|
if (caps)
|
2006-09-22 23:39:58 +02:00
|
|
|
return toupper((unsigned char) out);
|
2006-07-13 06:15:25 +02:00
|
|
|
else
|
|
|
|
return out;
|
|
|
|
}
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_2comp(buf, len) */
|
|
|
|
|
|
|
|
static void
|
|
|
|
s_2comp(unsigned char *buf, int len)
|
|
|
|
{
|
|
|
|
int i;
|
|
|
|
unsigned short s = 1;
|
|
|
|
|
|
|
|
for (i = len - 1; i >= 0; --i)
|
|
|
|
{
|
|
|
|
unsigned char c = ~buf[i];
|
|
|
|
|
|
|
|
s = c + s;
|
|
|
|
c = s & UCHAR_MAX;
|
|
|
|
s >>= CHAR_BIT;
|
|
|
|
|
|
|
|
buf[i] = c;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* last carry out is ignored */
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_tobin(z, buf, *limpos) */
|
|
|
|
|
|
|
|
static mp_result
|
|
|
|
s_tobin(mp_int z, unsigned char *buf, int *limpos, int pad)
|
|
|
|
{
|
|
|
|
mp_size uz;
|
|
|
|
mp_digit *dz;
|
|
|
|
int pos = 0,
|
|
|
|
limit = *limpos;
|
|
|
|
|
|
|
|
uz = MP_USED(z);
|
|
|
|
dz = MP_DIGITS(z);
|
|
|
|
while (uz > 0 && pos < limit)
|
|
|
|
{
|
|
|
|
mp_digit d = *dz++;
|
|
|
|
int i;
|
|
|
|
|
|
|
|
for (i = sizeof(mp_digit); i > 0 && pos < limit; --i)
|
|
|
|
{
|
|
|
|
buf[pos++] = (unsigned char) d;
|
|
|
|
d >>= CHAR_BIT;
|
|
|
|
|
|
|
|
/* Don't write leading zeroes */
|
|
|
|
if (d == 0 && uz == 1)
|
|
|
|
i = 0; /* exit loop without signaling truncation */
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Detect truncation (loop exited with pos >= limit) */
|
|
|
|
if (i > 0)
|
|
|
|
break;
|
|
|
|
|
|
|
|
--uz;
|
|
|
|
}
|
|
|
|
|
|
|
|
if (pad != 0 && (buf[pos - 1] >> (CHAR_BIT - 1)))
|
|
|
|
{
|
|
|
|
if (pos < limit)
|
|
|
|
buf[pos++] = 0;
|
|
|
|
else
|
|
|
|
uz = 1;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* Digits are in reverse order, fix that */
|
|
|
|
REV(unsigned char, buf, pos);
|
|
|
|
|
|
|
|
/* Return the number of bytes actually written */
|
|
|
|
*limpos = pos;
|
|
|
|
|
|
|
|
return (uz == 0) ? MP_OK : MP_TRUNC;
|
|
|
|
}
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* {{{ s_print(tag, z) */
|
|
|
|
|
|
|
|
#if 0
|
|
|
|
void
|
|
|
|
s_print(char *tag, mp_int z)
|
|
|
|
{
|
|
|
|
int i;
|
|
|
|
|
|
|
|
fprintf(stderr, "%s: %c ", tag,
|
|
|
|
(MP_SIGN(z) == MP_NEG) ? '-' : '+');
|
|
|
|
|
|
|
|
for (i = MP_USED(z) - 1; i >= 0; --i)
|
|
|
|
fprintf(stderr, "%0*X", (int) (MP_DIGIT_BIT / 4), z->digits[i]);
|
|
|
|
|
|
|
|
fputc('\n', stderr);
|
|
|
|
|
|
|
|
}
|
|
|
|
|
|
|
|
void
|
|
|
|
s_print_buf(char *tag, mp_digit *buf, mp_size num)
|
|
|
|
{
|
|
|
|
int i;
|
|
|
|
|
|
|
|
fprintf(stderr, "%s: ", tag);
|
|
|
|
|
|
|
|
for (i = num - 1; i >= 0; --i)
|
|
|
|
fprintf(stderr, "%0*X", (int) (MP_DIGIT_BIT / 4), buf[i]);
|
|
|
|
|
|
|
|
fputc('\n', stderr);
|
|
|
|
}
|
|
|
|
#endif
|
|
|
|
|
|
|
|
/* }}} */
|
|
|
|
|
|
|
|
/* HERE THERE BE DRAGONS */
|