[笔记]中国剩余定理(CRT) & 扩展中国剩余定理(exCRT)

中国剩余定理(CRT)

P1495 【模板】中国剩余定理(CRT)/ 曹冲养猪

对于线性同余方程组:

\[\begin{cases} x\equiv a_1\pmod{b_1}\\ x\equiv a_2\pmod{b_2}\\ \dots\\ x\equiv a_n\pmod{b_n}\\ \end{cases}\]

其中\(b_i\)两两互质。

下文中,令\(N=\prod\limits_{i=1}^n b_i\)


考虑如何求解其中一个特解\(x\)

根据同余运算的可加性,我们可以将原方程拆成\(n\)个:

\[\begin{cases} x_1\equiv a_1\pmod{b_1}\\ x_1\equiv 0\pmod{b_2}\\ \dots\\ x_1\equiv 0\pmod{b_n}\\ \end{cases} \begin{cases} x_2\equiv 0\pmod{b_1}\\ x_2\equiv a_2\pmod{b_2}\\ \dots\\ x_2\equiv 0\pmod{b_n}\\ \end{cases} \dots \begin{cases} x_n\equiv 0\pmod{b_1}\\ x_n\equiv 0\pmod{b_2}\\ \dots\\ x_n\equiv a_n\pmod{b_n}\\ \end{cases}\]

\(x=\sum_{i=1}^n x_i\)

对于第一个方程组,我们令\(x_1=y_1\times a_1\),则有:

\[\begin{cases} y_1\equiv 1\pmod{b_1}\\ y_1\equiv 0\pmod{b_2}\\ \dots\\ y_1\equiv 0\pmod{b_n}\\ \end{cases}\]

上面的转化过程是CRT的核心思想。

由后\(n-1\)个式子可知\(y_1\equiv 0\pmod{\text{lcm}_{i=2}^n\ b_i}\),即\(y_1\equiv 0\pmod{\frac{N}{b_1}}\),令\(c_i=\frac{N}{b_i}\),则\(y_1=c_1\times k_1(k_1\in \mathbb{Z})\)

由第\(1\)个式子可知\(c_1\times k_1\equiv 1\pmod{b_1}\),于是\(k_1\equiv c_1^{-1}\pmod{b_1}\),即\(k_1\)\(c_1\)在模\(b_1\)意义下的乘法逆元。

类比可得:

\[x_i=y_i\times a_i=c_i\times k_i\times a_i\\ x=\sum\limits_{i=1}^n c_i\times k_i\times a_i\]

至此完成\(x\)的求解。


再考虑如何用特解\(x\)表示所有的通解。

考虑\(x,y\)同为该方程组的解,则根据同余的传递性,有:

\[\begin{cases} x\equiv y\pmod{b_1}\\ x\equiv y\pmod{b_2}\\ \dots\\ x\equiv y\pmod{b_n} \end{cases}\]

也即\(x\equiv y\pmod{\text{lcm}_{i=1}^n b_i}\),即\(x\equiv y\pmod N\)

所以通解可以表示为\(kN+x(k\in \mathbb{Z})\)

这也告诉我们,对于\(k\in \mathbb{Z}\),该方程组的解在\([kN,(k+1)N)\)上是唯一的。


由于\(b_i\)不一定是质数,所以乘法逆元采用exGCD

点击查看代码
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=15;
int n,fac=1,a[N],b[N];
__int128 ans;
int exgcd(int a,int b,int& x,int& y){
	int d;
	if(b==0) x=1,y=0,d=a;
	else d=exgcd(b,a%b,y,x),y-=a/b*x;
	return d;
}
signed main(){
	ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
	cin>>n;
	for(int i=1;i<=n;i++){
		cin>>b[i]>>a[i];
		fac*=b[i];
	}
	for(int i=1;i<=n;i++){
		int c=fac/b[i],k,t;
		exgcd(c,b[i],k,t);
		(ans+=(__int128)c*k*a[i])%=fac;
	}
	cout<<(long long)((ans+fac)%fac)<<"\n";
	return 0;
}

扩展中国剩余定理(exCRT)

CRT不能解决\(b_i\)不两两互质的线性同余方程组。

这是因为模\(p\)意义下\(a\)的逆元存在\(\iff \gcd(a,p)=1\),即\(a,p\)互质。

而如果存在\(j\neq i\)使得\(\gcd(b_i,b_j)\neq 1\),则容易知道\(\gcd(b_i,c_i)\neq 1\),自然就不存在\(k_i\)

这一缺陷来源于CRT的核心思想:构造\(x_i\)使得\(\begin{aligned}\begin{cases} x_i\equiv 1\pmod{b_k} & k=i\\ x_i\equiv 0\pmod{b_k} & k\neq i \end{cases}\end{aligned}\),这就注定了我们要使用乘法逆元。

exCRT则可以处理\(b_i\)不两两互质下的问题。

所谓exCRT,其实与CRT本身没有什么关系,它采用的是另一种思想——将同余式进行合并。

我们考虑\(2\)个同余式构成的方程组:

\[\begin{cases} x\equiv a_1\pmod{b_1}\\ x\equiv a_2\pmod{b_2} \end{cases}\]

要将其合并为\(x\equiv a\pmod b\)的形式。

我们将同余式重新表示为:\(x=a_1+k_1\times b_1=a_2+k_2\times b_2(k_1,k_2\in \mathbb{Z})\)

变形可得:\(k_1\times b_1-k_2\times b_2=a_2-a_1\)

这是一个典型的不定方程\(ax+by=c\)的形式,未知量是\((k_1,k_2)\),可以用exGCD求解。

如何解$ax+by=c$

具体来说,该方程有解\(\iff \gcd(a,b)\mid c\),必要性显然,充分性由裴蜀定理可得。

然后我们用exGCD求出\(ax+by=\gcd(a,b)\)的解,再将\(x,y\)同时乘\(\frac{c}{\gcd(a,b)}\)即可。


若存在\(ax+by=c\),则可以根据特解\(x,y\)求出任意通解\(x',y'\)
\(\begin{cases} x'=x+k\times \frac{b}{\gcd(a,b)}\\ y'=y-k\times \frac{a}{\gcd(a,b)} \end{cases}(k\in \mathbb{Z})\)

若解得特解为\((k_1',k_2')\),通解\((k_1,k_2)\)即为\((k_1'+k\times\frac{b_2}{\gcd(b_1,b_2)},k_2'+k\times\frac{b_1}{\gcd(b_1,b_2)})\),其中\(k\in \mathbb{Z}\)

选择其一带入:

\[\begin{aligned} x&=a_1+k_1\times b_1\\ &=a_1+(k_1'+k\times\frac{b_2}{\gcd(b_1,b_2)})\times b_1\\ &=a_1+k_1'\times b_1+k\times\text{lcm}(b_1,b_2) \end{aligned}\]

即:

\[x\equiv a_1+k_1'\times b_1\pmod{\text{lcm}(b_1,b_2)} \]

然后重复上述步骤,直到所有式子合并为\(1\)个为止。

点击查看代码
#include<bits/stdc++.h>
#define int __int128
using namespace std;
const int N=1e5+10;
long long tmp;
int n,a[N],b[N];
int exgcd(int a,int b,int& x,int& y){
	int d;
	if(b==0) x=1,y=0,d=a;
	else d=exgcd(b,a%b,y,x),y-=a/b*x;
	return d;
}
signed main(){
	ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
	cin>>tmp,n=tmp;
	for(int i=1;i<=n;i++) cin>>tmp,b[i]=tmp,cin>>tmp,a[i]=tmp;
	int pb=b[1],pa=a[1],k1,k2;
	for(int i=2;i<=n;i++){
		int g=exgcd(pb,b[i],k1,k2);//实际上应该是-b[i],为了避免负数所以这样写,注意此时k2要取相反数
		//if((a[i]-pa)%g) cout<<"no solution\n";//保证有解
		k1=(a[i]-pa)/g*k1;
		pa=pa+pb*k1;
		pb=pb/g*b[i];
		pa%=pb;
	}
	tmp=(pa+pb)%pb;
	cout<<tmp<<"\n";
	return 0;
}

(实际上,也有在CRT的基础上动“小手术”来解决exCRT的方法,见此文,我没有学就不放了)

例题

会不断添加。

P3868 [TJOI2009] 猜数字

变形一下:

\[b_i\mid (n-a_i)\\ n-a_i\equiv 0\pmod{b_i}\\ n\equiv a_i\pmod{b_i}\]

然后是CRT的板子了。

点击查看代码
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=15;
int n,fac=1,a[N],b[N];
__int128 ans;
int exgcd(int a,int b,int& x,int& y){
	int d;
	if(b==0) x=1,y=0,d=a;
	else d=exgcd(b,a%b,y,x),y-=a/b*x;
	return d;
}
signed main(){
	ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
	cin>>n;
	for(int i=1;i<=n;i++) cin>>a[i];
	for(int i=1;i<=n;i++) cin>>b[i],fac*=b[i];
	for(int i=1;i<=n;i++){
		int c=fac/b[i],k,t;
		exgcd(c,b[i],k,t);
		(ans+=(__int128)c*k*a[i])%=fac;
	}
	cout<<(long long)((ans+fac)%fac)<<"\n";
	return 0;
}

P4774 [NOI2018] 屠龙勇士

我们发现,用于对付每条龙的剑是固定的,因此我们预处理出对付第\(i\)条龙的剑的攻击力,记为\(b_i\)

\(i\)条龙能被击杀,当且仅当\(p_i\mid a_i-b_ix\),也即\(b_ix\equiv a_i\pmod{p_i}\)

所以我们实际上是要求解下面的线性同余方程组:

\[\begin{cases} b_1x\equiv a_1\pmod{p_1}\\ b_2x\equiv a_2\pmod{p_2}\\ \dots\\ b_nx\equiv a_n\pmod{p_n}\\ \end{cases}\]

我们发现,相比之前的方程,此题在\(x\)前面加了一个系数\(b_i\)

由于\(p\)不互质,所以采用exCRT。

首先求\(b_i\)在模\(p_i\)意义下的逆元是不可行的,因为\(b_i,p_i\)未必互质。

让我们看看能否沿用之前的思路来分析。

合并成带参的形式

仍然取两个方程,试图将它们合并。

\[\begin{cases} b_1x\equiv a_1\pmod{p_1}\\ b_2x\equiv a_2\pmod{p_2} \end{cases}\]

\[\begin{cases} b_1x=a_1+p_1k_1\\ b_2x=a_2+p_2k_2 \end{cases}(k_1,k_2\in \mathbb{Z})\]

发现\(x\)前系数不同,无法直接合并,因此不妨让系数变成\(\text{lcm}(b_1,b_2)=L\)

\[\begin{cases} Lx=\frac{L}{b_1}(a_1+p_1k_1)\\ Lx=\frac{L}{b_2}(a_2+p_2k_2) \end{cases}\]

统一系数后,就和exCRT基本相同了,可以类比上面的分析。

将两式联立并整理,得到:

\[\frac{L}{b_1}(a_1+p_1k_1)=\frac{L}{b_2}(a_2+p_2k_2)\\ k_1\frac{L}{b_1}p_1-k_2\frac{L}{b_2}p_2=\frac{L}{b_2}a_2-\frac{L}{b_1}a_1\]

解线性方程组得到特解\((k_1',k_2')\),根据上面的结论,通解\((k_1,k_2)\)可以表示为\((k\in \mathbb{Z})\)

\[\large(k_1'+k\frac{\frac{L}{b_2}p_2}{\gcd(\frac{L}{b_1}p_1,\frac{L}{b_2}p_2)},k_2'+k\frac{\frac{L}{b_1}p_1}{\gcd(\frac{L}{b_1}p_1,\frac{L}{b_2}p_2)}) \]

任选其一带入:

\[\begin{aligned} \large Lx&=\frac{L}{b_1}(a_1+p_1k_1)\\ &=\dots\\ &=\frac{L}{b_1}(a_1+p_1k_1)+k\times \text{lcm}(\frac{L}{b_1}p_1,\frac{L}{b_2}p_2) \end{aligned}\]

变成同余式的形式:

\[Lx\equiv \frac{L}{b_1}(a_1+p_1k_1) \pmod{\text{lcm}(\frac{L}{b_1}p_1,\frac{L}{b_2}p_2)} \]

至此完成合并。

最终得到的是一个线性同余方程\(bx\equiv a\pmod p\),理论上应该可以化成不定方程用exGCD来解。不过代码到现在还没写出来,有一些奇怪的问题(^^;

合并成不带参的形式

我们当然更希望最终合并成的同余方程式中,\(x\)系数是\(1\)

因此假设之前的方程式已经合并成了该形式,现在向其中添加一个带系数的同余方程式:

\[\begin{aligned}\begin{cases} x&\equiv a\pmod p\\ b_ix&\equiv a_i\pmod{p_i} \end{cases}\end{aligned}\]

\[\begin{aligned}\begin{cases} x&=a+k_1p&(1)\\ b_ix&=a_i+k_2p_i&(2) \end{cases}\end{aligned}\]

\((1)\)式乘\(b_i\),两式联立并整理,得到:

\[k_1b_ip-k_2p_i=a_i-ab_i \]

exGCD解不定方程得到\((k_1',k_2')\),通解\((k_1,k_2)\)表示为\((k\in\mathbb{Z})\)

\[(k_1'+k\frac{p_i}{\gcd(b_ip,p_i)},k_2+k\frac{b_ip}{\gcd(b_ip,p_i)}) \]

由于我们想要\(x\)的系数为\(1\),所以选择\(k_1\)带入\((1)\)式:

\[x=a+k_1p+k\frac{p_ip}{\gcd(b_ip,p_i)} \]

变成同余式的形式:

\[x\equiv a+k_1p\pmod{\frac{p_ip}{\gcd(b_ip,p_i)}} \]

至此完成合并。

点击查看代码
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int N=1e5+10;
int t,n,m,a[N],p[N],z[N],b[N];
multiset<int> se;
int exgcd(int a,int b,int& x,int& y){
	int d;
	if(b==0) x=1,y=0,d=a;
	else d=exgcd(b,a%b,y,x),y-=a/b*x;
	return d;
}
int solve(){
	se.clear();
	cin>>n>>m;
	for(int i=1;i<=n;i++) cin>>a[i];
	for(int i=1;i<=n;i++) cin>>p[i];
	for(int i=1;i<=n;i++) cin>>z[i];
	for(int i=1,x;i<=m;i++) cin>>x,se.insert(x);
	int pp=1,pa=0,C,g,k1,k2,mx=0,tmp;
	for(int i=1;i<=n;i++){
		auto it=se.upper_bound(a[i]);
		if(it!=se.begin()) --it;
		b[i]=*it,se.erase(it),se.insert(z[i]);
		mx=max(mx,(a[i]-1)/b[i]+1);
	}
	for(int i=1;i<=n;i++){
		g=exgcd((__int128)b[i]*pp%p[i],p[i],k1,k2);
		C=(a[i]-b[i]*pa)%p[i];
		if(C%g) return -1;
		pa+=(__int128)k1*(C/g)%(p[i]/g)*pp%(tmp=p[i]/g*pp);
		pa%=(pp=tmp);
	}
	if(pa<mx) pa+=((mx-pa-1)/pp+1)*pp;
	return pa;
}
signed main(){
	ios::sync_with_stdio(0),cin.tie(0),cout.tie(0);
	cin>>t;
	while(t--) cout<<solve()<<"\n";
	return 0;
}

\[\mathbb{\mid TO\ BE\ CONTINUED\mid} \]

posted @ 2025-08-03 15:36  Sinktank  阅读(67)  评论(0)    收藏  举报
★CLICK FOR MORE INFO★ TOP-BOTTOM-THEME
Enable/Disable Transition
Copyright © 2023 ~ 2025 Sinktank - 1328312655@qq.com
Illustration from 稲葉曇『リレイアウター/Relayouter/中继输出者』,by ぬくぬくにぎりめし.