1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149
|
#ifdef PLAN9
#pragma src "/sys/src/libmp"
#pragma lib "libmp.a"
#endif
#define _MPINT 1
/*
* the code assumes mpdigit to be at least an int
* mpdigit must be an atomic type. mpdigit is defined
* in the architecture specific u.h
*/
typedef struct mpint mpint;
struct mpint
{
int sign; /* +1 or -1 */
int size; /* allocated digits */
int top; /* significant digits */
mpdigit *p;
char flags;
};
enum
{
MPstatic= 0x01,
Dbytes= sizeof(mpdigit), /* bytes per digit */
Dbits= Dbytes*8 /* bits per digit */
};
/* allocation */
void mpsetminbits(int n); /* newly created mpint's get at least n bits */
mpint* mpnew(int n); /* create a new mpint with at least n bits */
void mpfree(mpint *b);
void mpbits(mpint *b, int n); /* ensure that b has at least n bits */
void mpnorm(mpint *b); /* dump leading zeros */
mpint* mpcopy(mpint *b);
void mpassign(mpint *old, mpint *new);
/* random bits */
mpint* mprand(int bits, void (*gen)(uchar*, int), mpint *b);
/* conversion */
mpint* strtomp(char*, char**, int, mpint*); /* ascii */
int mpfmt(Fmt*);
char* mptoa(mpint*, int, char*, int);
mpint* letomp(uchar*, uint, mpint*); /* byte array, little-endian */
int mptole(mpint*, uchar*, uint, uchar**);
mpint* betomp(uchar*, uint, mpint*); /* byte array, little-endian */
int mptobe(mpint*, uchar*, uint, uchar**);
uint mptoui(mpint*); /* unsigned int */
mpint* uitomp(uint, mpint*);
int mptoi(mpint*); /* int */
mpint* itomp(int, mpint*);
uvlong mptouv(mpint*); /* unsigned vlong */
mpint* uvtomp(uvlong, mpint*);
vlong mptov(mpint*); /* vlong */
mpint* vtomp(vlong, mpint*);
/* divide 2 digits by one */
void mpdigdiv(mpdigit *dividend, mpdigit divisor, mpdigit *quotient);
/* in the following, the result mpint may be */
/* the same as one of the inputs. */
void mpadd(mpint *b1, mpint *b2, mpint *sum); /* sum = b1+b2 */
void mpsub(mpint *b1, mpint *b2, mpint *diff); /* diff = b1-b2 */
void mpleft(mpint *b, int shift, mpint *res); /* res = b<<shift */
void mpright(mpint *b, int shift, mpint *res); /* res = b>>shift */
void mpmul(mpint *b1, mpint *b2, mpint *prod); /* prod = b1*b2 */
void mpexp(mpint *b, mpint *e, mpint *m, mpint *res); /* res = b**e mod m */
void mpmod(mpint *b, mpint *m, mpint *remainder); /* remainder = b mod m */
/* quotient = dividend/divisor, remainder = dividend % divisor */
void mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder);
/* return neg, 0, pos as b1-b2 is neg, 0, pos */
int mpcmp(mpint *b1, mpint *b2);
/* extended gcd return d, x, and y, s.t. d = gcd(a,b) and ax+by = d */
void mpextendedgcd(mpint *a, mpint *b, mpint *d, mpint *x, mpint *y);
/* res = b**-1 mod m */
void mpinvert(mpint *b, mpint *m, mpint *res);
/* bit counting */
int mpsignif(mpint*); /* number of sigificant bits in mantissa */
int mplowbits0(mpint*); /* k, where n = 2**k * q for odd q */
/* well known constants */
extern mpint *mpzero, *mpone, *mptwo;
/* sum[0:alen] = a[0:alen-1] + b[0:blen-1] */
/* prereq: alen >= blen, sum has room for alen+1 digits */
void mpvecadd(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *sum);
/* diff[0:alen-1] = a[0:alen-1] - b[0:blen-1] */
/* prereq: alen >= blen, diff has room for alen digits */
void mpvecsub(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *diff);
/* p[0:n] += m * b[0:n-1] */
/* prereq: p has room for n+1 digits */
void mpvecdigmuladd(mpdigit *b, int n, mpdigit m, mpdigit *p);
/* p[0:n] -= m * b[0:n-1] */
/* prereq: p has room for n+1 digits */
int mpvecdigmulsub(mpdigit *b, int n, mpdigit m, mpdigit *p);
/* p[0:alen*blen-1] = a[0:alen-1] * b[0:blen-1] */
/* prereq: alen >= blen, p has room for m*n digits */
void mpvecmul(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p);
/* sign of a - b or zero if the same */
int mpveccmp(mpdigit *a, int alen, mpdigit *b, int blen);
/* divide the 2 digit dividend by the one digit divisor and stick in quotient */
/* we assume that the result is one digit - overflow is all 1's */
void mpdigdiv(mpdigit *dividend, mpdigit divisor, mpdigit *quotient);
/* playing with magnitudes */
int mpmagcmp(mpint *b1, mpint *b2);
void mpmagadd(mpint *b1, mpint *b2, mpint *sum); /* sum = b1+b2 */
void mpmagsub(mpint *b1, mpint *b2, mpint *sum); /* sum = b1+b2 */
/* chinese remainder theorem */
typedef struct CRTpre CRTpre; /* precomputed values for converting */
/* twixt residues and mpint */
typedef struct CRTres CRTres; /* residue form of an mpint */
#ifdef PLAN9
#pragma incomplete CRTpre
#endif
struct CRTres
{
int n; /* number of residues */
mpint *r[1]; /* residues */
};
CRTpre* crtpre(int, mpint**); /* precompute conversion values */
CRTres* crtin(CRTpre*, mpint*); /* convert mpint to residues */
void crtout(CRTpre*, CRTres*, mpint*); /* convert residues to mpint */
void crtprefree(CRTpre*);
void crtresfree(CRTres*);
#ifdef VARARGCK
#pragma varargck type "B" mpint*
#endif
|