zoukankan      html  css  js  c++  java
  • 中国剩余定理及其扩展学习笔记

    Chinese Remainder Theorem(CRT)

    又称孙子定理,是用来解决一类同余方程组问题的方法。

    这个要用到欧几里得算法相关知识,不会走这里GO


    我们从一个简单的例子入手:

    现有如下式子(p1,p2为质数):

    {xa1(mod p1)xa2(mod p2)

    按照同余式,上式可以写成:

    {x=a1+k1p1x=a2+k2p2其中k1,k2Z

    通过联立起来,将其写成a1+k1p1=a2+k2p2,k1p1k2p2=a2a1

    又因为gcd(p1,p2)=1,所以此方程肯定有解。

    对应ax+by=c式子,那么{a=p1b=p2c=a2a1,所以根据这种形式,我们用exgcd可以求出一组解(k^1,k^2),然后我们将其带回原方程组可得x0,然后通过exgcd的知识我们可以由这组(k^1,k^2)特解,推导出其它的解:{k1=p2t+k^1k2=p1t+k^2,k1,k2Z

    带入原方程求解得:x=a1+k1p1=a1+p1p2t+k^1=x0+p1p2t等价于xx0(mod p1p2)
    我们就得到了方程的答案,同时将模数合并了起来。


    所以对于多个方程构成的方程组,我们每次选择两个合并出一个新的,然后再用新的去和另一个合并(注意:模数相乘后小心爆范围long long)。

    那么对于如下的模方程组:

    {xa1 (mod m1)xa2 (mod m2)xa3 (mod m3)xan (mod mn)

    所有的mi均互质,所以用exgcd求出特解即可解决。(LRJ的蓝书上也有讲解)。

    模板代码:

    \[TJOI2009猜数字]
    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #define ll long long
    using namespace std;
    const int M=21;
    int n;
    ll a[M],b[M];
    void exgcd(ll a,ll b,ll &x,ll &y){
      if(!b){x=1;y=0;}
      else {exgcd(b,a%b,y,x);y-=x*(a/b);}
    }
    ll mul(ll a,ll b,ll mod){
      ll c=0;
      for(;b>0;b>>=1,a=(a+a)%mod)if(b&1) c=(c+a)%mod;
      return c;
    }
    ll CRT(ll *a,ll *b,int n){
      ll M=1,d=0,x,y;
      for(int i=1;i<=n;i++) M*=b[i];
        for(int i=1;i<=n;i++){
          ll w=M/b[i];
          exgcd(b[i],w,x,y);y=(y%b[i]+b[i])%b[i];//exgcd求解
          d=(d+mul(mul(w,y,M),(a[i]%b[i]+b[i])%b[i],M))%M;
        }
        return (d+M)%M;
    }
    int main(){
      scanf("%d",&n);
      for(int i=1;i<=n;i++)scanf("%lld",&a[i]);
      for(int i=1;i<=n;i++)scanf("%lld",&b[i]);
      ll ans=CRT(a,b,n);
      printf("%lld
    ",ans);
      return 0;
    }

    扩展中国剩余定理

    我们现在解决的一类模方程的求法,但是这是在模数互质的条件下。大多数时候模数都不会互质,那么应该如何求解呢?
    我们仍然从简单的入手

    看如下一个方程组:

    {xa1 (mod m1)xa2 (mod m2)

    同样的我们将其变形得到a1+k1m1=a2+k2m2,如果该方程要有整数解,那么根据裴蜀定理,要满足gcd(m1,m2)|(a2a1)

    这里我们先假设用exgcd求出了一组解特解(k1,k2),那么通解有是什么呢?
    接下来令g=gcd(m1,m2),又因为g|(a2a1)(因为我们先假设它有解),所以原方程可以写成

    k1m1gk2m2g=a2a1g
    ,这样我们就可得知gcd(m1g,m2g)=1,系数就转换为之前的互质的情况了,互质后我们再通过之前的方法得到一组特解(k^1,k^2),通解便为:
    {k1=m2gt+k^1k2=m2gt+k^2tZ

    回带入原方程,我们就可以得到:
    x=a1+k1m1

    =x0+m1m2gt

    =x0+lcm(m1,m2)t

    lcm为最小公倍数。

    实质上就等于解:xx0 (mod lcm(m1,m2))这个方程。

    举个例子,解如下方程:

    {x2 (mod 3)x3 (mod 5)x2 (mod 7)

    先自己试着解一下再看下面解答。

    先联立x 8(mod 15),再联立得到最终答案x23 (mod 105),那么答案就为x=23,反带回去发现确实成立。

    【模板IN

    模板代码:

    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #define ll long long
    using namespace std;
    const int M=201;
    ll xx(ll a,ll b,ll mod){
        if(a>b) swap(a,b);ll c=0;
        for(ll i=0;(1ll<<i)<=b;i++){if(b&(1ll<<i)){c=(c+(1ll<<i)*a)%mod;}}
        return c;
    }
    void exgcd(ll a,ll b,ll &d,ll &x,ll &y){
        if(!b){d=a;x=1;y=0;}
        else {exgcd(b,a%b,d,y,x);y-=x*(a/b);}
    }
    ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
    ll inv(ll a,ll b){
        ll g,x,y;
        exgcd(a,b,g,x,y);
        if(g!=1) return -1;
        return (x%b+b)%b;
    }
    bool merge(ll a1,ll b1,ll a2,ll b2,ll &a,ll &b){
        ll g=gcd(b1,b2);
        ll c=a2-a1;
        if(c%g) return 0;//判断有无解
        c=(c%b2+b2)%b2;
        b1/=g;b2/=g;
        c/=g;c=xx(c,inv(b1,b2),b2);
        b=b1*b2*g;
        c=xx(c,xx(b1,g,b),b);
        c+=a1;
        a=(c%b+b)%b;
        return 1;
    }
    ll CRT(ll *a,ll *b,int n){
        ll aa=a[1],bb=b[1];
        for(int i=2;i<=n;i++){
            ll da=a[i],db=b[i];
            ll nexa,nexb;
            //合并
            if(!merge(aa,bb,da,db,nexa,nexb)) return -1;
            aa=nexa;bb=nexb;
        }
        return (aa%bb+bb)%bb;
    }
    ll a[M],b[M];
    int n;
    int main(){
        scanf("%d",&n);
        for(int i=1;i<=n;i++){scanf("%lld%lld",&b[i],&a[i]);}
        ll ans=CRT(a,b,n);
        printf("%lld
    ",ans);
        return 0;
    }
    //模板2
    #include<cstdio>
    #include<cstring>
    #include<iostream>
    #include<algorithm>
    #define LL long long
    using namespace std;
    const int N=210;
    int n;
    LL x,y,lcm,equ[N][2];
    LL multi(LL a,LL b,LL p)
    {
        a=(a%p+p)%p;b=(b%p+p)%p;LL ans=0;
        for(;a;a>>=1,b=(b*2)%p)if(a&1)ans=(ans+b)%p;
        return ans;
    }
    LL exgcd(LL a,LL b,LL &x,LL &y)
    {
        if(!b)
        {
            x=1,y=0;
            return a;
        }
        LL val=exgcd(b,a%b,x,y);
        LL t=x;x=y;y=t-a/b*y;return val;
    }
    int main()
    {
        scanf("%d",&n);
        for(int i=1;i<=n;i++)
         scanf("%lld%lld",&equ[i][1],&equ[i][0]);
        for(int i=1;i<n;i++)
        {
            LL val=exgcd(equ[i][1],equ[i+1][1],x,y);
            lcm=equ[i][1]/val*equ[i+1][1];
            if((equ[i+1][0]-equ[i][0])%val)
            {
                printf("No Answer
    ");
                return 0;
            }
            val=multi(y,(equ[i+1][0]-equ[i][0])/val,lcm);
            equ[i+1][0]=(multi(equ[i+1][1],-val,lcm)+equ[i+1][0]+lcm)%lcm;
            equ[i+1][1]=lcm;
        }
        printf("%lld
    ",(equ[n][0]%equ[n][1]+equ[n][1])%equ[n][1]);
        return 0;
    }

    一个额外的小技巧:
    例如NOI2018的屠龙勇士。
    当同余方程组有系数,如下式子:

    {a1xb1 (mod m1)a2xb2 (mod m2)anxbn(mod mn)

    如果a,b都互质,我们用扩展欧几里得定理就可以求得逆元,然后把系数除过去就可以计算了。

    但是当a,b不一定互质怎么办?

    我们观察式子,axb(mod m),可以将其写成ax=km+b,我们记g=gcd(a,m)(gcd表示求最大公约数),那么左右可以写成(gw1)x=k(gw2)+b,由于左边为g的倍数,右边那么同理也应该为g的倍数,所以显然b也应该为g的倍数,若不为则该方程组无整数解,有的话将a,b,m同时除以g,就得到了系数互质的方程,用原来的方法解就可以了。


    End

    Thanks For Reading!

  • 相关阅读:
    Comparable与Comparator
    【源码】String类、jdk1.6subString内存泄漏、字符串拼接几种区别、
    JAVA整型包装类的缓存策略
    通过tomcat把项目http请求转为https请求
    git rebase总结及git使用规范
    记一次对象序列化不打印value值为null的属性问题
    Layui-Tables+PHP分页
    Python操作字符串-截取IP地址
    命令行启动VMware虚拟机
    bat批处理备份桌面所有文档
  • 原文地址:https://www.cnblogs.com/VictoryCzt/p/10053429.html
Copyright © 2011-2022 走看看