【模板】扩展中国剩余定理(EXCRT)

算法介绍

孙子定理是中国古代求解一次同余式方程组的方法。是数论中一个重要定理。又称中国余数定理。一元线性同余方程组问题最早可见于中国南北朝时期的数学著作《孙子算经》卷下第二十六题,叫做“物不知数”问题,原文如下:

有物不知其数,三三数之剩二,五五数之剩三,七七数之剩二。问物几何?即,一个整数除以三余二,除以五余三,除以七余二,求这个整数。《孙子算经》中首次提到了同余方程组问题,以及以上具体问题的解法,因此在中文数学文献中也会将中国剩余定理称为孙子定理。

用现代数学的语言来分析这个问题,中国剩余定理给出了以下的一元线性同余方程组:

[left{begin{matrix} xequiv a_1pmod {m_1}\xequiv a_2pmod {m_2} \vdots \xequiv a_kpmod {m_k} end{matrix}right.]

在中国剩余定理中给出了 (m_1,m_2,dots,m_k) 两两互质的条件,但是在扩展中国剩余定理中并没有这个条件,相较于前者,后者更难解决。

问题简述

问题就是,给定一个 (k) 个方程的线性同余方程组:

[left{begin{matrix} xequiv a_1pmod {m_1}\xequiv a_2pmod {m_2} \vdots \xequiv a_kpmod {m_k} end{matrix}right.]

其中 (m_1,m_2,dots,m_k) 不一定两两互质。

解题思路

我们的大致解题思路为将 (2) 个方程合并为一个新的方程,以此类推,最终我们会得到一个 (xequiv ypmod z) 的一个方程,易见上面的方程组的最小正整数解就是 (y)

正确性证明

接下来我们来解决合并方程的问题,我们考虑如下两个方程:

[left{begin{matrix} xequiv a_1pmod {m_1}\xequiv a_2pmod {m_2} end{matrix}right.]

我们根据第一个式子可以写出 (x) 的通解 (x=a_1+m_1times k) 其中 (k) 为任意整数,我们将这个通解带入第二个式子就可以得到 (a_1+m_1times kequiv a_2pmod {m_2}) 我们移一下项就可以得到 (m_1times kequiv a_2-a_1pmod {m_2}),这就是上面的方程组合并后的结果。

而这个方程有解的充要条件是 (gcd(m_1,m_2)mid a_2-a_1),这个其实就是裴蜀定理,这里不再概述。

我们继续讲,我们得到这个充要条件后我们可以判断这个方程是否有解,如果有解我们就继续进行接下来的操作。

我们设 (d=gcd(m_1,m_2)),然后将我们合并的方程变换一下就是:

[frac{m_1times k}{d}equiv frac{a_2-a_1}{d}pmod {frac{m_2}{d}} ]

然后,我们设 (m_1'=frac{m_1}{d},c=frac{a_2-a_1}{d},m_2'=frac{m_2}{d}) 于是我们就有:

[m_1'times kequiv cpmod {m_2'} ]

注意到此时 (m_1',m_2') 互质,所以 (m_1') 在模 (m_2') 的意义下存在乘法逆元,我们可以使用扩展欧几里得算法来求出逆元,即求出整数 (inv) 使得 (m_1'times invequiv 1pmod {m_2'}),所以我们继续将这个方程变换就变成了:

[kequiv ctimes invpmod {m_2'} ]

如果我们记 (k_0=ctimes inv)(k) 的通解为 (k_0+m_2'times t) 其中 (t) 为任意整数。

然后我们将这个 (k) 带回一开始的式子就可以得出:

[begin{aligned} x&=a_1+m_1times(k_0+m_2'times t)\ &=(a_1+m_1times k_0)+(m_1times m_2')times t\ &=(a_1+m_1times k_0)+frac{m_1times m_2}{d}times t\ &=(a_1+m_1times k_0)+mathrm{lcm}(m_1,m_2)times tend{aligned}]

我们设 (x_0=a_1+m_1times k_0,L=mathrm{lcm}(m_1,m_2)) 所以我们就愉快地得出了:

[left{begin{matrix} xequiv a_1pmod {m_1}\xequiv a_2pmod {m_2} end{matrix}right.Longleftrightarrow xequiv x_0pmod L]

于是,我们完成了合并方程的使命!

最后其实就是一个递推的过程我们一次合并前 (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; } 

发表评论

评论已关闭。

相关文章