算法介绍
孙子定理是中国古代求解一次同余式方程组的方法。是数论中一个重要定理。又称中国余数定理。一元线性同余方程组问题最早可见于中国南北朝时期的数学著作《孙子算经》卷下第二十六题,叫做“物不知数”问题,原文如下:
有物不知其数,三三数之剩二,五五数之剩三,七七数之剩二。问物几何?即,一个整数除以三余二,除以五余三,除以七余二,求这个整数。《孙子算经》中首次提到了同余方程组问题,以及以上具体问题的解法,因此在中文数学文献中也会将中国剩余定理称为孙子定理。
用现代数学的语言来分析这个问题,中国剩余定理给出了以下的一元线性同余方程组:
在中国剩余定理中给出了 \(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}{d},c=\frac{a_2-a_1}{d},m_2'=\frac{m_2}{d}\) 于是我们就有:
注意到此时 \(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;
}