poj1061-扩展欧几里得算法

来源:互联网 发布:js获取文件大小 编辑:程序博客网 时间:2024/06/14 21:59

首先声明:本题分析来自他人,不属本人创作。

题目分析:

     根据题意可知,题目要求的是求(x+k*m)modl=(y+k*n)mod l这个模方程的k的最小整数解;跟据模的性质对这个方程变形得:k*(m-n)+t*l=y-x(t为引进的一整数);

令a=(m-n);b=y-x;得:

    k*a+t*l=b; 0)

看到这个式子是不是感觉很眼熟,对,它就是扩展欧几里得法能解的标准方程a*x+b*y=d(此a,b非彼a,b)的标准形;唯一的不同是扩展欧几里得d是前面两系数a,b的最大公约数;而0)式中的b就不能保证了。

不过没关系,先用扩展欧几里得得出k*a+t*l=gcd(a,l)的一个解k0,t0;再观察两个式子:

                       k*a+t*l=b 0)

                      k0*a+t0*l=d 1)(令d=gcd(a,l));

将0)式两边都除以d,得:

         k*(a/d)+t*(l/d)=b/d;

因为d为a,b的最大公约数,所以a/d,l/d为整数,所以b/d必定为整数,否则方程无解;

将1)式两边都乘以b/d与0)式比较:

        k*a+t*l=b  0)

     (k0*b/d)*a+t0*(b/d)*l=b 1)

得k=k0*b/d;t=t0*b/d;

接下来就是由一个解扩展成一个解系:

       1)除以b得: k0*(a/d)+t0*(l/d)=1;

                          因为(a/d)和(l/d)互质,将方程写成:

                                      (k0+u*(l/d))*(a/d)+(t0+v*(a/d))*(l/d)=1;

                           只要u*(l/d)*(a/d)+v*(a/d)*(l/d)=0成立就行;

所以k,t的解系可以写成k=(k0+u*(l/d))*(b/d);t=(t0+v*(a/d))*b/d;

接下来就是求k的最小整数解了;

因为k=(k0 mod (l/d)+(u+k0/(l/d))*(l/d))*(b/d);

所以kmin=(k0*(b/d)) mod (l/d);

 

大功基本高成了,还差最后一步,那就是如何用扩展欧几里得解k*a+t*l=gcd(a,l):

对于a*x+b*y=d;
递归调用时
 令a=b;b=a%b;
将其变为形式2) b*x+a%b*y=d;
           变形:b*x+a*y-(a/b)*b*y=d;
           再变:a*y+b(x-a/b*y)=d; 3)

 与2)式比较:b*x+a%b*y=d;   2)
得:
   当a=b;b=a%b时:x=y;y=x-a/b*y;

调用过程中的x,y就是对应的a,b的解;
当回到顶层时,a,b就是最初的a,b所以此时的x,y就是所求解;

算法:

int extended_euclidean(int n, int m, int &x, int &y) {  

   if (m == 0) {  

   x = 1; y = 0; return n;  

1.      }  

2.     int g = extended_euclidean(m, n % m, x, y);  

3.     int t = x - n / m * y;  

4.        x = y;  

5.       y = t;  

6.      return g;  

7.  }     

以上均来自他人,本人只是想学习他人的见解。

本题代码如下:

#include <stdio.h>typedef long long INT;INT k1,k;INT gcd(INT a,INT b)               //求最大公约数。 {    return a%b?gcd(b,a%b):b;}void exgcd(INT a,INT b)          //扩展欧几里得函数。 {    if(!b)     { k=1,k1=0;return;}     exgcd(b,a%b);     INT t=k;     k=k1;     k1=t-a/b*k1;}int main(){    INT x,y,m,n,l,a,b,c,r;    while(~scanf("%lld%lld%lld%lld%lld",&x,&y,&m,&n,&l))    {       a=n-m,c=x-y;                    //步差,位差。        b=l;                            //圆环的长度。        r=gcd(a,b);                     //步差与圆环长度的最大公约数r。        if(c%r){                        //r最大公约数不是位差的约数,答案无解。         printf("Impossible\n");        continue;        }        a/=r,b/=r,c/=r;                //把等式除以公约数,使gcd(a,b)=1;         exgcd(a,b);                    //求出方程的特解。                            k=c*k-c*k/b*b;                 //求出最小解。(看分析)就知道了。         if(k<0)                        //判断最小解是否小于零,如果小于零,需要加上一周。          k+=b;        printf("%lld\n",k);         }      return 0;}                       



原创粉丝点击