2022杭电多校7 Connectivity of Erdős-Rényi Graph

题目描述

https://www.luogu.com.cn/problem/P4841

简要题意：现在有一个 $n$ 个点的无向图，每条边存在的概率都为 $p$，求图的连通块的期望个数

$n\le 5\times 10^5$

Solution

我们令 $f_n$ 表示 $n$ 个点联通的概率，$g_n$ 表示 $n$ 个点不连通的概率，那么我们知道最后的答案就是 $\sum_{i=1}^n\binom{n}{i}f_i(1-p)^{i(n-i)}$，我们首先考虑如何求 $f_n$

对于一个不连通的图，我们枚举 $1$ 号点所在连通块的大小，可以得到 $g_n=\sum_{i=1}^{n-1}\binom{n-1}{i-1}f_i(1-p)^{i(n-i)}$，我们知道 $g_n=1-f_n$，那么这个东西很像是一个分治 $NTT$ 的式子，但 $(1-p)^{i(n-i)}$ 无法拆开，我们考虑其组合意义，容易得到 $i(n-i)=\binom{n}{2}-\binom{i}{2}-\binom{n-i}{2}$，那么就可以拆成 $\frac{(1-p)^{\binom{n}{2}}}{(1-p)^{\binom{i}{2}}(1-p)^{\binom{n-i}{2}}}$，这样我们就可以做分治 $NTT$ 了，时间复杂度 $O(n\log^2 n)$，我们考虑继续化简

$$\begin{aligned} g_n&=\sum_{i=1}^{n-1}\binom{n-1}{i-1}f_i(1-p)^{i(n-i)}\\ \frac{1-f_n}{(n-1)!(1-p)^{\binom{n}{2}}}&=\sum_{i=1}^{n-1}\frac{f_i}{(i-1)!(1-p)^{\binom{i}{2}}}\frac{1}{(n-i)!(1-p)^{\binom{n-i}{2}}} \end{aligned}$$

我们令 $F(x)=\sum_{i=1}^n\frac{f_i}{(i-1)!(1-p)^{\binom{i}{2}}}x^i,G(x)=\sum_{i=1}^n\frac{1}{(i-1)!(1-p)^{\binom{i}{2}}}x^i,H(x)=\sum_{i=0}^n\frac{1}{i!(1-p)^{\binom{i}{2}}}x^i$，则可以得到 $\frac{G(x)}{H(x)+1}=F(x)$，需要注意 $H(x)=0$，那么这个式子可以用多项式求逆解决，时间复杂度 $O(n\log n)$

我们考虑指数公式定理，令 $g_n$ 表示 $n$ 个点无向图的概率，$f_n$ 表示 $n$ 个点联通的概率，显然有 $g_n=1$，同时我们令 $G(x)$ 和 $F(x)$ 分别为 $g_n$ 和 $f_n$ 的 $EGF$，那么我们应该有 $exp(F(x))=G(x)$，但实际上这并不成立，因为 $f_n$ 和 $f_m$ 合并贡献到 $g_{n+m}$ 时还需要乘上 $(1-p)^{nm}$ 这个系数，所以我们需要把 $\frac{1}{(1-p)^{\binom{n}{2}}}$ 的系数乘到 $f_n$ 和 $g_n$ 上，所以 $G(x)$ 和 $F(x)$ 应该为 $\frac{f_n}{(1-p)^{\binom{n}{2}}}$ 和 $\frac{g_n}{(1-p)^{\binom{n}{2}}}$ 的 $EGF$，那么我们只需要做一个 $\ln$ 即可，时间复杂度 $O(n\log n)$

下面给出求逆的代码

#include <iostream>
#include <vector>
#include <random>
#include <chrono>
#include <map>
#include <algorithm>
#define maxn 500010
#define ll long long
using namespace std;

const int p = 998244353;
inline int add(int x, int y) { return (x += y) >= p ? x - p : x; }
inline int mul(int x, int y) { return 1ll * x * y % p; }
inline int add(initializer_list<int> lst) { int s = 0; for (auto t : lst) s = add(s, t); return s; }
inline int mul(initializer_list<int> lst) { int s = 1; for (auto t : lst) s = mul(s, t); return s; }
ll pow_mod(ll x, ll n) {
    ll s = 1;
    for (; n; n >>= 1, x = x * x % p)
        if (n & 1) s = s * x % p;
    return s; 
}

#define Poly vector<int>
#define len(A) ((int) A.size())
namespace Pol {
    inline int add(int a, int b) { return (a += b) >= p ? a -= p : a; }
    inline int mul(int a, int b) { return 1ll * a * b % p; } 
    Poly operator - (const int &v, const Poly &a) {
        Poly res(a); 
        for (int i = 0; i < len(res); ++i) res[i] = p - res[i]; 
        res[0] = add(res[0], v); return res;
    }
    Poly operator - (const Poly &a, const int &v) {
        Poly res(a); res[0] = add(res[0], p - v); return res;
    }
    Poly operator * (const Poly &a, const int &v) {
        Poly res(a);
        for (int i = 0; i < len(res) ; ++i) res[i] = mul(res[i], v); 
        return res;
    }
 
    const int N = 4200000;
    const int G = 3;
    
    int P[N], inv[N], fac[N];
    void init_P(int n) {
        int l = 0; while ((1 << l) < n) ++l;
        for (int i = 0; i < n; ++i) P[i] = (P[i >> 1] >> 1) | ((i & 1) << l - 1);
    }
    void init_C() {
        if (fac[0]) return ;
        fac[0] = 1; for (int i = 1; i < N; ++i) fac[i] = mul(fac[i - 1], i);
        inv[N - 1] = pow_mod(fac[N - 1], p - 2); for (int i = N - 2; ~i; --i) inv[i] = mul(inv[i + 1], i + 1);
    }
    vector<int> init_W(int n) {
        vector<int> w(n); w[1] = 1;
        for (int i = 2; i < n; i <<= 1) {
            auto w0 = w.begin() + i / 2, w1 = w.begin() + i;
            int wn = pow_mod(G, (p - 1) / (i << 1));
            for (int j = 0; j < i; j += 2)
                w1[j] = w0[j >> 1], w1[j + 1] = mul(w1[j], wn);
        }
        return w; 
    } auto w = init_W(1 << 21);
    void DIT(Poly &a) {
        int n = len(a);
        for (int k = n >> 1; k; k >>= 1)
            for (int i = 0; i < n; i += k << 1)
                for (int j = 0; j < k; ++j) {
                    int x = a[i + j], y = a[i + j + k];
                    a[i + j + k] = mul(add(x, p - y), w[k + j]), a[i + j] = add(x, y);
                }
    }
    void DIF(Poly &a) {
        int n = len(a); 
        for (int k = 1; k < n; k <<= 1)
            for (int i = 0; i < n; i += k << 1)
                for (int j = 0; j < k; ++j) {
                    int x = a[i + j], y = mul(a[i + j + k], w[k + j]);
                    a[i + j + k] = add(x, p - y), a[i + j] = add(x, y);
                }
        int inv = pow_mod(n, p - 2);
        for (int i = 0; i < n; ++i) a[i] = mul(a[i], inv);
        reverse(a.begin() + 1, a.end());
    }
    Poly operator * (const Poly &A, const Poly &B) {
        int n = 1, n1 = len(A), n2 = len(B); while (n < n1 + n2 - 1) n <<= 1; init_P(n);
        Poly a(n), b(n);
        for (int i = 0; i < n1; ++i) a[i] = add(A[i], p);
        for (int i = 0; i < n2; ++i) b[i] = add(B[i], p);
        DIT(a); DIT(b);
        for (int i = 0; i < n; ++i) a[i] = mul(a[i], b[i]);
        DIF(a); a.resize(n1 + n2 - 1); return a; 
    }
    Poly MMul(const Poly &A, const Poly &B) { // 差卷积, 默认 A 和 B 的长度相同
        int n = 1, L = len(A); while (n < 2 * L - 1) n <<= 1; init_P(n);
        Poly a(n), b(n);
        for (int i = 0; i < L; ++i) a[i] = add(A[i], p);
        for (int i = 0; i < L; ++i) b[i] = add(B[i], p);
        reverse(b.begin(), b.begin() + L);
        DIT(a); DIT(b);
        for (int i = 0; i < n; ++i) a[i] = mul(a[i], b[i]);
        DIF(a); a.resize(L); reverse(a.begin(), a.end()); return a; 
    }
    Poly Der(const Poly &a) {
        Poly res(a);
        for (int i = 0; i < len(a) - 1; ++i) res[i] = mul(i + 1, res[i + 1]);
        res[len(a) - 1] = 0; return res;
    }
    Poly Int(const Poly &a) {
        static int inv[N];
        if (!inv[1]) {
            inv[1] = 1; 
            for (int i = 2; i < N; ++i) inv[i] = mul(p - p / i, inv[p % i]);
        }
        Poly res(a); res.resize(len(a) + 1);
        for (int i = len(a); i; --i) res[i] = mul(res[i - 1], inv[i]);
        res[0] = 0; return res;
    }
    Poly Inv(const Poly &a) {
        Poly res(1, pow_mod(a[0], p - 2)); 
        int n = 1; while (n < len(a)) n <<= 1; 
        for (int k = 2; k <= n; k <<= 1) {
            int L = 2 * k; init_P(L); Poly t(L);
            copy_n(a.begin(), min(k, len(a)), t.begin()); 
            t.resize(L); res.resize(L);
            DIT(res); DIT(t);
            for (int i = 0; i < L; ++i) res[i] = mul(res[i], add(2, p - mul(t[i], res[i])));
            DIF(res); res.resize(k);
        } res.resize(len(a)); return res;
    }
    Poly Offset(const Poly &a, int c) {
        int n = len(a); init_C();
        Poly t1(n), t2(n);
        for (int i = 0; i < n; ++i) t1[i] = mul(pow_mod(c, i), inv[i]);
        for (int i = 0; i < n; ++i) t2[i] = mul(a[i], fac[i]);
        t1 = MMul(t1, t2);
        for (int i = 0; i < n; ++i) t1[i] = mul(t1[i], inv[i]);
        return t1; 
    }
    pair<Poly, Poly> Divide(const Poly &a, const Poly &b) {
        int n = len(a), m = len(b); 
        Poly t1(a.rbegin(), a.rbegin() + n - m + 1), t2(b.rbegin(), b.rend()); t2.resize(n - m + 1);
        Poly Q = Inv(t2) * t1; Q.resize(n - m + 1); reverse(Q.begin(), Q.end());
        Poly R = Q * b; R.resize(m - 1); for (int i = 0; i < len(R); ++i) R[i] = add(a[i], p - R[i]);
        return make_pair(Q, R); 
    }
    Poly Ln(const Poly &a) {
        Poly res = Int(Der(a) * Inv(a));
        res.resize(len(a)); return res;
    }
    Poly Exp(const Poly &a) {
        Poly res(1, 1);
        int n = 1; while (n < len(a)) n <<= 1;
        for (int k = 2; k <= n; k <<= 1) {
            Poly t(res.begin(), res.end()); t.resize(k); t = Ln(t);
            for (int i = 0; i < min(len(a), k); ++i) t[i] = add(a[i], p - t[i]); t[0] = add(t[0], 1);
            res = res * t; res.resize(k);
        } res.resize(len(a)); return res;
    }
    Poly Sqrt(const Poly &a) { // a[0] = 1
        Poly res(1, 1); ll inv2 = pow_mod(2, p - 2); 
        int n = 1; while (n < len(a)) n <<= 1;
        for (int k = 2; k <= n; k <<= 1) {
            Poly t(res.begin(), res.end()), ta(a.begin(), a.begin() + min(len(a), k));
            t.resize(k); t = Inv(t) * ta; 
            res.resize(k); for (int i = 0; i < k; ++i) res[i] = mul(add(res[i], t[i]), inv2);
        } res.resize(len(a)); return res;
    }
    Poly Pow(const Poly &a, int k) { // a[0] = 1
        return Exp(Ln(a) * k);
    }
    Poly Pow(const Poly &a, int k, int kk) { 
        int n = len(a), t = n, m, v, inv, powv; Poly res(n);
        for (int i = n - 1; ~i; --i) if (a[i]) t = i, v = a[i];
        if (k && t >= (n + k - 1) / k) return res;
        if (t == n) { if (!k) res[0] = 1; return res; }
        m = n - t * k; res.resize(m);
        inv = pow_mod(v, p - 2); powv = pow_mod(v, kk); 
        for (int i = 0; i < m; ++i) res[i] = mul(a[i + (k > 0) * t], inv);
        res = Exp(Ln(res) * k); res.resize(n); 
        for (int i = m - 1; ~i; --i) {
            int tmp = mul(res[i], powv);
            res[i] = 0, res[i + t * k] = tmp;
        }
        return res;
    }
}  // namespace Pol

int fac[maxn], inv[maxn];
void init_C(int n) {
    fac[0] = 1; for (int i = 1; i <= n; ++i) fac[i] = mul(fac[i - 1], i);
    inv[n] = pow_mod(fac[n], p - 2); for (int i = n - 1; ~i; --i) inv[i] = mul(inv[i + 1], i + 1);
}

int a, b, n, m, q[maxn];

int pr, rp, pp[maxn], invpp[maxn];

void work() {
    cin >> m >> a >> b; pr = a * pow_mod(b, p - 2) % p; rp = add(1, p - pr);
    for (int i = 1; i <= m; ++i) cin >> q[i]; n = *max_element(q + 1, q + m + 1) + 1;
    for (int i = 0; i <= n; ++i) {
        pp[i] = pow_mod(rp, 1ll * i * (i - 1) / 2);
        invpp[i] = pow_mod(pp[i], p - 2);
    }
    Poly A(n), B(n);
    for (int i = 0; i < n; ++i) A[i] = mul(invpp[i], inv[i]); A[0] = 1; //A[0] = add(A[0], 1);
    for (int i = 1; i < n; ++i) B[i] = mul(invpp[i], inv[i - 1]);
    A = Pol::operator*(Pol::Inv(A), B);
    for (int i = 0; i < n; ++i) A[i] = mul({ A[i], fac[i - 1], inv[i] });
    for (int i = 0; i < n; ++i) B[i] = mul(invpp[i], inv[i]);
    Poly res = Pol::operator*(A, B);
    for (int i = 1; i <= m; ++i) cout << mul({ res[q[i]], fac[q[i]], pp[q[i]] }) << " \n"[i == m];
}

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

    int T; cin >> T; init_C(500000); 
    while (T--) work(); 
    return 0;
}