每种堆法（理解成名次序列，举例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。

也感谢来自其它大佬的恶意嘲讽