老师说NOIP会考,那还是复习&总结一下吧
求方程的最小正整数解:
\[\begin{cases} x\equiv a_1 \pmod {b_1}\\ x\equiv a_2 \pmod {b_2}\\ \cdots\\ x\equiv a_n \pmod {b_n}\\ \end{cases} \]假设现在已经求出了前 \(i-1\) 个方程的最小正整数解 \(res\) ,记 \(lcm=\operatorname{lcm}(b_1,b_2,\cdots,b_{i-1})\)
当加入第 \(i\) 个方程的时候,我们需要求出一个整数 \(x\) 满足 \(res+x\times lcm \equiv a_i \pmod {b_i}\)
移项可得 \(x\times lcm + y \times b_i = a_i - res\)
这是一个扩展欧几里得的形式,方程组无解当且仅当存在一个扩展欧几里得方程无解,即 \((a_i - res) \% \gcd(lcm, b_i) \not = 0\)
否则根据小奥知识,\(res' = res + x \times lcm \% lcm'\) ,\(lcm'\) 是指加入 \(b_i\) 之后的 \(lcm\) 值
但是直接拿 \(x\) 去乘 \(lcm\) 很可能爆 long long
。
我们其实可以把 \(x\) 通过取模限制在 \(\dfrac{b_i}{\gcd(b_i, lcm)}\) 以内的,这样只要出题人保证所有数的 \(lcm\) 最后不爆 long long
就行了
具体方法如下:
首先先看一个化简版的(变量与上面无关)
\[xy \% yz = y \times (x \% z) \]这个其实很好理解,就是余式两边同时乘上 \(y\) 结果不变
那么回到本题
\[x \times lcm \% lcm' = x \times lcm \% \dfrac{lcm \times b_i}{\gcd(lcm, b_i)}=lcm \times (x \% \dfrac{b_i}{\gcd(lcm, b_i)}) \]直接取模就行了
当然在实现的时候,由 \(exgcd\) 特性我们会先解出 \(x_0,y_0\) 满足 \(x_0 \times lcm + y_0 \times b_i = \gcd(lcm, b_i)\) ,我们的 \(x\) 应该等于 \(x_0 \times \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\n",excrt());
return 0;
}