思路与算法

k位存储意味着每次加C后按M=2^k回绕。执行t次后的值为(A+tC)%M，要等于B，即Ct与(B-A+M)%M模M同余。使用扩展欧几里得判断可解性并求最小非负解。无解输出FOREVER，A=B时不进入循环所以次数0；C=0也按同一方程处理。

正确性说明

循环第t次判断之前的变量值满足该模递推。退出恰好等价于模方程成立。设g=gcd(C,M)，若g不整除位置差则任何t都不能退出；否则约去g后的逆元给出唯一剩余类，规范化最小非负代表即首次退出次数。C=0且A=B时所有t可行，最小0；C=0且A不同于B则无解。终止四个0在构造2^k之前排除。

复杂度

欧几里得步骤O(log M)=O(k)，栈空间O(k)，其余O(1)，每组初始化O(1)。用1LL<<k构造M，k=32时仍正确。逆元和约分差都小于2^32，乘积最大(2^32-1)^2，超过signed64但小于unsigned64上界，所以先转unsigned long long乘法再取模。答案最大2^32-1，用long long输出，不截断为int。

易错点

不能写1<<32，不能用signed64直接乘两个32位归一化值。A=B必须输出0，不是一个周期。C=0、偶数C无解和k=32全覆盖。FOREVER全大写，四个0是结束标记。独立答案用Python任意精度，并对600个小位宽状态环直接模拟。

参考代码（C++11）

#include <iostream>
#include <string>
using namespace std;
long long exgcd(long long a, long long b, long long &x, long long &y) {
    if (!b) {
        x = 1;
        y = 0;
        return a;
    }
    long long u, v;
    long long g = exgcd(b, a % b, u, v);
    x = v;
    y = u - (a / b) * v;
    return g;
}
long long solve(long long c, long long delta, long long mod) {
    long long x, y, g = exgcd(c, mod, x, y);
    if (delta % g)
        return -1;
    long long reduced = mod / g;
    if (reduced == 1)
        return 0;
    x = (x % reduced + reduced) % reduced;
    unsigned long long product = (unsigned long long)x * (unsigned long long)(delta / g);
    return product % (unsigned long long)reduced;
}
int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    long long a, b, c, k;
    while (cin >> a >> b >> c >> k) {
        if (a == 0 && b == 0 && c == 0 && k == 0)
            break;
        long long mod = 1LL << k;
        long long answer = solve(c, (b - a + mod) % mod, mod);
        if (answer < 0)
            cout << "FOREVER\n";
        else
            cout << answer << '\n';
    }
}