扩展欧几里得 (exgcd)

这是 CRT 的前置知识，用于解决如下问题：

题意

给定方程 ax + by = c，其中 $a,b,c$ 为已知的正整数，求一组正整数解 $x,y$？

Follow up:

  1. 如果解不存在，输出 $-1$。
  2. 如果解存在，但没有正整数解（正整数解表示 $x,y$ 均 $>0$）求出所有整数解中，$x,y$ 的最小正整数值。
  3. 如果解存在，并且有正整数解，求出正整数解的数量，并求出 $x,y$ 在所有正整数解中，可能的最大最小值。

推导过程如下：

首先判断是否有解，只要知道 $c$ 是否为 $\gcd(a,b)$ 的倍数即可。

以下，我们设 $c=k * \gcd(a,b)$，那么我们只要先求出 $ax+by = \gcd(a,b)$ 的一组可行解 $x,y$，求出原方程的 $x,y$ 只要乘上 $\frac{c}{\gcd(a,b)}$ 即可。

$$ax_1 + by_1 = \gcd(a,b)$$

$$bx_2 + (a\% b)y_2 = \gcd(a,b)$$

注意到 $a \% b = a - \frac{a}{b} * b$，所以化简有：

$$ay_2 + b(x_2 - \frac{a}{b}y_2) = \gcd(a,b)$$

对比第一个式子的系数 $ax_1 + by_1 = \gcd(a,b)$ 可得

$$x_1 = y_2, y_1 = x_2 - \frac{a}{b} y_2$$

递归求解即可，base case为 $x = 1, y = 0$。

ll exgcd(ll a, ll b, ll& x, ll& y) {
    if (!b) {
        x = 1, y = 0;
        return a;  // 这是 gcd(a,b)
    }
    ll g = exgcd(b, a%b, x, y);
    ll x2 = x, y2 = y;
    x = y2, y = x2 - a/b * y2;
    return g;  // 这是 gcd(a,b)
}

int main() {
    ll a,b; cin >> a >> b;
    ll x,y;  // 无需定义，exgcd结束后 x,y 即为所求的 ax+by = gcd(a,b) 中的x,y值
    exgcd(a, b, x, y);
}
Follow Up
  1. 判断是否有解，只要判断 $c$ 是否为 $\gcd(a,b)$ 的倍数即可。

在考虑剩下两个问题之前，我们先将 $a,b,c$ 同除以 $\gcd(a,b)$，这样保证 $\gcd(a,b) = 1$，以下均遵循这个假设（主要是为了保证 $a,b$ 尽可能小）

  1. 要求 $x,y$ 的最小正整数解 $x_{min}, y_{min}$，注意到： $$a(x\pm kb) + b(y \mp ka) = c$$ 其中 $k$ 是任意整数。 并且 $x$ 越小，$y$ 越大，反之亦然。

    所以要求 $x_{min}$，只要保证 $k$ 的值使得 $(x+kb) > 0$ 并且尽可能小即可，此时 $y$ 的值取到所有正整数解（若存在）中的最大值，如果这个最大值仍然 $\leq 0$ 说明对于 $x_{min}$ 正整数解不存在。对于 $y_{min}$ 也要check一下，如果两种情况下都不存在正整数解，我们即可得出结论原方程无正整数解。

    • 注意到上面让 $a = \frac{a}{\gcd(a,b)}, b = \frac{b}{\gcd(a,b)}$ 使得 $a,b$ 尽可能小，保证了 $x_{min}, y_{min}$ 的正确性。

  2. 如果存在正整数解，那么根据上面的描述，解的数量和最大最小值就很容易得出了。

中国剩余定理 (crt)

定义

给定方程组： $$\begin{cases} x\equiv a_1 (\text{mod }m_1) \\
x\equiv a_2 (\text{mod }m_2) \\
… \\
x\equiv a_k (\text{mod }m_k) \\
\end{cases}$$

其中 $a_i \geq 0, m_i > 0, a_i,m_i \in \mathbb{Z}$, 且 $m_i$ 之间两两互质。

求满足条件的最小非负整数解 $x$ ？

结论

令 $M = m_1 \cdot m_2 ~… ~m_k$, $M_i = \frac{M}{m_i}$

令 $M_i^{-1}$ 为 $M_i$ 在 $\color{red} \text{mod }m_i$ 意义下的逆元

则，答案为: $x = \sum\limits_{i=1}^{k}a_iM_iM_i^{-1} ~ (\text{mod }M)$

证明

因为 $\forall i \neq j, M_i \equiv 0 ~(\text{mod } m_j)$

所以 $\forall i \neq j, a_iM_iM_i^{-1} \equiv 0~ (\text{mod } m_j)$

所以 $\forall i, x \equiv \sum\limits_{i=1}^{k}a_iM_iM_i^{-1} \equiv a_iM_iM_i^{-1} ~ (\bmod ~m_i)$

又因为 $\forall i, M_iM_i^{-1} \equiv 1 ~(\text{mod } m_i)$

所以 $\forall i, x \equiv a_iM_iM_i^{-1} \equiv a_i ~ (\text{mod } m_i)$

证明推导过程

证明推导过程

这里直接放参考链接一中的片段：

image

由于 $M = 3\times5\times7 = 105$

所以最后的解就是 $x = (n_1 + n_2 + n_3) ~\bmod ~ 105 = 233 ~ \bmod 105 = 23$

例题

  1. https://www.luogu.com.cn/problem/P3868

快速乘

在题目中，有可能出现 $\bmod$ 接近 $1e18$ 的情况，有可能会爆 long long，所以在计算大数乘法时，要用快速乘（原理类似于快速幂）：

ll qmul(ll a, ll b) {
    ll ans = 0;
    while (b) {
        if (b & 1) ans = (ans + a) % mod;
        b >>= 1;
        a = (a + a) % mod;
    }
    return ans;
}

模版

luogu-P3868-AC代码

题目链接: https://www.luogu.com.cn/problem/P3868

#include <bits/stdc++.h>
 
#define fastio ios::sync_with_stdio(false); cin.tie(0);
#define ll long long

ll M = 1;

ll qmul(ll a, ll b) {
    ll ans = 0;
    while (b) {
        if (b & 1) ans = (ans + a) % M;
        b >>= 1;
        a = (a + a) % M;
    }
    return ans;
}

ll exgcd(ll a, ll b, ll& x, ll& y) {
    if (!b) {
        x = 1, y = 0;
        return a;
    }
    ll g = exgcd(b, a%b, x, y);
    ll curx = y;
    ll cury = x - (a/b) * y;
    x = curx; y = cury;
    return g;
}

ll a[15], b[15];
int k;

ll crt() {
    ll ans = 0;
    for (int i = 1; i <= k; i++) {
        ll m = M / b[i];
        ll x,y;
        exgcd(m, b[i], x, y);
        if (x < 0) x += b[i];

        ans = (ans + qmul(qmul(a[i], m), x)) % M;
    }
    return ans;
}

int main() {
    fastio;

    cin >> k;
    for (int i = 1; i <= k; i++) cin >> a[i];
    for (int i = 1; i <= k; i++) cin >> b[i], M *= b[i];

    for (int i = 1; i <= k; i++) {
        a[i] -= (a[i]/b[i]) * b[i];  //将a[i]变成正数
        a[i] += b[i];
        a[i] %= b[i];
    }
    cout << crt() << endl;
}

参考链接

  1. https://zhuanlan.zhihu.com/p/103394468

拓展中国剩余定理 (excrt)

定义

用于 $m_1,m_2,…,m_k$ 并不互质的情况

给定方程组： $$\begin{cases} x\equiv a_1 (\text{mod }m_1) \\
x\equiv a_2 (\text{mod }m_2) \\
… \\
x\equiv a_k (\text{mod }m_k) \\
\end{cases}$$

其中 $a_i \geq 0, m_i > 0, a_i,m_i \in \mathbb{Z}$，（$m_i$ 之间无特殊关联）

求满足条件的最小非负整数解 $x$ ？

Intuition

先看一些简单的例子：

$$\begin{cases} x\equiv 2 ~(\text{mod }4) \\
\\
x\equiv 4 ~(\text{mod }6) \\
\end{cases}\Longrightarrow x \equiv 10 ~(\text{mod } 12)$$

$$\begin{cases} x\equiv 4 ~(\text{mod }6) \\
\\
x\equiv 3 ~(\text{mod }5) \\
\end{cases}\Longrightarrow ~\varnothing$$

可以看出几个特点：

  1. 答案的 $\text{mod}$ 为 $\text{lcm}(m_1, m_2)$
  2. 有可能无解 （若$m_1,m_2$互质则必然有解）

推导

考虑 $x$ 满足：

$$\begin{cases} x\equiv a_1 ~(\text{mod }m_1) \\
\\
x\equiv a_2 ~(\text{mod }m_2) \\
\end{cases}$$

则有 $x = k_1m_1 + a_1 = k_2m_2 + a_2$

即：$k_1m_1 - k_2m_2 = a_2 - a_1$

这个方程有解 $\iff \gcd(m_1,m_2) | (a_2 - a_1)$ （即 $(a_2 - a_1)$ 是 $\gcd(m_1,m_2)$ 的倍数）。

如果无解，就直接退出。

如果有解，就用 $exgcd(m_1, m_2, x’, y’)$ 解出 $k_1$ 的值。

设 $\gcd(m_1, m_2) = g$，则 $k_1 = x’ \cdot\frac{a_2 - a_1}{g}$

然后将 $k_1$ 代入 $x = k_1m_1 + a_1$，得到 $x$ 的值，答案就是 $x \equiv (k_1m_1 + a_1) (\text{mod ~lcm}(m_1,m_2))$

这样，$2$个同余方程就转化为了$1$个，用同样的方法继续进行下去，即可化为一个同余方程。

最终的 $\text{mod}$ 数就是 $M = \text{lcm}(m_1, m_2, m_3 … , m_n)$

模版

luogu-P4777-AC代码

题目链接: https://www.luogu.com.cn/problem/P4777

using namespace std;
#include <bits/stdc++.h>

#define ll long long
const int maxn = 1e5+5;

ll exgcd(ll a, ll b, ll& x, ll& y) {
    if (!b) {
        x = 1, y = 0;
        return a;
    }
    ll g = exgcd(b, a%b, x, y);
    ll x2 = x, y2 = y;
    x = y2, y = x2 - a/b * y2;
    return g;
}

ll mul(ll a, ll b, ll mod) {
    ll res = 0;
    while (b) {
        if (b&1) (res += a) %= mod;
        (a += a) %= mod;
        b >>= 1;
    }
    return res;
}

int n;
ll a[maxn], m[maxn];

int main() {
    fastio;

    cin >> n;

    for (int i = 1; i <= n; i++) cin >> m[i] >> a[i];
    ll cura = a[1], curm = m[1];

    for (int i = 2; i <= n; i++) {
        /* m1 = curm, m2 = m[i], a1 = cura, a2 = a[i] */
        ll x,y;
        ll g = exgcd(curm, m[i], x, y);  // g = gcd(curm, m[i])
        ll c = a[i] - cura;  // c = a2 - a1
        c = ((c % m[i]) + m[i]) % m[i];  // 处理负数，注意是 mod m2 (即 m[i])

        if (c % g) {  // 无解
            cout << -1 << endl;
            return 0;
        }

        ll M = (m[i] * (curm / g));  // M = lcm(m1, m2)
        ll k = mul(x, c / g, M);  // k1 = x' * (a2-a1) / g

        cura = (mul(k, curm, M) + cura) % M;  // x = k1m1 + a1
        curm = M;
    }
    cout << (cura % curm + curm) % curm << endl;
}

注意事项

  1. 我们需要处理可能出现的负数！（比如在快速乘的时候，如果 $b$ 为负数就会炸）

    1. 在 $exgcd(a,b,x’,y’)$ 中，$x'$ 有可能为负数。
    2. $c = a_2 - a_1$，则 $c$ 也有可能为负数。

    那么如何处理负数？我们需要找到这个负数对应的 $\text{mod}$为多少，我们可以在等式 $k_1m_1 - k_2m_2 = a_2 - a_1$ 中找到！

    因为我们要求的是 $k_1$ 的值，实际上 $k_1m_1 = (a_2-a_1) + k_2m_2$，即 $k_1m_1 \equiv (a_2-a_1) \text{ (mod } m_2)$

    所以无论是 $k_1 = x’ \cdot\frac{a_2 - a_1}{g}$ 还是 $(a_2-a_1)$，都取一个 $\text{mod } m_2$ 即可！

    代码段：

     ll c = a[i] - cura;
 c = ((c % m[i]) + m[i]) % m[i];  // 处理负数，注意是 mod m2 (即 m[i])

 ll k = mul(x, c / g, M);
 // k = ((k % m[i]) + m[i]) % m[i];  可以处理，但是没有必要，因为快速乘中的 a 可以为负数

  2. 注意快速乘过程中的 $\text{mod}$ ！因为快速乘求的实际上是 $x$ 的值，所以 $\text{mod } M$，其中 M = lcm(curm, m[i])

  3. 我们用 cura, curm 来维护当前的 同余方程，之后记得更新一下 curm = lcm(curm, m[i])

参考链接

  1. https://www.luogu.com.cn/blog/blue/kuo-zhan-zhong-guo-sheng-yu-ding-li