中国剩余定理(CRT)学习笔记

中国剩余定理(CRT)是用来求这样的一个不定方程组的解的：

$\begin{cases}x\equiv a_1\pmod{b_1}\\x\equiv a_2\pmod{b_2}\\x\equiv a_3\pmod{b_3}\\\qquad\dots\dots\\x\equiv a_n\pmod{b_n}\end{cases}$

其中 $b_i$ 两两互质。

我们可以先考虑一个简单点的问题：

$\begin{cases}x\equiv 2\pmod{3}\\x\equiv 3\pmod{5}\\x\equiv 2\pmod{7}\end{cases}$

我们可以：

• 令 $y_1$ 除以 $3$ 余 $2$，除以 $5$ 余 $0$，除以 $7$ 余 $0$

• 令 $y_2$ 除以 $3$ 余 $0$，除以 $5$ 余 $3$，除以 $7$ 余 $0$

• 令 $y_3$ 除以 $3$ 余 $0$，除以 $5$ 余 $0$，除以 $7$ 余 $2$

那么这个方程组的一个解就是 $x=y_1+y_2+y_3$。

继续对 $y_1$ 分解（其实对 $y_2$ 和 $y_3$ 也同理），令 $z$ 除以 $3$ 余 $1$，除以 $5$ 余 $0$，除以 $7$ 余 $0$，那么显然有 $y_1=z*2$。

令 $M=\prod\limits_{i=1}^nb_i=105$，那么 $z$ 肯定能被 $\dfrac{M}{b_1}=\dfrac{105}{3}=35$ 整除，不妨令 $z=35k$。那么 $35k\equiv 1\pmod{b_1}$ 即 $35k\equiv 1\pmod{3}$。很明显，$k$ 就是 $35$ 在模 $3$ 意义下的逆元。

所以 $y_1=a_1\cdot\dfrac{M}{b_1}\cdot\operatorname{inv}(\dfrac{M}{b_1},b_1)$ 即 $y_1=2*35*2=140$。

归纳一下，$y_i=a_i\cdot\dfrac{M}{b_i}\cdot\operatorname{inv}(\dfrac{M}{b_i},b_i)$，那么 $x=\sum\limits_{i=1}^n a_i\cdot\dfrac{M}{b_i}\cdot\operatorname{inv}(\dfrac{M}{b_i},b_i)$。

但由于这时候求出的 $x$ 只是任意解，最小解 $x_{min}$ 即为 $x\mod{\operatorname{lcm}(b_1,b_2,b_3,\dots,b_n)}$。

所以方程组

$\begin{cases}x\equiv a_1\pmod{b_1}\\x\equiv a_2\pmod{b_2}\\x\equiv a_3\pmod{b_3}\\\qquad\dots\dots\\x\equiv a_n\pmod{b_n}\end{cases}$

的最小解为：

$x_{min}=\sum\limits_{i=1}^n a_i\cdot\dfrac{M}{b_i}\cdot\operatorname{inv}(\dfrac{M}{b_i},b_i)\mod{\operatorname{lcm}(b_1,b_2,b_3,\dots,b_n)}$

其中 $M=\prod\limits_{i=1}^nb_i$。

可以发现 $\dfrac{M}{b_i}$ 其实能用前缀和维护。

模板题代码：

#include <iostream>
#include <cstdio>

using namespace std;

int n;
long long a[15],b[15];
long long sum[2][15];

long long exgcd(long long a,long long b,long long &x,long long &y)
{
	if(a<b)
	{
		return exgcd(b,a,y,x);
	}
	if(b==0)
	{
		x=1;
		y=0;
		return a;
	}
	long long tmpx,tmpy;
	long long res=exgcd(b,a%b,tmpx,tmpy);
	x=tmpy;
	y=tmpx-a/b*tmpy;
	return res;
}

inline long long inv(long long a,long long p) // 扩欧求逆元 
{
	long long x,y;
	if(exgcd(a,p,x,y)!=1)
	{
		return -1;
	}
	return (x%p+p)%p;
}

inline long long lcm(long long a,long long b)
{
	long long tmpa=a,tmpb=b;
	long long t=a%b;
	while(t!=0)
	{
		a=b;
		b=t;
		t=a%b;
	}
	return tmpa*tmpb/b;
}

int main()
{
	scanf("%d",&n);
	long long llcm=1;
	for(int i=1;i<=n;i++)
	{
		scanf("%lld%lld",&b[i],&a[i]);
		llcm=lcm(llcm,b[i]);
	}
	// 前缀和维护 M/b_i 
	sum[0][0]=sum[1][n+1]=1;
	for(int i=1;i<=n;i++)
	{
		sum[0][i]=sum[0][i-1]*b[i];
	}
	for(int i=n;i>=1;i--)
	{
		sum[1][i]=sum[1][i+1]*b[i];
	}
	long long ans=0;
	for(int i=1;i<=n;i++)
	{
		ans+=a[i]*sum[0][i-1]*sum[1][i+1]*inv(sum[0][i-1]*sum[1][i+1],b[i]);
	}
	printf("%lld\n",ans%llcm); // 求最小解 
	return 0;
}