CRT在算法竞赛中算是一个比较重要的模块,他的基本形式如下:
给出n个式子:
x≡a1 mod p1 x ≡ a 1 m o d p 1
x≡a2 mod p2 x ≡ a 2 m o d p 2
x≡a3 mod p3 x ≡ a 3 m o d p 3
… …
求x的最小正整数解。(pi互质)
那么如何求解x呢?
我们先考虑求解x的一个可行解 x0 x 0
因为同余的可加性,我们可以把 x0 x 0 拆成 x0=x1+x2+x3+…+xn x 0 = x 1 + x 2 + x 3 + … + x n ,每个 xi x i 对应的式子为:
xi≡ai mod pi x i ≡ a i m o d p i
对于不等于i的值j
xi≡0 mod pj x i ≡ 0 m o d p j
再由于同余的可乘性,我们可以把x拆成 x=a1x1+a2x2+…+anxn x = a 1 x 1 + a 2 x 2 + … + a n x n ,把式子改为:
xi≡1 mod pi x i ≡ 1 m o d p i
然后我们观察每个 xi x i :对于每个 xi x i 和不等于i的所有j,由于
xi≡0 mod pj x i ≡ 0 m o d p j
设所有p的乘积为P,那么 xi x i 一定是 Ppi P p i 的倍数,那么我们可以把 xi x i 改为 g(Ppi)≡1 mod pi g ( P p i ) ≡ 1 m o d p i ,通过这个式子可以看出g其实是 Ppi P p i 在模 pi p i 意义下的逆元, xi x i 就可以求出来, x0 x 0 也就求出来了, x0 x 0 的公式为:
x0=a1∗P1∗G1+a2∗P2∗G2+⋯+an∗Pn∗Gn x 0 = a 1 ∗ P 1 ∗ G 1 + a 2 ∗ P 2 ∗ G 2 + ⋯ + a n ∗ P n ∗ G n
其中 Pi P i 为 Ppi P p i , Gi G i 为 Pi P i 在模 pi p i 意义下的逆元
x的所有解为 x=x0+tP(t∈Z) x = x 0 + t P ( t ∈ Z )
需要注意的是这种方法只能用于p均为质数的情况下,如果p不是质数的话需要用其他方法,这不是CRT的内容,这里就不讲了,有兴趣的朋友可以去看看我前面写的一篇基础数论入门,那里面有所提及。
来道题练练手?
bzoj 1951-古代猪文
题目大意:
给出g,n,求 g∑d|nCdn g ∑ d | n C n d %999911659
Solution:
解决这个问题我们需要求出 ∑d|nCdn ∑ d | n C n d ,那么不妨设 x=∑d|nCdn x = ∑ d | n C n d
通过欧拉定理我们知道 gp−1≡1 mod p g p − 1 ≡ 1 m o d p
所以说 gx≡gx%(p−1) mod p g x ≡ g x % ( p − 1 ) m o d p
我们只需要求出 x%(p−1) x % ( p − 1 ) 即可,由于p-1不是一个质数,但是它可以分解成多个质数
( 999911658=2∗3∗4679∗35617 999911658 = 2 ∗ 3 ∗ 4679 ∗ 35617 )
那么这道题就可以用CRT了。
代码:
//999911658=2*3*4679*35617
#include<cstdio>
#include<cmath>
#include<iostream>
using namespace std;
int n,g;
int md[4]={
2,3,4679,35617};
int tot=999911658;
int ni[4][40000];
int mi[4][40000];
int nn=0;
int f[4];
int fast_pow(int x,int a,int mod)
{
if (x==mod) return 0;//特判g^0=0(因为指数是取模取出来的)
int ans=1;
for (;a;a>>=1,x=1ll*x*x%mod)
{
if (a&1) ans=1ll*ans*x%mod;
}
return ans;
}
int C(int n,int m,int mod)
{
if (n>m) return 0;
if (n<md[mod]&&m<md[mod]) return 1ll*mi[mod][m]*ni[mod][n]*ni[mod][m-n]%md[mod];
return C(n/md[mod],m/md[mod],mod)*C(n%md[mod],m%md[mod],mod)%md[mod];
}
int main()
{
scanf("%d%d",&n,&g);
for (int j=0;j<4;j++)
{
mi[j][0]=1;
for (int i=1;i<md[j];i++)
mi[j][i]=(mi[j][i-1]*i)%md[j];
}
for (int j=0;j<4;j++)
{
ni[j][md[j]-1]=fast_pow(mi[j][md[j]-1],md[j]-2,md[j]);
for (int i=md[j]-2;i>=1;i--)
{
ni[j][i]=ni[j][i+1]*(i+1)%md[j];
}
ni[j][0]=1;
}
//cout<<C(2,4,0)<<endl;
int qn=sqrt(n);
for (int i=1;i<=qn;i++)
{
if (n%i==0)
{
for (int j=0;j<4;j++)
{
f[j]=(f[j]+C(i,n,j))%md[j];
if (n/i!=i)
f[j]=f[j]+C(n/i,n,j)%md[j];
}
}
}
int ans=0;
for (int j=0;j<4;j++)
{
//cout<<f[j]<<endl;
ans=(ans+(1ll*f[j]*(tot/md[j]))%tot*fast_pow(tot/md[j],md[j]-2,md[j])%tot)%tot;
}
printf("%d",fast_pow(g,ans,tot+1));
}