拉格朗日反演（暂时鸽）与CF1349F2（xtq F2）

多项式真难

拉格朗日反演

先鸽着。

接下来看个小 例 题

CF1349F Slime and Sequences

这就是世界顶尖的计数水平么

先来看简单版。

序列不好计数，我们考虑把它转化成排列。

考虑对于一个长度为(n)排列(a)，我们记(s_i)表示满足(1le j<i,a_j<a_{j+1})的(j)的数量。

那么，如果我们在序列的第(a_i)位填上(s_i+1)，我们就得到了一个合法的序列。

然后，我们对于一个长度为(n)的合法序列(p)，我们可以将所有(i)以(p_i)从小到大为第一关键字，(i)从大到小为第二关键字排序，就得到了一个排列。

于是我们就将排列和合法的序列一一对应了。

然后我们考虑统计答案，记(ans_x)表示(x+1)在所有合法序列中出现的次数和。

那么，我们有

[ans_x=sum_{i=max(1,x)}^nf(i,x)inom{n}{i}(n-i)! ]

其中(f(i,j))表示长度为(i)的排列中有(j)个满足(a_k<a_{k+1})的(k)的方案数。

也即我们枚举这个数在排列的(s_i)中的出现位置(i)（也即出现在对应序列的(a_i)位置），然后我们要求(s_i=x)，然后从(n)个数中选出(i)个填入前(i)位，后面(n-i)位随便填。

然后我们考虑(f(x,y))怎么求。

如果是简单版本，xjb算一下就能做到(O(n^2))了。

对于难的版本，我们考虑优化。

我们记(g(x,y))表示长度为(x)的排列，至少有(y)个满足条件位置的方案数。

那么我们有

[g(x,y)=sum_{jge y}^xinom{j}{y}f(x,j) ]

二项式反演得

[f(x,y)=sum_{jge y}^x(-1)^{j-y}inom{j}{y}g(x,j) ]

然后我们考虑(g)怎么求，如果我们确定位置(i)有(a_i<a_{i+1})，我们就在(i)，(i+1)之间连一条边。那么由于我们有(y)条边，那么有(x-y)个连通块。由于我们要将(x)个数分到(x-y)个连通块中，而连通块内部的顺序是确定的，因此考虑使用EGF

这个EGF应该形如

[T(z)=sum_{ige1}frac{z^i}{i!}=e^z-1 ]

也即

[g(x,y)=x![z^x]T^{x-y}(z) ]

那么我们现在考虑答案长啥样。

由于(ans_0)比较特殊，我们可以把它去掉，这样求和符号下面的(max)就没有了。

[egin{split} ans_x&=sum_{i=x}^nf(i,x)inom{n}{i}(n-i)! \ &=sum_{i=x}^ninom{n}{i}(n-i)!sum_{j=x}^{i-1}(-1)^{j-x}inom{j}{x}g(i,j) \ &=sum_{i=x}^ninom{n}{i}(n-i)!sum_{j=x}^{i-1}(-1)^{j-x}inom{j}{x}i![z^i]T^{i-j}(z) \ &=n!sum_{i=x}^nsum_{j=x}^i(-1)^{j-x}inom{j}{x}[z^i]T^{i-j}(z) \ &=n!sum_{j=x}^n(-1)^{j-x}frac{j!}{x!(j-x)!}sum_{i=j}^n[z^i]T^{i-j}(z) end{split} ]

我们记(V_{x}=sum_{i=x}^n[z^i]T^{i-x}(z))那么上式可以化简为

[ans_x=frac{n!}{x!}sum_{i=x}^nfrac{(-1)^{i-x}}{(i-x)!} imes i!V_i ]

这是一个类似卷积的形式，可以通过一些翻转变成可以卷积的形式。

那么我们现在来求(V_x)。

[egin{split} V_x&=sum_{i=x}^n[z^i]T^{i-x}(z) \ &=sum_{i=x}^n[z^i](e^z-1)^{i-x} \ &=sum_{i=x}^n[z^x](frac{e^z-1}{z})^{i-x} \ &=sum_{i=0}^{n-x}[z^x](frac{e^z-1}{z})^i end{split} ]

记(F(z)=frac{e^z-1}{z})，那么这个式子形如等比数列求和

[egin{split} V_x&=sum_{i=0}^{n-x}[z^x]F^i(z) \ &=[z^x]frac{1-F^{n-x+1}(z)}{1-F(z)} \ &=[z^x]frac{1}{1-F(z)}+[z^x]frac{F^{n-x+1}(z)}{F(z)-1} end{split} ]

前半部分就是求个逆，不过由于分母的常数项为(0)，应该将分母的幂次降低(1)后求(z^{x-1})项系数。

再来考虑后半部分，设(S_x=[z^x]frac{F^{n-x+1}(z)}{F(z)-1})。

由于不太好求了，我们考虑引入一个新的未知数(w)

[egin{split} S_x&=[z^x]frac{F^{n-x+1}(z)}{F(z)-1} \ &=[z^{n+1}]frac{(zF(z))^{n-x+1}}{F(z)-1} \ &=[z^{n+1},w^{n-x+1}]frac{(wzF(z))^{n-x+1}}{F(z)-1} \ &=[z^{n+1},w^{n-x+1}]frac{1}{F(z)-1}sum_{ige0}(wzF(z))^i \ &=[z^{n+1},w^{n-x+1}]frac{1}{F(z)-1} imesfrac{1}{1-wzF(z)} end{split} ]

考虑拉格朗日反演。我们构造(G(z),H(z),P(z))，满足(H(G(z))=frac{1}{F(z)-1} imesfrac{1}{1-wzF(z)})，(P(G(z))=z)。

首先令(G(z)=zF(z))，我们再构造(T(z))满足(T(G(z))=F(z))，显然(T(z)=frac{z}{P(z)})。

于是(H(z)=frac{1}{T(z)-1} imesfrac{1}{1-wz})。

再来考虑(P(z))，由于(G(z)=zF(z)=e^z-1)，所以(P(e^z-1)=z)，(P(z)=ln(z+1))。

于是(T(z)=frac{z}{ln(z+1)})。

那么现在我们做拉格朗日反演（注意之后的求导将(w)作为常数）

[[z^n]H(G(z))=frac{1}{n}[z^{n-1}]H'(z)(frac{z}{P(z)})^n ]

鉴于我们要(z^{n+1})项，我们把式子里的(n)加上(1)

[[z^{n+1}]H(G(z))=frac{1}{n+1}[z^n]H'(z)(frac{z}{P(z)})^{n+1} ]

然后考虑式子里的(H'(z))

[egin{split} H'(z)&=frac{1}{1-wz}(frac{1}{T(z)-1})'+frac{1}{T(z)-1}(frac{1}{1-wz})' \ &=frac{1}{1-wz} imes-frac{T'(z)}{(T(z)-1)^2}+frac{1}{T(z)-1} imesfrac{w}{(1-wz)^2} \ &=-frac{T'(z)}{(1-wz)(T(z)-1)^2}+frac{w}{(1-wz)^2(T(z)-1)} end{split} ]

代入（公式可能会超出页面...有点难看）

[egin{split} &[z^{n+1}]H(G(z)) \ &=frac{1}{n+1}[z^n](-frac{T'(z)}{(1-wz)(T(z)-1)^2}+frac{w}{(1-wz)^2(T(z)-1)})(frac{z}{P(z)})^{n+1} \ &=frac{1}{n+1}[z^n](-frac{T'(z)}{(T(z)-1)^2} imesfrac{1}{1-wz}+frac{1}{(T(z)-1)} imesfrac{w}{(1-wz)^2})(frac{z}{P(z)})^{n+1} \ &=frac{1}{n+1}[z^n,w^m](-frac{T'(z)}{(T(z)-1)^2} imessum_{ige0}(wz)^i+frac{1}{(T(z)-1)} imes wsum_{ige0}(wz)^i imessum_{ige0}(wz)^i)(frac{z}{P(z)})^{n+1} \ &=frac{1}{n+1}[z^n,w^m](-frac{T'(z)}{(T(z)-1)^2} imessum_{ige0}(wz)^i+frac{1}{(T(z)-1)} imessum_{ige0}(i+1)z^iw^{i+1})(frac{z}{P(z)})^{n+1} \ &=frac{1}{n+1}[z^n](-frac{T'(z)}{(T(z)-1)^2} imes z^m+frac{1}{(T(z)-1)} imes mz^{m-1})(frac{z}{P(z)})^{n+1} end{split} ]

注意到(P(z)=frac{z}{T(z)})

[egin{split} [z^{n+1}]H(G(z))&=frac{1}{n+1}[z^n](-frac{T'(z)}{(T(z)-1)^2} imes z^m+frac{1}{(T(z)-1)} imes mz^{m-1})(frac{z}{P(z)})^{n+1} \ &=frac{1}{n+1}[z^n](-frac{T'(z)}{(T(z)-1)^2} imes z^m+frac{1}{(T(z)-1)} imes mz^{m-1})T^{n+1}(z) \ &=frac{1}{n+1}[z^n](-frac{T'(z)T^{n+1}(z)}{(T(z)-1)^2} imes z^m+frac{T^{n+1}(z)}{(T(z)-1)} imes mz^{m-1}) \ &=frac{1}{n+1}([z^{n-m}](-frac{T'(z)T^{n+1}(z)}{(T(z)-1)^2})+[z^{n-m+1}](frac{mT^{n+1}(z)}{(T(z)-1)})) end{split} ]

突然发现(T(z)-1)的常数项也是(0)（悲）

于是

[[z^{n+1}]H(G(z))=frac{1}{n+1}([z^{n-m+2}](-frac{T'(z)T^{n+1}(z)}{(z^{-1}(T(z)-1))^2})+[z^{n-m+2}](frac{mT^{n+1}(z)}{z^{-1}(T(z)-1)})) ]

然后照着做就是了。

真是简单呢（口区）

code:

注：如果你的实现不太优秀可能会被卡常（

#include<bits/stdc++.h>
#define vi vector<int>
using namespace std;
namespace poly{
	#define vi vector<int>
	#define ci const int&
	#define Red(x) (x+=(x>>31)&mod)
	const int LM=(1<<22),mod=998244353;
	int lm,lg[LM+10],rev[LM+10],rt[LM+10][2],iv[LM+10],*p,*q;
	int POW(int x,int y){
		int ret=1;
		while(y)y&1?ret=1ll*ret*x%mod:0,x=1ll*x*x%mod,y>>=1;
		return ret;
	}
	void NTT(vi&f,ci op){
		int tn=f.size(),l=lg[tn],r,t1,t2;
		long long nr;
		for(int i=0;i<tn;++i)rev[i]=(rev[i>>1]>>1)+(i&1)*(1<<l-1),rev[i]<i?swap(f[rev[i]],f[i]),0:0;
		for(int i=2;i<=tn;i<<=1){
			r=rt[i][op];
			for(int j=0;j<tn;j+=i){
				nr=1,p=&f[j],q=&f[j+(i>>1)];
				for(int k=j;k<j+(i>>1);++k,nr=nr*r%mod,++p,++q)t1=*p,t2=nr*(*q)%mod,(*p)-=(((*p)=t1+t2)>=mod?mod:0),(*q)=t1-t2,Red((*q));
			}
		}
		if(op)for(int i=0;i<tn;++i)f[i]=1ll*f[i]*iv[tn]%mod;
	}
	vi Poly(ci x){
		vi ret;
		return ret.push_back(x),ret;
	}
	vi Plus(vi x,vi y){
		int sz=max(x.size(),y.size());
		x.resize(sz),y.resize(sz);
		for(int i=0;i<sz;++i)(x[i]+=y[i])>=mod?x[i]-=mod:0;
		return x;
	}
	vi Minus(vi x,vi y){
		int sz=max(x.size(),y.size());
		x.resize(sz),y.resize(sz);
		for(int i=0;i<sz;++i)x[i]-=y[i],Red(x[i]);
		return x;
	}
	vi Mul(vi x,ci y){
		for(int i=0;i<x.size();++i)x[i]=1ll*x[i]*y%mod;
		return x;
	}
	vi Mul(vi x,vi y,ci sz){
		int tl=x.size()+y.size()-1,lth=1;
		while(lth<tl)lth<<=1;
		x.resize(lth),y.resize(lth),NTT(x,0),NTT(y,0);
		for(int i=0;i<lth;++i)x[i]=1ll*x[i]*y[i]%mod;
		NTT(x,1),x.resize(sz);
		return x;
	}
	vi Inv(vi x){
		if(x.size()==1)return x[0]=POW(x[0],mod-2),x;
		vi tmp=x;
		int ts=x.size(),sz=(ts+1>>1);
		tmp.resize(sz),tmp=Inv(tmp);
		int tl=ts+tmp.size()+tmp.size()-2,lth=1;
		while(lth<tl)lth<<=1;
		x.resize(lth),tmp.resize(lth),NTT(x,0),NTT(tmp,0);
		for(int i=0;i<lth;++i)tmp[i]=(2-1ll*x[i]*tmp[i])%mod*tmp[i]%mod,Red(tmp[i]);
		NTT(tmp,1);
		return tmp.resize(ts),tmp;
	}
	vi Ln(vi x){
		vi tmp=x;
		for(int i=1;i<tmp.size();++i)tmp[i-1]=1ll*i*tmp[i]%mod;
		tmp.pop_back(),tmp=Mul(tmp,Inv(x),x.size());
		for(int i=x.size()-1;i>0;--i)tmp[i]=1ll*tmp[i-1]*iv[i]%mod;
		tmp[0]=0;
		return tmp;
	}
	vi Exp(vi x){
		if(x.size()==1)return x[0]=1,x;
		int sz=(x.size()+1>>1);
		vi tmp=x,t2;
		tmp.resize(sz),tmp=Exp(tmp),t2=tmp;
		t2.resize(x.size());
		return Mul(tmp,Plus(Minus(Poly(1),Ln(t2)),x),x.size());
	}
	vi POW(vi x,ci y,ci yc){
		return Exp(Mul(Ln(x),y));
	}
	void init(ci x){
		lm=1;
		while(lm<x)lm<<=1;
		for(int i=2;i<=lm;++i)lg[i]=lg[i>>1]+1,lg[i]!=lg[i-1]?rt[i][0]=POW(3,(mod-1)/i),rt[i][1]=POW(rt[i][0],mod-2):0;
		iv[1]=1;
		for(int i=2;i<=lm;++i)iv[i]=(-1ll*(mod/i)*iv[mod%i])%mod+mod;
	}
}
using namespace poly;
int n,tmp,s[100010],V[100010],fac[100010],invf[100010],a0;
vi m,m1,md,t1,t2,val1,val2,f,ff,v,v2,tv,mp,im;
int main(){
	scanf("%d",&n),poly::init((n+10)*3),t1.push_back(0),t1.push_back(1),t2=t1,t2[0]=1,t2.resize(n+5);
	fac[0]=1;
	for(int i=1;i<=n+5;++i)fac[i]=1ll*fac[i-1]*i%mod;
	invf[n+5]=POW(fac[n+5],mod-2);
	for(int i=n+4;i>=0;--i)invf[i]=1ll*invf[i+1]*(i+1)%mod;
	m=Ln(t2);
	for(int i=0;i<m.size()-1;++i)m[i]=m[i+1];
	m.pop_back(),m=Inv(m),md.resize(m.size()-1),m1.resize(m.size()-1);
	for(int i=0;i<m.size()-1;++i)md[i]=1ll*m[i+1]*(i+1)%mod,m1[i]=m[i+1];
	mp=POW(m,n+1,n+1),im=Inv(m1);
	val1=Mul(Mul(md,mp,n+2),Mul(im,im,n+2),n+2);
	val2=Mul(mp,im,n+2),tmp=POW(n+1,mod-2);
	for(int i=0;i<=n+1;++i)val1[i]=1ll*val1[i]*tmp%mod,val2[i]=1ll*val2[i]*tmp%mod;
	for(int i=0;i<=n;++i)s[i]=(val1[i+1]+(mod-1ll)*(n-i+1)%mod*val2[i+1])%mod;
	f.resize(n+2);
	for(int i=0;i<n+2;++i)f[i]=mod-invf[i+2];
	f=Inv(f),v.resize(n+1),v2.resize(n+1);
	for(int i=0;i<=n;++i)(v[i]=f[i+1]-s[i]-(!i))>=mod?v[i]-=mod:(v[i]<0?v[i]+=mod:0),v[i]=1ll*v[i]*fac[i]%mod;
	for(int i=0;i<=n-i;++i)swap(v[i],v[n-i]);
	for(int i=0;i<=n;++i)v2[i]=((i&1?mod-1ll:1ll)*invf[i])%mod;
	tv=Mul(v,v2,n+1);
	for(int i=0;i<=n-i;++i)swap(tv[i],tv[n-i]);
	for(int i=1;i<=n;++i)a0=(a0+1ll*fac[n]*invf[i])%mod;
	printf("%d",a0);
	for(int i=1;i<n;++i)printf(" %lld",1ll*tv[i]*fac[n]%mod*invf[i]%mod);
	return 0;
}