【模板】擴展中國剩余定理(EXCRT)
算法介紹
孫子定理是中國古代求解一次同余式方程組的方法。是數論中一個重要定理。又稱中國余數定理。一元線性同余方程組問題最早可見于中國南北朝時期的數學著作《孫子算經》卷下第二十六題,叫做“物不知數”問題,原文如下:
有物不知其數,三三數之剩二,五五數之剩三,七七數之剩二。問物幾何?即,一個整數除以三余二,除以五余三,除以七余二,求這個整數?!秾O子算經》中首次提到了同余方程組問題,以及以上具體問題的解法,因此在中文數學文獻中也會將中國剩余定理稱為孫子定理。
用現代數學的語言來分析這個問題,中國剩余定理給出了以下的一元線性同余方程組:
在中國剩余定理中給出了 \(m_1,m_2,\dots,m_k\) 兩兩互質的條件,但是在擴展中國剩余定理中并沒有這個條件,相較于前者,后者更難解決。
問題簡述
問題就是,給定一個 \(k\) 個方程的線性同余方程組:
其中 \(m_1,m_2,\dots,m_k\) 不一定兩兩互質。
解題思路
我們的大致解題思路為將 \(2\) 個方程合并為一個新的方程,以此類推,最終我們會得到一個 \(x\equiv y\pmod z\) 的一個方程,易見上面的方程組的最小正整數解就是 \(y\)。
正確性證明
接下來我們來解決合并方程的問題,我們考慮如下兩個方程:
我們根據第一個式子可以寫出 \(x\) 的通解 \(x=a_1+m_1\times k\) 其中 \(k\) 為任意整數,我們將這個通解帶入第二個式子就可以得到 \(a_1+m_1\times k\equiv a_2\pmod {m_2}\) 我們移一下項就可以得到 \(m_1\times k\equiv a_2-a_1\pmod {m_2}\),這就是上面的方程組合并后的結果。
而這個方程有解的充要條件是 \(\gcd(m_1,m_2)\mid a_2-a_1\),這個其實就是裴蜀定理,這里不再概述。
我們繼續講,我們得到這個充要條件后我們可以判斷這個方程是否有解,如果有解我們就繼續進行接下來的操作。
我們設 \(d=\gcd(m_1,m_2)\),然后將我們合并的方程變換一下就是:
然后,我們設 \(m_1'=\frac{m_1}w0obha2h00,c=\frac{a_2-a_1}w0obha2h00,m_2'=\frac{m_2}w0obha2h00\) 于是我們就有:
注意到此時 \(m_1',m_2'\) 互質,所以 \(m_1'\) 在模 \(m_2'\) 的意義下存在乘法逆元,我們可以使用擴展歐幾里得算法來求出逆元,即求出整數 \(inv\) 使得 \(m_1'\times inv\equiv 1\pmod {m_2'}\),所以我們繼續將這個方程變換就變成了:
如果我們記 \(k_0=c\times inv\) 則 \(k\) 的通解為 \(k_0+m_2'\times t\) 其中 \(t\) 為任意整數。
然后我們將這個 \(k\) 帶回一開始的式子就可以得出:
我們設 \(x_0=a_1+m_1\times k_0,L=\mathrm{lcm}(m_1,m_2)\) 所以我們就愉快地得出了:
于是,我們完成了合并方程的使命!
最后其實就是一個遞推的過程我們一次合并前 \(2\) 個方程,最后就能得到答案!
代碼實現
#include<bits/stdc++.h>
#define LL __int128
#define R register
using namespace std;
namespace fastIO{char *p1,*p2,buf[100000];
#define nc() (p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++)
inline void read(LL&n){LL x=0,f=1;char ch=nc();while(ch<48||ch>57){if(ch=='-'){f=-1;}ch=nc();}while(ch>=48&&ch<=57){x=(x<<3)+(x<<1)+(ch^48),ch=nc();}n=x*f;}inline void read(string&s){s="";char ch=nc();while(ch==' '||ch=='\n'){ch=nc();}while(ch!=' '&&ch!='\n'){s+=ch,ch=nc();}}inline void read(char&ch){ch=nc();while(ch==' '||ch=='\n'){ch=nc();}}inline void write(LL x){if(x<0){putchar('-'),x=-x;}if(x>9){write(x/10);}putchar(x%10+'0');return;}inline void write(const string&s){for(R LL i=0;i<(int)s.size();i++){putchar(s[i]);}}inline void write(const char&c){putchar(c);}
}using namespace fastIO;
inline LL mul(LL a,LL b,const LL&mod){
a=(a%mod+mod)%mod;
b=(b%mod+mod)%mod;
LL res=0;
while(b){
if(b&1)res=(res+a)%mod;
a=(a+a)%mod;
b>>=1;
}
return res;
}
void exgcd(LL a,LL b,LL&x,LL&y){
if(b==0){
x=1;
y=0;
}
else{
exgcd(b,a%b,y,x);
y-=x*(a/b);
}
}
LL inv_mod(LL a,LL m){
LL x,y;
exgcd(a,m,x,y);
return (x%m+m)%m;
}
LL gcd(LL a,LL b){
return b?gcd(b,a%b):a;
}
LL n,a[100005],b[100005];
signed main(){
read(n);
for(int i=0;i<n;i++){
read(a[i]);
read(b[i]);
}
LL a0=a[0];
LL b0=(b[0]%a0+a0)%a0;
for(int i=1;i<n;i++){
LL ai=a[i];
LL bi=(b[i]%ai+ai)%ai;
LL d=gcd(a0,ai);
LL dif=bi-b0;
LL a0_=a0/d;
LL ai_=ai/d;
LL dif_=dif/d;
LL c=(dif_%ai_+ai_)%ai_;
LL inv=inv_mod(a0_,ai_);
LL t0=mul(inv,c,ai_);
LL a0__=(a0/d)*ai;
LL mod__=a0__;
LL p=mul(a0,t0,mod__);
LL b0__=(b0+p)%mod__;
a0=mod__;
b0=b0__;
}
write(b0);
return 0;
}

浙公網安備 33010602011771號