zoukankan      html  css  js  c++  java
  • 扩展欧几里得(exgcd)与同余详解

    exgcd入门以及同余基础

    gcd,欧几里得的智慧结晶,信息竞赛的重要算法,数论的...(编不下去了

    讲exgcd之前,我们先普及一下同余的性质:

    1. left ( a+b 
ight )mod p=left (amod p+ bmod p
ight )mod p
    2. left ( a-b 
ight )mod p=left (amod p- bmod p
ight )mod p
    3. left ( a	imes b 
ight )mod p=left (amod p)	imes (bmod p
ight )mod p
    4. aequiv b(mod p),那么a^{n}equiv b^{n}(mod p)
    5. amod p1=xamod p2=x,且p1,p2互质,amod (p	imes q)=x

     有了这三个式子,就不用怕在计算时溢出了。

    下面我会用gcdleft ( a,b 
ight )lcmleft ( a,b 
ight )分别表示a与b的最大公约数与最小公倍数。

    首先会来学扩欧的同学肯定都会欧几里得算法(即辗转相除法)了吧

    而通过观察发现:lcmleft ( a,b 
ight )=frac{a	imes b}{gcdleft ( a,b 
ight )}=frac{a}{gcdleft ( a,b 
ight )}	imes b,先除后乘防溢出。

    所以gcdleft ( a,b 
ight )lcmleft ( a,b 
ight )的代码如下:

    1 inline int gcd(int a,int b)
    2 {return (b==0)?a:gcd(b,a%b);}
    3 inline int lcm(int a,int b)
    4 {return a/gcd(a,b)*b;

    讲exgcd之前先引入一种方程——不定方程

     所谓不定方程,是指未知数的个数多于方程个数,且未知数受到某些限制(如要求是有理数、整数或正整数等等)的方程或方程组。

                                                                                             ————百度百科 

    就是形如ax+by=c的方程,其中a,b,c已知。

    1.判断是否有解

    如果cmod gcdleft ( a,b 
ight )
eq 0,那么方程无解。

    2.转化

    方程可转化为{a}'x+{b}'y={c}'

    其中{a}'=frac{a}{gcdleft (a,b 
ight )}{b}'=frac{b}{gcdleft (a,b 
ight )}{c}'=frac{c}{gcdleft (a,b 
ight )}

    3.求一组特解

    接着就用到了exgcd。

    我们知道gcd有一个性质gcd(a,b) = gcd(b, amod b)

    如果,一直循环下去,b将等于0,那么x将等于c/a,y=0。

    1 inline void exgcd(int a,int b,int c)
    2 {
    3     if(!b)
    4     {x=c/a;y=0;return;}
    5     exgcd(b,a%b,c);
    6     x=y;
    7     y=(c-a*x)/b;
    8     return;
    9 }

    这就求出了一组特解。

    exgcd的模板我也在这摆出来

    1 inline void exgcd(int a,int b)
    2 {
    3     if(!b)
    4     {x=1;y=0;return;}
    5     exgcd(b,a%b);
    6     k=x;x=y;
    7     y=k-a/b*y;
    8     return;
    9 }

    这是求ax+by 
ight =gcdleft ( a,b 
ight )时用的,后面讲同余方程会讲。

    4.构造通解

    我们假设x1,y1是我们求出的一组特解,那么

    left{egin{matrix} & &x=ccdot x1+ccdot bcdot k \ & & y=ccdot y1-ccdot acdot k end{matrix}
ight.     left ( kepsilon z 
ight )


    同余类问题

    1.单个同余方程

    求的是axequiv b( mod p)关于x的解

    转化一下,就成了ax+py=b,就可以直接套exgcd模板。

    2.同余方程组

    left{egin{matrix} & &xequiv b1left (mod p1 
ight )\ & & xequiv b2left (mod p2 
ight ) end{matrix}
ight.

    1.有解的充要条件left ( b1-b2 
ight )mod gcd(p1,p2)=0

    2.left{egin{matrix} & &x-p1	imes y1=b1\ & & x-p2	imes y2=b2end{matrix}
ight.

    下式减上式得p1	imes y1-p2	imes y2=b2-b1

    再用exgcd求出y1和y2,{x}'=b1+p1	imes y1=b2+p2	imes y2

    3.关于通解

    所有的x mod lcm(p1,p2)有唯一解,这样就可以通过特解,求通解了。

    4.至于式子更多的同余方程组,就先联立两个,就可以得出新的方程xequiv {x}'left ( mod lcm(p1,p2) 
ight )

    再联立下一个。


    exgcd用处及题目讲解

    1.求同余方程的解

    例如这道题P1082

    这是一道裸的扩欧模板题,变形之后就是求ax+by=1

    套模板即可。

     1 inline void exgcd(int a,int b)
     2 {
     3     if(b==0)
     4     {x=1;y=0;return;}
     5     exgcd(b,a%b);
     6     k=x;x=y;
     7     y=k-a/b*y;
     8     return;
     9 }
    10 int main()
    11 {
    12     int n,m;
    13     read(n),read(m);
    14     exgcd(n,m);
    15     printf("%d",(x+m)%m);
    16 }

    还有一道模板P1516

    仔细观察,推一下后我们发现,这在就是在求x+k	imes mequiv y+k	imes n(mod l)的解。

    进而可以推出x+k	imes m-y-k	imes n=l	imes p

    合并同类项后(x-y)+k	imes(m-n)-l	imes p=0

    把一些东西移到左边来后(x-y)=k	imes(n-m)+l	imes p

    把(x-y),(n-m)各看成一个整体后,问题就成了解一个不定方程。

     1 inline int exgcd(long long a,long long b)
     2 {
     3     if(b==0)
     4     {x=1;y=0;return a;}
     5     ans=exgcd(b,a%b);
     6     k=x;x=y;
     7     y=k-a/b*y;
     8     return ans;
     9 }
    10 int main()
    11 {
    12     long long x1,y1,m,n,l;
    13     read(x1),read(y1),read(m),read(n),read(l);
    14     if(n-m<0)swap(x1,y1);
    15     exgcd(std::abs(n-m),l);
    16     if((x1-y1)%ans!=0)
    17     printf("Impossible");
    18     else
    19     printf("%lld",((x*((x1-y1)/ans))%(l/ans)+(l/ans))%(l/ans));
    20 }

    还有一道也是模板P4777,涉及同余方程组求解,上面已详细的讲了,近期我也会发一篇中国剩余定理的博客

     1 inline long long mul(long long a,long long b,long long mod)
     2 {
     3     long long res=0;
     4     while(b>0)
     5     {
     6         if(b&1) res=(res+a)%mod;
     7         a=(a+a)%mod;
     8         b>>=1;
     9     }
    10     return res;
    11 }
    12 long long exgcd(long long a,long long b,long long &x,long long &y)
    13 {
    14     if(!b)
    15     {x=1;y=0;return a;}
    16     long long gcd=exgcd(b,a%b,x,y);
    17     k1=x;x=y;
    18     y=k1-a/b*y;
    19     return gcd;
    20 }
    21 int main()
    22 {
    23     io::begin();
    24     io::read(n);
    25     for(register int i=1;i<=n;i++)
    26     io::read(b1[i]),io::read(a1[i]);
    27     long long x,y,k;
    28     long long m=b1[1],ans=a1[1];
    29     for(int i=2;i<=n;i++)
    30     {
    31         long long a=m,b=b1[i],c=(a1[i]-ans%b+b)%b;
    32         long long gcd=exgcd(a,b,x,y);
    33         long long p=b/gcd;
    34         x=mul(x,c/gcd,p);
    35         ans+=x*m;
    36         m*=p;
    37         ans=(ans%m+m)%m;
    38     }
    39     printf("%lld",(ans%m+m)%m);
    40 }

    2.扩欧求逆元

    这是一种很重要的算法,至于逆元怎么跟扩欧扯上关系,大家可以点这里乘法逆元及两道模板题详解

    这里就不多赘述了,大家可以用扩欧a一下P3811,P2613。


    我要讲的讲完了,如果觉得讲的还好,请关注我的blog,谢谢

  • 相关阅读:
    (转)Docker network命令
    (转)Navicat连接MySQL8.0亲测有效
    (转)Docker 网络
    Docker问题方案收集
    (转)docker run的--rm选项详解
    (转)docker-compose安装
    (转)教你分分钟搞定Docker私有仓库Registry
    (转)Docker入门——Dockerfile详解
    (转)Windows下安装Docker, GitBash环境配置
    (转)教你分分钟搞定Docker私有仓库Registry
  • 原文地址:https://www.cnblogs.com/mashiro-/p/9526323.html
Copyright © 2011-2022 走看看