雖然中國剩餘定理和拓展中國剩餘定理只差兩個字,但他倆的解法相差十萬八千里,因此會不會CRT無所謂c++
求相似$$\begin{cases}x \equiv b_{1}\pmod{a_{1}} \\x \equiv b_{2}\pmod{a_{2}} \\...\\x \equiv b_{n}\pmod{a_{n}} \\ \end{cases}$$的線性同餘方程組的解git
假設如今咱們只有兩個同餘方程$$x \equiv b_{1}\pmod{a_{1}}$$ $$x \equiv b_{2}\pmod{a_{2}}$$ide
改寫它們得$$x=a_{1}k_{1}+b_{1} ①$$ $$x=a_{2}k_{2}+b_{2}$$ui
聯立$$a_{1}k_{1}+b_{1}=a_{2}k_{2}+b_{2}$$spa
因此$$a_{1}k_{1}=b_{2}-b_{1}+a_{2}k_{2}$$code
記$$c=gcd(a_{1},a_{2})$$ $$d_{1}=\frac{a_{1}}{c}$$ $$d_{2}=\frac{a_{2}}{c}$$blog
兩邊同除以$c$得$$d_{1}k_{1}=\frac{b_{2}-b_{1}}{c}+d_{2}k_{2}$$get
即$$d_{1}k_{1} \equiv \frac{b_{2}-b_{1}}{c}\%d_{2}\pmod{d_{2}}$$it
顯然,方程組有解的條件是$c|(b_{2}-b_{1})$io
爲何要在這個時候取模$\%d_{2}$?由於以後取模的時候模數和被模數已經不互質了
記$d_{1}$在模$d_{2}$意義下的逆元爲$p$,兩邊同乘$p$得$$k_{1} \equiv \frac{p(b_{2}-b_{1})}{c}\%d_{2} \pmod{d_{2}}$$
即$$k_{1} = \frac{p(b_{2}-b_{1})}{c}\%d_{2}+d_{2}*y$$
帶入①式得$$x=\frac{p(b_{2}-b_{1})}{c}\%d_{2}*a_{1}+d_{2}a_{1}y+b_{1}$$
即$$x \equiv \frac{p(b_{2}-b_{1})}{c}\%d_{2}*a_{1}+b_{1}\pmod{a_{1}d_{2}}$$
注意,上式中的$a_{1}$萬萬不可乘到$\%d_{2}$前面
因此咱們就弄出了一個新的同餘方程。這個方程也能夠寫成$$x \equiv b \pmod a$$其中$$b=(inv(\frac{a_{1}}{(a_{1},a_{2})},\frac{a_{2}}{(a_{1},a_{2})})*\frac{b_{2}-b_{1}}{(a_{1},a_{2})}\%\frac{a_{2}}{(a_{1},a_{2})}*a_{1}+b_{1})$$ $$a=\frac{a_{1}a_{2}}{(a_{1},a_{2})}$$
其中$inv(m,n)$表示$m$在模$n$意義下的逆元,$(m,n)$表示$m$和$n$的最大公約數
而後把新求出的方程和後面的方程聯立,迭代求解,直到只剩一個方程
毫無疑問的是,在求解方程組的時候很容易爆$long long$,那該怎麼辦呢?
__int128 快(龜)速(速)乘(乘)是個不錯的主意。
何爲快速乘?
就是一種和快速冪相似的東西,利用二進制在log時間內求出兩個數相乘對另外一個數取模的結果而且不會爆long long
但要注意的是,千萬不要讓負數出如今快速乘裏面,尤爲是被拆分的那個乘數,會死循環的
1 int mul(int x,int k,int mod) { 2 ll ans=0; 3 while(k) { 4 if(k&1)(ans+=x)%=mod; 5 k>>=1; 6 (x<<=1)%=mod; 7 } 8 return ans; 9 }
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define INF 0x7fffffff 4 #define ME 0x7f 5 #define FO(s) freopen(s".in","r",stdin);freopen(s".out","w",stdout) 6 #define fui(i,a,b,c) for(int i=(a);i<=(b);i+=(c)) 7 #define fdi(i,a,b,c) for(int i=(a);i>=(b);i-=(c)) 8 #define fel(i,a) for(register int i=hd[a];i;i=dg[i].nxt) 9 #define ll long long 10 #define int long long 11 #define MEM(a,b) memset(a,b,sizeof(a)) 12 #define maxn 100010 13 ll n; 14 ll a[maxn],b[maxn],a0,b0; 15 template<class T> 16 inline T read(T &n){ 17 n=0;int t=1;double x=10;char ch; 18 for(ch=getchar();!isdigit(ch)&&ch!='-';ch=getchar());(ch=='-')?t=-1:n=ch-'0'; 19 for(ch=getchar();isdigit(ch);ch=getchar()) n=n*10+ch-'0'; 20 if(ch=='.') for(ch=getchar();isdigit(ch);ch=getchar()) n+=(ch-'0')/x,x*=10; 21 return (n*=t); 22 } 23 template<class T> 24 T write(T n){ 25 if(n<0) putchar('-'),n=-n; 26 if(n>=10) write(n/10);putchar(n%10+'0');return n; 27 } 28 template<class T> 29 T writeln(T n){write(n);putchar('\n');return n;} 30 int exgcd(int a,int b,int &x,int &y){ 31 if(!b)return a+((x=1)*(y=0)); 32 int xx,yx,z;z=exgcd(b,a%b,xx,yx); 33 x=yx,y=xx-a/b*yx;return z; 34 } 35 int gcd(int a,int b){if(!b)return a;return gcd(b,a%b);} 36 int mul(int x,int k,int mod){ 37 ll ans=0; 38 while(k){ 39 if(k&1)(ans+=x)%=mod; 40 k>>=1; 41 (x<<=1)%=mod; 42 } 43 return ans; 44 } 45 signed main(){ 46 read(n); 47 fui(i,1,n,1)read(a[i]),read(b[i]); 48 a0=a[1],b0=b[1]; 49 fui(i,2,n,1){ 50 int c=gcd(a0,a[i]),d1=a0/c,d2=a[i]/c,p,x,y; 51 int a1=a0,a2=a[i],b1=b0,b2=b[i]; 52 if(((b2-b1)/c*c)^(b2-b1))return 0*writeln(-1); 53 exgcd(d1,d2,p,y);((p%=d2)+=d2)%=d2; 54 x=mul(p,(((b2-b1)/c)%d2+d2)%d2,d2)*a1+b1; 55 y=a1*d2;b0=x,a0=y;((b0%=a0)+=a0)%=a0; 56 } 57 int x,y;exgcd(1,a0,x,y); 58 ((x%=a0)+=a0)%=a0; 59 x=mul(x,b0,a0); 60 writeln(x); 61 return 0; 62 }
好像還有NOI2018屠龍勇士?23333先等我去A掉它