mpeuclid.c (1326B)
1 #include "os.h" 2 #include <mp.h> 3 4 /* extended euclid */ 5 /* */ 6 /* For a and b it solves, d = gcd(a,b) and finds x and y s.t. */ 7 /* ax + by = d */ 8 /* */ 9 /* Handbook of Applied Cryptography, Menezes et al, 1997, pg 67 */ 10 11 void 12 mpeuclid(mpint *a, mpint *b, mpint *d, mpint *x, mpint *y) 13 { 14 mpint *tmp, *x0, *x1, *x2, *y0, *y1, *y2, *q, *r; 15 16 if(a->sign<0 || b->sign<0) 17 sysfatal("mpeuclid: negative arg"); 18 19 if(mpcmp(a, b) < 0){ 20 tmp = a; 21 a = b; 22 b = tmp; 23 tmp = x; 24 x = y; 25 y = tmp; 26 } 27 28 if(b->top == 0){ 29 mpassign(a, d); 30 mpassign(mpone, x); 31 mpassign(mpzero, y); 32 return; 33 } 34 35 a = mpcopy(a); 36 b = mpcopy(b); 37 x0 = mpnew(0); 38 x1 = mpcopy(mpzero); 39 x2 = mpcopy(mpone); 40 y0 = mpnew(0); 41 y1 = mpcopy(mpone); 42 y2 = mpcopy(mpzero); 43 q = mpnew(0); 44 r = mpnew(0); 45 46 while(b->top != 0 && b->sign > 0){ 47 /* q = a/b */ 48 /* r = a mod b */ 49 mpdiv(a, b, q, r); 50 /* x0 = x2 - qx1 */ 51 mpmul(q, x1, x0); 52 mpsub(x2, x0, x0); 53 /* y0 = y2 - qy1 */ 54 mpmul(q, y1, y0); 55 mpsub(y2, y0, y0); 56 /* rotate values */ 57 tmp = a; 58 a = b; 59 b = r; 60 r = tmp; 61 tmp = x2; 62 x2 = x1; 63 x1 = x0; 64 x0 = tmp; 65 tmp = y2; 66 y2 = y1; 67 y1 = y0; 68 y0 = tmp; 69 } 70 71 mpassign(a, d); 72 mpassign(x2, x); 73 mpassign(y2, y); 74 75 mpfree(x0); 76 mpfree(x1); 77 mpfree(x2); 78 mpfree(y0); 79 mpfree(y1); 80 mpfree(y2); 81 mpfree(q); 82 mpfree(r); 83 mpfree(a); 84 mpfree(b); 85 }