Use gcc -ansi -pedantic in 9c. Fix many non-C89-isms.
This commit is contained in:
rsc
parent
226d80b821
commit
cbeb0b26e4
492 changed files with 3218 additions and 3167 deletions
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// convert a big-endian byte array (most significant byte first) to an mpint
|
||||
/* convert a big-endian byte array (most significant byte first) to an mpint */
|
||||
mpint*
|
||||
betomp(uchar *p, uint n, mpint *b)
|
||||
{
|
||||
|
|
@ -12,18 +12,18 @@ betomp(uchar *p, uint n, mpint *b)
|
|||
if(b == nil)
|
||||
b = mpnew(0);
|
||||
|
||||
// dump leading zeros
|
||||
/* dump leading zeros */
|
||||
while(*p == 0 && n > 1){
|
||||
p++;
|
||||
n--;
|
||||
}
|
||||
|
||||
// get the space
|
||||
/* get the space */
|
||||
mpbits(b, n*8);
|
||||
b->top = DIGITS(n*8);
|
||||
m = b->top-1;
|
||||
|
||||
// first digit might not be Dbytes long
|
||||
/* first digit might not be Dbytes long */
|
||||
s = ((n-1)*8)%Dbits;
|
||||
x = 0;
|
||||
for(; n > 0; n--){
|
||||
|
|
|
|||
|
|
@ -1,20 +1,20 @@
|
|||
#include "os.h"
|
||||
#include <mp.h>
|
||||
|
||||
// chinese remainder theorem
|
||||
//
|
||||
// handbook of applied cryptography, menezes et al, 1997, pp 610 - 613
|
||||
/* chinese remainder theorem */
|
||||
/* */
|
||||
/* handbook of applied cryptography, menezes et al, 1997, pp 610 - 613 */
|
||||
|
||||
struct CRTpre
|
||||
{
|
||||
int n; // number of moduli
|
||||
mpint **m; // pointer to moduli
|
||||
mpint **c; // precomputed coefficients
|
||||
mpint **p; // precomputed products
|
||||
mpint *a[1]; // local storage
|
||||
int n; /* number of moduli */
|
||||
mpint **m; /* pointer to moduli */
|
||||
mpint **c; /* precomputed coefficients */
|
||||
mpint **p; /* precomputed products */
|
||||
mpint *a[1]; /* local storage */
|
||||
};
|
||||
|
||||
// setup crt info, returns a newly created structure
|
||||
/* setup crt info, returns a newly created structure */
|
||||
CRTpre*
|
||||
crtpre(int n, mpint **m)
|
||||
{
|
||||
|
|
@ -30,18 +30,18 @@ crtpre(int n, mpint **m)
|
|||
crt->p = crt->c+n;
|
||||
crt->n = n;
|
||||
|
||||
// make a copy of the moduli
|
||||
/* make a copy of the moduli */
|
||||
for(i = 0; i < n; i++)
|
||||
crt->m[i] = mpcopy(m[i]);
|
||||
|
||||
// precompute the products
|
||||
/* precompute the products */
|
||||
u = mpcopy(mpone);
|
||||
for(i = 0; i < n; i++){
|
||||
mpmul(u, m[i], u);
|
||||
crt->p[i] = mpcopy(u);
|
||||
}
|
||||
|
||||
// precompute the coefficients
|
||||
/* precompute the coefficients */
|
||||
for(i = 1; i < n; i++){
|
||||
crt->c[i] = mpcopy(mpone);
|
||||
for(j = 0; j < i; j++){
|
||||
|
|
@ -70,7 +70,7 @@ crtprefree(CRTpre *crt)
|
|||
free(crt);
|
||||
}
|
||||
|
||||
// convert to residues, returns a newly created structure
|
||||
/* convert to residues, returns a newly created structure */
|
||||
CRTres*
|
||||
crtin(CRTpre *crt, mpint *x)
|
||||
{
|
||||
|
|
@ -88,7 +88,7 @@ crtin(CRTpre *crt, mpint *x)
|
|||
return res;
|
||||
}
|
||||
|
||||
// garners algorithm for converting residue form to linear
|
||||
/* garners algorithm for converting residue form to linear */
|
||||
void
|
||||
crtout(CRTpre *crt, CRTres *res, mpint *x)
|
||||
{
|
||||
|
|
@ -109,7 +109,7 @@ crtout(CRTpre *crt, CRTres *res, mpint *x)
|
|||
mpfree(u);
|
||||
}
|
||||
|
||||
// free the residue
|
||||
/* free the residue */
|
||||
void
|
||||
crtresfree(CRTres *res)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -11,19 +11,19 @@ testcrt(mpint **p)
|
|||
|
||||
fmtinstall('B', mpconv);
|
||||
|
||||
// get a modulus and a test number
|
||||
/* get a modulus and a test number */
|
||||
m = mpnew(1024+160);
|
||||
mpmul(p[0], p[1], m);
|
||||
x = mpnew(1024+160);
|
||||
mpadd(m, mpone, x);
|
||||
|
||||
// do the precomputation for crt conversion
|
||||
/* do the precomputation for crt conversion */
|
||||
crt = crtpre(2, p);
|
||||
|
||||
// convert x to residues
|
||||
/* convert x to residues */
|
||||
res = crtin(crt, x);
|
||||
|
||||
// convert back
|
||||
/* convert back */
|
||||
y = mpnew(1024+160);
|
||||
crtout(crt, res, y);
|
||||
print("x %B\ny %B\n", x, y);
|
||||
|
|
|
|||
|
|
@ -1,12 +1,12 @@
|
|||
#define mpdighi (mpdigit)((ulong)1<<(Dbits-1))
|
||||
#define DIGITS(x) ((Dbits - 1 + (x))/Dbits)
|
||||
|
||||
// for converting between int's and mpint's
|
||||
/* for converting between int's and mpint's */
|
||||
#define MAXUINT ((uint)-1)
|
||||
#define MAXINT (MAXUINT>>1)
|
||||
#define MININT (MAXINT+1)
|
||||
|
||||
// for converting between vlongs's and mpint's
|
||||
/* for converting between vlongs's and mpint's */
|
||||
#define MAXUVLONG (~0ULL)
|
||||
#define MAXVLONG (MAXUVLONG>>1)
|
||||
#define MINVLONG (MAXVLONG+1ULL)
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// convert a little endian byte array (least significant byte first) to an mpint
|
||||
/* convert a little endian byte array (least significant byte first) to an mpint */
|
||||
mpint*
|
||||
letomp(uchar *s, uint n, mpint *b)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -2,14 +2,14 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// sum = abs(b1) + abs(b2), i.e., add the magnitudes
|
||||
/* sum = abs(b1) + abs(b2), i.e., add the magnitudes */
|
||||
void
|
||||
mpmagadd(mpint *b1, mpint *b2, mpint *sum)
|
||||
{
|
||||
int m, n;
|
||||
mpint *t;
|
||||
|
||||
// get the sizes right
|
||||
/* get the sizes right */
|
||||
if(b2->top > b1->top){
|
||||
t = b1;
|
||||
b1 = b2;
|
||||
|
|
@ -34,7 +34,7 @@ mpmagadd(mpint *b1, mpint *b2, mpint *sum)
|
|||
mpnorm(sum);
|
||||
}
|
||||
|
||||
// sum = b1 + b2
|
||||
/* sum = b1 + b2 */
|
||||
void
|
||||
mpadd(mpint *b1, mpint *b2, mpint *sum)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -37,7 +37,7 @@ mpint *mpzero = &_mpzero;
|
|||
|
||||
static int mpmindigits = 33;
|
||||
|
||||
// set minimum digit allocation
|
||||
/* set minimum digit allocation */
|
||||
void
|
||||
mpsetminbits(int n)
|
||||
{
|
||||
|
|
@ -48,7 +48,7 @@ mpsetminbits(int n)
|
|||
mpmindigits = DIGITS(n);
|
||||
}
|
||||
|
||||
// allocate an n bit 0'd number
|
||||
/* allocate an n bit 0'd number */
|
||||
mpint*
|
||||
mpnew(int n)
|
||||
{
|
||||
|
|
@ -72,7 +72,7 @@ mpnew(int n)
|
|||
return b;
|
||||
}
|
||||
|
||||
// guarantee at least n significant bits
|
||||
/* guarantee at least n significant bits */
|
||||
void
|
||||
mpbits(mpint *b, int m)
|
||||
{
|
||||
|
|
@ -101,7 +101,7 @@ mpfree(mpint *b)
|
|||
return;
|
||||
if(b->flags & MPstatic)
|
||||
sysfatal("freeing mp constant");
|
||||
memset(b->p, 0, b->size*Dbytes); // information hiding
|
||||
memset(b->p, 0, b->size*Dbytes); /* information hiding */
|
||||
free(b->p);
|
||||
free(b);
|
||||
}
|
||||
|
|
@ -140,7 +140,7 @@ mpassign(mpint *old, mpint *new)
|
|||
memmove(new->p, old->p, Dbytes*old->top);
|
||||
}
|
||||
|
||||
// number of significant bits in mantissa
|
||||
/* number of significant bits in mantissa */
|
||||
int
|
||||
mpsignif(mpint *n)
|
||||
{
|
||||
|
|
@ -159,7 +159,7 @@ mpsignif(mpint *n)
|
|||
return 0;
|
||||
}
|
||||
|
||||
// k, where n = 2**k * q for odd q
|
||||
/* k, where n = 2**k * q for odd q */
|
||||
int
|
||||
mplowbits0(mpint *n)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// return neg, 0, pos as abs(b1)-abs(b2) is neg, 0, pos
|
||||
/* return neg, 0, pos as abs(b1)-abs(b2) is neg, 0, pos */
|
||||
int
|
||||
mpmagcmp(mpint *b1, mpint *b2)
|
||||
{
|
||||
|
|
@ -15,7 +15,7 @@ mpmagcmp(mpint *b1, mpint *b2)
|
|||
return mpveccmp(b1->p, b1->top, b2->p, b2->top);
|
||||
}
|
||||
|
||||
// return neg, 0, pos as b1-b2 is neg, 0, pos
|
||||
/* return neg, 0, pos as b1-b2 is neg, 0, pos */
|
||||
int
|
||||
mpcmp(mpint *b1, mpint *b2)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -2,9 +2,9 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
//
|
||||
// divide two digits by one and return quotient
|
||||
//
|
||||
/* */
|
||||
/* divide two digits by one and return quotient */
|
||||
/* */
|
||||
void
|
||||
mpdigdiv(mpdigit *dividend, mpdigit divisor, mpdigit *quotient)
|
||||
{
|
||||
|
|
@ -14,15 +14,15 @@ mpdigdiv(mpdigit *dividend, mpdigit divisor, mpdigit *quotient)
|
|||
hi = dividend[1];
|
||||
lo = dividend[0];
|
||||
|
||||
// return highest digit value if the result >= 2**32
|
||||
/* return highest digit value if the result >= 2**32 */
|
||||
if(hi >= divisor || divisor == 0){
|
||||
divisor = 0;
|
||||
*quotient = ~divisor;
|
||||
return;
|
||||
}
|
||||
|
||||
// at this point we know that hi < divisor
|
||||
// just shift and subtract till we're done
|
||||
/* at this point we know that hi < divisor */
|
||||
/* just shift and subtract till we're done */
|
||||
q = 0;
|
||||
x = divisor;
|
||||
for(i = Dbits-1; hi > 0 && i >= 0; i--){
|
||||
|
|
|
|||
|
|
@ -2,9 +2,9 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// division ala knuth, seminumerical algorithms, pp 237-238
|
||||
// the numbers are stored backwards to what knuth expects so j
|
||||
// counts down rather than up.
|
||||
/* division ala knuth, seminumerical algorithms, pp 237-238 */
|
||||
/* the numbers are stored backwards to what knuth expects so j */
|
||||
/* counts down rather than up. */
|
||||
|
||||
void
|
||||
mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
||||
|
|
@ -13,11 +13,11 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
mpdigit qd, *up, *vp, *qp;
|
||||
mpint *u, *v, *t;
|
||||
|
||||
// divide bv zero
|
||||
/* divide bv zero */
|
||||
if(divisor->top == 0)
|
||||
abort();
|
||||
|
||||
// quick check
|
||||
/* quick check */
|
||||
if(mpmagcmp(dividend, divisor) < 0){
|
||||
if(remainder != nil)
|
||||
mpassign(dividend, remainder);
|
||||
|
|
@ -26,8 +26,8 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
return;
|
||||
}
|
||||
|
||||
// D1: shift until divisor, v, has hi bit set (needed to make trial
|
||||
// divisor accurate)
|
||||
/* D1: shift until divisor, v, has hi bit set (needed to make trial */
|
||||
/* divisor accurate) */
|
||||
qd = divisor->p[divisor->top-1];
|
||||
for(s = 0; (qd & mpdighi) == 0; s++)
|
||||
qd <<= 1;
|
||||
|
|
@ -44,13 +44,13 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
vp = v->p+v->top-1;
|
||||
vn = v->top;
|
||||
|
||||
// D1a: make sure high digit of dividend is less than high digit of divisor
|
||||
/* D1a: make sure high digit of dividend is less than high digit of divisor */
|
||||
if(*up >= *vp){
|
||||
*++up = 0;
|
||||
u->top++;
|
||||
}
|
||||
|
||||
// storage for multiplies
|
||||
/* storage for multiplies */
|
||||
t = mpnew(4*Dbits);
|
||||
|
||||
qp = nil;
|
||||
|
|
@ -60,15 +60,15 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
qp = quotient->p+quotient->top-1;
|
||||
}
|
||||
|
||||
// D2, D7: loop on length of dividend
|
||||
/* D2, D7: loop on length of dividend */
|
||||
for(j = u->top; j > vn; j--){
|
||||
|
||||
// D3: calculate trial divisor
|
||||
/* D3: calculate trial divisor */
|
||||
mpdigdiv(up-1, *vp, &qd);
|
||||
|
||||
// D3a: rule out trial divisors 2 greater than real divisor
|
||||
/* D3a: rule out trial divisors 2 greater than real divisor */
|
||||
if(vn > 1) for(;;){
|
||||
memset(t->p, 0, 3*Dbytes); // mpvecdigmuladd adds to what's there
|
||||
memset(t->p, 0, 3*Dbytes); /* mpvecdigmuladd adds to what's there */
|
||||
mpvecdigmuladd(vp-1, 2, qd, t->p);
|
||||
if(mpveccmp(t->p, 3, up-2, 3) > 0)
|
||||
qd--;
|
||||
|
|
@ -76,21 +76,21 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
break;
|
||||
}
|
||||
|
||||
// D4: u -= v*qd << j*Dbits
|
||||
/* D4: u -= v*qd << j*Dbits */
|
||||
sign = mpvecdigmulsub(v->p, vn, qd, up-vn);
|
||||
if(sign < 0){
|
||||
|
||||
// D6: trial divisor was too high, add back borrowed
|
||||
// value and decrease divisor
|
||||
/* D6: trial divisor was too high, add back borrowed */
|
||||
/* value and decrease divisor */
|
||||
mpvecadd(up-vn, vn+1, v->p, vn, up-vn);
|
||||
qd--;
|
||||
}
|
||||
|
||||
// D5: save quotient digit
|
||||
/* D5: save quotient digit */
|
||||
if(qp != nil)
|
||||
*qp-- = qd;
|
||||
|
||||
// push top of u down one
|
||||
/* push top of u down one */
|
||||
u->top--;
|
||||
*up-- = 0;
|
||||
}
|
||||
|
|
@ -101,7 +101,7 @@ mpdiv(mpint *dividend, mpint *divisor, mpint *quotient, mpint *remainder)
|
|||
}
|
||||
|
||||
if(remainder != nil){
|
||||
mpright(u, s, remainder); // u is the remainder shifted
|
||||
mpright(u, s, remainder); /* u is the remainder shifted */
|
||||
remainder->sign = dividend->sign;
|
||||
}
|
||||
|
||||
|
|
|
|||
|
|
@ -1,12 +1,12 @@
|
|||
#include "os.h"
|
||||
#include <mp.h>
|
||||
|
||||
// extended euclid
|
||||
//
|
||||
// For a and b it solves, d = gcd(a,b) and finds x and y s.t.
|
||||
// ax + by = d
|
||||
//
|
||||
// Handbook of Applied Cryptography, Menezes et al, 1997, pg 67
|
||||
/* extended euclid */
|
||||
/* */
|
||||
/* For a and b it solves, d = gcd(a,b) and finds x and y s.t. */
|
||||
/* ax + by = d */
|
||||
/* */
|
||||
/* Handbook of Applied Cryptography, Menezes et al, 1997, pg 67 */
|
||||
|
||||
void
|
||||
mpeuclid(mpint *a, mpint *b, mpint *d, mpint *x, mpint *y)
|
||||
|
|
@ -44,16 +44,16 @@ mpeuclid(mpint *a, mpint *b, mpint *d, mpint *x, mpint *y)
|
|||
r = mpnew(0);
|
||||
|
||||
while(b->top != 0 && b->sign > 0){
|
||||
// q = a/b
|
||||
// r = a mod b
|
||||
/* q = a/b */
|
||||
/* r = a mod b */
|
||||
mpdiv(a, b, q, r);
|
||||
// x0 = x2 - qx1
|
||||
/* x0 = x2 - qx1 */
|
||||
mpmul(q, x1, x0);
|
||||
mpsub(x2, x0, x0);
|
||||
// y0 = y2 - qy1
|
||||
/* y0 = y2 - qy1 */
|
||||
mpmul(q, y1, y0);
|
||||
mpsub(y2, y0, y0);
|
||||
// rotate values
|
||||
/* rotate values */
|
||||
tmp = a;
|
||||
a = b;
|
||||
b = r;
|
||||
|
|
|
|||
|
|
@ -2,17 +2,17 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// res = b**e
|
||||
//
|
||||
// knuth, vol 2, pp 398-400
|
||||
/* res = b**e */
|
||||
/* */
|
||||
/* knuth, vol 2, pp 398-400 */
|
||||
|
||||
enum {
|
||||
Freeb= 0x1,
|
||||
Freee= 0x2,
|
||||
Freem= 0x4,
|
||||
Freem= 0x4
|
||||
};
|
||||
|
||||
//int expdebug;
|
||||
/*int expdebug; */
|
||||
|
||||
void
|
||||
mpexp(mpint *b, mpint *e, mpint *m, mpint *res)
|
||||
|
|
@ -47,7 +47,7 @@ mpexp(mpint *b, mpint *e, mpint *m, mpint *res)
|
|||
tofree |= Freem;
|
||||
}
|
||||
|
||||
// skip first bit
|
||||
/* skip first bit */
|
||||
i = e->top-1;
|
||||
d = e->p[i];
|
||||
for(bit = mpdighi; (bit & d) == 0; bit >>= 1)
|
||||
|
|
|
|||
|
|
@ -3,12 +3,12 @@
|
|||
|
||||
#define iseven(a) (((a)->p[0] & 1) == 0)
|
||||
|
||||
// extended binary gcd
|
||||
//
|
||||
// For a anv b it solves, v = gcd(a,b) and finds x and y s.t.
|
||||
// ax + by = v
|
||||
//
|
||||
// Handbook of Applied Cryptography, Menezes et al, 1997, pg 608.
|
||||
/* extended binary gcd */
|
||||
/* */
|
||||
/* For a anv b it solves, v = gcd(a,b) and finds x and y s.t. */
|
||||
/* ax + by = v */
|
||||
/* */
|
||||
/* Handbook of Applied Cryptography, Menezes et al, 1997, pg 608. */
|
||||
void
|
||||
mpextendedgcd(mpint *a, mpint *b, mpint *v, mpint *x, mpint *y)
|
||||
{
|
||||
|
|
@ -53,7 +53,7 @@ mpextendedgcd(mpint *a, mpint *b, mpint *v, mpint *x, mpint *y)
|
|||
D = mpcopy(mpone);
|
||||
|
||||
for(;;) {
|
||||
// print("%B %B %B %B %B %B\n", u, v, A, B, C, D);
|
||||
/* print("%B %B %B %B %B %B\n", u, v, A, B, C, D); */
|
||||
while(iseven(u)){
|
||||
mpright(u, 1, u);
|
||||
if(!iseven(A) || !iseven(B)) {
|
||||
|
|
@ -64,7 +64,7 @@ mpextendedgcd(mpint *a, mpint *b, mpint *v, mpint *x, mpint *y)
|
|||
mpright(B, 1, B);
|
||||
}
|
||||
|
||||
// print("%B %B %B %B %B %B\n", u, v, A, B, C, D);
|
||||
/* print("%B %B %B %B %B %B\n", u, v, A, B, C, D); */
|
||||
while(iseven(v)){
|
||||
mpright(v, 1, v);
|
||||
if(!iseven(C) || !iseven(D)) {
|
||||
|
|
@ -75,7 +75,7 @@ mpextendedgcd(mpint *a, mpint *b, mpint *v, mpint *x, mpint *y)
|
|||
mpright(D, 1, D);
|
||||
}
|
||||
|
||||
// print("%B %B %B %B %B %B\n", u, v, A, B, C, D);
|
||||
/* print("%B %B %B %B %B %B\n", u, v, A, B, C, D); */
|
||||
if(mpcmp(u, v) >= 0){
|
||||
mpsub(u, v, u);
|
||||
mpsub(A, C, A);
|
||||
|
|
|
|||
|
|
@ -23,7 +23,7 @@ to32(mpint *b, char *buf, int len)
|
|||
uchar *p;
|
||||
int n, rv;
|
||||
|
||||
// leave room for a multiple of 5 buffer size
|
||||
/* leave room for a multiple of 5 buffer size */
|
||||
n = b->top*Dbytes + 5;
|
||||
p = malloc(n);
|
||||
if(p == nil)
|
||||
|
|
@ -32,7 +32,7 @@ to32(mpint *b, char *buf, int len)
|
|||
if(n < 0)
|
||||
return -1;
|
||||
|
||||
// round up buffer size, enc32 only accepts a multiple of 5
|
||||
/* round up buffer size, enc32 only accepts a multiple of 5 */
|
||||
if(n%5)
|
||||
n += 5 - (n%5);
|
||||
rv = enc32(buf, len, p, n);
|
||||
|
|
|
|||
|
|
@ -3,12 +3,12 @@
|
|||
|
||||
#define iseven(a) (((a)->p[0] & 1) == 0)
|
||||
|
||||
// use extended gcd to find the multiplicative inverse
|
||||
// res = b**-1 mod m
|
||||
/* use extended gcd to find the multiplicative inverse */
|
||||
/* res = b**-1 mod m */
|
||||
void
|
||||
mpinvert(mpint *b, mpint *m, mpint *res)
|
||||
{
|
||||
mpint *dc1, *dc2; // don't care
|
||||
mpint *dc1, *dc2; /* don't care */
|
||||
|
||||
dc1 = mpnew(0);
|
||||
dc2 = mpnew(0);
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// res = b << shift
|
||||
/* res = b << shift */
|
||||
void
|
||||
mpleft(mpint *b, int shift, mpint *res)
|
||||
{
|
||||
|
|
@ -15,17 +15,17 @@ mpleft(mpint *b, int shift, mpint *res)
|
|||
return;
|
||||
}
|
||||
|
||||
// a negative left shift is a right shift
|
||||
/* a negative left shift is a right shift */
|
||||
if(shift < 0){
|
||||
mpright(b, -shift, res);
|
||||
return;
|
||||
}
|
||||
|
||||
// b and res may be the same so remember the old top
|
||||
/* b and res may be the same so remember the old top */
|
||||
otop = b->top;
|
||||
|
||||
// shift
|
||||
mpbits(res, otop*Dbits + shift); // overkill
|
||||
/* shift */
|
||||
mpbits(res, otop*Dbits + shift); /* overkill */
|
||||
res->top = DIGITS(otop*Dbits + shift);
|
||||
d = shift/Dbits;
|
||||
l = shift - d*Dbits;
|
||||
|
|
@ -46,7 +46,7 @@ mpleft(mpint *b, int shift, mpint *res)
|
|||
for(i = 0; i < d; i++)
|
||||
res->p[i] = 0;
|
||||
|
||||
// normalize
|
||||
/* normalize */
|
||||
while(res->top > 0 && res->p[res->top-1] == 0)
|
||||
res->top--;
|
||||
}
|
||||
|
|
|
|||
|
|
@ -2,9 +2,9 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// remainder = b mod m
|
||||
//
|
||||
// knuth, vol 2, pp 398-400
|
||||
/* remainder = b mod m */
|
||||
/* */
|
||||
/* knuth, vol 2, pp 398-400 */
|
||||
|
||||
void
|
||||
mpmod(mpint *b, mpint *m, mpint *remainder)
|
||||
|
|
|
|||
|
|
@ -2,18 +2,18 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
//
|
||||
// from knuth's 1969 seminumberical algorithms, pp 233-235 and pp 258-260
|
||||
//
|
||||
// mpvecmul is an assembly language routine that performs the inner
|
||||
// loop.
|
||||
//
|
||||
// the karatsuba trade off is set empiricly by measuring the algs on
|
||||
// a 400 MHz Pentium II.
|
||||
//
|
||||
/* */
|
||||
/* from knuth's 1969 seminumberical algorithms, pp 233-235 and pp 258-260 */
|
||||
/* */
|
||||
/* mpvecmul is an assembly language routine that performs the inner */
|
||||
/* loop. */
|
||||
/* */
|
||||
/* the karatsuba trade off is set empiricly by measuring the algs on */
|
||||
/* a 400 MHz Pentium II. */
|
||||
/* */
|
||||
|
||||
// karatsuba like (see knuth pg 258)
|
||||
// prereq: p is already zeroed
|
||||
/* karatsuba like (see knuth pg 258) */
|
||||
/* prereq: p is already zeroed */
|
||||
static void
|
||||
mpkaratsuba(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
||||
{
|
||||
|
|
@ -21,7 +21,7 @@ mpkaratsuba(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
int u0len, u1len, v0len, v1len, reslen;
|
||||
int sign, n;
|
||||
|
||||
// divide each piece in half
|
||||
/* divide each piece in half */
|
||||
n = alen/2;
|
||||
if(alen&1)
|
||||
n++;
|
||||
|
|
@ -39,7 +39,7 @@ mpkaratsuba(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
v0 = b;
|
||||
v1 = b + v0len;
|
||||
|
||||
// room for the partial products
|
||||
/* room for the partial products */
|
||||
t = mallocz(Dbytes*5*(2*n+1), 1);
|
||||
if(t == nil)
|
||||
sysfatal("mpkaratsuba: %r");
|
||||
|
|
@ -49,7 +49,7 @@ mpkaratsuba(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
res = t + 3*(2*n+1);
|
||||
reslen = 4*n+1;
|
||||
|
||||
// t[0] = (u1-u0)
|
||||
/* t[0] = (u1-u0) */
|
||||
sign = 1;
|
||||
if(mpveccmp(u1, u1len, u0, u0len) < 0){
|
||||
sign = -1;
|
||||
|
|
@ -57,35 +57,35 @@ mpkaratsuba(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
} else
|
||||
mpvecsub(u1, u1len, u0, u1len, u0v0);
|
||||
|
||||
// t[1] = (v0-v1)
|
||||
/* t[1] = (v0-v1) */
|
||||
if(mpveccmp(v0, v0len, v1, v1len) < 0){
|
||||
sign *= -1;
|
||||
mpvecsub(v1, v1len, v0, v1len, u1v1);
|
||||
} else
|
||||
mpvecsub(v0, v0len, v1, v1len, u1v1);
|
||||
|
||||
// t[4:5] = (u1-u0)*(v0-v1)
|
||||
/* t[4:5] = (u1-u0)*(v0-v1) */
|
||||
mpvecmul(u0v0, u0len, u1v1, v0len, diffprod);
|
||||
|
||||
// t[0:1] = u1*v1
|
||||
/* t[0:1] = u1*v1 */
|
||||
memset(t, 0, 2*(2*n+1)*Dbytes);
|
||||
if(v1len > 0)
|
||||
mpvecmul(u1, u1len, v1, v1len, u1v1);
|
||||
|
||||
// t[2:3] = u0v0
|
||||
/* t[2:3] = u0v0 */
|
||||
mpvecmul(u0, u0len, v0, v0len, u0v0);
|
||||
|
||||
// res = u0*v0<<n + u0*v0
|
||||
/* res = u0*v0<<n + u0*v0 */
|
||||
mpvecadd(res, reslen, u0v0, u0len+v0len, res);
|
||||
mpvecadd(res+n, reslen-n, u0v0, u0len+v0len, res+n);
|
||||
|
||||
// res += u1*v1<<n + u1*v1<<2*n
|
||||
/* res += u1*v1<<n + u1*v1<<2*n */
|
||||
if(v1len > 0){
|
||||
mpvecadd(res+n, reslen-n, u1v1, u1len+v1len, res+n);
|
||||
mpvecadd(res+2*n, reslen-2*n, u1v1, u1len+v1len, res+2*n);
|
||||
}
|
||||
|
||||
// res += (u1-u0)*(v0-v1)<<n
|
||||
/* res += (u1-u0)*(v0-v1)<<n */
|
||||
if(sign < 0)
|
||||
mpvecsub(res+n, reslen-n, diffprod, u0len+v0len, res+n);
|
||||
else
|
||||
|
|
@ -104,7 +104,7 @@ mpvecmul(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
mpdigit d;
|
||||
mpdigit *t;
|
||||
|
||||
// both mpvecdigmuladd and karatsuba are fastest when a is the longer vector
|
||||
/* both mpvecdigmuladd and karatsuba are fastest when a is the longer vector */
|
||||
if(alen < blen){
|
||||
i = alen;
|
||||
alen = blen;
|
||||
|
|
@ -119,10 +119,10 @@ mpvecmul(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *p)
|
|||
}
|
||||
|
||||
if(alen >= KARATSUBAMIN && blen > 1){
|
||||
// O(n^1.585)
|
||||
/* O(n^1.585) */
|
||||
mpkaratsuba(a, alen, b, blen, p);
|
||||
} else {
|
||||
// O(n^2)
|
||||
/* O(n^2) */
|
||||
for(i = 0; i < blen; i++){
|
||||
d = b[i];
|
||||
if(d != 0)
|
||||
|
|
|
|||
|
|
@ -22,7 +22,7 @@ mprand(int bits, void (*gen)(uchar*, int), mpint *b)
|
|||
betomp(p, n*Dbytes, b);
|
||||
free(p);
|
||||
|
||||
// make sure we don't give too many bits
|
||||
/* make sure we don't give too many bits */
|
||||
m = bits%Dbits;
|
||||
n--;
|
||||
if(m > 0){
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// res = b >> shift
|
||||
/* res = b >> shift */
|
||||
void
|
||||
mpright(mpint *b, int shift, mpint *res)
|
||||
{
|
||||
|
|
@ -15,7 +15,7 @@ mpright(mpint *b, int shift, mpint *res)
|
|||
return;
|
||||
}
|
||||
|
||||
// a negative right shift is a left shift
|
||||
/* a negative right shift is a left shift */
|
||||
if(shift < 0){
|
||||
mpleft(b, -shift, res);
|
||||
return;
|
||||
|
|
@ -27,13 +27,13 @@ mpright(mpint *b, int shift, mpint *res)
|
|||
r = shift - d*Dbits;
|
||||
l = Dbits - r;
|
||||
|
||||
// shift all the bits out == zero
|
||||
/* shift all the bits out == zero */
|
||||
if(d>=b->top){
|
||||
res->top = 0;
|
||||
return;
|
||||
}
|
||||
|
||||
// special case digit shifts
|
||||
/* special case digit shifts */
|
||||
if(r == 0){
|
||||
for(i = 0; i < b->top-d; i++)
|
||||
res->p[i] = b->p[i+d];
|
||||
|
|
|
|||
|
|
@ -2,14 +2,14 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// diff = abs(b1) - abs(b2), i.e., subtract the magnitudes
|
||||
/* diff = abs(b1) - abs(b2), i.e., subtract the magnitudes */
|
||||
void
|
||||
mpmagsub(mpint *b1, mpint *b2, mpint *diff)
|
||||
{
|
||||
int n, m, sign;
|
||||
mpint *t;
|
||||
|
||||
// get the sizes right
|
||||
/* get the sizes right */
|
||||
if(mpmagcmp(b1, b2) < 0){
|
||||
sign = -1;
|
||||
t = b1;
|
||||
|
|
@ -32,7 +32,7 @@ mpmagsub(mpint *b1, mpint *b2, mpint *diff)
|
|||
mpnorm(diff);
|
||||
}
|
||||
|
||||
// diff = b1 - b2
|
||||
/* diff = b1 - b2 */
|
||||
void
|
||||
mpsub(mpint *b1, mpint *b2, mpint *diff)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -2,9 +2,9 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// convert an mpint into a big endian byte array (most significant byte first)
|
||||
// return number of bytes converted
|
||||
// if p == nil, allocate and result array
|
||||
/* convert an mpint into a big endian byte array (most significant byte first) */
|
||||
/* return number of bytes converted */
|
||||
/* if p == nil, allocate and result array */
|
||||
int
|
||||
mptobe(mpint *b, uchar *p, uint n, uchar **pp)
|
||||
{
|
||||
|
|
@ -22,7 +22,7 @@ mptobe(mpint *b, uchar *p, uint n, uchar **pp)
|
|||
*pp = p;
|
||||
memset(p, 0, n);
|
||||
|
||||
// special case 0
|
||||
/* special case 0 */
|
||||
if(b->top == 0){
|
||||
if(n < 1)
|
||||
return -1;
|
||||
|
|
@ -46,7 +46,7 @@ mptobe(mpint *b, uchar *p, uint n, uchar **pp)
|
|||
}
|
||||
}
|
||||
|
||||
// guarantee at least one byte
|
||||
/* guarantee at least one byte */
|
||||
if(s == p){
|
||||
if(p >= e)
|
||||
return -1;
|
||||
|
|
|
|||
|
|
@ -2,10 +2,10 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// convert an mpint into a little endian byte array (least significant byte first)
|
||||
/* convert an mpint into a little endian byte array (least significant byte first) */
|
||||
|
||||
// return number of bytes converted
|
||||
// if p == nil, allocate and result array
|
||||
/* return number of bytes converted */
|
||||
/* if p == nil, allocate and result array */
|
||||
int
|
||||
mptole(mpint *b, uchar *p, uint n, uchar **pp)
|
||||
{
|
||||
|
|
@ -23,7 +23,7 @@ mptole(mpint *b, uchar *p, uint n, uchar **pp)
|
|||
return -1;
|
||||
memset(p, 0, n);
|
||||
|
||||
// special case 0
|
||||
/* special case 0 */
|
||||
if(b->top == 0){
|
||||
if(n < 1)
|
||||
return -1;
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// prereq: alen >= blen, sum has at least blen+1 digits
|
||||
/* prereq: alen >= blen, sum has at least blen+1 digits */
|
||||
void
|
||||
mpvecadd(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *sum)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -11,19 +11,19 @@ mpdigmul(mpdigit a, mpdigit b, mpdigit *p)
|
|||
mpdigit x, ah, al, bh, bl, p1, p2, p3, p4;
|
||||
int carry;
|
||||
|
||||
// half digits
|
||||
/* half digits */
|
||||
ah = HI(a);
|
||||
al = LO(a);
|
||||
bh = HI(b);
|
||||
bl = LO(b);
|
||||
|
||||
// partial products
|
||||
/* partial products */
|
||||
p1 = ah*bl;
|
||||
p2 = bh*al;
|
||||
p3 = bl*al;
|
||||
p4 = ah*bh;
|
||||
|
||||
// p = ((p1+p2)<<(Dbits/2)) + (p4<<Dbits) + p3
|
||||
/* p = ((p1+p2)<<(Dbits/2)) + (p4<<Dbits) + p3 */
|
||||
carry = 0;
|
||||
x = p1<<(Dbits/2);
|
||||
p3 += x;
|
||||
|
|
@ -33,12 +33,12 @@ mpdigmul(mpdigit a, mpdigit b, mpdigit *p)
|
|||
p3 += x;
|
||||
if(p3 < x)
|
||||
carry++;
|
||||
p4 += carry + HI(p1) + HI(p2); // can't carry out of the high digit
|
||||
p4 += carry + HI(p1) + HI(p2); /* can't carry out of the high digit */
|
||||
p[0] = p3;
|
||||
p[1] = p4;
|
||||
}
|
||||
|
||||
// prereq: p must have room for n+1 digits
|
||||
/* prereq: p must have room for n+1 digits */
|
||||
void
|
||||
mpvecdigmuladd(mpdigit *b, int n, mpdigit m, mpdigit *p)
|
||||
{
|
||||
|
|
@ -66,7 +66,7 @@ mpvecdigmuladd(mpdigit *b, int n, mpdigit m, mpdigit *p)
|
|||
*p = part[1] + carry;
|
||||
}
|
||||
|
||||
// prereq: p must have room for n+1 digits
|
||||
/* prereq: p must have room for n+1 digits */
|
||||
int
|
||||
mpvecdigmulsub(mpdigit *b, int n, mpdigit m, mpdigit *p)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -2,7 +2,7 @@
|
|||
#include <mp.h>
|
||||
#include "dat.h"
|
||||
|
||||
// prereq: a >= b, alen >= blen, diff has at least alen digits
|
||||
/* prereq: a >= b, alen >= blen, diff has at least alen digits */
|
||||
void
|
||||
mpvecsub(mpdigit *a, int alen, mpdigit *b, int blen, mpdigit *diff)
|
||||
{
|
||||
|
|
|
|||
|
|
@ -84,7 +84,7 @@ from10(char *a, mpint *b)
|
|||
|
||||
b->top = 0;
|
||||
for(;;){
|
||||
// do a billion at a time in native arithmetic
|
||||
/* do a billion at a time in native arithmetic */
|
||||
x = 0;
|
||||
for(i = 0; i < 9; i++){
|
||||
y = tab.t10[*(uchar*)a];
|
||||
|
|
@ -97,7 +97,7 @@ from10(char *a, mpint *b)
|
|||
if(i == 0)
|
||||
break;
|
||||
|
||||
// accumulate into mpint
|
||||
/* accumulate into mpint */
|
||||
uitomp(mppow10[i], pow);
|
||||
uitomp(x, r);
|
||||
mpmul(b, pow, b);
|
||||
|
|
@ -191,7 +191,7 @@ strtomp(char *a, char **pp, int base, mpint *b)
|
|||
break;
|
||||
}
|
||||
|
||||
// if no characters parsed, there wasn't a number to convert
|
||||
/* if no characters parsed, there wasn't a number to convert */
|
||||
if(e == a)
|
||||
return nil;
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue