中国剩余定理(CRT) & 扩展中国剩余定理(ExCRT)总结

中国剩余定理(CRT) & 扩展中国剩余定理(ExCRT)总结

标签：数学方法——数论
阅读体验：https://zybuluo.com/Junlier/note/1300035

前置浅讲

前置知识点：\(Exgcd\)
这两个东西都是用来解同余方程组的
形如
\[ \left\{ \begin{aligned} x\equiv B_1(mod\ W_1)\\ x\equiv B_2(mod\ W_2)\\ \cdots\\ x\equiv B_n(mod\ W_n)\\ \end{aligned} \right. \]给定\(B_i\)和\(W_i\)，求解唯一解\(x\)满足上述方程组
PS：个人认为自己的\(ExCRT\)讲得好一些

中国剩余定理(CRT)

心理准备：这里我觉得自己讲得不是很清晰，有点说不清的感觉
在上述方程中，存在一种特殊情况，即\(W_i\)全部互质
有什么用呢，用处就是可以用中国剩余定理（孙子定理）
首先，古人告诉我们：
解上面那个方程相当于对于每一个\(B_i\)：

把\(B_i\)变成\(1\)，其他的\(B\)变成\(0\)的解\(x\)
然后答案就是\(\sum(x*B_i)\)
也就是解每一个形如下面的方程组的解\(x\)在乘上\(B_i\)：
\[ \left\{ \begin{aligned} x\equiv 0(mod\ W_1)\\ x\equiv 0(mod\ W_2)\\ \cdots\\ x\equiv 1(mod\ W_i)\\ \cdots\\ x\equiv 0(mod\ W_n)\\ \end{aligned} \right. \]

别问我为什么，我不知道
那么考虑怎么求呢？
记\(M\)为\(\prod W_i\)
显然解\(x\)必须是\(M/W_i\)的倍数
那么方程变为
\[ (M/W_i)y\equiv 1(\mod W_i) \]
\(Exgcd\)求解即可，再加进答案

lst CRT()
{
    lst tot=1,Ans=0;
    for(int i=1;i<=n;++i)tot*=W[i];
    for(int i=1;i<=n;++i)
    {
        lst now=tot/W[i],x,y;
        Exgcd(now,W[i],x,y);//不需要我讲吧！
        x=(x%W[i]+W[i])%W[i];//这个取膜貌似很关键诶
        Ans=(Ans+Mult(Mult(x,now,tot),B[i],tot)+tot)%tot;//Mult是快速乘
    }return Ans>=0?Ans:Ans+tot;
}

扩展中国剩余定理(ExCRT)

再把方程组放一遍：
\[ \left\{ \begin{aligned} x\equiv B_1(mod\ W_1)\\ x\equiv B_2(mod\ W_2)\\ \cdots\\ x\equiv B_n(mod\ W_n)\\ \end{aligned} \right. \]这个时候我们的\(W_i\)不互质了
那么我们考虑对所有的方程一个一个求解
假设我们求解到了第\(i\)个方程
前面的方程组解出来答案是\(Ans\)
那我们是不是可以把之前求解的答案看做一个这样的同余方程：
\[ x\equiv Ans(\mod M) \]其中\(M\)是前\(i-1\)个\(W\)的最小公倍数，\(x\)使我们想得到的新解
如果看不出请补一下同余方程。。。
很显然这些都是已知量了吧
那我们为了求出前\(i\)个方程的解就相当于要解出下面这个方程组了对不对
\[ \left\{ \begin{aligned} x\equiv Ans(\mod M) \\ x\equiv B_i(\mod W_i)\\ \end{aligned} \right. \]
用\(A_1\)表示Ans，\(B_1\)表示M，\(A_2\)表示\(B_i\)，\(B_2\)表示\(W_i\)，\(Ans\)表示新答案
可以化为
\[ Ans=B_1x1+A_1=B_2x2+A_2 \]\[ \therefore B_1x1-B_2x2=A_2-A_1 \]为了使用扩展欧几里德我们设\[ B_1x-B_2y=Gcd(B_1,B_2) \]那么我们会解出一个\(x'\)，而真正的\(x1\)如下
\[x1=x'×\dfrac{A_2-A_1}{GCD(B_1,B_2)}\]我们再回代就得到了新解\(Ans\)\[ Ans=B_1x1+A_1 \]有点凌乱请认真看(trust me)！
模板题：洛谷P4777 【模板】扩展中国剩余定理（EXCRT）
代码：

#include<bits/stdc++.h>
#define lst long long
#define ldb double
#define N 100050
using namespace std;
const int Inf=1e9;
lst read()
{
    lst s=0,m=0;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-')m=1;ch=getchar();}
    while( isdigit(ch))s=(s<<3)+(s<<1)+(ch^48),ch=getchar();
    return m?-s:s;
}

lst n,Ans,x,y,M;
lst A[N],B[N];

lst qmul(lst x,lst y,lst p)
{
    lst ret=0;
    while(y)
    {
        if(y&1)ret=(ret+x)%p;
        x=(x+x)%p;y>>=1;
    }return ret;
}

lst Exgcd(lst a,lst b,lst &x,lst &y)
{
    if(!b){x=1,y=0;return a;}
    lst ss=Exgcd(b,a%b,x,y),t;
    t=x,x=y,y=t-(a/b)*y; return ss;
}

int main()
{
    n=read();
    for(int i=1;i<=n;++i)
        B[i]=read(),A[i]=read();
    Ans=A[1],M=B[1];
    //根据上面的详解一步对应一行
    for(int i=2;i<=n;++i)
    {
        lst Get=((A[i]-Ans)%B[i]+B[i])%B[i];
        lst GCD=Exgcd(M,B[i],x,y);
        x=qmul(x,Get/GCD,B[i]);//qmul是龟速乘
        Ans+=M*x;
        M*=B[i]/GCD;//这里是更新辅助下一次计算
        Ans=(Ans+M)%M;
    }printf("%lld\n",Ans);
    return 0;
}