WD与积木

转载

mb5ff980b461ced 2021-04-01 20:44:00

文章标签 i++ 组合数递推前缀和分治 文章分类 代码人生

XVI.WD与积木

本题有两种思路。

首先，两种思路共同的地方在于都将期望化成了\(\dfrac{\text{所有方案一共的层数}}{\text{总共的方案数}}\)。我们设其为\(\dfrac{f_n}{g_n}\)。

思路1：从DP开始

我们先考虑求出\(g_n\)。

我们有

\[g_n=\sum\limits_{i=1}^n\dbinom{n}{i}g_{n-i} \]

它的意义是枚举第一层里放了哪些积木。

拆开组合数，得到

\[\dfrac{g_n}{n!}=\sum\limits_{i=1}^n\dfrac{1}{i!}\times\dfrac{g_{n-i}}{(n-1)!} \]

我们考虑设\(G_i=\dfrac{g_n}{n!},H_i=\dfrac{1}{i!}\)，

则即有

\[G_n=\sum\limits_{i=1}^nH_iG_{n-i} \]

到这里已经可以分治FFT求出了。但是我们还有更好的方法。

我们计算\(GH\)，就得到了

\[(GH)_n=\sum\limits_{i=0}^nH_iG_{n-i} \]

带回上面的递推式，我们发现这实际上就是

\[(GH)_n=\sum\limits_{i=1}^nH_iG_{n-i}+H_0G_n=2G_n \]

但是这个递推式有一个特例，即\(n=0\)。此时有\((GH)_0=2G_0-1\)。

故最终有

\[GH=2G-1 \]

处理就得到了

\[G=\dfrac{1}{2-H} \]

求逆即可求出\(G\)。

我们再考虑求出\(f_n\)。

我们有

\[f_n=\sum\limits_{i=1}^n\dbinom{n}{i}(f_{n-i}+g_{n-i}) \]

其中，\(i\)还是枚举第一层里放了多少积木；\(f_{n-i}\)是剔除第一层后剩下部分的总共积木数；\(g_{n-i}\)因为对于剩下\(n-i\)块积木，它们的每一种方法都增加了\(1\)层。

还是仿照我们之前的样子，拆开组合数并设\(F_i=\dfrac{f_i}{i!}\)，就有

\[F_n=\sum\limits_{i=1}^nH_i(F_{n-i}+G_{n-i}) \]

计算\(H(F+G)\)，有

\[\Big(H(F+G)\Big)_n=\sum\limits_{i=0}^nH_i(F_{n-i}+G_{n-i})=2F_n+G_n \]

照例看一下\(f_0\)，发现此时仍有\(\Big(H(F+G)\Big)_0=2F_0+G_0\)（\(F_0=0,G_0=H_0=1\)）

故最终有

\[H(F+G)=2F+G \]

于是

\[F=\dfrac{G-GH}{H-2} \]

代入\(G=\dfrac{1}{2-H}\)，最终就得到

\[F=G(G-1) \]

然后最终答案就是\(\dfrac{f_n}{g_n}=\dfrac{F_n}{G_n}\)。

此种方案的代码：

#include<bits/stdc++.h>
using namespace std;
const int M=1e5;
const int N=1<<20;
const int mod=998244353;
const int all=18;
int ksm(int x,int y){
	int rt=1;
	for(;y;x=(1ll*x*x)%mod,y>>=1)if(y&1)rt=(1ll*rt*x)%mod;
	return rt;
}
namespace Poly{
	int rev[N];
	const int G=3;
	void NTT(int *a,int tp,int LG){
		int lim=(1<<LG),invlim=ksm(lim,mod-2);
		for(int i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(LG-1));
		for(int i=0;i<lim;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
		for(int md=1;md<lim;md<<=1){
			int rt=ksm(G,(mod-1)/(md<<1));
			if(tp==-1)rt=ksm(rt,mod-2);
			for(int stp=md<<1,pos=0;pos<lim;pos+=stp){
				int w=1;
				for(int i=0;i<md;i++,w=(1ll*w*rt)%mod){
					int x=a[pos+i],y=(1ll*w*a[pos+md+i])%mod;
					a[pos+i]=(x+y)%mod;
					a[pos+md+i]=(x-y+mod)%mod;
				}
			}
		}
		if(tp==-1)for(int i=0;i<lim;i++)a[i]=(1ll*a[i]*invlim)%mod;
	}
	int A[N],B[N],C[N];
	void mul(int *a,int *b,int *c,int LG){//using: Array A and B
		int lim=(1<<LG);
		for(int i=0;i<lim;i++)A[i]=B[i]=0;
		for(int i=0;i<(lim>>1);i++)A[i]=a[i],B[i]=b[i];
		NTT(A,1,LG),NTT(B,1,LG);
		for(int i=0;i<lim;i++)A[i]=1ll*A[i]*B[i]%mod;
		NTT(A,-1,LG);
		for(int i=0;i<lim;i++)c[i]=A[i];
	}
	void inv(int *a,int *b,int LG){//using: Array C
		b[0]=ksm(a[0],mod-2);
		for(int k=1;k<=LG+1;k++){
			mul(b,a,C,k);
			for(int i=0;i<(1<<k);i++)C[i]=(mod-C[i])%mod;
			(C[0]+=2)%=mod;
			mul(C,b,b,k);
		}
	}	
}
int a[N];
int f[N],g[N],h[N],res[N];
int fac[N],inv[N];
int T,n;
int main(){
	fac[0]=1;
	for(int i=1;i<=M;i++)fac[i]=1ll*fac[i-1]*i%mod;
	inv[M]=ksm(fac[M],mod-2);
	for(int i=M-1;i>=0;i--)inv[i]=1ll*inv[i+1]*(i+1)%mod;
	for(int i=0;i<=M;i++)h[i]=inv[i];
	for(int i=0;i<=M;i++)a[i]=mod-h[i];
	(a[0]+=2)%=mod;
	Poly::inv(a,g,all);
	for(int i=0;i<=M;i++)a[i]=g[i];
	(a[0]+=mod-1)%=mod;
	Poly::mul(a,g,f,all+1);
	for(int i=0;i<=M;i++)res[i]=1ll*f[i]*ksm(g[i],mod-2)%mod;
	scanf("%d",&T);
	while(T--)scanf("%d",&n),printf("%d\n",res[n]);
	return 0;
}

思路2：从生成函数开始。

仍然先考虑\(g_n\)。

我们考虑一个一层放多少积木的EGF，\(F(x)=\sum\limits_{i=1}^n\dfrac{1}{i!}=e^x-1\)

EGF的卷积就是将这些积木标号的过程；所以我们最终的答案就是

\[\hat{g}=\sum\limits_{i=0}^{\infty}F^i \]

它是\(F\)的一阶前缀和，则有

\[\hat{g}=\dfrac{1}{1-F}=\dfrac{1}{2-e^x} \]

通过复原\(e^x\)为\(\sum\limits_{i=0}^n\dfrac{1}{i!}\)可以直接通过多项式求逆求出\(\hat{g}\)。

我们再考虑\(f_n\)。

则我们只需要乘上一个层数即可；则有

\[\hat{f}=\sum\limits_{i=0}^{\infty}iF^i \]

它是\(F\)的二阶前缀和；于是就有