POJ 1061 青蛙的约会 (扩展欧几里德解不定方程)

题目链接：http://poj.org/problem?id=1061 题目大意：两只青蛙在一个首尾相接的轴上(1 - L)跳，并且其中一个起点在x，每步跳m，另一个起点在y，每步跳n，问他们经过某步后可不可以相遇，如果可以，找出步数. 思路： 经典的扩展欧几里德解线性不定方程的题： 显而易见的我们有方程 x + m*p = y + n*p (mod L),但是这个形式不太好解，未知数在两边，那么我们换个形式： (x + m*p) - (y + n*p) = t*L       ->      (n - m) * p + t * L = (x - y) 于是回到了我们熟悉的形式：解不定方程  ax + by = c. 我们用扩展欧几里德来解决这个问题： 1.扩展欧几里德算法解整数不定方程 ax + by = gcd(a, b)  （由定理可知其一定有解，证明略了.） （注：ext_gcd求出方程一个通解x₀, y₀， 其所有的整数解形式为：(x₀ + kb/g, y₀ - ka/g)   (g = gcd(a,b) ) ） 首先：根据定理 gcd(a, b) = gcd(b,  a mod b)我们可以得到等式ax₁ + by₁ = bx₂ + (a mod b)y₂ 整理一下：ax₁ + by₁ = bx₂ + (a - a/b*b)y₂    ->     ax₁ + by₁ = ay₂ + b(x₂ - a/b * y₂) 由恒等定理得： x₁ = y₂ y₁ = x₂ - a/b*y₂ 然后我们在欧几里德辗转相除法的基础上递归地利用gcd(b, a mod b)的解来得到gcd(a, b)的解：  

void ext_gcd(long long a, long long b, long long &x, long long &y){
    if (b == 0){
        x = 1;
        y = 0 ;
        return ;
    }
    ext_gcd(b, a%b, x, y);
    long long tmp = x;
    x = y;
    y = tmp - a/b * y;
    return ;
}

  2.解整数不定方程ax + by = c ★ ①首先由定理"所有的ax + by都一定是gcd(a,b)的倍数"可知方程有解的充要条件是gcd(a,b)必须能整除c. ②然后我们令等式两边同时除以gcd(a,b)得到等价方程   ->   a'x + b'y = c'        易知gcd(a', b') = 1.利用扩展欧几里德求出方程解(x₀, y₀)，则方程a'x + b'y = 1的通解为(x₀ + kb, y₀ - ka)，且(c' * x₀, c' * y₀)即为方程a'x + b'y = c'的解. ③注：题目要求满足方程的x为解中最小的整数，求最小整数过程中便很容易出现一个问题：如果x_#是方程a'x + b'y = 1的最小整数x解，我们不能说c' * x_#是方程a'x + b'y = c'的最小整数x解！比如方程413x + 1908y = 365的最小整数x解为(121, -26)，看出问题了吧！方程a'x + b'y = c'的所有整数解并不是(  c' * (x₀ + kb'),   c' * (y₀ - ka')   )！事实上我们应该这么做（虽然我也不会证为什么……）：我们求出了a'x + b'y = 1的解(x₀, y₀)后，先求出a'x + b'y = c'的解(c' * x₀, c' * y₀)，则方程a'x + b'y = c'的所有整数解为(c' * x₀ + kb', c' * y₀ - ka'). 代码:  

#include 
#include 
using namespace std;
long long gcd(long long a, long long b){
    return b ? gcd(b, a%b) : a;
}
void ext_gcd(long long a, long long b, long long &x, long long &y){
    if (b == 0){
        x = 1;
        y = 0 ;
        return ;
    }
    ext_gcd(b, a%b, x, y);
    long long tmp = x;
    x = y;
    y = tmp - a/b * y;
    return ;
}

int main(){
    long long x, y, m, n, L;
    scanf("%I64d%I64d%I64d%I64d%I64d", &x, &y, &m, &n, &L);
    long long a = n - m;
    long long b = L;
    long long g = x - y;
    if (g % gcd(a,b)){
        printf("Impossible\n");
        return 0;
    }
    long long p = gcd(a, b);
    a /= p;
    b /= p;
    g /= p;
    long long ansx, ansy;
    ext_gcd(a, b, ansx, ansy);
    ansx *= g;
    long long tmp = (long long)abs(double(b));
    ansx = (ansx % tmp + tmp) % tmp;
    printf("%I64d\n", ansx);
    return 0;
}