• 『转』扩展欧几里德算法求不定方程


     例题是 POJ 1061 青蛙的约会 

      题目大意是,一个周长为L的圆, A、B两只青蛙,分别位于 x 、y 处,每次分别能跳跃 m 、n ,问最少多少次能够相遇,如若不能输出 “ Impossible”

      此题其实就是扩展欧几里德算法-求解不定方程,线性同余方程。

      设过 k1 步后两青蛙相遇,则必满足以下等式:

        ( x + m*k1 ) - ( y + n*k1 ) = k2*L  ( k2 =0 , 1 , 2....)         //这里的k2:   存在一个整数k2, 使其满足上式

      稍微变一下形得:   ( m - n )*k1 - k2*L= y - x 

      令  a = m - n , b = -L , c = y - x. 即

        a*k1 + b*k2 = c

      转换成了线性同余方程,只要上式存在整数解,则两青蛙能相遇,否则不能。 

      欧几里德扩展原理 是解决线性同余方程的工具

      下面我们介绍  欧几里德定理的证明,以及如何推广到扩展欧几里德

     

      先来看下什么是欧几里德扩展原理:

      

      欧几里德算法又称辗转相除法,用于计算两个整数a,b的最大公约数。其计算原理依赖于下面的定理:

      定理:gcd(a,b) = gcd(b,a mod b)

      证明:

    复制代码
            令 d = gcd( a, b ) 

            则 d | a  ,d | b    (   x|y 表示  x能够整除y,即 y%x = 0 )

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

      前面已经知道  d | b
        
            又  d | (a%b) => d | r => d | ( a - bk ) 

            => ( d | a ) - ( d | bk )   

            所以 d | (a%b) 

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

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

     

     

      欧几里德算法就是根据这个原理来做的,其算法用C++语言描述为: 

    1 int Gcd(int a, int b)
    2 {
    3   if(b == 0)    return a;
    4          return Gcd(b, a % b);
    5 }
     

      当然你也可以写成迭代形式:

    复制代码
     1 int Gcd(int a, int b)
     2 {
     3   while(b != 0)
     4   {
     5        int r = b;
     6        b = a % b;
     7        a = r;
     8   }
     9   return a;
    10 }
    复制代码
     

     

     

      扩展欧几里德算法

     

      问题: 计算  a*x + b*y = Gcd( a , b )

    复制代码
    由 欧几里德定理

      得: Gcd( a , b ) = Gcd( b , a%b )

      又  Gcd( b, a%b ) = b*x0 + (a%b) * y0          // 这里( x0, y0 )指 存在一对(x0,y0) 使其满足 等式

      所以  a*x + b*y = b*x0 + (a%b) * y0  

      又 a % b = a - (a / b) * b                                      // 注:这里的 (a/b) 是程序设计语言中的除法

      带入式子得:

      a*x + b*y = b*x0 + [ a - (a / b) * b ] * y0

      将式子右边变换下得:

      a*x + b*y = a*y0 + b*( x0 - (a/b) * y0 ) 

      根据对称我们可以知道:

      x = y0 , y = x0 - (a/b) * y0

      这样我们就得到了求解 x,y 的方法:x,y 的值基于 x0,y0
      由此,我们可以通过类似  欧几里德定理 循环迭代地 求出 x , y 

      再考虑下边界条件: 因为 a*x + b*y = Gcd( a, b )   

      当 b = 0 时, Gcd( a, b ) = a

      所以 a*x + b*y = a

      所以 x = 1, y = 0

      由此,我们知道边界条件为 当 b = 0 时, x = 1, y = 0  
    复制代码
     

      有上面的推导过程,我们可以得出C++的实现:  

    复制代码
    int exGcd(int a, int b, int &x, int &y)  // 此 &x 表示引用,目的是像 C中 指针那样传递地址,达到修改值的目的
    {
      if(b == 0)
      {
          x = 1;
          y = 0;
          return a;
      }
      int r = exGcd(b, a % b, x, y);
      int t = x;  // 因为 x 本身值要被覆盖,之后又要使用,所以定义临时变量来存储信息
      x = y;
      y = t - a / b * y;
      return r;
    }
    复制代码
     

    不定方程 a * x + b * y = c 的求解过程:

    复制代码
      求 a * x + b * y = c 的整数解。

        设 d = Gcd(a,b) , 若 c % d != 0  则不定方程无整数解 , 由数论相关定理可以证明

        由上面的推导,我们知道,我们可以使用 ExGcd: a*x + b*y = Gcd( a, b )   求出一组解 ( x0 , y0 ) 使其满足 a * x0 + b * y0 = Gcd( a, b )

        令  A = a/d , B = b/d, C = c/d  带入原式

        A * x + B * y = C                        ------ 等式 1
        因为 Gcd( A, B ) = 1 , 由 扩展GCD 性质得

        A * xi + B * yi = Gcd( A, B ) = 1     ------ 等式2

        我们可以通过 扩展GCD 求出 等式2 的一组解 ( x0, y0 ) 
        使其满足 A * x0 + B * y0  = 1        ------- 等式3
        等式3 两边同乘以 C, 可得
        A * x0 * C+ B * y0 * C = C       ------- 等式4  
        通过观察,我们可以知道 
        二元组 ( x0*C, y0*C ) 是 等式1 的一组解 使其满足 A*x + B*y = C
        所以 ( x0*C, y0*C ) 是 a*x + b*y = c 的一组解
        

        另,我们将式子1 变形下:

        A * ( x + B )  + B * ( y + A ) = C

        我们可以知道 x , y 解有多组. 但都满足如下关系

        x = x0*C + k * B    ( 其中 k 为整数 )

        y = y0*C + k * A
    复制代码
       此时方程的所有解为:x=c*k1-b*t,x的最小的可能值是0,令x=0可求出当x最小时的t的取值,但由于x=0是可能的最小取值,实际上可能x根本取不到0,那么由计算机的取整除法可知:由 t=c*k1/b算出的t,代回x=c*k1-b*t中,求出的x可能会小于0,此时令t=t+1,求出的x必大于0;如果代回后x仍是大于等于0的,那么不需要再做修正。


  • 相关阅读:
    编写内核模块
    ubuntu安装虚拟磁带库mhvtl
    linux中断与异常
    jquery 强大表格插件 DataTables
    判断元素是否存时,使用isset会比in_array快得多
    MarkDown 语法
    接口测试、压力测试工具
    jquery 复制文本到剪切板插件(非 flash)
    fiddler抓包HTTPS请求
    php mongodb manager 查数据的各种姿势
  • 原文地址:https://www.cnblogs.com/CKboss/p/3350809.html
Copyright © 2020-2023  润新知