扩展欧几里得算法

来源:互联网 发布:电脑软件学校招生 编辑:程序博客网 时间:2024/04/19 17:46

根据此题 整理下 扩展欧几里得

扩展欧几里得是用来判断并求 ax + by = c 是否有解及其解的数学算法

首先列出定理

 1、ax + by = gcd(a,b)

 2、gcd(a,b) = gcd(b, a%b);


 本题 题意是 判断两个同向的青蛙 在初始位置(x,y),速度不同(m,n)的情况下能否在一个球形路线(周长为L)下相遇

 也就是可以列出等式   

    (x + m*t)%L = (y + n*t)%L;

==>     (x + m*t)%L -(y + n*t)%L =0

==>     (x + m*t) -(y + n*t)  + k*L =0

==>     (m - n )*t  + k*L = y-x

        a   *   x  + b*y = c;

判断是否有解以及解;

如果C对gcd(a,b)取余不是0那么就说明无论怎么走 都不会相遇 

如果是0

说明有解 带入扩展欧几里得算法求解

我们对 ax + by = c 化简 根据定理2

ax + by = gcd (a,b) = gcd (b,a%b)  =  bx +(a%b) *y = c  = ...

当化简到最后也就是 b = 0时候 原式可以表示成 a*x +0*y = gcd(a,0) = a

此时也就是说 ==> x=1,y = 0; 

求得一个关于a*x +0*y = gcd(a,0) = a的解

ax' + by'  

     =  bx + (a%b) *y

     = bx + (a - a/b*b ) *y

                     = bx + ay- a/b*by

                     = ay + b(x-a/b*y)

也就是说 当前ab下的解

通过下一层递归的x和y求出

x' = y

y' = (x-a/b*y)

于是可以写出扩展欧几里得算法

当递归回溯回来的时候 我们就得到了 ax' + by' = gcd(a,b)的解  

ll egcd(ll a,ll b,ll &x,ll &y){    if(b==0)     {        x = 1,y = 0;        return a;    }    ll ans = egcd(b,a%b,x,y);    ll t = x;    x = y;    y = t-a/b*y;    return ans;}

可是问题还没解决 我们要求的是等于C的解我们现在得到了gcdab的一组解,这就需要我们进一步转化答案 

设t = C/gcd(a,b)

可以将 ax + by  = c =   a * t *x + b * t *y = gcd(a,b) *t  

也就是 我们得到的解 的X值 也就是时间值 再乘以 t 就是对应C 下的x解

然后得到解后有可能解是负数 那么需要我们进一步找通解

通解公式

x = x0 + k * ( b / gcd(a,b) )

y = y0 + k * ( a / gcd(a,b) )

最终可以获得最小正整数解

例题  POJ 1061

code:

#include <iostream>#include <iomanip>#include <stdio.h>#include <stdlib.h>#include <algorithm>#include <functional>#include <vector>#include <cmath>#include <string>#include <stack>#include <queue>using namespace std;typedef long long ll;ll egcd(ll a,ll b,ll &x,ll &y){    if(b==0)    {        x = 1,y = 0;        return a;    }    ll ans = egcd(b,a%b,x,y);    ll t = x;    x = y;    y = t-a/b*y;    return ans;}int main(){    ll x,y,m,n,L;    cin>>x>>y>>m>>n>>L;    ll a,b,c;    a = m-n;    b = L;    c = y-x;    if(a<0)    {        a = -a;        c = -c;    }    ll gcd = egcd(a,b,x,y);    cout<<egcd(x,y,a,b)<<endl;    if(c%gcd!=0)puts("Impossible");    else    {        x = x*c/gcd;        ll t = b/gcd;        if(x>=0)            x %=t;        else            x = x%t+t;    }    return 0;}


==>     (x + m*t)%L -(y + n*t)%L =0