[P5162] WD與積木

每種堆法（理解成名次序列，舉例3,3,8,2和7,7,100,2都對應2,2,1,3這個名次序列）等概率出現；題目中“兩種堆法不同當且僅當某個積木在兩種堆法中處於不同的層中”可見這是個組合問題，於是設一個名次的生成函數F(x)=sum i=1->inf 1*x^i/(i!) 。因為F(x)的常數項為0，有F^inf(x)=0。

名次序列的方案數的生成函數：A(x) = sum i=0->inf F^i(x) = 1/(1-F(x))
名次序列的各種方案的名次個數之和的生成函數：B(x) = sum i=0->inf i*F^i(x) = F(x)/(1-F(x))^2 = A^2(x)-A(x) = A(x)*(A(x)-1)

預處理出A(x)，B(x)的前MAX_N項繫數，對於n，答案為B(x)[n]/A(x)[n]。

註意事項

關於從0開始求和

如果是從1開始的情形：

A(x)=sum i=1->inf F^i(x)
= F(x)*(1-F^inf(x))/(1-F(x))
= F(x)/(1-F(x))
= (1-(1-F(x)))/(1-F(x))
= 1/(1-F(x))-1

B(x)=sum i=1-> inf i*F^i(x)
= ...
= F(x)/(1-F(x))^2-1

拋開A(x)[0]和B(x)[0]，似乎關係並沒有什麼影響。

關於EGF的捲積運算

某個生成函數F(x)=sum i=0->inf a[i]*x^i/(i!)

對於EGF意義下的捲積
F(x)*F(x)=sum i=0->inf (sum j=0->i C(i,j)*a[j]*a[i-j]) x^i/(i!)
= sum i=0->inf (sum j=0->i (i!)/(j!)/((i-j)!)*a[j]*a[i-j]) x^i/(i!)
= sum i=0->inf (sum j=0->i a[j]/(j!)*a[i-j]/((i-j)!)) x^i
與在OGF意義下的捲積
F(x)=sum i=0->inf b[i]*x^i ,b[i]=a[i]/(i!)
F(x)*F(x)=sum i=0->inf (sum j=0->i b[j]*b[i-j]) x^i
完全一致.

這告訴我們在使用EGF的具體計算時，可以將1/(i!)放到繫數中去用一般多項式演算法計算得到一個OGF結構，
所得到的函數的各繫數再乘上i!就是正確的EGF結構下的標誌序列了；顯然，該性質也適用於EGF的連續捲積！

換句話說,EGF的標誌序列儘管除以一個i!,然後與OGF無異地各種操作,最終結果再讓標誌序列乘上i!就好了。

代碼實現中因為除法關係沒有乘上i!。

參考實現

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

const int MAX_N=1e5+10;
const int P=998244353,G=3;

int a[MAX_N<<2];
int b[MAX_N<<2];
int c[MAX_N<<2];
int f[MAX_N<<2];
int r[MAX_N<<2];
int frr[MAX_N];

int qpow(long long x,int y) {
    long long c=1;
    for(; y; y>>=1, x=x*x%P)
        if(y&1) c=c*x%P;
    return c;
}
void ntt(int *a,int lmt,int tp) {
    for(int i=0; i<lmt; ++i) if(i<r[i]) swap(a[i],a[r[i]]);
    for(int m=1; m<lmt; m<<=1) {
        int wm=qpow(G,(P-1)/(m<<1));
        if(tp==-1) wm=qpow(wm,P-2);
        for(int i=0; i<lmt; i+=(m<<1)) {
            long long w=1, tmp;
            for(int j=0; j<m; ++j, w=w*wm%P) {
                tmp=w*a[i+j+m]%P;
                a[i+j+m]=(a[i+j]-tmp+P)%P;
                a[i+j]=(a[i+j]+tmp)%P;
            }
        }
    }
    if(tp==-1) {
        long long tmp=qpow(lmt,P-2);
        for(int i=0; i<lmt; ++i) a[i]=tmp*a[i]%P;
    }
}
void polyrev(int deg,int *a,int *b) {
    if(deg==1) {
        b[0]=qpow(a[0],P-2);
        return;
    }
    polyrev((deg+1)>>1,a,b);
    int lmt=1, l=0;
    while(lmt<(deg<<1)) lmt<<=1, ++l;
    for(int i=0; i<lmt; ++i) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
    for(int i=0; i<deg; ++i) c[i]=a[i];
    fill(c+deg,c+lmt,0);
    ntt(c,lmt,1), ntt(b,lmt,1);
    for(int i=0; i<lmt; ++i) b[i]=(2LL-1LL*c[i]*b[i]%P+P)%P*b[i]%P; 
    ntt(b,lmt,-1);
    fill(b+deg,b+lmt,0);
}

int main() {
    long long fac=1;
    for(int i=1; i<MAX_N; ++i) fac=fac*i%P;
    frr[MAX_N-1]=qpow(fac,P-2);
    for(int i=MAX_N-1; i; --i) frr[i-1]=1LL*i*frr[i]%P;
    for(int i=0; i<MAX_N; ++i) a[i]=P-frr[i];
    a[0]=(a[0]+2)%P;
    polyrev(MAX_N,a,b);
    for(int i=0; i<MAX_N; ++i) a[i]=f[i]=b[i];
    int lmt=1;
    while(lmt<=MAX_N) lmt<<=1; lmt<<=1;
    a[0]=(a[0]+P-1)%P;
    ntt(f,lmt,1), ntt(a,lmt,1);
    for(int i=0; i<lmt; ++i) f[i]=1LL*f[i]*a[i]%P;
    ntt(f,lmt,-1);
    
    int T,n; 
    scanf("%d",&T);
    while(T--) {
        scanf("%d",&n);
        printf("%lld\n",1LL*f[n]*qpow(b[n],P-2)%P);
    }
    return 0;
}

感謝@Weng_Weiji的答疑，思路也是來自於他（她？）的題解；代碼參考了@Owen_codeisking。
也感謝來自其它大佬的惡意嘲諷