tree
对于 (n) 个点带标号的无根森林，计算所有森林的树的个数的 (k) 次方，对 (998244353) 取模。
自闭，错了一堆关于长度的问题，这里以后一定要注意
比如需要 (n) 次多项式的点值，但是却使用了乘法后的 (2 imes n) 多项式的前 (n) 项点值；对 (2 imes n) 多项式按长度 (2 imes n) DFT 直接什么的
以及还是容易写错 exp 板子的问题...
由 (prufer) 序列，我们知道树的生成函数是 (A(x)=sum_{ige 0}frac{i^{i-2}}{i!}x^i) ，注意这里我们要算的话需要定义 (0) 个点的树个数是 (0)
设森林个数的 (k) 次方的生成函数是 (B_k(x)) ，我们有
[B_k(x)=sum_{ige 0}i^kfrac{A^i(x)}{i!}
]
我们可以知道 (B_0=e^A)
因为 (k) 比较小，并且我们知道 (B_0) ，所以考虑建立 (B_k) 与 (B_{k-1}) 之间的递推关系式
因为关键在于 (i^k) 与 (i^{k-1}) 不一样，所以考虑求导去凑
[egin{aligned}
B_{k-1}'(x)&=sum_{ige 0}frac{i^{k-1}}{i!}(A^i(x))'\
&=sum_{ige 0}frac{i^k}{i!}A^{i-1}(x)A'(x)\
&=frac{A'(x)}{A(x)}sum_{ige 0}frac{i^k}{i!}A^i(x)\
&=frac{A'(x)}{A(x)}B_k(x)
end{aligned}
]
因此
[B_k(x)=frac{A(x)B'_{k-1}(x)}{A'(x)}
]
然后套上多项式板子即可
复杂度 (O(knlog n))
Code:
#include <cstdio>
#include <cctype>
#include <algorithm>
using std::min;
const int SIZE=1<<21;
char ibuf[SIZE],*iS,*iT;
//#define gc() (iS==iT?(iT=(iS=ibuf)+fread(ibuf,1,SIZE,stdin),iS==iT?EOF:*iS++):*iS++)
#define gc() getchar()
template <class T>
void read(T &x)
{
int f=0;x=0;char c=gc();
while(!isdigit(c)) f|=c=='-',c=gc();
while(isdigit(c)) x=x*10+c-'0',c=gc();
if(f) x=-x;
}
const int N=(1<<17)+10;
const int mod=998244353,G=3,Gi=332748118;
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
#define mul(x,y) (1ll*(x)*(y)%mod)
int qp(int d,int k){int f=1;while(k){if(k&1)f=mul(f,d);d=mul(d,d),k>>=1;}return f;}
int turn[N],fac[N],inv[N];
void NTT(int *a,int lim,int typ)
{
int L=-1;for(int i=1;i<lim;i<<=1) ++L;
for(int i=0;i<lim;i++)
{
turn[i]=turn[i>>1]>>1|(i&1)<<L;
if(i<turn[i]) std::swap(a[i],a[turn[i]]);
}
for(int le=1;le<lim;le<<=1)
{
int wn=qp(typ?G:Gi,(mod-1)/(le<<1));
for(int i=0;i<lim;i+=le<<1)
{
int w=1;
for(int j=i;j<i+le;j++,w=mul(w,wn))
{
int x=a[j],y=mul(w,a[j+le]);
a[j]=add(x,y),a[j+le]=add(x,mod-y);
}
}
}
if(!typ)
{
int inv=qp(lim,mod-2);
for(int i=0;i<lim;i++) a[i]=mul(a[i],inv);
}
}
int iva[N],ivb[N];
void polyinv(int *a,int *b,int lim)
{
if(lim==1){b[0]=qp(a[0],mod-2);return;}
polyinv(a,b,lim>>1);
for(int i=0;i<lim<<1;i++) iva[i]=ivb[i]=0;
for(int i=0;i<lim;i++) iva[i]=a[i],ivb[i]=b[i];
NTT(iva,lim<<1,1),NTT(ivb,lim<<1,1);
for(int i=0;i<lim<<1;i++) iva[i]=mul(ivb[i],add(2,mod-mul(iva[i],ivb[i])));
NTT(iva,lim<<1,0);
for(int i=0;i<lim;i++) b[i]=iva[i];
}
void polyqd(int *a,int lim)
{
for(int i=0;i<lim-1;i++) a[i]=mul(a[i+1],i+1);
a[lim-1]=0;
}
void polyint(int *a,int lim)
{
for(int i=lim-1;i;i--)
a[i]=mul(a[i-1],mul(fac[i-1],inv[i]));
a[0]=0;
}
int lna[N],lnb[N];
void polyln(int *a,int lim)
{
for(int i=0;i<lim<<1;i++) lna[i]=lnb[i]=0;
for(int i=0;i<lim;i++) lna[i]=a[i];
polyinv(lna,lnb,lim);
polyqd(lna,lim);
NTT(lna,lim<<1,1),NTT(lnb,lim<<1,1);
for(int i=0;i<lim<<1;i++) lna[i]=mul(lna[i],lnb[i]);
NTT(lna,lim<<1,0);
polyint(lna,lim);
for(int i=0;i<lim;i++) a[i]=lna[i];
}
int exa[N],exb[N];
void polyexp(int *a,int *b,int lim)
{
if(lim==1){b[0]=1;return;}
polyexp(a,b,lim>>1);
for(int i=0;i<lim<<1;i++) exa[i]=exb[i]=0;
for(int i=0;i<lim;i++) exa[i]=exb[i]=b[i];
polyln(exb,lim);
for(int i=0;i<lim;i++) exb[i]=add(a[i]+(i==0),mod-exb[i]);
NTT(exa,lim<<1,1),NTT(exb,lim<<1,1);
for(int i=0;i<lim<<1;i++) exa[i]=mul(exa[i],exb[i]);
NTT(exa,lim<<1,0);
for(int i=0;i<lim;i++) b[i]=exa[i];
}
int A[N],B[21][N],C[N],D[N];
int main()
{
int n,k;
read(n),read(k);
int lim=1;
while(lim<=n) lim<<=1;
fac[0]=1;for(int i=1;i<=lim<<1;i++) fac[i]=mul(fac[i-1],i);
inv[lim<<1]=qp(fac[lim<<1],mod-2);
for(int i=(lim<<1)-1;~i;i--) inv[i]=mul(inv[i+1],i+1);
A[0]=C[0]=0;A[1]=C[1]=1;
for(int i=2;i<lim;i++) A[i]=C[i]=mul(qp(i,i-2),inv[i]);
polyexp(A,B[0],lim);
polyqd(C,lim);
polyinv(C,D,lim);
NTT(A,lim<<1,1),NTT(D,lim<<1,1);
for(int i=0;i<lim<<1;i++) A[i]=mul(A[i],D[i]);
NTT(A,lim<<1,0);//warning
for(int i=lim;i<lim<<1;i++) A[i]=0;
NTT(A,lim<<1,1);
for(int i=1;i<=k;i++)
{
polyqd(B[i-1],lim);
for(int j=lim;j<lim<<1;j++) B[i-1][j]=0;//warning
NTT(B[i-1],lim<<1,1);
for(int j=0;j<lim<<1;j++) B[i][j]=mul(B[i-1][j],A[j]);
NTT(B[i],lim<<1,0);
}
printf("%lld
",mul(B[k][n],fac[n]));
return 0;
}
2019.6.24