-
Notifications
You must be signed in to change notification settings - Fork 522
/
Euclid.java
143 lines (124 loc) · 4.46 KB
/
Euclid.java
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
package numbertheory;
import java.math.BigInteger;
import java.util.Arrays;
import java.util.Random;
public class Euclid {
public static int gcd(int a, int b) {
return b == 0 ? Math.abs(a) : gcd(b, a % b);
}
public static int gcd2(int a, int b) {
while (b != 0) {
int t = b;
b = a % b;
a = t;
}
return Math.abs(a);
}
public static int lcm(int a, int b) {
return Math.abs(a / gcd(a, b) * b);
}
// returns { gcd(a,b), x, y } such that gcd(a,b) = a*x + b*y
public static long[] euclid(long a, long b) {
long x = 1, y = 0, x1 = 0, y1 = 1;
// invariant: a=a_orig*x+b_orig*y, b=a_orig*x1+b_orig*y1
while (b != 0) {
long q = a / b;
long _x1 = x1;
long _y1 = y1;
long _b = b;
x1 = x - q * x1;
y1 = y - q * y1;
b = a - q * b;
x = _x1;
y = _y1;
a = _b;
}
return a > 0 ? new long[] {a, x, y} : new long[] {-a, -x, -y};
}
public static long[] euclid2(long a, long b) {
if (b == 0)
return a > 0 ? new long[] {a, 1, 0} : new long[] {-a, -1, 0};
long[] r = euclid2(b, a % b);
return new long[] {r[0], r[2], r[1] - a / b * r[2]};
}
public static int mod(int a, int m) {
a %= m;
return a >= 0 ? a : a + m;
}
public static int mod(long a, int m) {
a %= m;
return (int) (a >= 0 ? a : a + m);
}
// precondition: m > 1 && gcd(a, m) = 1
public static int modInverse(int a, int m) {
a = mod(a, m);
return a == 0 ? 0 : mod((int) ((1 - (long) modInverse(m, a) * m) / a), m);
}
// precondition: m > 0 && gcd(a, m) = 1
public static int modInverse2(int a, int m) {
return mod((int) euclid(a, m)[1], m);
}
// precondition: p is prime
public static int[] generateInverses(int p) {
int[] res = new int[p];
res[1] = 1;
for (int i = 2; i < p; ++i) res[i] = (int) ((long) (p - p / i) * res[p % i] % p);
return res;
}
// returns x ≡ a[i] (mod p[i]), where gcd(p[i], p[j]) == 1
public static BigInteger garnerRestore(int[] a, int[] p) {
int[] x = a.clone();
for (int i = 0; i < x.length; ++i)
for (int j = 0; j < i; ++j)
x[i] = mod(
BigInteger.valueOf(p[j]).modInverse(BigInteger.valueOf(p[i])).longValue() * (x[i] - x[j]), p[i]);
BigInteger res = BigInteger.valueOf(x[0]);
BigInteger m = BigInteger.ONE;
for (int i = 1; i < x.length; i++) {
m = m.multiply(BigInteger.valueOf(p[i - 1]));
res = res.add(m.multiply(BigInteger.valueOf(x[i])));
}
return res;
}
// returns x ≡ a[i] (mod p[i]), where gcd(p[i], p[j]) == 1
public static int simpleRestore(int[] a, int[] p) {
int res = 0;
for (int i = 0, m = 1; i < a.length; i++, m *= p[i])
while (res % p[i] != a[i]) res += m;
return res;
}
// Usage example
public static void main(String[] args) {
Random rnd = new Random(1);
for (int steps = 0; steps < 10000; steps++) {
int a = rnd.nextInt(20) - 10;
int b = rnd.nextInt(20) - 10;
BigInteger xa = BigInteger.valueOf(a);
BigInteger xb = BigInteger.valueOf(b);
long gcd1 = gcd(a, b);
long gcd2 = gcd2(a, b);
long gcd = xa.gcd(xb).longValue();
long[] euclid1 = euclid(a, b);
long[] euclid2 = euclid2(a, b);
int inv1 = 0;
int inv2 = 0;
int inv = 0;
if (gcd == 1 && b > 0) {
inv1 = modInverse(a, b);
inv2 = modInverse2(a, b);
inv = xa.modInverse(xb).intValue();
}
if (gcd1 != gcd || gcd2 != gcd || !Arrays.equals(euclid1, euclid2) || euclid1[0] != gcd || inv1 != inv
|| inv2 != inv) {
System.err.println(a + " " + b);
}
}
long a = 6;
long b = 9;
long[] res = euclid(a, b);
System.out.println(res[1] + " * (" + a + ") "
+ " + " + res[2] + " * (" + b + ") = gcd(" + a + "," + b + ") = " + res[0]);
System.out.println(Arrays.toString(generateInverses(7)));
System.out.println(garnerRestore(new int[] {200_000_125, 300_000_333}, new int[] {1000_000_007, 1000_000_009}));
}
}