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