中国剩余定理及其拓展
\(\text{CRT & exCRT}\)
\(\text{CRT}\)
中国剩余定理(\(Chinese~Remainder~Theorem, \text{CRT}\))是数论中的一个关于一元线性同余方程组的定理,说明了一元线性同余方程组有解的准则以及求解方法。
一元线性同余方程组的一般形式:
令 \(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\) 意义下的逆元。
为什么它成立呢?
模板(洛谷 \(\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\) 个方程,现在要合并的是:
整理一下,有:
移项,得:
再写成一元线性同余方程的形式:
很容易就可以用 \(\text{exgcd}\) 求出 \(k_1\) 的特解,进而更新此时 \(x\) 的特解 \(res\)(即前 \(k\) 个一元线性同余方程的特解)。
而前 \(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\)。
那么原方程组的解可写成:
化成一元线性同余方程的形式:
以此类推,一直合并,直到剩下最后一个一元线性同余方程,此时的 \(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;
}