結果
問題 | No.950 行列累乗 |
ユーザー | 37zigen |
提出日時 | 2019-12-14 18:08:19 |
言語 | Java21 (openjdk 21) |
結果 |
WA
|
実行時間 | - |
コード長 | 20,137 bytes |
コンパイル時間 | 3,250 ms |
コンパイル使用メモリ | 94,908 KB |
実行使用メモリ | 94,336 KB |
最終ジャッジ日時 | 2024-06-28 03:18:51 |
合計ジャッジ時間 | 33,147 ms |
ジャッジサーバーID (参考情報) |
judge3 / judge5 |
(要ログイン)
テストケース
テストケース表示入力 | 結果 | 実行時間 実行使用メモリ |
---|---|---|
testcase_00 | WA | - |
testcase_01 | AC | 258 ms
48,728 KB |
testcase_02 | WA | - |
testcase_03 | AC | 165 ms
42,608 KB |
testcase_04 | AC | 135 ms
41,604 KB |
testcase_05 | WA | - |
testcase_06 | AC | 220 ms
43,948 KB |
testcase_07 | WA | - |
testcase_08 | WA | - |
testcase_09 | AC | 216 ms
43,544 KB |
testcase_10 | WA | - |
testcase_11 | WA | - |
testcase_12 | AC | 172 ms
42,528 KB |
testcase_13 | WA | - |
testcase_14 | WA | - |
testcase_15 | AC | 172 ms
42,380 KB |
testcase_16 | WA | - |
testcase_17 | WA | - |
testcase_18 | AC | 173 ms
42,396 KB |
testcase_19 | WA | - |
testcase_20 | WA | - |
testcase_21 | AC | 214 ms
43,632 KB |
testcase_22 | WA | - |
testcase_23 | WA | - |
testcase_24 | AC | 213 ms
43,708 KB |
testcase_25 | WA | - |
testcase_26 | WA | - |
testcase_27 | AC | 172 ms
42,692 KB |
testcase_28 | WA | - |
testcase_29 | WA | - |
testcase_30 | AC | 218 ms
43,796 KB |
testcase_31 | WA | - |
testcase_32 | WA | - |
testcase_33 | AC | 164 ms
42,808 KB |
testcase_34 | WA | - |
testcase_35 | WA | - |
testcase_36 | AC | 171 ms
42,784 KB |
testcase_37 | WA | - |
testcase_38 | WA | - |
testcase_39 | AC | 171 ms
42,612 KB |
testcase_40 | WA | - |
testcase_41 | WA | - |
testcase_42 | AC | 215 ms
43,880 KB |
testcase_43 | TLE | - |
testcase_44 | TLE | - |
testcase_45 | AC | 216 ms
43,836 KB |
testcase_46 | TLE | - |
testcase_47 | AC | 216 ms
43,948 KB |
testcase_48 | TLE | - |
testcase_49 | TLE | - |
testcase_50 | WA | - |
testcase_51 | WA | - |
testcase_52 | WA | - |
testcase_53 | WA | - |
testcase_54 | AC | 170 ms
42,876 KB |
testcase_55 | AC | 136 ms
41,200 KB |
testcase_56 | WA | - |
testcase_57 | AC | 135 ms
41,632 KB |
testcase_58 | AC | 133 ms
41,060 KB |
testcase_59 | AC | 168 ms
42,720 KB |
testcase_60 | AC | 167 ms
42,940 KB |
ソースコード
import java.io.*; import java.util.*; class Main { public static void main(String[] args) { new Main().run(); } long MOD; class K {// quadratic field long a, b, base; public K(long a, long b, long base) { this.a = (a % MOD + MOD) % MOD; this.b = (b % MOD + MOD) % MOD; this.base = (base % MOD + MOD) % MOD; } K add(K o) { return new K((a + o.a) % MOD, (b + o.b) % MOD, base); } K sub(K o) { return new K((a - o.a + MOD) % MOD, (b - o.b + MOD) % MOD, base); } K mul(K o) { return new K((a * o.a % MOD + b * o.b % MOD * base % MOD) % MOD, (a * o.b % MOD + o.a * b % MOD) % MOD, base); } // (a+b√w)/(c+d√w)=(a+b√w)(c-d√w)/(cc-ddw) K div(K o) { long det = o.det(); K ret = this.mul(new K(o.a, (MOD - o.b) % MOD, base)); ret.a = ret.a * inv(det, MOD) % MOD; ret.b = ret.b * inv(det, MOD) % MOD; return ret; } K rev() { return new K(1, 0, base).div(this); } long det() { return (a * a % MOD - b * b % MOD * base % MOD + MOD) % MOD; } boolean equiv(K o) { return a == o.a && b == o.b; } List<Long> toList() { List<Long> list = new ArrayList<Long>(); list.add(a); list.add(b); return list; } } long[][] pow(long[][] a, long n) { long[][] ret = new long[2][2]; ret[0][0] = ret[1][1] = 1; for (; n > 0; n >>= 1, a = mul(a, a)) { if (n % 2 == 1) ret = mul(ret, a); } return ret; } K[][] mul(K coe, K[][] a) { K[][] ret = new K[a.length][a[0].length]; for (int i = 0; i < a.length; ++i) for (int j = 0; j < a[i].length; ++j) ret[i][j] = a[i][j].mul(coe); return ret; } long[][] mul(long[][] a, long[][] b) { long[][] ret = new long[a.length][b[0].length]; for (int i = 0; i < a.length; ++i) { for (int j = 0; j < b[i].length; ++j) { for (int k = 0; k < a[i].length; ++k) { ret[i][j] += a[i][k] * b[k][j] % MOD; ret[i][j] = (ret[i][j] % MOD + MOD) % MOD; } } } return ret; } K[][] mul(K[][] a, K[][] b) { K[][] ret = new K[a.length][b[0].length]; for (int i = 0; i < ret.length; ++i) { for (int j = 0; j < ret[i].length; ++j) { ret[i][j] = new K(0, 0, a[0][0].base); } } for (int i = 0; i < a.length; ++i) { for (int j = 0; j < b[i].length; ++j) { for (int k = 0; k < a[i].length; ++k) { ret[i][j] = ret[i][j].add(a[i][k].mul(b[k][j])); } } } return ret; } long det(long[][] a) { return (a[0][0] * a[1][1] % MOD - a[0][1] * a[1][0] % MOD + MOD) % MOD; } K det(K[][] a) { return (a[0][0].mul(a[1][1])).sub(a[0][1].mul(a[1][0])); } long[][] rndmat(long p) { long[][] ret = new long[2][2]; Random rnd = new Random(); for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) ret[i][j] = rnd.nextInt((int) 5) % MOD; return ret; } void run() { // long p = (long) 5; // MOD = p; // for (int i = 0; i < 10000; ++i) { // long[][] a = rndmat(p); // long[][] b = rndmat(p); // long ans0 = exact(a, b, p); // long ans1 = solve(a, b, p); // if (ans0 != ans1) { // tr(a, b); // System.err.println("exact=" + ans0); // System.err.println("fast =" + ans1); // throw new AssertionError(); // } // } Scanner sc = new Scanner(System.in); long p = sc.nextLong(); long[][] a = new long[2][2]; long[][] b = new long[2][2]; for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) a[i][j] = sc.nextLong(); for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) b[i][j] = sc.nextLong(); // System.out.println(exact(a, b, p)); System.out.println(solve(a, b, p)); } long solve2(long[][] a, long[][] b, long p) { if (equiv(a, b)) return 1; if (p == 2) return exact(a, b, p); Scanner sc = new Scanner(System.in); MOD = p; long w = (a[0][0] * a[0][0] % MOD - 2 * a[0][0] * a[1][1] % MOD + a[1][1] * a[1][1] % MOD + 4 * a[0][1] * a[1][0] % MOD) % MOD; K eigen1 = new K(inv(2, MOD) * (a[0][0] + a[1][1]) % MOD, MOD - inv(2, MOD), w); K eigen2 = new K(inv(2, MOD) * (a[0][0] + a[1][1]) % MOD, inv(2, MOD), w); K[][] S = new K[2][2]; K[][] J = new K[2][2]; K[][] eigenvec1 = new K[2][1]; K[][] eigenvec2 = new K[2][1]; K[][] ak = new K[][] { { new K(a[0][0], 0, w), new K(a[0][1], 0, w) }, { new K(a[1][0], 0, w), new K(a[1][1], 0, w) } }; K[][] bk = new K[][] { { new K(b[0][0], 0, w), new K(b[0][1], 0, w) }, { new K(b[1][0], 0, w), new K(b[1][1], 0, w) } }; System.err.println("eigen1=" + eigen1.a + "+" + eigen1.b + "sqrt" + eigen1.base); System.err.println("eigen2=" + eigen2.b + "+" + eigen2.b + "sqrt" + eigen2.base); if (!eigen1.equiv(eigen2)) { eigenvec1 = eigenvec(ak, eigen1); eigenvec2 = eigenvec(ak, eigen2); S = new K[][] { { eigenvec1[0][0], eigenvec2[0][0] }, { eigenvec1[1][0], eigenvec2[1][0] } }; J = new K[][] { { eigen1, new K(0, 0, w) }, { new K(0, 0, w), eigen2 } }; } else { throw new AssertionError(); } bk = mul(invmat(S), mul(bk, S)); System.err.println("S^(-1)bS=" + bk[0][0].a + "+" + bk[0][0].b + " " + bk[0][1].a + "+" + bk[0][1].b); System.err.println(" " + bk[1][0].a + "+" + bk[1][0].b + " " + bk[1][1].a + "+" + bk[1][1].b); System.err.println("J=" + J[0][0].a + "+" + J[0][0].b + " " + J[0][1].a + "+" + J[0][1].b); System.err.println(" " + J[1][0].a + "+" + J[1][0].b + " " + J[1][1].a + "+" + J[1][1].b); System.err.println("SAS^(-1)=" + mul(mul(bk, S), invmat(S))[0][0].a + "+" + mul(mul(bk, S), invmat(S))[0][0].b + " " + mul(mul(bk, S), invmat(S))[0][1].a + "+" + J[0][1].b); System.err.println(" " + mul(mul(bk, S), invmat(S))[1][0].a + "+" + mul(mul(bk, S), invmat(S))[1][0].b + " " + mul(mul(bk, S), invmat(S))[1][1].a + "+" + mul(mul(bk, S), invmat(S))[1][1].b); if (!bk[0][1].equiv(new K(0, 0, w)) || !bk[1][0].equiv(new K(0, 0, w))) { return -1; } long sol1 = discretelog(J[0][0], bk[0][0]); long sol2 = discretelog(J[1][1], bk[1][1]); System.err.println("sol1=" + sol1 + " sol2" + sol2); if (sol1 == -1 || sol2 == -1) return -1; long ord1 = ord(J[0][0], MOD); long ord2 = ord(J[1][1], MOD); long ans1 = sol1; long ans2 = sol2; if (Math.abs(ans1 - ans2) % gcd(ord1, ord2) != 0) return -1; while (ans1 != ans2) { if (ans1 < ans2) { ans1 += (ans2 - ans1 + ord1 - 1) / ord1 * ord1; } else { ans2 += (ans1 - ans2 + ord2 - 1) / ord2 * ord2; } } while (!equiv(pow(J, ans1), bk)) { ans1 += lcm(ord1, ord2); ans2 += lcm(ord1, ord2); } return ans1; } long solve(long[][] a, long[][] b, long p) { MOD = p; if (equiv(a, b)) return 1; else if (a[0][0] == 0 && a[0][1] == 0 && a[1][0] == 0 && a[1][1] == 0) return -1; else if (a[0][0] == a[1][1] && a[0][0] != 0 && a[0][1] == 0 && a[1][0] == 0) { if (b[0][0] == b[1][1] && b[0][0] != 0 && b[0][1] == 0 && b[1][0] == 0) { return discretelog(a[0][0], b[0][0]); } else { return -1; } } if (p == 2) return exact(a, b, p); Scanner sc = new Scanner(System.in); System.err.println("p=" + p); System.err.println("r=" + rank(a)); System.err.println("a=" + a[0][0] + " " + a[0][1]); System.err.println(" " + a[1][0] + " " + a[1][1]); long w = (a[0][0] * a[0][0] % MOD - 2 * a[0][0] % MOD * a[1][1] % MOD + a[1][1] * a[1][1] % MOD + 4 * a[0][1] % MOD * a[1][0] % MOD) % MOD; if (!isQuadraticResidue(w)) { System.err.println("sqrt[w] is not quadratic residue."); return solve2(a, b, p); } System.err.println("w=" + w); w = sqrt((w + MOD) % MOD, p); long det0 = det(a); long det1 = det(b); if (det0 == 0 && det1 != 0) return -1; long eigen1 = inv(2, MOD) * (MOD - w + a[0][0] + a[1][1]) % MOD; long eigen2 = inv(2, MOD) * (MOD + w + a[0][0] + a[1][1]) % MOD; System.err.println("eigen1=" + eigen1 + " eigen2=" + eigen2); if (eigen1 == 0 && eigen2 == 0) return exact(a, b, p); long[][] S = new long[2][2]; long[][] J = new long[2][2]; long[][] eigenvec1 = new long[2][1]; long[][] eigenvec2 = new long[2][1]; if (eigen1 != eigen2) { eigenvec1 = eigenvec(a, eigen1); eigenvec2 = eigenvec(a, eigen2); S = new long[][] { { eigenvec1[0][0], eigenvec2[0][0] }, { eigenvec1[1][0], eigenvec2[1][0] } }; J = new long[][] { { eigen1, 0 }, { 0, eigen2 } }; } else { // diagonalization is impossible. System.err.println("diagonalization is impossible."); eigenvec1 = eigenvec(a, eigen1); long[][] tmp = new long[][] { { (a[0][0] - eigen1 + MOD) % MOD, a[0][1] }, { a[1][0], (a[1][1] - eigen1 + MOD) % MOD } }; eigenvec2 = solveLinearEq(tmp, eigenvec1); S = new long[][] { { eigenvec1[0][0], eigenvec2[0][0] }, { eigenvec1[1][0], eigenvec2[1][0] } }; J = new long[][] { { eigen1, 1 }, { 0, eigen2 } }; } System.err.println("eigenvec1=" + eigenvec1[0][0] + " " + eigenvec1[1][0]); System.err.println("eigenvec2=" + eigenvec2[0][0] + " " + eigenvec2[1][0]); System.err.println("J=" + J[0][0] + " " + J[0][1]); System.err.println(" " + J[1][0] + " " + J[1][1]); b = mul(invmat(S), mul(b, S)); System.err.println("S^(-1)BS=" + b[0][0] + " " + b[0][1]); System.err.println(" " + b[1][0] + " " + b[1][1]); if (b[1][0] != 0 && J[1][0] == 0) { return -1; } if (eigen1 != eigen2 && eigen1 != 0 && eigen2 != 0) { if (b[0][1] != 0 || b[1][0] != 0) return -1; } if ((b[0][0] != 0 && J[0][0] == 0) || (b[1][1] != 0 && J[1][1] == 0)) return -1; if (J[0][1] == 0 && b[0][1] != 0) return -1; long sol1 = discretelog(J[0][0], b[0][0]); long sol2 = discretelog(J[1][1], b[1][1]); System.err.println("sol1=" + sol1 + " sol2=" + sol2); if (sol1 == -1 || sol2 == -1) return -1; long ord1 = ord(J[0][0], MOD); System.err.println("ord1=" + ord1); long ord2 = ord(J[1][1], MOD); System.err.println("ord2=" + ord2); long ans1 = sol1; long ans2 = sol2; if (Math.abs(ans1 - ans2) % gcd(ord1, ord2) != 0) return -1; while (ans1 != ans2) { if (ans1 < ans2) { ans1 += (ans2 - ans1 + ord1 - 1) / ord1 * ord1; } else { ans2 += (ans1 - ans2 + ord2 - 1) / ord2 * ord2; } } System.err.println("J^ans1=" + pow(J, ans1)[0][0] + " " + pow(J, ans1)[0][1]); System.err.println(" " + pow(J, ans1)[1][0] + " " + pow(J, ans1)[1][1]); long[][] powJ = pow(J, ans1); if (powJ[0][1] != 0) { long N = (b[0][1] * inv(pow(J[0][0], ans1 - 1), MOD) % MOD - ans1) % MOD * inv(lcm(ord1, ord2), MOD) % MOD; N = (N + MOD) % MOD; return ans1 + N * lcm(ord1, ord2); } else { return ans1; } } boolean equiv(long[][] a, long[][] b) { boolean ret = true; if (a[0].length != b[0].length || a[1].length != b[1].length) throw new AssertionError(); for (int i = 0; i < a.length; ++i) for (int j = 0; j < a[0].length; ++j) ret &= a[i][j] == b[i][j]; return ret; } long[][] eigenvec(long[][] a, long eigenval) { long[][] tmp = new long[][] { { (a[0][0] - eigenval + MOD) % MOD, a[0][1] }, { a[1][0], (a[1][1] - eigenval + MOD) % MOD } }; long[][] ret = new long[2][1]; if (tmp[0][0] != 0 || tmp[0][1] != 0) { ret = new long[][] { { (MOD - tmp[0][1]) % MOD }, { tmp[0][0] } }; } else { ret = new long[][] { { (MOD - tmp[1][1]) % MOD }, { tmp[1][0] } }; } if (!equiv(mul(a, ret), mul(new long[][] { { eigenval, 0 }, { 0, eigenval } }, ret))) { throw new AssertionError(); } return ret; } K[][] eigenvec(K[][] a, K eigenval) { long base = a[0][0].base; K[][] tmp = new K[][] { { a[0][0].sub(eigenval), a[0][1] }, { a[1][0], a[1][1].sub(eigenval) } }; if (!tmp[0][0].equiv(new K(0, 0, base)) || !tmp[0][1].equiv(new K(0, 0, base))) { return new K[][] { { new K(0, 0, base).sub(tmp[0][1]) }, { tmp[0][0] } }; } else { return new K[][] { { new K(0, 0, base).sub(tmp[1][1]) }, { tmp[1][0] } }; } } boolean equiv(K[][] a, K[][] b) { boolean ret = true; if (a[0].length != b[0].length || a[1].length != b[1].length) throw new AssertionError(); for (int i = 0; i < a.length; ++i) for (int j = 0; j < a[0].length; ++j) { ret &= a[i][j].equiv(b[i][j]); } return ret; } long exact(long[][] a, long[][] b, long p) { MOD = p; for (int i = 1; i < p * p; ++i) { long[][] pw_a = pow(a, i); boolean equiv = true; for (int j = 0; j < 2; ++j) for (int k = 0; k < 2; ++k) equiv &= pw_a[j][k] == b[j][k]; if (equiv) return i; } return -1; } long inv(long a, long mod) { return pow(a, mod - 2); } long ord(long a, long p) { if (a == 0 || a == 1) return 1; long ret = p - 1; for (long div = 2; div * div <= p - 1; ++div) { if ((p - 1) % div != 0) continue; if (pow(a, div) == 1) ret = Math.min(ret, div); else if (pow(a, (p - 1) / div) == 1) ret = Math.min(ret, (p - 1) / div); } return ret; } long ord(K a, long p) { long base = a.base; if (a.equiv(new K(0, 0, base)) || a.equiv(new K(1, 0, base))) return 1; long ret = 2 * p - 1;// (a+bw)^2p=a+bw for (long div = 2; div * div <= 2 * p - 1; ++div) { if ((p - 1) % div != 0) continue; if (pow(a, div).equiv(new K(1, 0, base))) ret = Math.min(ret, div); else if (pow(a, (p - 1) / div).equiv(new K(1, 0, base))) ret = Math.min(ret, (p - 1) / div); } return ret; } long pow(long a, long n) { long ret = 1; for (; n > 0; n >>= 1, a = a * a % MOD) { if (n % 2 == 1) ret = ret * a % MOD; } return ret; } K pow(K a, long n) { K ret = new K(1, 0, a.base); for (; n > 0; n >>= 1, a = a.mul(a)) { if (n % 2 == 1) ret = ret.mul(a); } return ret; } K[][] pow(K[][] a, long n) { K[][] ret = new K[2][2]; for (int i = 0; i < 2; ++i) { for (int j = 0; j < 2; ++j) { ret[i][j] = new K(i == j ? 1 : 0, 0, a[0][0].base); } } for (; n > 0; n >>= 1, a = mul(a, a)) { if (n % 2 == 1) ret = mul(ret, a); } return ret; } // return x s.t. a^x = b && x>0 long discretelog(long a, long b) { if (a == 1) { if (b == 1) return 1; else return -1; } else if (a == 0) { if (b == 0) return 1; else return -1; } // a^(um+v) = b // a^v = b a^(-m)^u System.err.println("solve " + a + "^x=" + b); int m = (int) (Math.sqrt(MOD) + 1); long pw = pow(a, m); HashMap<Long, Integer> map = new HashMap<>(); for (int v = m; v >= 1; --v) { map.put(pw, v); pw = pw * inv(a, MOD) % MOD; } long ima = pow(inv(a, MOD), m); long ipw = 1; for (int i = 0; i <= m; ++i) { boolean one = b * ipw % MOD == 1; if (map.containsKey(b * ipw % MOD)) { long ret = i * m + (one && i == 0 ? map.get(b * ipw % MOD) : 0); if (ret != 0) return ret; } ipw = ipw * ima % MOD; } return -1; } long exact_discretelog(K a, K b) { for (int i = 1; i <= 2 * MOD * MOD - 1; ++i) { if (pow(a, i).equiv(b)) { return i; } } return -1; } // return x s.t. a^x = b && x>0 long discretelog(K a, K b) { System.err .println("discretelog " + a.a + "+" + a.b + "sqrt" + b.base + " " + b.a + "+" + b.b + "sqrt" + b.base); // System.err.println("exact_discretlog " + exact_discretelog(a, b)); long base = a.base; if (a.equiv(new K(1, 0, base))) { if (b.equiv(new K(1, 0, base))) return 1; else return -1; } else if (a.equiv(new K(0, 0, base))) { if (b.equiv(new K(0, 0, base))) return 1; else return -1; } System.err.println("det(a)=" + a.det() + " det(b)=" + b.det()); long e = discretelog(a.det(), b.det()); long ord = ord(a.det(), MOD); if (pow(a, e).equiv(b)) return e; System.err.println("ord(det(a))=" + ord); System.err.println("det(a)^" + e + "=det(b)"); int m = (int) (Math.sqrt(2 * MOD + 2)); HashMap<List<Long>, Integer> map = new HashMap<>(); K oa = pow(a, ord).rev(); K pw = pow(a, ord * (m - 1)); for (int v = m - 1; v >= 0; --v) { map.put(pw.toList(), v); pw = pw.mul(oa); } K ima = pow(a.rev(), ord * m); K ipw = new K(1, 0, base); for (int i = 0; i <= m; ++i) { if (map.containsKey(b.mul(ipw).mul(pow(a.rev(), e)).toList())) { long get = map.get(b.mul(ipw).mul(pow(a.rev(), e)).toList()); long ret = ord * (i * m + get) + e; if (ret != 0) { System.err.println("fast_discretlog" + ret); return ret; } } ipw = ipw.mul(ima); } return -1; } long[][] invmat(long[][] a) { if (det(a) == 0) throw new AssertionError(); long[][] ret = new long[2][2]; ret[0][0] = a[1][1]; ret[1][1] = a[0][0]; ret[0][1] = (MOD - a[0][1]) % MOD; ret[1][0] = (MOD - a[1][0]) % MOD; for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) ret[i][j] = ret[i][j] * inv(det(a), MOD) % MOD; return ret; } K[][] invmat(K[][] a) { long base = a[0][0].base; if (det(a).equiv(new K(0, 0, base))) throw new AssertionError(); K[][] ret = new K[2][2]; ret[0][0] = new K(a[1][1].a, a[1][1].b, base); ret[1][1] = new K(a[0][0].a, a[0][0].b, base); ret[0][1] = new K(-a[0][1].a, -a[0][1].b, base); ret[1][0] = new K(-a[1][0].a, -a[1][0].b, base); for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) ret[i][j] = ret[i][j].mul(det(a).rev()); return ret; } long gcd(long a, long b) { if (a > b) return gcd(b, a); if (a == 0) return b; return gcd(a, b % a); } long lcm(long a, long b) { return a / gcd(a, b) * b; } long rank(long[][] a) { long[][] tmp = new long[2][2]; for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) tmp[i][j] = a[i][j]; int p = 0; for (int i = 0; i < 2; ++i) { boolean find = false; for (int j = p; j < 2; ++j) { if (tmp[j][i] != 0) { if (p != j) { for (int k = 0; k < 2; ++k) { tmp[p][k] ^= tmp[j][k]; tmp[j][k] ^= tmp[p][k]; tmp[p][k] ^= tmp[j][k]; } } find = true; break; } } if (find) { for (int j = 0; j < 2; ++j) { if (p == j) continue; long pre = tmp[j][i]; for (int k = 0; k < 2; ++k) { tmp[j][k] -= tmp[p][k] * inv(tmp[p][i], MOD) % MOD * pre % MOD; tmp[j][k] = (tmp[j][k] + MOD) % MOD; } } ++p; } } return p; } long[][] solveLinearEq(long[][] A, long[][] v) { long[][] x = new long[2][1]; long[][] A_ = new long[2][2]; long[][] v_ = new long[2][1]; for (int i = 0; i < 2; ++i) for (int j = 0; j < 2; ++j) A_[i][j] = A[i][j]; for (int i = 0; i < 2; ++i) v_[i][0] = v[i][0]; int p = 0; for (int i = 0; i < 2; ++i) { boolean find = false; for (int j = p; j < 2; ++j) { if (A_[j][i] != 0) { if (p != j) { for (int k = 0; k < 2; ++k) { A_[p][k] ^= A_[j][k]; A_[j][k] ^= A_[p][k]; A_[p][k] ^= A_[j][k]; } v_[p][0] ^= v_[j][0]; v_[j][0] ^= v_[p][0]; v_[p][0] ^= v_[j][0]; } find = true; break; } } if (find) { for (int j = 0; j < 2; ++j) { if (p == j) continue; long coe = A_[j][i]; long normalize = inv(A_[p][i], MOD); for (int k = 0; k < 2; ++k) A_[p][k] = A_[p][k] * normalize % MOD; v_[p][0] = v_[p][0] * normalize % MOD; for (int k = 0; k < 2; ++k) { A_[j][k] -= A_[p][k] * coe % MOD; A_[j][k] = (A_[j][k] + MOD) % MOD; } v_[j][0] = (v_[j][0] - v_[p][0] * coe % MOD + MOD) % MOD; } x[i][0] = v_[p][0]; ++p; } } if (!equiv(mul(A, x), v)) { tr(A); tr(x); tr(v); throw new AssertionError(); } return x; } boolean isQuadraticResidue(long a) { return a == 0 || pow(a, (MOD - 1) / 2) == 1; } long sqrt(long a, long p) { if (a == 0) return 0; int b = 0; while (pow((b * b % p - a + p) % p, (p - 1) / 2) != p - 1) ++b; long[] d = { 1, 0 }; long[] m = { b, 1 }; long n = (p + 1) / 2; for (; n > 0; n >>= 1, m = poly_mul(m, m, b, a, p)) { if (n % 2 == 1) d = poly_mul(d, m, b, a, p); } return d[0]; } long[] poly_mul(long[] u, long[] v, long b, long a, long p) { long[] ret = new long[3]; for (int i = 0; i < 2; ++i) { for (int j = 0; j < 2; ++j) { ret[i + j] += u[i] * v[j]; ret[i + j] %= p; } } ret[0] += ret[2] * (b * b - a); ret[0] %= p; for (int i = 0; i < ret.length; ++i) { while (ret[i] < 0) ret[i] += p; } return Arrays.copyOf(ret, 2); } void tr(Object... objects) { System.out.println(Arrays.deepToString(objects)); } }