(扩展)中国剩余定理
中国剩余定理
概念
中国剩余定理(\(\texttt{Chinese Remainder Theorem}\), \(\tt CRT\))可以求解关于 \(x\) 的线性同余方程组
其中 \(p_1, p_2, ..., p_k\) 两两互质。
思想
设 \(M = \prod\limits_{i = 1}^k p_i\),\(m_i = \frac{M}{p_i}\)
另设 \(m_i^{-1}\) 为 \(m_i\) 在模 \(p_i\) 意义下的逆元,即 \(m_im_i^{-1} \equiv 1(\bmod p_i)\)
则方程的唯一解为 \(x = \sum\limits_{i = 1}^k a_im_im_i^{-1}\)
证明
首先证明方程解的正确性。
然后证明方程解在模 \(M\) 意义下的唯一性。
考虑采用反证法。
设原方程组在模 \(M\) 意义下有两个不同的整数解 \(x, y\),
则 \(\forall 1 \leq i \leq k\),有 \(x \equiv (x \bmod m) \equiv y(\bmod p_i)\)
所以 \(x, y\) 在模 \(M\) 意义下是同一个解。
模板
#include <cstdio>
using namespace std;
typedef long long ll;
const int maxn = 15;
int n;
ll mul;
ll a[maxn], p[maxn], inv[maxn], m[maxn];
ll exgcd(ll a, ll b, ll &x, ll &y)
{
ll d = a;
if (b == 0)
x = 1, y = 0;
else
{
d = exgcd(b, a % b, y, x);
y -= (a / b * x);
}
return d;
}
int main()
{
ll x, y, ans = 0;
mul = 1;
scanf("%d", &n);
for (int i = 1; i <= n; i++)
{
scanf("%lld%lld", &p[i], &a[i]);
mul *= p[i];
}
for (int i = 1; i <= n; i++)
{
m[i] = mul / p[i];
exgcd(m[i], p[i], inv[i], y);
ans = (ans + a[i] * m[i] * (inv[i] < 0 ? inv[i] + p[i] : inv[i]));
}
printf("%lld\n", ans % mul);
return 0;
}
扩展中国剩余定理
概念
扩展中国剩余定理(\(\tt EXCRT\))可以求解关于 \(x\) 的线性同余方程组:
其中 \(p_1, p_2, ..., p_k\) 不一定 两两互质。
思路
如果按照 \(\tt CRT\) 的思路求解,会出现 \(m_i\) 的逆元根本不存在的情况,因此考虑合并同余方程。
假设前 \(i - 1\) 个同余方程的一个特解为 \(x^{\prime}\) 且 \(\operatorname{lcm}_{j = 1}^{i - 1} p_j = M\)
则前 \(i - 1\) 个同余方程的通解形式为 \(x^{\prime} + tM, t \in \mathbb{Z}\)
合并第 \(i\) 个同余方程实际上等价于求满足要求的 \(t\),使得 \(x^{\prime} + tM \equiv a_i(\bmod p_i)\)
即 \(tM \equiv a_i - x^{\prime}(\bmod p_i)\)
原方程可进一步化为 \(tM + kp_i = a_i - x^{\prime}, k \in \mathbb{Z}\)
当该方程无解时,原方程组无解。
反之,考虑用扩展欧几里得求出 \(t, k\) 的一组特解,则前 \(i\) 个同余方程的一个特解为 \(x^{\prime} + tM\)
显然联立方程 \(x \equiv 0(\bmod 1)\) 与原方程组对解没有影响
所以在代码中可以把 \(\operatorname{lcm}_{j = 1}^{i - 1}p_j\) 的初始值赋值为 \(1\),将特解赋值为 \(0\),直接从原方程组中的第 \(1\) 个方程开始合并
代码
#include <cstdio>
using namespace std;
typedef long long ll;
const int maxn = 1e5 + 5;
int n;
ll a[maxn], p[maxn];
ll fmul(ll a, ll b, ll mod)
{
ll res = 0;
while (b)
{
if (b & 1)
res = (res + a) % mod;
a = (a + a) % mod;
b /= 2;
}
return res;
}
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);
ll t = x;
x = y;
y = t - a / b * y;
return g;
}
ll excrt()
{
ll ans = 0, mul = 1, x, y;
for (int i = 1; i <= n; i++)
{
ll A = mul, B = p[i], C = (a[i] - ans % B + B) % B;
ll g = exgcd(A, B, x, y), bg = B / g;
if (C % g != 0)
return -1;
x = fmul(x, C / g, bg);
ans = ans + x * mul;
mul = mul * bg;
ans = (ans % mul + mul) % mul;
}
return ans;
}
int main()
{
scanf("%d", &n);
for (int i = 1; i <= n; i++)
scanf("%lld%lld", &p[i], &a[i]);
printf("%lld\n", excrt());
return 0;
}