Initial check in of big integer library, v2010.04.30
[pdfium.git] / third_party / bigint / BigIntegerAlgorithms.cc
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 invariant.
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 invariant.
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 }