SRM 529 div1 med: MinskyMystery

問題

TopCoder Statistics - Problem Statement

解法

よくわからないので, とりあえずいろいろ実験するためのコードを書きます。

ll bag[5];

int main() {
    cin.tie(0);
    ios::sync_with_stdio(false);
    ll N;
    cin >> N;
    bag[0] = N;
    bag[1]++;
    ll cnt = 1;
    while (1) {
        cout << bag[0] << "  " << bag[1] << endl;
        bag[1]++;
        cnt++;
        bag[4] = 0;
        while (bag[0] != 0) {
            ll mini = min(bag[0], bag[1]);
            cnt += 2*mini;
            bag[0] -= mini;
            bag[1] -= mini;
            bag[2] += mini;
            bag[3] += mini;
            cnt++;
            bag[4]++;
            if (bag[0] == 0 && bag[1] == 0) {
                bag[4] += bag[3];
                cnt += bag[3];
                bag[3] = 0;
                cout << cnt << endl;
                return 0;
            }
            bag[1] += bag[3];
            cnt += bag[3];
            bag[3] = 0;
        }
        bag[0] += bag[2];
        cnt += bag[2];
        bag[2] = 0;
    }
    cout << cnt << endl;
    return 0;
}

while (1) の直後に
cout << bag[0] << " " << bag[1] << endl;
と書くと気づきますが, bag[0] は常に N という値を取り, bag[1] は, 1 ずつ増えていきます。で, bag[1] の値が bag[0] の値の約数になると, 操作は終了します。

ここまで気づくとやることはわかります。コードにすると, 以下のようなことをやればよいです。

ll bag[5];

int main() {
    cin.tie(0);
    ios::sync_with_stdio(false);
    ll N;
    cin >> N;
    bag[0] = N;
    bag[1]++;
    ll cnt = 1;
    while (1) {
        bag[1]++;
        cnt++;
        bag[4] = 0;
        if (bag[0]%bag[1]) {
            // while 前の bag[2], bag[3] の分
            cnt += 3*bag[0];
            // bag[4] の分
            cnt += bag[0]/bag[1]+1;
            // while 後の bag[2] の分
            cnt += bag[0];
        } else {
            cnt += bag[0]/bag[1] * (3*bag[1]+1);
            cout << cnt << endl;
            break;
        }
    }
    cout << cnt << endl;
    return 0;
}

これをこのままやると, N が 10^9 以上の素数の場合, bag[1] の値がかなりでかくなってしまうので間に合いません。
そこで, bag[0]/bag[1] の値の候補がそれほど多くないことを利用します。
実際, 1, 2, ..., までは普通に計算して, それより後は N/2+1 〜 N-1, N/3+1 〜 N/2, ... などと計算すると候補は O() 個だけしかないので, ここは二分探索すれば間に合います(といいつつ下のコードは 1.9 秒くらいかかってるので危ない)。

const ll MOD = 1e9+9;

class MinskyMystery {
public:
    int computeAnswer(long long N) {
        if (N < 2) return -1;
        vector<ll> div;
        for (ll d = 1; d*d <= N; d++) {
            if (N%d == 0) {
                div.push_back(d);
                if (d*d != N) div.push_back(N/d);
            }
        }
        sort(div.begin(), div.end());
        ll d;
        for (ll el : div) {
            if (el > 1) {
                d = el;
                break;
            }
        }
        cout << d << endl;
        ll ret = d+(N%MOD) * ((d-2) % MOD);
        ret += ((N/d) % MOD) * ((3*d+1)%MOD) % MOD;
        ret += (3*N%MOD) * ((d-2)%MOD) % MOD;
        ret %= MOD;
        for (ll p = 2; p < d; ) {
            ll low = p, high = d;
            while (high - low > 1) {
                ll med = (high+low)/2;
                if (N/p == N/med) low = med;
                high = med;
            }
            low = min(low, d-1);
            ll num = low-p+1;
            ret += (num%MOD) * ((N/p+1)%MOD) % MOD;
            ret %= MOD;
            p = low+1;
        }
        return ret;
    }
};