Submission #61843132

---

Source Code Expand

#include <bits/stdc++.h>
using namespace std;
using ll=long long;
#define rep(i,a,b) for (int i=(int)(a);i<(int)(b);i++)
#define all(p) p.begin(),p.end()
#include <atcoder/modint>
using mint = atcoder::modint;

namespace po167{
template<class T>
struct Binomial{
    std::vector<T> fact_vec, fact_inv_vec;
    void extend(int m = -1){
        int n = fact_vec.size();
        if (m == -1) m = n * 2;
        if (n >= m) return;
        fact_vec.resize(m);
        fact_inv_vec.resize(m);
        for (int i = n; i < m; i++){
            fact_vec[i] = fact_vec[i - 1] * T(i);
        }
        fact_inv_vec[m - 1] = T(1) / fact_vec[m - 1];
        for (int i = m - 1; i > n; i--){
            fact_inv_vec[i - 1] = fact_inv_vec[i] * T(i);
        }
    }
    Binomial(int MAX = 0){
        fact_vec.resize(1, T(1));
        fact_inv_vec.resize(1, T(1));
        extend(MAX + 1);
    }

    T fact(int i){
        if (i < 0) return 0;
        while (int(fact_vec.size()) <= i) extend();
        return fact_vec[i];
    }
    T invfact(int i){
        if (i < 0) return 0;
        while (int(fact_inv_vec.size()) <= i) extend();
        return fact_inv_vec[i];
    }
    T C(int a, int b){
        if (a < b || b < 0) return 0;
        return fact(a) * invfact(b) * invfact(a - b);
    }
};
}

// https://ferin-tech.hatenablog.com/entry/2019/08/11/ラグランジュ補間
// x座標が相異なるn+1点(x_i,y_i)を通るn次以下の多項式f(x)を返す
// O(n^2)
#define REP(i, a) for (int i = 0; i < a; i++)
#define FOR(i, a, b) rep(i, a, b)
vector<mint> lagrange_interpolation(vector<mint> x, vector<mint> y) {
    const ll n = x.size() - 1;
    REP(i, n+1) {
        mint t = 1;
        REP(j, n+1) if(i != j) t *= x[i]-x[j];
        y[i] /= t;
    }
    ll cur = 0, nxt = 1;
    vector<vector<mint>> dp(2, vector<mint>(n+2));
    dp[0][0] = -1 * x[0], dp[0][1] = 1;
    FOR(i, 1, n+1) {
        REP(j, n+2) {
            dp[nxt][j] = dp[cur][j] * -1 * x[i];
            if(j >= 1) dp[nxt][j] += dp[cur][j-1];
        }
        swap(nxt, cur);
    }
    vector<mint> xinv(n+1);
    REP(i, n+1) xinv[i] = x[i].inv();
    vector<mint> ret(n+1);  // f(x)
    REP(i, n+1) {
        if(y[i]==0) continue;
        // 0割り対策の場合分け
        if(x[i] == 0) {
            REP(j, n+1) ret[j] += dp[cur][j+1] * y[i];
        } else {
            ret[0] -= dp[cur][0] * xinv[i] * y[i];
            mint pre = -1 * dp[cur][0] * xinv[i];
            FOR(j, 1, n+1) {
                ret[j] -= (dp[cur][j] - pre) * xinv[i] * y[i];
                pre = -1 * (dp[cur][j] - pre) * xinv[i];
            }
        }
    }
    return ret;
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int N, P;
    cin >> N >> P;
    mint::set_mod(P);
    int H = N / 2 ;
    int M = N * (N - 1) / 2;
    vector<mint> ans(M + 1);
    po167::Binomial<mint> table(N);
    
    rep(rp, 0, M + 1) {
        vector dp(N + 1, vector(H + 1, vector<mint>(H + 1)));
        dp[1][1][1] = 1;
        vector<mint> p2(N * N + 1, 1);
        rep(i, 0, N * N) p2[i + 1] = p2[i] * (rp + 2);
        vector p3(N, vector<mint>(N));
        rep(i, 0, N){
            p3[i][0] = 1;
            rep(j, 1, N){
                p3[i][j] = p3[i][j - 1] * (p2[i] - 1);
            }
        }
        rep(sum, 1, N) rep(fr, 1, min(H, sum) + 1) rep(a, fr, H + 1) if(dp[sum][fr][a].val()){
            rep(k, 1, min(N - sum, H) + 1) {
                int b = sum - a + k;
                if (b > H) break;
                mint tmp = dp[sum][fr][a] * table.C(N - sum, k);
                tmp *= p3[fr][k];
                tmp *= p2[k * (k - 1) / 2];
                dp[sum + k][k][b] += tmp;
            }
        }
        rep(i, 0, H + 1) ans[rp] += dp[N][i][H];
    }
    
    vector<mint> X(M + 1);
    rep(i, 0, M + 1) X[i] = i + 1;
    auto res = lagrange_interpolation(X, ans);
    rep(i, N - 1, M + 1){
        cout << res[i].val() << (i == M ? "\n" : " ");
    }
}

Submission Info

Judge Result