| OLD | NEW |
| (Empty) |
| 1 #include "BigIntegerAlgorithms.hh" | |
| 2 | |
| 3 BigUnsigned gcd(BigUnsigned a, BigUnsigned b) { | |
| 4 BigUnsigned trash; | |
| 5 // Neat in-place alternating technique. | |
| 6 for (;;) { | |
| 7 if (b.isZero()) | |
| 8 return a; | |
| 9 a.divideWithRemainder(b, trash); | |
| 10 if (a.isZero()) | |
| 11 return b; | |
| 12 b.divideWithRemainder(a, trash); | |
| 13 } | |
| 14 } | |
| 15 | |
| 16 void extendedEuclidean(BigInteger m, BigInteger n, | |
| 17 BigInteger &g, BigInteger &r, BigInteger &s) { | |
| 18 if (&g == &r || &g == &s || &r == &s) | |
| 19 throw "BigInteger extendedEuclidean: Outputs are aliased"; | |
| 20 BigInteger r1(1), s1(0), r2(0), s2(1), q; | |
| 21 /* Invariants: | |
| 22 * r1*m(orig) + s1*n(orig) == m(current) | |
| 23 * r2*m(orig) + s2*n(orig) == n(current) */ | |
| 24 for (;;) { | |
| 25 if (n.isZero()) { | |
| 26 r = r1; s = s1; g = m; | |
| 27 return; | |
| 28 } | |
| 29 // Subtract q times the second invariant from the first invarian
t. | |
| 30 m.divideWithRemainder(n, q); | |
| 31 r1 -= q*r2; s1 -= q*s2; | |
| 32 | |
| 33 if (m.isZero()) { | |
| 34 r = r2; s = s2; g = n; | |
| 35 return; | |
| 36 } | |
| 37 // Subtract q times the first invariant from the second invarian
t. | |
| 38 n.divideWithRemainder(m, q); | |
| 39 r2 -= q*r1; s2 -= q*s1; | |
| 40 } | |
| 41 } | |
| 42 | |
| 43 BigUnsigned modinv(const BigInteger &x, const BigUnsigned &n) { | |
| 44 BigInteger g, r, s; | |
| 45 extendedEuclidean(x, n, g, r, s); | |
| 46 if (g == 1) | |
| 47 // r*x + s*n == 1, so r*x === 1 (mod n), so r is the answer. | |
| 48 return (r % n).getMagnitude(); // (r % n) will be nonnegative | |
| 49 else | |
| 50 throw "BigInteger modinv: x and n have a common factor"; | |
| 51 } | |
| 52 | |
| 53 BigUnsigned modexp(const BigInteger &base, const BigUnsigned &exponent, | |
| 54 const BigUnsigned &modulus) { | |
| 55 BigUnsigned ans = 1, base2 = (base % modulus).getMagnitude(); | |
| 56 BigUnsigned::Index i = exponent.bitLength(); | |
| 57 // For each bit of the exponent, most to least significant... | |
| 58 while (i > 0) { | |
| 59 i--; | |
| 60 // Square. | |
| 61 ans *= ans; | |
| 62 ans %= modulus; | |
| 63 // And multiply if the bit is a 1. | |
| 64 if (exponent.getBit(i)) { | |
| 65 ans *= base2; | |
| 66 ans %= modulus; | |
| 67 } | |
| 68 } | |
| 69 return ans; | |
| 70 } | |
| OLD | NEW |