[数论]扩展欧几里得算法

http://poj.org/problem?id=2115

void extgcd(ll a, ll b, ll &g, ll &x, ll &y)
{
    if (b == 0) g = a, x = 1, y = 0;
    else{
        extgcd(b, a % b, g, y, x);
        y -= (a / b) * x;
    }
}

此算法能计算出ax + by = gcd(a, b)的整数解。
可以解决ax + by = c的问题，如果c % gcd(a, b) == 0，一定有整数解，否则一定无整数解。
结果乘上c / gcd(a, b)就可以了。

int main()
{
    ll a, b, c, d, x, y, g;
    while (cin >> a >> b >> c >> d){
        if (a == 0 && b == 0 && c == 0 && d == 0) break;
        if (b - a == 0){
            cout << "0" << endl;
            continue;
        }
        extgcd(c, (ll)1 << d, g, x, y);
        if ((b - a) % g) cout << "FOREVER" << endl;
        else{
            ll t = ((ll)1 << d) / g;
            x = (((x * (b - a) / g) % t) + t) % t;
            cout << x << endl;
        }
    }
    return 0;
}

这题可以得到这个方程：(a + c * x) % 2^d = b
化简得：c * x – 2^d * y = b – a;
因为gcd(a, b) = gcd(|a|, |b|)
于是就可以使用扩展欧几里得算法
判断一下是否有解
然后再同乘倍数：x = (((x * (b – a) / g) % t) + t) % t
为了防止出现负数或者出现t，先把结果%t再加t再%t。
关键是t是什么东西，因为要在最早等于b的时候就退出循环，因为g已经整除(b – a)了，所以在模2^d的一个循环里就必定有答案，x每次加1，c * x就每次加c，模2^{d的循环周期就是2}d / gcd(c, 2^d)，所以x对这个取余就是正确答案了。