#127922 (C++11) No.215 素数サイコロと合成数サイコロ (3-Hard)

提出ソース
結果

問題	No.215 素数サイコロと合成数サイコロ (3-Hard)
ユーザー	pekempey
提出日時	2016-10-31 20:33:45
言語	C++11 (gcc 11.4.0)
結果	AC
実行時間	3,385 ms / 4,000 ms
コード長	11,457 bytes
コンパイル時間	2,055 ms
コンパイル使用メモリ	179,452 KB
実行使用メモリ	24,416 KB
最終ジャッジ日時	2024-11-25 00:14:21
合計ジャッジ時間	12,498 ms
ジャッジサーバーID （参考情報）	judge1 / judge2
このコードへのチャレンジ
（要ログイン）
テストケース

テストケース表示
入力	結果	実行時間実行使用メモリ
testcase_00	AC	3,251 ms 24,100 KB
testcase_01	AC	3,385 ms 24,416 KB
権限があれば一括ダウンロードができます
ソースコード

raw source code
// #pragma GCC optimize ("O3")
// #pragma GCC target ("avx")
#include <bits/stdc++.h>
typedef __int128_t int128;

// worst: 1152921504606839475 300 300 => 610505353
//        1000000000000000000 300 300 => 817328043

// 今回は不使用
class NumberTheoreticTransform {
public:
    int mod;
    int root;

public:
    NumberTheoreticTransform(int mod, int root) : mod(mod), root(root) { }

    std::vector<int> mul(std::vector<int> a, std::vector<int> b) {
        int s = a.size() + b.size() - 1;
        int t = 1;
        while (t < s) t *= 2;
        for (int i = 0; i < a.size(); i++) a[i] %= mod;
        for (int i = 0; i < b.size(); i++) b[i] %= mod;
        a.resize(t);
        b.resize(t);
        ntt(a);
        ntt(b);
        for (int i = 0; i < a.size(); i++) a[i] = mul(a[i], b[i]);
        ntt(a, true);
        a.resize(s);
        return a;
    }

private:
    int mul(int x, int y) {
        return int64_t(x) * y % mod;
    }

    int add(int x, int y) {
        return (x += y) >= mod ? x - mod : x;
    }

    int pow(int x, int y) {
        int res = 1;
        for (; y > 0; x = mul(x, x), y >>= 1) if (y & 1) res = mul(res, x);
        return res;
    }

    int inv(int x) {
        return pow(x, mod - 2);
    }

    void ntt(std::vector<int> &a, bool rev = false) {
        int n = a.size();
        int h = 0;
        for (int i = 0; 1 << i < n; i++) h++;
        for (int i = 0; i < n; i++) {
            int j = 0;
            for (int k = 0; k < h; k++) j |= (i >> k & 1) << (h - 1 - k);
            if (i < j) std::swap(a[i], a[j]);
        }
        for (int i = 1; i < n; i *= 2) {
            int w = pow(root, (mod - 1) / (i * 2));
            if (rev) w = inv(w);
            for (int j = 0; j < n; j += i * 2) {
                int wn = 1;
                for (int k = 0; k < i; k++) {
                    int s = a[j + k + 0];
                    int t = mul(a[j + k + i], wn);
                    a[j + k + 0] = add(s, t);
                    a[j + k + i] = add(s, mod - t);
                    wn = mul(wn, w);
                }
            }
        }
        int v = inv(n);
        if (rev) for (int i = 0; i < n; i++) a[i] = mul(a[i], v);
    }
};

// 今回は不使用
// mod を取りつつ Karatsuba 法
class Karatsuba {
private:
    int mod;
    
public:
    Karatsuba(int mod) : mod(mod) {}

    std::vector<int> mul(std::vector<int> a, std::vector<int> b) {
        int s = std::max<int>(a.size(), b.size());
        int t = 1;
        int u = a.size() + b.size() - 1;
        while (t < s) t *= 2;
        a.resize(t);
        b.resize(t);
        std::vector<int> c(t * 6);
        mul(a.data(), b.data(), t, c.data());
        c.resize(u);
        return c;
    }

private:
    void mul(int a[], int b[], int n, int res[]) {
        int *a0 = a, *a1 = a + n / 2;
        int *b0 = b, *b1 = b + n / 2;
        int *c0 = res + n * 5, *c1 = res + n * 5 + n / 2;
        int *x0 = res, *x1 = res + n, *x2 = res + n * 2;
        if (n <= 8) {
            for (int i = 0; i < n * 2; i++) res[i] = 0;
            for (int i = 0; i < n; i++) for (int j = 0; j < n; j++) upd(res[i + j], int64_t(a[i]) * b[j] % mod);
            return;
        }
        for (int i = 0; i < n / 2; i++) {
            c0[i] = add(a0[i], a1[i]);
            c1[i] = add(b0[i], b1[i]);
        }
        mul(a0, b0, n / 2, x0);
        mul(a1, b1, n / 2, x1);
        mul(c0, c1, n / 2, x2);
        for (int i = 0; i < n; i++) upd(x2[i], mod - add(x0[i], x1[i]));
        for (int i = 0; i < n; i++) upd(res[i + n / 2], x2[i]);
    }

    int add(int x, int y) {
        return (x += y) >= mod ? x - mod : x;
    }

    void upd(int &x, int y) {
        x = add(x, y);
    }
};

// int128 を使って mod を省いた Karatsuba 法
class Karatsuba128 {
private:
    int mod;
    
public:
    Karatsuba128(int mod) : mod(mod) {}

    std::vector<int> mul(std::vector<int> a, std::vector<int> b, int pmod = -1) {
        if (pmod == -1) pmod = a.size() + b.size() - 1;
        else {
            if (a.size() > pmod) a.resize(pmod);
            if (b.size() > pmod) b.resize(pmod);
        }
        int s = std::max<int>(a.size(), b.size());
        int t = 1;
        while (t < s) t *= 2;
        std::vector<int128> A(t), B(t), C(t * 6);
        for (int i = 0; i < a.size(); i++) A[i] = a[i];
        for (int i = 0; i < b.size(); i++) B[i] = b[i];
        mul(A.data(), B.data(), t, C.data());
        std::vector<int> c(pmod);
        for (int i = 0; i < std::min<int>(pmod, C.size()); i++) {
            c[i] = C[i] % mod;
            if (c[i] < 0) c[i] += mod;
        }
        return c;
    }

private:
    void mul(int128 a[], int128 b[], int n, int128 res[]) {
        if (n <= 8) {
            memset(res, 0, sizeof(res[0]) * n * 2);
            for (int i = 0; i < n; i++) for (int j = 0; j < n; j++) {
                res[i + j] += a[i] * b[j];
            }
            return;
        }
        int128 *a0 = a, *a1 = a + n / 2;
        int128 *b0 = b, *b1 = b + n / 2;
        int128 *c0 = res + n * 5, *c1 = res + n * 5 + n / 2;
        int128 *x0 = res, *x1 = res + n, *x2 = res + n * 2;
        for (int i = 0; i < n / 2; i++) {
            c0[i] = a0[i] + a1[i];
            c1[i] = b0[i] + b1[i];
        }
        mul(a0, b0, n / 2, x0);
        mul(a1, b1, n / 2, x1);
        mul(c0, c1, n / 2, x2);
        for (int i = 0; i < n; i++) x2[i] -= x0[i] + x1[i];
        for (int i = 0; i < n; i++) res[i + n / 2] += x2[i];
    }
};

// 今回は不使用
// 3種の mod を使って CRT で復元する畳み込み
class AnyModConvolution {
private:
    std::vector<NumberTheoreticTransform> ntt;
    std::vector<int> mods;
    std::vector<int> roots;
    int mod;

public:
    AnyModConvolution(int mod) : mod(mod) {
        mods.push_back(469762049);
        mods.push_back(167772161);
        mods.push_back(754974721);
        roots.push_back(3);
        roots.push_back(3);
        roots.push_back(11);
        for (int i = 0; i < 3; i++) {
            ntt.emplace_back(mods[i], roots[i]);
        }
    }

    std::vector<int> mul(std::vector<int> a, std::vector<int> b) {
        std::vector<std::vector<int>> res(ntt.size());
        for (int i = 0; i < ntt.size(); i++) {
            res[i] = ntt[i].mul(a, b);
        }
        for (int i = 0; i < res[0].size(); i++) {
            std::vector<int> x;
            for (int j = 0; j < ntt.size(); j++) {
                x.push_back(res[j][i]);
            }
            res[0][i] = linear_congruence(x, mods, mod);
        }
        return res[0];
    }

    int pow(int64_t x, int64_t y, int mod) {
        int64_t res = 1;
        for (; y > 0; x = x * x % mod, y >>= 1) if (y & 1) res = res * x % mod;
        return res;
    }

    int inv(int x, int mod) {
        return pow(x, mod - 2, mod);
    }

    int linear_congruence(std::vector<int> x, std::vector<int> m, int mod) {
        m.push_back(mod);
        std::vector<int> u(x.size());
        for (int i = 0; i <= x.size(); i++) {
            int64_t s = 0, v = 1;
            for (int j = 0; j < i; j++) {
                (s += u[j] * v) %= m[i];
                (v *= m[j]) %= m[i];
            }
            if (i == x.size()) return s;
            u[i] = (x[i] + m[i] - s) * inv(v, m[i]) % m[i];
        }
        abort();
        return -1;
    }
};

class AnyModPolynomial {
private:
    int mod;
    Karatsuba128 kara128;

    std::vector<int> t;

public:
    AnyModPolynomial(int mod) : mod(mod), kara128(mod) {}

    std::vector<int> mul(const std::vector<int> &a, const std::vector<int> &b) {
        return kara128.mul(a, b);
    }

    std::vector<int> add(const std::vector<int> &a, const std::vector<int> &b) {
        std::vector<int> res(std::max<int>(a.size(), b.size()));
        for (int i = 0; i < res.size(); i++) {
            res[i] = add(i < a.size() ? a[i] : 0, i < b.size() ? b[i] : 0);
        }
        return res;
    }

    std::vector<int> sub(const std::vector<int> &a, const std::vector<int> &b) {
        std::vector<int> res(std::max<int>(a.size(), b.size()));
        for (int i = 0; i < res.size(); i++) {
            res[i] = add(i < a.size() ? a[i] : 0, i < b.size() ? (mod - b[i]) : 0);
        }
        return res;
    }

    std::vector<int> quot(std::vector<int> a, std::vector<int> b) {
        if (a.size() < b.size()) return{ 0 };
        while (b.back() == 0) b.pop_back();
        int s = a.size() - b.size() + 1;
        while (!is_pow2(a.size() - b.size() + 1)) a.push_back(0);
        int n = a.size() - b.size() + 1;
        if (t.empty()) t = inv(rev(b), n);
        std::vector<int> q = rev(mul(rev(a), t, n));
        q.resize(s);
        return q;
    }

    std::vector<int> rem(std::vector<int> a, std::vector<int> b) {
        while (b.back() == 0) b.pop_back();
        std::vector<int> r = sub(a, mul(quot(a, b), b, b.size() - 1));
        r.resize(b.size() - 1);
        return r;
    }

    std::vector<int> remmul(std::vector<int> a, std::vector<int> b, std::vector<int> c) {
        return rem(mul(a, b), c);
    }

private:
    int pow(int64_t x, int64_t y, int mod) {
        int64_t res = 1;
        for (; y > 0; x = x * x % mod, y >>= 1) if (y & 1) res = res * x % mod;
        return res;
    }

    int inv(int x, int mod) {
        return pow(x, mod - 2, mod);
    }

    int is_pow2(int n) {
        return (n & (n - 1)) == 0;
    }

    int mul(int x, int y) {
        return int64_t(x) * y % mod;
    }

    int add(int x, int y) {
        return (x += y) >= mod ? x - mod : x;
    }

    void upd(int &x, int y) {
        x = add(x, y);
    }

    int pow(int x, int y) {
        int res = 1;
        for (; y > 0; x = mul(x, x), y >>= 1) if (y & 1) res = mul(res, x);
        return res;
    }

    int inv(int x) {
        return pow(x, mod - 2);
    }

    std::vector<int> mul(std::vector<int> a, std::vector<int> b, int n) {
        return kara128.mul(a, b, n);
    }

    std::vector<int> inv(std::vector<int> a, int n) {
        assert(is_pow2(n));
        std::vector<int> b(n);
        b[0] = inv(a[0]);
        for (int i = 0; (1 << i) < a.size(); i++) {
            b = sub(add(b, b), mul(mul(b, b, n), a, n));
        }
        return b;
    }

    std::vector<int> rev(std::vector<int> a) {
        reverse(a.begin(), a.end());
        return a;
    }
};

const int mod = 1e9 + 7;

int64_t dp[2][301][4000], way[10000];

void calc(int64_t dp[301][4000], int64_t n, std::vector<int> dice) {
    dp[0][0] = 1;
    for (int k = 0; k < dice.size(); k++) {
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < 3910; j++) {
                (dp[i + 1][j + dice[k]] += dp[i][j]) %= mod;
            }
        }
    }
}

int main() {
    using namespace std;
    int64_t N, P, C;
    cin >> N >> P >> C;
    calc(dp[0], P, { 2, 3, 5, 7, 11, 13 });
    calc(dp[1], C, { 4, 6, 8, 9, 10, 12 });
    for (int i = 0; i <= 3900; i++) {
        for (int j = 0; j <= 3900; j++) {
            (way[i + j] += dp[0][P][i] * dp[1][C][j]) %= mod;
        }
    }
    int64_t K = P * 13 + C * 12 + 1;
    vector<int> a(K + 1);
    for (int i = 0; i < K; i++) a[i] = mod - way[K - i];
    a[K] = 1;
    int64_t M = N + K - 1;
    AnyModPolynomial poly(mod);
    vector<int> f = { 1 }, g = { 0, 1 };
    for (; M > 0; g = poly.remmul(g, g, a), M >>= 1) if (M & 1) f = poly.remmul(f, g, a);
    int64_t ans = 0;
    for (int i = 0; i < K; i++) ans += int64_t(f[i]);
    cout << ans % mod << endl;
}
yukicoder

結果

テストケース

ソースコード
