CRT&EXCRT学习笔记
学\(EXCRT\)的时间挺久了,有点忘了.
写一篇博客记录一下.
CRT
首先,我们要知道中国剩余定理是什么
它是用来求解这样一个同余方程的
其中,\(p_1,p_2...p_n\)两两互质.
我们先规定一些数.\(M=\prod_{i=1}^np_i,L=lcm(p_1,p_2...p_n)\)
\(CRT\)的流程是这样的.
假设现在执行到\(i\)
令\(m=\frac{M}{p_i}\)
找到一个\(x\),使得\(mx\equiv 1(mod\ p_i)\)
找\(x\)可以用逆元的方法实现.
那么,显然地,\(a_imx\equiv a_i(mod\ p_i)\)
我们令\(t_i=a_imx\)
那么,答案就是\((\sum_{i=1}^nt_i)\% L\)
这个证明显然吧
EXCRT
有的时候,出题人并不保证所有的\(p_i\)互质.
那么我们该怎么办呢?
其实也不难.
假设有\(2\)个方程\(x\equiv c_1(mod\ p_1),x\equiv c_2(mod\ p_2)\)
我们可以列出相应的方程.
\(x=c_1+k_1p_1,x=c_2+k_2p_2\)
那么\(c_1+k_1p_1=c_2+k_2p_2\)
\(k_1p_1=c_2-c_1+k_2p_2\)
令\(g=gcd(p_1,p_2)\),那么\(\frac{p_1}{g}k_1=\frac{c_2-c_1}{g}+\frac{p_2}{g}k_2\).
因此只有\(g\)能整除\(c_2-c_1\)时才可以合并,否则无解.
此时,可以化简方程为\(\frac{p_1}{g}k_1\equiv\frac{c_2-c_1}{g}(mod\ \frac{p_2}{g})\).
由于\((\frac{p_1}{g},\frac{p_2}{g})=1\),我们可以用\(exgcd\)求出\(\frac{p_1}{g}\)在模\(\frac{p_2}{g}\)意义下的逆元,假设为\(inv\).
那么\(k_1\equiv\frac{c_2-c_1}{g}*inv(mod\ \frac{p_2}{g})\).
我们把这个方程回代到第一个方程中,得到
\(x\equiv c_1+(\frac{c_2-c_1}{g}*inv\% \frac{p_2}{g})(mod\ \frac{p_1*p_2}{g})\)
但是,这个方程中所有的量都是已知的,因此我们就一直合并下去即可.
代码如下
#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>
#include<vector>
#define N (100010)
#define inf (0x7f7f7f7f)
#define rg register int
#define Label puts("NAIVE")
#define spa print(' ')
#define ent print('\n')
#define rand() (((rand())<<(15))^(rand()))
typedef long double ld;
typedef long long LL;
typedef unsigned long long ull;
using namespace std;
inline char read(){
static const int IN_LEN=1000000;
static char buf[IN_LEN],*s,*t;
return (s==t?t=(s=buf)+fread(buf,1,IN_LEN,stdin),(s==t?-1:*s++):*s++);
}
template<class T>
inline void read(T &x){
static bool iosig;
static char c;
for(iosig=false,c=read();!isdigit(c);c=read()){
if(c=='-')iosig=true;
if(c==-1)return;
}
for(x=0;isdigit(c);c=read())x=((x+(x<<2))<<1)+(c^'0');
if(iosig)x=-x;
}
inline char readchar(){
static char c;
for(c=read();!isalpha(c);c=read())
if(c==-1)return 0;
return c;
}
const int OUT_LEN = 10000000;
char obuf[OUT_LEN],*ooh=obuf;
inline void print(char c) {
if(ooh==obuf+OUT_LEN)fwrite(obuf,1,OUT_LEN,stdout),ooh=obuf;
*ooh++=c;
}
template<class T>
inline void print(T x){
static int buf[30],cnt;
if(x==0)print('0');
else{
if(x<0)print('-'),x=-x;
for(cnt=0;x;x/=10)buf[++cnt]=x%10+48;
while(cnt)print((char)buf[cnt--]);
}
}
inline void flush(){fwrite(obuf,1,ooh-obuf,stdout);}
int n; LL c,m,c2,mo;
LL mult(LL a,LL b,LL mod){//a*b
LL res=0,fu=1;
if(a<0)fu=-fu,a=-a;
if(b<0)fu=-fu,b=-b;
while(b){
if(b&1)res=(res+a)%mod;
a=(a+a)%mod,b>>=1;
}
res*=fu;
if(res<0)(res+=((-res-1)/mod+1)*mod);
return res;
}
void exgcd(LL a,LL b,LL &x,LL &y){
if(!b){x=1,y=0;return;}
exgcd(b,a%b,y,x),y-=(a/b)*x;
}
LL gcd(LL a,LL b){return b?gcd(b,a%b):a;}
LL inv(LL a,LL b){
LL x,y;
exgcd(a,b,x,y);
return (x<=0)?(x+b):x;
}
int main(){
read(n),read(m),read(c);
for(int i=1;i<n;i++){
read(mo),read(c2);
LL t=gcd(m,mo),s=inv(m/t,mo/t),tc=c,tm=m;
m=(mo/t*tm),c=(tc+mult(tm,mult(s,(c2-c)/t,(mo/t)),m))%m;
}
if(c<0)c+=((-c-1)/m+1)*m;
printf("%lld\n",c%m);
}
