zoukankan      html  css  js  c++  java
  • 扩展欧几里得算法

    【转载】http://blog.csdn.net/qq_34494458/article/details/52637193

    一:欧几里得算法(辗转相除法)

                           基本算法:设a=qb+r,其中a,b,q,r都是整数,则gcd(a,b)=gcd(b,r),即gcd(a,b)=gcd(b,a%b)。

    证明:

           a可以表示成a = kb + r,则r = a mod b

      假设d是a,b的一个公约数,则有

      d|a, d|b,而r = a - kb,因此d|r

      因此d是(b,a mod b)的公约数

      假设d 是(b,a mod b)的公约数,则

      d | b , d |r ,但是a = kb +r

      因此d也是(a,b)的公约数

      因此(a,b)和(b,a mod b)的公约数是一样的,其最大公约数也必然相等,得证

    算法实现:

    int gcd(int a,int b) {             //递归算法  
        return b ? gcd(b, a%b) : a;  
    }  
      
    int gcd(int a, int b) {      //迭代算法  
        while(b) {  
            int c = a%b;  
            a = b;  
            b = c;  
        }  
        return a;  
    }  

    二   扩展欧几里得算法:

                   基本算法:对于不完全为 0 的非负整数 a,b,gcd(a,b)表示 a,b 的最大公约数,必然存在整数对 x,y ,使得 gcd(a,b)=ax+by。

       证明:设 a>b。

      1,显然当 b=0,gcd(a,b)=a。此时 x=1,y=0;

      2,ab!=0 时

      设 ax1+by1=gcd(a,b);

      bx2+(a mod b)y2=gcd(b,a mod b);

      根据朴素的欧几里德原理有 gcd(a,b)=gcd(b,a mod b);

      则:ax1+by1=bx2+(a mod b)y2;

      即:ax1+by1=bx2+(a-(a/b)*b)y2=ay2+bx2-(a/b)*by2;

      根据恒等定理得:x1=y2; y1=x2-(a/b)*y2;

          这样我们就得到了求解 x1,y1 的方法:x1,y1 的值基于 x2,y2.

       上面的思想是以递归定义的,因为 gcd 不断的递归求解一定会有个时候 b=0,所以递归可以结束。

    //递归代码  
    int exgcd(int a, int b, int&x, int&y) {  
       if (b == 0) {  
           x = 1;  
           y = 0;  
           return a;  
       }  
       int r = exgcd(b, a%b, y, x);  
       int t = x;  
       y = y - a/b*t;  
       return r;  
      }
       // 非递归算法  
    int exgcd(int m, int n, int &x, int &y) {  
       int x1, y1, x0, y0;  
       x0 = 1; y0 = 0;  
       x1 = 0; y1 = 1;  
       x = 0; y = 1;  
       int r = m%n;  
       int q = (m-r)/n;  
       while(r) {  
           x = x0 - q*x1;  
           y = y0 - q*y1;  
           x0 = x1; y0 = y1;  
           x1 = x; y1 = y;  
           m = n; n = r; r = m%n;  
           q = (m-r)/n;  
       }  
       return n;  
    }

    刘汝佳的:

    void gcd(int a, int b, int& d, int& x, int& y) {  
        if (!b) {d = a; x = 1; y = 0; }  
        else {gcd(b, a%b, d, y, x); y -= x*(a/b); }  
    }  

    上面求出的 x(当a>b时)都是最接近0的正整数。(不知道为什么)

    扩展欧几里德算法的应用主要有以下两个方面:

    (1)求解不定方程;

    (2)求解模线性方程(线性同余方程)与逆元;

    (1)求解不定方程;

      1.对于不定整数方程ax+by=m,若 m mod Gcd(a, b)=0,则该方程存在整数解,否则不存在整数解。

         证明:

       首先扩展欧几里德主要是用来与求解线性方程相关的问题,所以我们从一个线性方程开始分析。现在假设这个线性方程为a*x+b*y=m,如果这个线性方程有解,那么一定有gcd(a,b) | m,即a,b的最大公约数能够整除m(m%gcd(a,b)==0)。证明很简单,由于a%gcd(a,b)==b%gcd(a,b)==0,所以a*x+b*y肯定能够整除gcd(a,b),如果线性方程成立,那么就可以用m代替a*x+b*y,从而得到上面的结论,利用上面的结论就可以用来判断一个线性方程是否有解。
      2.ax+by=Gcd(a, b) 两边同时乘以 m/Gcd(a, b)得a(x*m/Gcd(a, b))+b(y*m/Gcd(a, b))=m;

    上面已经列出找一个整数解的方法,在找到 a*x+  b*y = Gcd(a, b)的一组解x0,y0后,不定方程ax+by=m的一组解满足 x = m(x0)/gcd , y = m*(y0)/gcd;
    所以通解为  x = m*(x0)/gcd + k*lcm/a , y = m*(y0)/gcd + k*lcm/b (其中k为任意整数);
    证明:

          令a1=a/gcd(a,b),b1=b/gcd(a,b),m1=m/gcd(a,b)。那么x,y的一组解就是x0*m1,y0*m1,但是由于满足方程的解无穷多个,在实际的解题中一般都会去求解x或是y的最小正数的值。以求x为例,又该如何求解呢?还是从方程入手,现在的x,y已经满足a*x+b*y=m,那么a*(x+n*b)+b*(y-n*a)=m显然也是成立的。可以得出x+n*b(n=…,-2,-1,0,1,2,…)就是方程的所有x解的集合,由于每一个x都肯定有一个y和其对应,所以在求解x的时候可以不考虑y的取值。取k使得x+k*b>0,x的最小正数值就应该是(x+k*b)%b,但是这个值真的是最小的吗??如果我们将方程最有两边同时除以gcd(a,b),则方程变为a1*x+b1*y=m1,同上面的分析可知,此时的最小值应该为(x+k*b1)%b1,由于b1<=b,所以这个值一定会小于等于之前的值。在实际的求解过程中一般都是用while(x<0)x+=b1来使得为正的条件满足,为了更快的退出循环,可以将b1改为b(b是b1的倍数),并将b乘以一个倍数后再加到x上。

    代码:

        //不定方程  
    void linear_equation(int a, int b, int c, int &x, int &y) {  
        int d = exgcd(a, b, x, y);  
        if (c%d)  
            return ;  
        int k = c/d;  
        x *= k; y *= k; // 一组解  
        return ;  
    }  

    相关题目:poj2115 poj2142 poj1061。

    (2)求解模线性方程(线性同余方程)与逆元;

            同余方程 ax≡b (mod n)对于未知数 x 有解,当且仅当 gcd(a,n) | b。且方程有解时,方程有 gcd(a,n) 个解。

               求解方程 ax≡b (mod n) 相当于求解方程 ax+ (-ny)= b, (x, y为整数)。

       //ax = b (mod n) 即 (ax) mod n = b mod n

     算法导论上有两个定理:

      定理一:设d = gcd(a, n), 假定对整数x', y', 有d = ax' + ny', 如果d | b, 则方程ax = b(mod n)有一个解的值为x0, 满足:、

          x0 = x'(b / d)(mod n)

      定理二:假设方程ax = b(mod n)有解, x0是方程的任意一个解, 则方程对模n恰有d个不同的解, 分别为: xi = x0 + i * (n / d), 其中 i = 1,2,3......d - 1

      有了这两个定理, 解方程就不难了。

    代码:

       // 模线性方程  
    bool modular_linear_equation(int a, int b, int n) {  
        int x, y, x0, i;  
        int d = exgcd(a, n, x, y);  
        if (b%d)  
            return false;  
        x0 = x*(b/d)%n; //特解  
        for(i = 0; i < d; i++)  
            printf("%d
    ", (x0 + i*(n/d))%n);  
        return true;  
    }  

    相关题目  hdu1576.

    同余方程ax≡b (mod n),如果 gcd(a,n)== 1,则方程只有唯一解。
          在这种情况下,如果 b== 1,同余方程就是 ax=1 (mod n ),gcd(a,n)= 1。
          这时称求出的 x 为 a 的对模 n 乘法的逆元。
          对于同余方程 ax= 1(mod n ), gcd(a,n)= 1 的求解就是求解方程
          ax+ ny= 1,x, y 为整数。这个可用扩展欧几里德算法求出,原同余方程的唯一解就是用扩展欧几里德算法得出的 x 。

    代码:

     
    ll inv(ll a, ll n) {  
        ll d, x, y;  
        d = exgcd(a, n, x, y);  
        return (d == 1) ? (x%n + n)%n : -1;  
    }  

    相关题目: hdu2115.

  • 相关阅读:
    C# Dictionary几种遍历方式
    Android 代码中文字在手机上显示乱码问题解决方法
    【Android学习5】Clean 之后R文件丢失
    【android学习4】Eclipse中Clean作用
    Java基础语法(练习)
    Java基础语法(自定义类、ArrayList集合)
    Java基础语法(方法)
    Java基础语法(数组)
    Java基础(Scanner、Random、流程控制语句)
    Java基础(变量、运算符)
  • 原文地址:https://www.cnblogs.com/pach/p/5911163.html
Copyright © 2011-2022 走看看