【算法学习笔记】37：扩展中国剩余定理（EXCRT）求解任意线性同余方程组

中国剩余定理仅适用于模数$m_1...m_n$两两互质的情况，如果没有这个限制条件则不能使用。

下面进行一般情况的解法推导，注意虽然叫扩展中国剩余定理（EXCRT），实际上只是扩宽了中国剩余定理（CRT）的求解限制，解法和原理上和CRT完全没啥关系。

1 推导过程

对于一般的线性同余方程组：
$$\begin{gathered} xmod{m}_1\equiv a_1 \\ xmod{m}_2\equiv a_2 \\ ... \\ xmod{m}_n\equiv a_n \end{gathered}$$

1.1 对于任意两个方程求得一组解$k_1$、$k_2$

先考虑其中两个方程组成的方程组：
$$\begin{gathered} xmod{m}_1\equiv a_1 \\ xmod{m}_2\equiv a_2 \end{gathered}$$

说明存在$k_1$、$k_2$分别是$x$除以$m_1$、$m_2$的解，并留下余数是$a_1$和$a_2$，也即：
$$\begin{gathered} x=k_1\cdot m_1+a_1 \\ x=k_2\cdot m_2+a_2 \end{gathered}$$

从而可以用$x$把两式关联起来：
$$k_1\cdot m_1+a_1=k_2\cdot m_2+a_2$$

整理即得一个不定方程：
$$k_1\cdot m_1-k_2\cdot m_2=a_2-a_1$$

其中$m_1$、$m_2$、$a_2$、$a_1$都是给定的，只要求出$k_1$和$k_2$就能算出$x$了，只要用扩展欧几里得算法就能计算$k_1$和$-k_2$，仅当$a_2-a_1$能被$m_1$和$m_2$的最大公约数整除时候有解，即$gcd(m_1,m_2)|a_2-a_1$。

1.2 从$k_1$、$k_2$得到$x$的所有解的表示形式

从前面的推导可以得出：
$$\begin{gathered} x=k_1\cdot m_1+a_1 \\ x=k_2\cdot m_2+a_2 \end{gathered}$$

只要得到了一组$k_1$、$k_2$，已经满足了$k_1\cdot m_1-k_2\cdot m_2=a_2-a_1$，对于最大公约数$d=gcd(m_1,m_2)$，只要找任何一个倍数$k$，让$k_1$加上$k\cdot \frac{m_2}{d}$，让$k_2$加上$k\cdot \frac{m_1}{d}$，该等式还是成立的。

所以引入整数$k$来表示$x$的所有解，应该形如（只看第一个式子的话）：
$$x=(k_1+k\frac{m_2}{d})m_1+a_1$$

整理一下，得到：
$$x=m_1{k}_1+a_1+k\frac{m_1{m}_2}{d}$$

其中$d=gcd(m_1,m_2)$，所以$\frac{m_1{m}_2}{d}$也即$lcm(m_1,m_2)$，因此：
$$x=m_1{k}_1+a_1+k\cdot lcm(m_1,m_2)$$

1.3 任意两个不定方程的合并

其中$m_1{k}_1+a_1$这一项在求出$k_1$之后就是固定的了，记为$A$，后面的$lcm(m_1,m_2)$记为$M$，则$x$的所有解的形式就可以表示为：
$$x=k\cdot M+A$$

这个方程不仅和前面的不定方程形式一样，而且其中$M$和$A$都是确定的，$k$是变化的。所以总是可以通过这样的方式把两个不定方程合并成一个，也即对于前两个方程
$$\begin{gathered} x=k_1\cdot m_1+a_1 \\ x=k_2\cdot m_2+a_2 \end{gathered}$$

可以合并为$x=kM+A$，其中$M=lcm(m_1,m_2)$，$A=m_1{k}_1+a_1$，然后可以拿着新方程和第三个方程继续合并，即合并：
$$\begin{gathered} x=k\cdot M+A \\ x=k_3\cdot m_3+a_3 \end{gathered}$$

一共有$n$个方程，只要这样合并$n-1$次就只剩一个方程$x=k\cdot M+A$了，这时要找$x$的最小正整数取值就很简单了，一定是$A$模$M$的正的余数。

2 例题：AcWing 204. 表达整数的奇怪方式

注意这题里$a$和$m$的含义和前面的推导是相反的。

这题里限制了方程的数量$n$是不超过25个的，因为每次合并两个方程，得到的新方程的$M$都是前两个方程里$m_i$和$m_j$的最小公倍数，所以最后的$M$就是$m_1...m_n$的最小公倍数，如果方程太多了这个最小公倍数就太大了，可能超出整数范围，所以数据里做了这样的限制。

为了防止得到的新的方程里$k$越来越大爆掉，每次都取一个最小的正整数$k$。

小技巧：要取任何一个整数p对于模q的数学模（最小正整数），只要(p % q + q) % q即可

#include <bits/stdc++.h>

using namespace std;

typedef long long LL;

LL exgcd(LL a, LL b, LL& x, LL& y) {
    if (!b) {
        x = 1, y = 0;
        return a;
    }
    LL d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}

int main() {
    int n; cin >> n;
    
    // 读取第一组方程
    LL m1, a1; cin >> m1 >> a1;
    for (int i = 1; i < n; i ++ ) {
        // 每次合并完读取下一组方程
        LL m2, a2; cin >> m2 >> a2;
        // 尝试合并两个方程
        // x = k1 m1 + a1
        // x = k2 m2 + a2
        // k1 m1 - k2 m2 = a2 - a1
        // 用exgcd计算k1和-k2(丢弃)
        LL k1, k2;
        LL d = exgcd(m1, m2, k1, k2);
        // 如果有解，d必须是a2 - a1的倍数
        if ((a2 - a1) % d) {
            puts("-1");
            return 0;
        }
        // 把d变成a2 - a1时候，k1也跟着倍增
        k1 *= (a2 - a1) / d;
        // 由于k1可以是任意的k1 + k * (m2 / d)
        // 这里为了防止太大，每次把k1变成满足条件的最小正值
        LL t = m2 / d;
        k1 = (k1 % t + t) % t;
        // 接下来有了k1，合并好的新方程是
        // x = k1 M + A，其中
        // A = m1 k1 + a1，M = lca(m2, m1)
        // 这里M就是下一轮的m1，A就是下一轮的a1
        a1 = m1 * k1 + a1;
        m1 = m1 * m2 / d;
    }
    
    // 最后的解是上一轮的A = m1 k1 + a1，放在了a1里
    // 为了调整k1让A是最小正整数，即是对m1取模时的最小正数
    cout << (a1 % m1 + m1) % m1 << endl;
    
    return 0;
}