题解 luogu P5162 WD与积木

传送门

---

【分析】

令 (f_{n, m}) 表示 (n) 个积木划分成 (m) 层的方案数

由于层之间的有序性，不难写出 (displaystyle f_{n, m}=sum_{i=1}^n dbinom n i f_{n-i, m-1}=sum_{i+j=n}dbinom n i f_{n-i, m-1}cdot [i>0])

将第二维作为生成函数，则 (displaystyle hat F_m[n]=sum_{i+j=n}hat F_{m-1}[j]cdot(e^x-1)[i])

故 (displaystyle hat F_m(x)=hat F_{m-1}(x)cdot (e^x-1)=hat F_1(x)cdot (e^x-1)^{m-1}=(e^x-1)^{m-1}cdot sum_{n=1}^infty {1over n!}x^n=(e^x-1)^m)

所以层数的期望可以等价为所有方案的层数，乘上所有的方案数的逆元

而 (hat F_m[n]) 表示 (n) 个积木划分成 (m) 层的方案数

则 (displaystyle left(sum_{m=0}^infty hat F_m(x) ight)[n]cdot n!) 表示 (n) 个积木的方案数

且 (displaystyle left(sum_{m=0}^infty mhat F_m(x) ight)[n]cdot n!) 表示 (n) 个积木所有方案的层数和

前者为 (displaystyle sum_{m=0}^infty f^m(x)={1over 1-f(x)}) 的形式，答案代入得到 (displaystyle {1over 2-e^x})

后者为 (displaystyle sum_{m=0}^infty mf^m(x)=f(x)cdot { ext dover ext df(x)}({1over 1-f(x)})={f(x)over (1-f(x))^2}) ，答案代入得到 ((displaystyle {1over (2-e^x)^2}-{1over 2-e^x}))

因此我们记 (displaystyle F(x)={1over 2-e^x}, G(x)={1over (2-e^x)^2})

则答案应为 (displaystyle {(G[n]-F[n])cdot n!over F[n]cdot n!}={G[n]over F[n]}-1)

因此直接求逆、卷积，最后 (O(nlog n)) 预处理答案即可

---

【代码】

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef pair<ll, ll> pii;
typedef double db;
#define fi first
#define se second

const int LimBit=18, M=1<<LimBit<<1, P=998244353;
inline ll kpow(ll a,ll x) { ll ans=1; for(;x;x>>=1,a=a*a%P) if(x&1) ans=ans*a%P; return ans; }
int a[M], b[M], c[M];
struct NTT{
    static const int G=3;
    int N_, N, na, nb, invN;
    int Stk[M], curStk;
    NTT() {
        curStk=0;
    }

    int w[2][M], rev[M];
    inline void work() {
        if(N==N_) return ;
        N_=N;
        int d=__builtin_ctz(N);
        w[0][0]=w[1][0]=1;
        for(int i=1, x=kpow(G, (P-1)/N), y=kpow(x, P-2); i<N; ++i) {
            rev[i] = (rev[i>>1]>>1)|((i&1)<<d-1);
            w[0][i]=(ll)x*w[0][i-1]%P, w[1][i]=(ll)y*w[1][i-1]%P;
        }
        invN=kpow(N, P-2);
    }
    inline void FFT(int *a, int f) {
        for(int i=0; i<N; ++i) if(i<rev[i]) swap(a[i], a[rev[i]]);
        for(int i=1; i<N; i<<=1)
            for(int j=0, t=N/(i<<1); j<N; j+=i<<1)
                for(int k=0, l=0, x, y; k<i; ++k, l+=t)
                    x=(ll)w[f][l]*a[j+k+i]%P, y=a[j+k], a[j+k+i]=(y-x+P)%P, a[j+k]=(y+x)%P;
        if(f) for(int i=0; i<N; ++i) a[i]=(ll)a[i]*invN%P;
    }
    inline void doit(int *a, int *b, int na, int nb) {
        for(N=1; N<na+nb-1; N<<=1);
        for(int i=na; i<N; ++i) a[i]=0;
        for(int i=0; i<nb; ++i) c[i]=b[i];
        for(int i=nb; i<N; ++i) c[i]=0;
        work(); FFT(a, 0); FFT(c, 0);
        for(int i=0; i<N; ++i) a[i]=(ll)a[i]*c[i]%P;
        FFT(a, 1);
    }
    inline void get_inv(int *f, int *g, int n) {
        int pos=curStk;
        for(int i=n; i>1; i=i+1>>1) Stk[++curStk]=i;
        g[0]=kpow(f[0], P-2);
        N=2;
        for(int l=Stk[curStk]; curStk>pos; l=Stk[--curStk]) {
            N<<=(N<l+l-1);
            for(int i=0; i<l; ++i) a[i]=f[i];
            for(int i=l; i<N; ++i) a[i]=g[i]=0;
            work(); FFT(a, 0); FFT(g, 0);
            for(int i=0; i<N; ++i) g[i]=g[i]*(2-(ll)g[i]*a[i]%P+P)%P;
            FFT(g, 1);
            for(int i=l; i<N; ++i) g[i]=0;
        }
    }
}ntt;

int f[M], g[M], lim=1e5+1;;
inline void init() {
    f[0]=P-1;
    for(int i=1; i<lim; ++i)
        f[i]=(ll)i*f[i-1]%P;
    f[lim-1]=kpow(f[lim-1], P-2);
    for(int i=lim-1; i; --i)
        f[i-1]=(ll)i*f[i]%P;
    f[0]=1;

    ntt.get_inv(f, g, lim);
    for(int i=0; i<lim; ++i) f[i]=g[i];
    ntt.doit(g, f, lim, lim);

    for(int i=0; i<lim; ++i)
        f[i]=( g[i]*kpow(f[i], P-2)+P-1 )%P;
}

int main(){
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);
    init();
    int T, n; cin>>T;
    while(T--&&cin>>n) cout<<f[n]<<"
";
    cout.flush();
    return 0;
}

---

可以考虑到 (displaystyle { ext dover ext dx}({1over 2-e^x})={2over (2-e^x)^2}-{1over 2-e^x})

即 (F'(x)=2G(x)-F(x))

对比系数得到 (displaystyle {g[n]-f[n]over f[n]}={f[n+1](n+1)+f[n]-2f[n]over 2f[n]}={f[n+1](n+1)over 2f[n]}-{1over 2})

省去一次多项式乘法

int f[M], g[M], lim=1e5+2;
inline void init() {
    f[0]=P-1;
    for(int i=1; i<lim; ++i)
        f[i]=(ll)i*f[i-1]%P;
    f[lim-1]=kpow(f[lim-1], P-2);
    for(int i=lim-1; i; --i)
        f[i-1]=(ll)i*f[i]%P;
    f[0]=1;

    ntt.get_inv(f, g, lim);
    for(int i=0; i<lim; ++i)
        f[i]=( g[i+1]*(i+1ll)%P*kpow(g[i]<<1, P-2)%P+(P>>1) )%P;
}