扩展中国剩余定理(EXCRT)
老师说NOIP会考,那还是复习&总结一下吧
求方程的最小正整数解:
假设现在已经求出了前 (i-1) 个方程的最小正整数解 (res) ,记 (lcm=operatorname{lcm}(b_1,b_2,cdots,b_{i-1}))
当加入第 (i) 个方程的时候,我们需要求出一个整数 (x) 满足 (res+x imes lcm equiv a_i pmod {b_i})
移项可得 (x imes lcm + y imes b_i = a_i - res)
这是一个扩展欧几里得的形式,方程组无解当且仅当存在一个扩展欧几里得方程无解,即 ((a_i - res) \% gcd(lcm, b_i) ot = 0)
否则根据小奥知识,(res' = res + x imes lcm \% lcm') ,(lcm') 是指加入 (b_i) 之后的 (lcm) 值
但是直接拿 (x) 去乘 (lcm) 很可能爆 long long
。
我们其实可以把 (x) 通过取模限制在 (dfrac{b_i}{gcd(b_i, lcm)}) 以内的,这样只要出题人保证所有数的 (lcm) 最后不爆 long long
就行了
具体方法如下:
首先先看一个化简版的(变量与上面无关)
这个其实很好理解,就是余式两边同时乘上 (y) 结果不变
那么回到本题
直接取模就行了
当然在实现的时候,由 (exgcd) 特性我们会先解出 (x_0,y_0) 满足 (x_0 imes lcm + y_0 imes b_i = gcd(lcm, b_i)) ,我们的 (x) 应该等于 (x_0 imes dfrac{a_i - res}{gcd(lcm, b_i)}) ,这一步要龟速乘或者光速乘,可以在这一步取模。
但是由于 (exgcd) 的缘故 (EXCRT) 必然带 (log) ,用光速乘意义不大。
另外,这东西严格吊打 CRT
啊,到现在为止还没见过 CRT
能做 EXCRT
不能做的事。
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
typedef double db;
#define x first
#define y second
#define sz(v) (int)v.size()
#define pb(x) push_back(x)
#define mkp(x,y) make_pair(x,y)
#define int long long//------------------Warning!!!-------------------------------------------
//#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
//char buf[1<<21],*p1=buf,*p2=buf;
inline int read(){
int x=0,f=1;char c=getchar();
while(!isdigit(c)){if(c=='-')f=0;c=getchar();}
while(isdigit(c))x=x*10+c-'0',c=getchar();
return f?x:-x;
}
#define N 100005
int n,a[N],b[N];
void exgcd(int&x,int&y,int&d,int a,int b){
if(!b)return x=1,y=0,d=a,void();
exgcd(y,x,d,b,a%b),y-=a/b*x;
}
int mul(int x,int y,int mod){
// int res=x*y-(long long)((long double)x/mod*y+0.5)*mod;
// return res<0?res+mod:res;//光速乘
int res=0;
while(y){
if(y&1)(res+=x)%=mod;
y>>=1,(x<<=1)%=mod;
}
return res;
}//龟速乘
int excrt(){
int res=a[1],lcm=b[1];
for(int i=2;i<=n;++i){
int x,y,d,t=((a[i]-res)%b[i]+b[i])%b[i];
exgcd(x,y,d,lcm,b[i]);
if(t%d)return -1;
res+=lcm*mul(x,t/d,b[i]/d);
lcm=lcm/d*b[i];
res=(res+lcm)%lcm;
}
return res;
}
signed main(){
n=read();
for(int i=1;i<=n;++i)b[i]=read(),a[i]=read();
printf("%lld
",excrt());
return 0;
}