Luogu5748 集合划分计数

发布于 22 天前  41 次阅读


Link

Sol

贝尔数就是一行第二类斯特林数的和,关于贝尔数和第二类斯特林数之间一些转换关系,此处不再赘述

首先给出贝尔数的递推公式

B_n = \sum_{i = 1}^{n} B_{n - i} \binom{n - 1}{i - 1}

组合意义显然意见,把集合划分成 n 个非空子集的方案数就可以先枚举某一特定元素所在集合的大小,比如说考虑 1 这个元素所在的集合是 \binom{n - 1}{i - 1} 种,那么剩下的元素划分的方案数就是 B_{n - i}

i - 1 不好搞,把枚举 i 换成枚举 i - 1 ,就变成了

B_n = \sum_{i = 0}^{n - 1} B_{n - i - 1} \binom{n - 1}{i}

把组合数拆开,并且把两边同时乘上 x^{n - 1} , 得到:

B_n x^{n - 1} = \sum_{i = 0}^{n - 1} \frac{(n-1)! x^{n-1}}{i! (n - i - 1)!} B_{n - i - 1}

这玩意儿你一看十分不好做,但是如果观察看见 x^{i} \times x^{n - i - 1} = x^{n - 1} ,事情就简单了不少

B_n x^{n - 1} = \sum_{i = 0} ^{n - 1} \frac{x^{i}}{i!} \cdot \frac{B_{n - i - 1}}{(n - i - 1)!} x^{n - i - 1} \cdot (n - 1)!

(n - 1)! 项除过去,现在已经变成了:

\frac{B_n}{(n - 1)!} x^{n - 1} = \sum_{i = 0}^{n - 1} \frac{x^{i}}{i!} \cdot \frac{B_{n - i - 1}}{(n - i - 1)!} x^{n - i - 1}

对其求导,得:

\frac{\mathrm{d}}{\mathrm{d}x} (\frac{B_n}{n!} x ^{n}) = \sum_{i + j = n - 1} \frac{x^i}{i!} x^i \cdot \frac{B_j}{j!} x^j

设贝尔数的EGF为F(x) , 上面的式子已经是一个 EGF 的形式了,所以重写变成

\frac{\mathrm{d} }{\mathrm{d} x} F'(x) = F(x) \cdot e^x

u = F(x) , 两边同时乘以算子 \mathrm{d} x

\mathrm{d} u = u \cdot e^x \cdot \mathrm{d} u

把右边的 u 除过去

\frac{\mathrm{d} u}{u} = e^x \cdot \mathrm{d} u

积分得

\int \frac{\mathrm{d} u}{u} = \int e^x \mathrm{d} x

\ln u = e^x + C

x = 1 带入 EGF 的时候可以算出一个 C = -1 ,所以

F(x) = e^{e^x - 1}

因为 e^x - 1 = \sum_{n \geq 1} \frac{1}{n!} x^n

所以把阶乘逆元作为初始多项式再做 exp 就好了

Code

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

const int SIZE = 2e6 + 5, siz = 1e5;
const int G = 3, Gi = 332748118, mod = 998244353;

int n;
int f[SIZE], g[SIZE], c[SIZE], lne[SIZE], iv[SIZE], cal[SIZE], pos[SIZE], W[SIZE], WN[SIZE], fac[SIZE], inv[SIZE];

namespace GTR {
    const int bufl = 1 << 15;
    char buf[bufl], *s = buf, *t = buf;
    inline int fetch() {
        if (s == t) { t = (s = buf) + fread(buf, 1, bufl, stdin); if (s == t) return EOF; }
        return *s++;
    }
    inline int read() {
        int a = 0, b = 1, c = fetch();
        while (c < 48 || c > 57) b ^= c == '-', c = fetch();
        while (c >= 48 && c <= 57) a = (a << 1) + (a << 3) + c - 48, c = fetch();
        return b ? a : -a;
    }
}
using GTR::read;

int qPow(int a, int b) {
    int ans = 1;
    for ( ; b; b >>= 1, a = a * a % mod) {
        if (b & 1) ans = ans * a % mod;
    }
    return ans % mod;
}

namespace poly {
    void NTT(int *a, int n, int opt) {
        int w, wn, t, *a0, *a1;
        for (int i = 0; i < n; ++ i) if (i < pos[i]) std::swap(a[i], a[pos[i]]);
        for (int i = 1, j, k; i < n; i <<= 1) {
            if (opt == 1) W[i << 1] ? wn = W[i << 1] : W[i << 1] = wn = qPow(G, (mod - 1) / (i << 1));
            else if (opt == -1) WN[i << 1] ? wn = WN[i << 1] : WN[i << 1] = wn = qPow(Gi, (mod - 1) / (i << 1));
            for (j = 0; j < n; j += (i << 1))
                for (w = 1, k = 0, a1 = i + (a0 = a + j); k < i; ++ k, ++ a0, ++ a1, (w *= wn) %= mod)
                    t = *a1 * w % mod, *a1 = (*a0 - t + mod) % mod, (*a0 += t) %= mod;
        }
        if (opt == 1) return;
        int inv = qPow(n, mod - 2);
        for (int i = 0; i < n; ++ i) (a[i] *= inv) %= mod;
    }

    void derive(int n, int *a, int *b) {
        for (int i = 1; i < n; ++ i) b[i - 1] = i * a[i] % mod; b[n - 1] = 0;
    }

    void calculus(int n, int *a, int *b) {
        for (int i = 1; i < n; ++ i) b[i] = a[i - 1] * qPow(i, mod - 2) % mod; b[0] = 0;
    }

    void inv(int ind, int *a, int *b) {
        if (ind == 1) { return b[0] = qPow(a[0], mod - 2), void(); }
        inv((ind + 1) >> 1, a, b);
        int N = 1, M = ind << 1, len = 0;
        for ( ; N <= M; N <<= 1, ++ len);
        for (int i = 0; i < N; ++ i)
            pos[i] = (pos[i >> 1] >> 1 | ((i & 1) << (len - 1))), c[i] = (i >= ind ? 0 : a[i]);
        NTT(c, N, 1); NTT(b, N, 1);
        for (int i = 0; i < N; ++ i) b[i] = (2ll - (b[i] * c[i] % mod) + mod) % mod * b[i] % mod;
        NTT(b, N, -1);
        for (int i = ind; i < N; ++ i) b[i] = 0;
    }

    void ln(int bas, int *a, int *b) {
        derive(bas, a, cal); inv(bas, a, iv);
        int N = 1, M = bas << 1, len = 0;
        for ( ; N <= M; N <<= 1, ++ len);
        for (int i = 0; i < N; ++ i) pos[i] = (pos[i >> 1] >> 1 | ((i & 1) << (len - 1)));
        NTT(cal, N, 1); NTT(iv, N, 1);
        for (int i = 0; i < N; ++ i) cal[i] = cal[i] * iv[i] % mod;
        NTT(cal, N, -1);
        calculus(N, cal, b);
    }

    void exp(int ind, int *a, int *b) {
        if (ind == 1) { return b[0] = 1, void(); }
        exp((ind + 1) >> 1, a, b);
        int N = 1, M = ind << 1, len = 0;
        for ( ; N <= M; N <<= 1, ++ len);
        for (int i = 0; i < (M << 1); ++ i) iv[i] = lne[i] = 0;
        ln(ind, b, lne);
        for (int i = 0; i < N; ++ i) pos[i] = (pos[i >> 1] >> 1 | ((i & 1) << (len - 1)));
        for (int i = 0; i < N; ++ i) c[i] = (i >= ind ? 0 : a[i]);
        NTT(b, N, 1); NTT(lne, N, 1); NTT(c, N, 1);
        for (int i = 0; i < N; ++ i) b[i] = (1ll - (lne[i] - c[i] + mod) % mod + mod) % mod * b[i] % mod;
        NTT(b, N, -1);
        for (int i = ind; i < N; ++ i) b[i] = 0;
    }
}

signed main() {
    fac[0] = 1ll, f[0] = 0;
    for (int i = 1; i <= siz; ++ i) fac[i] = fac[i - 1] * i % mod;
    inv[siz] = qPow(fac[siz], mod - 2);
    for (int i = siz - 1; i; -- i) inv[i] = inv[i + 1] * (i + 1) % mod;
    for (int i = 1; i <= siz; ++ i) f[i] = inv[i];
    poly::exp(siz + 1, f, g);   
    int tn = read();
    while (tn -- ) {
        n = read();
        printf("%lld\n", g[n] * fac[n] % mod);
    }
    return 0;
}

生存 关心