結果

問題 No.2829 GCD Divination
ユーザー Tatsu_mrTatsu_mr
提出日時 2024-08-02 22:42:13
言語 C++17
(gcc 12.3.0 + boost 1.83.0)
結果
AC  
実行時間 109 ms / 2,000 ms
コード長 4,535 bytes
コンパイル時間 2,898 ms
コンパイル使用メモリ 217,580 KB
実行使用メモリ 237,340 KB
最終ジャッジ日時 2024-08-02 22:42:19
合計ジャッジ時間 5,232 ms
ジャッジサーバーID
（参考情報） 		judge1 / judge2
このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 AC 2 ms
6,816 KB
testcase_01 AC 2 ms
6,940 KB
testcase_02 AC 72 ms
198,528 KB
testcase_03 AC 2 ms
6,940 KB
testcase_04 AC 90 ms
237,340 KB
testcase_05 AC 109 ms
205,708 KB
testcase_06 AC 89 ms
171,672 KB
testcase_07 AC 82 ms
154,548 KB
testcase_08 AC 59 ms
104,456 KB
testcase_09 AC 47 ms
87,356 KB
testcase_10 AC 41 ms
123,068 KB
testcase_11 AC 5 ms
7,920 KB
testcase_12 AC 17 ms
41,496 KB
testcase_13 AC 61 ms
182,544 KB
testcase_14 AC 44 ms
125,732 KB
testcase_15 AC 40 ms
99,764 KB
testcase_16 AC 18 ms
46,880 KB
testcase_17 AC 12 ms
26,608 KB
testcase_18 AC 11 ms
20,796 KB
testcase_19 AC 13 ms
33,212 KB
testcase_20 AC 2 ms
6,940 KB
testcase_21 AC 39 ms
87,440 KB
testcase_22 AC 30 ms
88,472 KB
testcase_23 AC 50 ms
124,772 KB
testcase_24 AC 56 ms
163,132 KB
testcase_25 AC 65 ms
195,388 KB
testcase_26 AC 30 ms
85,152 KB
testcase_27 AC 58 ms
168,400 KB
testcase_28 AC 14 ms
33,468 KB
testcase_29 AC 36 ms
100,652 KB
testcase_30 AC 69 ms
201,456 KB
testcase_31 AC 12 ms
25,428 KB
testcase_32 AC 40 ms
112,524 KB
testcase_33 AC 30 ms
81,900 KB
testcase_34 AC 56 ms
159,556 KB
権限があれば一括ダウンロードができます

ソースコード

diff # 

#include <bits/stdc++.h>
using namespace std;

using ld = long double;

template <class T>
T modpow(T a, T b, T mod) {
    T cur = a % mod, res = 1 % mod;
    while (b) {
        if (b & 1) {
            res = (res * cur) % mod;
        }
        cur = (cur * cur) % mod;
        b >>= 1;
    }
    return res;
}

bool MillerRabin(long long n) {
    if (n <= 1) {
        return false;
    }
    if (n == 2 || n == 7 || n == 61) {
        return true;
    }
    if (n % 2 == 0) {
        return false;
    }
    
    vector<long long> A;
    if (n < 4759123141) {
        A = {2, 7, 61};
    } else {
        A = {2, 325, 9375, 28178, 450775, 9780504, 1795265022};
    }
    long long s = 0, d = n - 1;
    while (d % 2 == 0) {
        s++;
        d >>= 1;
    }
    for (auto a : A) {
        if (a % n == 0) {
            return true;
        }
        long long x = modpow<__int128_t>(a, d, n);
        if (x == 1) {
            continue;
        }
        bool ok = false;
        for (int i = 0; i < s; i++) {
            if (x == n - 1) {
                ok = true;
                break;
            }
            x = (__int128_t)x * x % n;
        }
        if (!ok) {
            return false;
        }
    }
    return true;
}

long long gcd(long long x, long long y) {
    if (y == 0) {
        return x;
    }
    return gcd(y, x % y);
}

unsigned int xorshift() {
    static unsigned int x = 123456789, y = 362436069, z = 521288629, w = 88675123;
    unsigned int t = (x ^ (x << 11));
    x = y;
    y = z;
    z = w;
    return (w = (w ^ (w >> 19)) ^ (t ^ (t >> 8)));
}

long long Pollard(long long n) {
    if (n % 2 == 0) {
        return 2LL;
    }
    if (MillerRabin(n)) {
        return n;
    }
    
    long long i = 0;
    while (true) {
        i++;
        long long r = xorshift();
        auto f = [&](long long x) {
            return (__int128_t(x) * x + r) % n;
        };
        long long x = i, y = f(x);
        while (true) {
            long long p = gcd(abs(y - x + n), n);
            if (p == 0 || p == n) {
                break;
            }
            if (p != 1) {
                return p;
            }
            x = f(x);
            y = f(f(y));
        }
    }
}

vector<long long> prime_factorize(long long n) {
    if (n == 1) {
        return {};
    }
    long long p = Pollard(n);
    if (p == n) {
        return {p};
    }
    vector<long long> l = prime_factorize(p);
    vector<long long> r = prime_factorize(n / p);
    for (auto x : r) {
        l.emplace_back(x);
    }
    sort(l.begin(), l.end());
    return l;
}

vector<long long> divisors(long long n) {
    if (n == 1) {
        return {1LL};
    }
    auto divisor_dfs = [&](auto divisor_dfs, vector<pair<long long, long long>> &p, long long t, int cur, vector<long long> &res) -> void {
        if (cur == p.size()) {
            res.push_back(t);
            return;
        }
        divisor_dfs(divisor_dfs, p, t, cur + 1, res);
        for (int i = 0; i < p[cur].second; i++) {
            t *= p[cur].first;
            divisor_dfs(divisor_dfs, p, t, cur + 1, res);
        }
    };
    
    vector<long long> res, pf = prime_factorize(n);
    
    vector<pair<long long, long long>> p;
    long long cnt = 1, now = pf[0];
    for (int i = 1; i < (int)pf.size(); i++) {
        if (pf[i] == now) {
            cnt++;
        } else {
            p.push_back({now, cnt});
            now = pf[i];
            cnt = 1;
        }
    }
    p.push_back({now, cnt});
    
    divisor_dfs(divisor_dfs, p, 1, 0, res);
    sort(res.begin(), res.end());
    return res;
}

int main() {
    long long n;
    cin >> n;
    vector<ld> dp(n + 1, -1);
    dp[1] = 0.0;
    auto f = [&](auto f, long long x) -> ld {
        if (dp[x] != -1) {
            return dp[x];
        }
        ld res = (ld)x;
        vector<long long> divs = divisors(x);
        vector<int> cnt(divs.back() + 1);
        for (int i = (int)divs.size() - 1; i >= 0; i--) {
            int d = divs[i];
            if (i == (int)divs.size() - 1) {
                cnt[d] = 1;
                continue;
            }
            int num = x / d;
            for (int j = i + 1; j < (int)divs.size(); j++) {
                if (divs[j] % divs[i] == 0) {
                    num -= cnt[divs[j]];
                }
            }
            cnt[divs[i]] = num;
            res += f(f, divs[i]) * (ld)num;
        }
        res /= (ld)(x - 1);
        dp[x] = res;
        return res;
    };
    cout << fixed << setprecision(15) << f(f, n) << endl;
}
0