题目传送门：洛谷 P5345。
很荣幸为 X Round 1 贡献了自己的一题。
题意简述：
给定 (n) 组 (k_i,g_i,r_i)（(0le k_i,r_i<g_ile 10^7)）。
求关于 (x) 的方程组 (left{egin{matrix}k_1^xequiv r_1pmod{g_1}\k_2^xequiv r_2pmod{g_2}\vdots\k_n^xequiv r_npmod{g_n}end{matrix} ight.)（定义 (0^0=1)）在 ([0,10^9]) 内的最小整数解，或判断在这个范围内无解。
题解：
分为两个部分解决。
第一部分：分别求出每个方程的解。
首先观察 (k^jmod g) 的规律，以 (g=495616)，(k=124) 为例：
定义 (f(i)=kimod g)，则把每个在 ([0,g)) 之间的整数看作一个节点后，形成了一个基环内向森林。
从 (1mod g) 开始走唯一的出边，必然形成一个 ( ho) 形结构，让我们模拟这个过程：
可以看到，出现了长度为 (5) 的循环。
-
如果 (r=245760)，则 (xequiv 7pmod{5})，且 (xge 7)。
-
如果 (r=12544)，则 (x=4)。
-
如果 (r=2)，则无解。
以上就是每个方程的 (3) 种解，无限解，唯一解或无解。
那么，问题就是，如何区分这 (3) 种解，以及如何求出解。
首先，需要将 ( ho) 形的“尾巴”和循环分开考虑，不难证明“尾巴”的长度不超过 (log_2g)。
而且可以发现，若在“尾巴”上找到了解，那么就只有一组解。所以需要判断一个值是否在尾巴上。
其实这是很简单的，如果 (gcd(x,m) eqgcd(f(x),m))，那么 (x) 就在尾巴上，反之就在循环内。
据此，我们区分了在“尾巴”上的解，也就是只有一组解的情况，接下来考虑不在“尾巴”上的情况。
因为解不在“尾巴”上，则如果解不在循环内就是无解了，所以首先判断解是否在循环内。
假设第一个在循环上的数为 (c)，尾巴长度为 (o)，则有 (cequiv k^opmod{g})。
令 (d=gcd(c,g))，则原方程 (k^xequiv rpmod{g}) 可化为 (ccdot k^{x-o}equiv rpmod{g}) 且 (xge o)，这是因为我们假定解一定在循环内。
可以进一步化为 (k^{x-o}equivfrac{r}{d}left(frac{c}{d} ight)^{-1}pmod{frac{g}{d}})，若 (r) 不是 (d) 的倍数则无解。
这是因为循环内的值均是 (d) 的倍数，可以让方程同除 (d)，又因为 (frac{c}{d}perpfrac{g}{d})，可以取逆元。
令 (a=k,b=frac{r}{d}left(frac{c}{d} ight)^{-1},m=frac{g}{d})，则有 (aperp m)。
使用 BSGS 求出 (a^yequiv bpmod{m}) 的最小自然数解 (y) 后，有原方程的最小自然数解为 (y+o)。
再使用 BSGS 求出 (a^zequiv 1pmod{m}) 的最小正整数解 (z)，即 (a) 对 (m) 的阶。
则原方程的解为 (xequiv y+opmod{z}) 且 (xge y+o)。
至此，第一部分解决。
第二部分：合并每个方程的解。
首先，若前面的方程出现了无解的情况，则方程组也无解。
若前面的方程出现了只有唯一解的情况，只需要检查此唯一解是否满足所有方程即可。
接下来讨论以上每个方程的解均形如 (xequiv rpmod{q}) 且 (xge r) 的形式。
分开考虑前半部分和后半部分，对于前半部分，显然是 ExCRT 的形式。
对于后半部分，可以化为 (xgemax r_i)，令 (max r_i=x_0)，放到最后考虑。
考虑使用 ExCRT 解决前半部分，令前 (i) 个方程组的解为 (xequiv P_ipmod{Q_i})，有 (P_0=0,Q_0=1)。
但是这里出现一个问题，那就是 (P_i,Q_i) 可能很大，但是这里其实没必要使用高精度，要怎么处理这种情况呢？
因为以上方程的解 (xequiv rpmod{q}) 中，(q) 均小于 (10^7)，所以不需要担心这部分。
考虑这样一种情况：(Q_{i-1}>10^9)，而再合并进 (xequiv r_ipmod{q_i}) 可能会让 (Q_i) 超出 long long 能够表示的范围。
这时其实不需要再合并了，只需要判断 (P_{i-1}) 是否满足第 (i) 个方程即可。因为若不满足，合并后 (x) 必然会至少增加一倍的 (Q_{i-1})，这就超出了 (10^9) 的范围，从而可以直接输出无解。
最后，当方程成功合并完之后，(Q_n) 可能不是真实值，但是当 (Q_nle 10^9) 时必然是真实值。
此时可以再考虑 (x_0)，真实的解应该为 (x=leftlceilfrac{x_0-P_n}{Q_n} ight ceilcdot Q_n+P_n)。
若 (Q_n> 10^9)，则必然要满足 (P_nge x_0)，否则无解。
据此写出代码，复杂度 (mathcal{O}left(sum_{i=1}^{n}left(sqrt{g_i}+log g_i ight) ight))：
#include <cstdio>
#include <algorithm>
#include <cstring>
#include <cmath>
#define Fail return puts("Impossible"), 0
#define mp std::make_pair
typedef long long LL;
typedef std::pair<int, int> pii;
const int MG = 10000005, MS = 3175;
const int U = 1e9;
int gcd(int a, int b) { return b ? gcd(b, a % b) : a; }
template<typename T>
T exgcd(T a, T b, T &x, T &y) {
if (!b) return x = 1, y = 0, a;
int d = exgcd(b, a % b, y, x);
return y -= a / b * x, d;
}
int buk[MG], stk[MS];
inline int BSGS(int a, int b, int m) {
int S = sqrt(m - 1) + 1;
int A = 1, f = -1, t = 0;
for (int i = 0; i < S; ++i) {
buk[stk[++t] = (LL)b * A % m] = i;
A = (LL)A * a % m;
}
int C = 1;
for (int i = 1; !~f && i <= S; ++i)
if (~buk[C = (LL)C * A % m])
f = i * S - buk[C];
while (t) buk[stk[t--]] = -1;
return f;
}
inline pii ExBSGS(int a, int b, int m) {
int o = 0, A = 1 % m, d = 1, nd, x, y;
while (1) {
if (d == (nd = gcd((LL)A * a % m, m))) break;
if (A == b) return mp(o, -1);
++o, A = (LL)A * a % m, d = nd;
}
if (b % d) return mp(-1, -1);
m /= d, b /= d, A /= d;
exgcd(A, m, x, y);
b = (LL)b * (x + m) % m;
x = BSGS(a, b, m);
if (!~x) return mp(-1, -1);
y = BSGS(a, 1 % m, m);
return mp(x % y + o, y);
}
inline bool Combine(LL &a1, LL &m1, LL a2, LL m2) {
LL k1, k2, g = exgcd(m1, m2, k1, k2);
if ((a2 - a1) % g) return 0;
a1 += (k1 * ((a2 - a1) / g) % m2 + m2) * m1 % (m1 / g * m2);
return a1 %= m1 *= m2 / g, 1;
}
const int MN = 1005;
int N;
int k[MN], r[MN], g[MN];
int x[MN], q[MN], X, MaxX;
int main() {
memset(buk, -1, sizeof buk), X = -1;
scanf("%d", &N);
for (int i = 1; i <= N; ++i) {
scanf("%d%d%d", &k[i], &g[i], &r[i]);
pii ans = ExBSGS(k[i] % g[i], r[i] % g[i], g[i]);
x[i] = ans.first, q[i] = ans.second;
if (!~x[i]) Fail;
if (!~q[i]) X = x[i];
MaxX = std::max(MaxX, x[i]);
}
if (~X) {
for (int i = 1; i <= N; ++i) {
if (~q[i] && (X < x[i] || (X - x[i]) % q[i])) Fail;
if (!~q[i] && X != x[i]) Fail;
}
return printf("%d
", X), 0;
}
LL P = 0, Q = 1;
for (int i = 1; i <= N; ++i) {
if (Q > U && (P - x[i]) % q[i]) Fail;
if (Q <= U && !Combine(P, Q, x[i] % q[i], q[i])) Fail;
}
if (P < MaxX) P += ((MaxX - P - 1) / Q + 1) * Q;
if (P > U) Fail;
printf("%lld
", P);
return 0;
}