[总结] （扩展）中国剩余定理

转载

mb5fed6f751964c 2021-08-12 16:54:00

文章标签 数学数论 i++ #include 中国剩余定理 文章分类 代码人生

中国剩余定理

有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二。问物几何？

解线性方程组：

\(\left(\begin{matrix} x\equiv b_1 \pmod{a_1}\\x \equiv b_2 \pmod{a_2}\\ \cdots \\ x\equiv b_n \pmod{a_n}\end{matrix} \right)\)

解法

\(int \ lcj=a_1*a_2\cdots*a_n,m_i=lcj/a[i],m_i^{-1}=m_i\)关于\(a[i]\)的逆元

\(ans=\sum_{i=1}^nb[i]*m_i*m_i^{-1}\)

如果是最小的非负整数解，\(ans=(ans+lcj)mod\ lcj\)

P1495 【模板】中国剩余定理(CRT)/曹冲养猪

#include <iostream>
#include <cstdio>
#include <cstring>
#define int long long
using namespace std;
const int maxn = 15;
int a[maxn],b[maxn],lcj=1;
void exgcd(int a,int b,int &x,int &y){
if(!b){
x=1;y=0;
return ;
  }
exgcd(b,a%b,x,y);
int tmp=x;x=y;
y=tmp-(a/b)*y;
}
int n,ans=0;
signed main(){
scanf("%lld",&n);
for(int i=1;i<=n;i++)scanf("%lld%lld",a+i,b+i),lcj=1LL*lcj*a[i];
for(int i=1;i<=n;i++){
int mi=lcj/a[i];
int mi_=0,y0=0;
exgcd(mi,a[i],mi_,y0);
    (ans+=1LL*b[i]*mi*mi_)%=lcj;
  }
printf("%lld\n",(ans+lcj)%lcj);
return 0;
}

扩展中国剩余定理

P4777 【模板】扩展中国剩余定理（EXCRT）

与中国剩余定理的区别

扩展中国剩余定理用来解决模数不互质的情况。

解法

采用合并答案的思想。

如果我们解决了前 \(k-1\) 个同余方程，令 \(M=\Pi_{i=1}^{k-1}a[i]\) ，那么前 k-1 个方程的通解为 \(x=x+t*M,(t\in Z)\)。

现在考虑第 \(k\) 个方程。

\(x\equiv b_k \pmod{a_k}\)

等价于：

\(x+t* M\equiv b_k \pmod{a_k}\) 有解。

即 \(t* M \equiv b_k-x \pmod{a_k}\) 有解，可以解得 \(t\)。

现在考虑合并答案。

\(x=x+t* M,M=lcm(M,a[k])\)。

实现

需要用到快 \((gui)\) 速乘，可以防止溢出

int mul(int a,int b,int P){
int res=0;
while(b){
if(b&1)res=(res+a)%P;
a=(a+a)%P;
b>>=1;
  }
return res;
}

#include <iostream>
#include <cstdio>
#include <cstring>
#define int long long
using namespace std;
const int maxn = 1e5 + 10;
int n,ai[maxn],bi[maxn];
int mul(int a,int b,int P){
int res=0;
while(b){
if(b&1)res=(res+a)%P;
a=(a+a)%P;
b>>=1;
  }
return res;
}
int exgcd(int a,int b,int &x,int &y){
if(!b){
x=1;y=0;return a;
  }
int gcd=exgcd(b,a%b,x,y);
int tmp=x;x=y;
y=tmp-(a/b)*y;
return gcd;
}
int exCRT(){
int x,y;
int ans=bi[1],M=ai[1];
for(int i=2;i<=n;i++){
int a=M,b=ai[i],c=(bi[i]-ans%b+b)%b;
int gcd=exgcd(a,b,x,y),bg=b/gcd;
if(c%gcd!=0)return -1;

x=mul(x,c/gcd,bg);
ans+=x*M;
M*=bg;
ans=(ans%M+M)%M;
  }
return (ans%M+M)%M;
}
signed main(){
scanf("%lld",&n);
for(int i=1;i<=n;i++)scanf("%lld%lld",ai+i,bi+i);
printf("%lld",exCRT());
}

例题

P4774 【NOI2018】屠龙勇士

前言

中国剩余定理起的是合并答案的作用。

题解

解同余方程组：

\(\left(\begin{matrix}b_1* x\equiv a_1 \pmod{p_1}\\b_2*x \equiv a_2 \pmod{p_2}\\ \cdots \\ b_n*x\equiv a_n \pmod{p_n}\end{matrix} \right)\)

有系数了怎么办？？

为了合并答案，我们需要把前面的系数 \(b_i\) 去掉

无非就是解 \(n\) 个不定方程然后合并答案，这亦然是扩展中国剩余定理的根本

因此在处理中国剩余定理的系数时，可以尝试解完不定方程再合并答案。

因此问题变成了。

\(\left(\begin{matrix} x\equiv x_0 \pmod{P_1}\\x \equiv x_0 \pmod{P_2}\\ \cdots \\ x\equiv x_n \pmod{P_n}\end{matrix} \right)\)

其中 \(P\) 为每一个不定方程的剩余系通解的模数，即 \(b/gcd\) , \(x0\) 为解得的不定方程的最小非负整数解。

细节部分

根据中国剩余定理可知，最后的通解为：

\(x=x+k* M,k\in Z\)

这里求的就不是最小非负整数解了

每攻击一个巨龙都有一定的次数，所以最终次数一定是 \(\geq\) 最大次数的。

因此在剩余系里找一个合法的最小的解就 ok 了

这里用 \(multiset\) 处理前后缀，千万别用 \(set\) ！！

#include <iostream>
#include <cstdio>
#include <set>
#include <cstring>
#include <algorithm>
using namespace std;
#define LL long long
const int maxn = 1e5 + 10; 
LL read()
{
LL f=1,x=0;
char ss=getchar();
while(ss<'0'||ss>'9'){if(ss=='-')f=-1;ss=getchar();}
while(ss>='0'&&ss<='9'){x=x*10+ss-'0';ss=getchar();}
return f*x;
}
int T;
LL n,m,a[maxn],p[maxn],atk[maxn],mi[maxn],maxx=0;
LL mod[maxn],x0[maxn];
multiset<LL> s;
LL exgcd(LL a,LL b,LL &x,LL &y){
if(!b){
x=1;y=0;
return a;
    }
LL g=exgcd(b,a%b,x,y);
LL tmp=x;x=y;
y=tmp-(a/b)*y;
return g;
}
LL mul(LL a,LL b,LL P){
LL res=0;
while(b){
if(b&1)res=(res+a)%P;
a=(a+a)%P;
b>>=1LL;
    }
return res;
}
LL ans=0,M=0;
void clear(){
s.clear();
ans=0;M=0;
maxx=0;
memset(atk,0,sizeof atk);
memset(mod,0,sizeof mod);
memset(x0,0,sizeof x0);
return ;
}
void init(){
clear();
n=read();m=read();
LL input;
for(int i=1;i<=n;i++)a[i]=read();
for(int i=1;i<=n;i++)p[i]=read();
for(int i=1;i<=n;i++)mi[i]=read();
for(int i=1;i<=m;i++){
input=read();s.insert(input);
    }
for(int i=1;i<=n;i++){
multiset<LL> :: iterator it=s.upper_bound(a[i]);
if(it!=s.begin())it--;
atk[i]=(*it);
s.erase(it);
s.insert(mi[i]);
    }
//for(int i=1;i<=n;i++)cerr<<atk[i]<<" ";
return ;
}
LL exCRT(LL *ai,LL *bi){
ans=bi[1];
M=ai[1];
LL x,y;
for(int i=2;i<=n;i++){
LL a=M,b=ai[i],c=(bi[i]-ans%b+b)%b;
LL gcd=exgcd(a,b,x,y),bg=b/gcd;
if(c%gcd)return -1;
//
x=mul(x,c/gcd,bg)%bg;
ans+=x*M;
M*=bg;
ans=(ans%M+M)%M;
    }
if(ans>=maxx)return ans;
else return ans+((maxx-ans)/M+((maxx-ans)%M ? 1 : 0))*M;
}
LL solve(){
LL x,y;
for(int i=1;i<=n;i++){
LL gcd=exgcd(atk[i],p[i],x,y);
LL c=a[i];
if(c%gcd!=0)return -1;
//
mod[i]=p[i]/gcd;
x=(mul(x,c/gcd,mod[i])+mod[i])%mod[i];
x0[i]=x;
maxx=max(maxx,a[i]/atk[i]+(a[i]%atk[i]?1:0));
    }
LL ans=exCRT(mod,x0);
return ans;
}
int main(){
//freopen("1.in","r",stdin);
//freopen("1.out","w",stdout); 
scanf("%d",&T);
while(T--){
init();
ans=solve();
printf("%lld\n",ans);
    }
return 0;
}