Loading

中国剩余定理及其拓展

\(\text{CRT & exCRT}\)

\(\text{CRT}\)

中国剩余定理(\(Chinese~Remainder~Theorem, \text{CRT}\))是数论中的一个关于一元线性同余方程组的定理,说明了一元线性同余方程组有解的准则以及求解方法。

一元线性同余方程组的一般形式:

\[\left\{ \begin{aligned} &x \equiv a_1 \pmod{m_1} \\ &x \equiv a_2 \pmod{m_2} \\ &\dots \\ &x \equiv a_n \pmod{m_n} \end{aligned} \right. \]

\(M = \prod\limits_{i = 1}^n m_i\), 当自然数 \(m\) 两两互质时,上述方程组在模 \(M\) 意义下有唯一解。

这个解为 \(x \equiv \sum\limits_{i = 1}^{n} (a_i \times p_i^{-1} \times p_i)\),其中 \(p_i = \frac{M}{m_i}, p_i^{-1}\)\(p_i\) 在模 \(m_i\) 意义下的逆元。

为什么它成立呢?

\[x \bmod m_t = \sum_{i = 1}^n (a_i \times p_i^{-1} \times p_i \bmod m_t) \bmod m_t\\ \because \forall k \in [1, t - 1] \cup [t + 1, n], m_t|p_k \\ \therefore \forall k \in [1, t - 1] \cup [t + 1, n], a_k \times p_k^{-1} \times p_k \bmod m_t = 0 \\ \therefore x \bmod m_t = a_t \times p_t^{-1} \times p_t \bmod m_t = a_t \bmod m_t \\ \text{即} x \equiv a_t \pmod{m_t} \]

模板(洛谷 \(\text{P1495}\)

\(Code\)

#include <bits/stdc++.h>

#define MAXN 100100

using namespace std;

typedef long long ll;

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

template<typename _T>
void read(_T &_x) {
    _x = 0;
    _T _f = 1;
    char _ch = getchar();
    while (_ch < '0' || '9' < _ch) {
        if (_ch == '-') _f = -1;
        _ch = getchar();
    }
    while ('0' <= _ch && _ch <= '9') {
        _x = (_x << 3) + (_x << 1) + (_ch & 15);
        _ch = getchar();
    }
    _x *= _f;
}

template<typename _T>
void write(_T _x) {
    if (_x < 0) {
        putchar('-');
        _x = -_x;
    }
    if (_x > 9) write(_x / 10);
    putchar('0' + _x % 10);
}

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

ll CRT() {
    ll M = 1;
    for (int i = 1; i <= n; i++) M *= m[i];
    ll x = 0;
    for (int i = 1; i <= n; i++) {
        ll Mi = M / m[i], y, z;
        exgcd(Mi, m[i], y, z);
        x = (x + a[i] * Mi * y) % M;
    }
    return (x + M) % M;
}

int main() {
    read(n);
    for (int i = 1; i <= n; i++) {
        read(m[i]), read(a[i]);
    }
    write(CRT());
    return 0;
}

\(\text{exCRT}\)

\(\text{CRT}\) 只能求解 \(m\) 两两互质时的一元线性同余方程组(如果 \(m\) 不两两互质,\(p_i^{-1}\) 不存在),那么当 \(m\) 不一定两两互质时,还能解出一元线性同余方程组吗?

答案是肯定的——\(\text{exCRT}\)

\(\text{exCRT}\) 的核心思想就是一元线性同余方程组的合并。假设我们现在已经合并了前 \(k - 1\) 个方程,现在要合并的是:

\[\left\{ \begin{aligned} &x \equiv res \pmod{p} \\ &x \equiv a_k \pmod{m_k} \end{aligned} \right. \]

整理一下,有:

\[x = k_1p + res = k_2m_k + a \Rightarrow k_1p + res = k_2m_k + a \]

移项,得:

\[k_1p + k_2m_k = a - res \]

再写成一元线性同余方程的形式:

\[k_1p \equiv a - res \pmod{m_k} \]

很容易就可以用 \(\text{exgcd}\) 求出 \(k_1\) 的特解,进而更新此时 \(x\) 的特解 \(res\)(即前 \(k\) 个一元线性同余方程的特解)。

\[res = res + k_1p \]

而前 \(k\) 个一元线性同余方程组的通解一定是 \(x = res + k \times \text{lcm}(p, m_k), k \in \Z\),因为当前的 \(res\) 已经满足上述两个同余方程了,其后加上的常数 \(c\) 一定要满足 \(c \bmod p = c \bmod m_k = 0\) 才能使得 \(res + c\) 也是上述两个同余方程的解,解得 \(c = k \times \text{lcm}(p, m_k), k \in \Z\)

那么原方程组的解可写成:

\[x = res + k \times \text{lcm}(p, m_k) \]

化成一元线性同余方程的形式:

\[x \equiv res \pmod{\text{lcm}(p, m_k)} \]

以此类推,一直合并,直到剩下最后一个一元线性同余方程,此时的 \(res\) 就是通解啦~

模板(洛谷 \(\text{P4777}\)

\(Code\)

#include <bits/stdc++.h>

#define MAXN 100100

using namespace std;

typedef long long ll;

int n;
ll a[MAXN], b[MAXN];

template<typename _T>
void read(_T &_x) {
    _x = 0;
    _T _f = 1;
    char _ch = getchar();
    while (_ch < '0' || '9' < _ch) {
        if (_ch == '-') _f = -1;
        _ch = getchar();
    }
    while ('0' <= _ch && _ch <= '9') {
        _x = (_x << 3) + (_x << 1) + (_ch & 15);
        _ch = getchar();
    }
    _x *= _f;
}

template<typename _T>
void write(_T _x) {
    if (_x < 0) {
        putchar('-');
        _x = -_x;
    }
    if (_x > 9) write(_x / 10);
    putchar('0' + _x % 10);
}

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

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

ll exCRT() {
    ll p = a[1], res = b[1];
    for (int i = 2; i <= n; i++) {
        ll w = (b[i] - res % a[i] + a[i]) % a[i];
        ll x, y;
        ll gcd = exgcd(p, a[i], x, y);
        if (w % gcd) return -1;
        x = sm(x, w / gcd, a[i]);
        res += x * p;
        p *= a[i] / gcd;
        res = (res + p) % p;
    }
    return (res % p + p) % p;
}

int main() {
    read(n);
    for (int i = 1; i <= n; i++) read(a[i]), read(b[i]);
    write(exCRT());
    return 0;
}
posted @ 2022-09-27 10:04  Chy12321  阅读(67)  评论(0编辑  收藏  举报