生成函数Euler变换学习笔记（无标号有根树计数）

众所周知，对于有标号计数的指数型生成函数 $f(x)$，将其任意地进行无顺序的组合，得到的生成函数是$exp(f(x))$。

而对于无标号计数的这样的组合，我们就需要引入一个叫 $\text{Eular}$ 变换的东西 $\mathcal{E}(f(x))$

（这个符号打法：\mathcal(E)）

为了方便，本文中记$f_i$表示$[x^i]f(x)$。两种表述可能混用。

我们先引入一个实际问题：$n$个点无标号有根树的计数。两棵树被认为相同，当且仅当存在一种重标号方案，使得它们的根标号相同且边集相等。

对于一棵树，我们去掉它的根，然后就分成若干棵独立的子树。所以对无标号有根树的$\text{OGF}$ $F(x)$，有

$$F(x)=x\cdot \mathcal{E}(F(x))$$

先思考这个问题的主要矛盾在什么地方：对于分出来的若干子树，如果有一些子树完全一致，那么它们只会贡献一次方案。但如果只是大小相同而本质不同，仍然算作不同的方案。

换句话说，最终的方案只和每种本质不同的子树使用的个数有关。考虑背包。

对每种本质不同的大小为$k$的树构造完全背包的生成函数，即

$$\frac{1}{1-x^k}$$

那么我们可以得到一个简单粗暴的定义式：

$$\mathcal{E}(F(x))=\prod _{i=1}^{\infty }\frac{1}{(1-x^i{)}^{F_i}}$$

$F_i$是因为有$F_i$种本质不同大小为$i$的物品

然后开始愉快地推式子

$$G(x)=\mathcal{E}(F(x))=\prod _{i=1}^{\infty }\frac{1}{(1-x^i{)}^{F_i}}$$

既然有连乘，先求个$\ln$

$$\ln (G(x))=-\sum _{i=1}^{\infty }{F}_i\ln (1-x^i)$$

$\ln$还是看不顺眼，求个导

$$\frac{G^{\prime }(x)}{G(x)}=\sum _{i=1}^{\infty }{F}_i\frac{i{x}^{i-1}}{1-x^i}$$

又看见了喜闻乐见的完全背包，直接展开

$$\frac{G^{\prime }(x)}{G(x)}=\sum _{i=1}^{\infty }{F}_i\cdot i{x}^{i-1}\sum _{j=0}^{\infty }{x}^{ij}$$

$$\frac{G^{\prime }(x)}{G(x)}=\sum _{i=1}^{\infty }{F}_i\cdot i\sum _{j=0}^{\infty }{x}^{i(j+1)-1}$$

$$\frac{G^{\prime }(x)}{G(x)}=\sum _{i=1}^{\infty }{F}_i\cdot i\sum _{j=1}^{\infty }{x}^{ij-1}$$

再积分回来

$$\ln (G(x))=\sum _{i=1}^{\infty }{F}_i\cdot i\sum _{j=1}^{\infty }\frac{x^{ij}}{ij}$$

整理一下

$$\ln (G(x))=\sum _{i=1}^{\infty }{F}_i\sum _{j=1}^{\infty }\frac{x^{ij}}{j}$$

$$\ln (G(x))=\sum _{j=1}^{\infty }\frac{1}{j}\sum _{i=1}^{\infty }{F}_i{x}^{ij}$$

$$\ln (G(x))=\sum _{j=1}^{\infty }\frac{F(x^j)}{j}$$

$$\ln (G(x))=\sum _{i=1}^{\infty }\frac{F(x^i)}{i}$$

得到 $\text{Eular}$ 变换的计算式（雾）

$$\mathcal{E}(F(x))=\exp (\sum _{i=1}^{\infty }\frac{F(x^i)}{i})$$

然后回到原问题

$$F(x)=x\cdot \mathcal{E}(F(x))$$

$$F(x)=x\cdot \exp (\sum _{i=1}^{\infty }\frac{F(x^i)}{i})$$

两边求导(如果高数不好建议自己多写几步)

$$F^{\prime }(x)=\exp (\sum _{i=1}^{\infty }\frac{F(x^i)}{i})+x\cdot \exp (\sum _{i=1}^{\infty }\frac{F(x^i)}{i})\cdot \sum _{i=1}^{\infty }\frac{i{x}^{i-1}\cdot F^{\prime }(x^i)}{i}$$

$$F^{\prime }(x)=\exp (\sum _{i=1}^{\infty }\frac{F(x^i)}{i})+F(x)\sum _{i=1}^{\infty }{x}^{i-1}\cdot F^{\prime }(x^i)$$

两边同时乘上$x$

$$x\cdot F^{\prime }(x)=F(x)+F(x)\sum _{i=1}^{\infty }{x}^i\cdot F^{\prime }(x^i)$$

#undef G

记

$$G(x)=\sum _{i=1}^{\infty }{x}^i\cdot F^{\prime }(x^i)$$

$$x\cdot F^{\prime }(x)=F(x)+F(x)G(x)$$

考虑求这东西的第$x^n$项

$$n{F}_n=F_n+\sum _{k=1}^{n-1}{F}_k{G}_{n-k}$$

$$F_n=\frac{1}{n-1}\sum _{k=1}^{n-1}{F}_k{G}_{n-k}$$

如果我们知道$G$，就可以直接分治 NTT 算出$F$

现在考虑这个$G$

$$G(x)=\sum _{i=1}^{\infty }{x}^i\cdot F^{\prime }(x^i)$$

$$G_n=\sum _{i=1}^{\infty }[x^{n-i}]F^{\prime }(x^i)$$

经过大胆的观察和仔细的猜想

$$G_n=\sum _{d\mid n}(\frac{n-d}{d}+1)F_{\frac{n-d}{d}+1}$$

$$G_n=\sum _{d\mid n}(\frac{n}{d})F_{\frac{n}{d}}$$

$$G_n=\sum _{d\mid n}d{F}_d$$

然后每求出一个$F$，就可以调和级数复杂度求出之前的$G$，就可以分治 NTT 了。

但是分治时如果$L=1$，$G$可能不够用。解决办法是先不计算不够的地方的贡献，在后面再把它加上。

具体的讲，如果我们要计算$F_{1\sim mid}$对$F_{mid+1\sim R}$的贡献，我们只计算$i$不超过$mid$的$G_i$的贡献，这样少算了$G_{mid+1,R-1}$的贡献。而这部分贡献可以通过$R$之前交换$F$和$G$做一次类似的卷积计算出来。详见代码。

复杂度$O(n{\log }^{2}n)$

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cctype>
#define MAXN 524288
using namespace std;
const int MOD=998244353;
typedef long long ll;
inline int add(const int& x,const int& y){return x+y>=MOD? x+y-MOD:x+y;}
inline int dec(const int& x,const int& y){return x<y? x-y+MOD:x-y;}
inline int qpow(int a,int p)
{
	int ans=1;
	while (p)
	{
		if (p&1) ans=(ll)ans*a%MOD;
		a=(ll)a*a%MOD;p>>=1;
	}
	return ans;
}
#define inv(x) qpow(x,MOD-2)
int rt[2][24],r[MAXN],l,lim;
inline void init(){lim=1<<l;for (int i=0;i<lim;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));}
inline void NTT(int* a,int type)
{
	for (int i=0;i<lim;i++) if (i<r[i]) swap(a[i],a[r[i]]);
	for (int L=0;L<l;L++)
	{
		int mid=1<<L,len=mid<<1;
		ll Wn=rt[type][L+1];
		for (int s=0;s<lim;s+=len)
			for (int k=0,w=1;k<mid;k++,w=w*Wn%MOD)
			{
				int x=a[s+k],y=(ll)w*a[s+mid+k]%MOD;
				a[s+k]=add(x,y),a[s+mid+k]=dec(x,y);
			}
	}
	if (type)
	{
		int t=inv(lim);
		for (int i=0;i<lim;i++) a[i]=(ll)a[i]*t%MOD;
	}
}
int f[MAXN],g[MAXN],n;
int a[MAXN],b[MAXN];
void solve(int L,int R)
{
	if (L==R)
	{
		if (L>1) f[L]=(ll)f[L]*inv(L-1)%MOD;
		int t=(ll)L*f[L]%MOD;
		for (int i=L;i<=n;i+=L) g[i]=add(g[i],t);
		return;
	}
	int mid=(L+R)>>1;
	solve(L,mid);
	int len=R-L;
	for (l=0;(1<<l)<=(len<<1);l++);
	init();
	if (L==1)
	{
		for (int i=0;i<lim;i++) a[i]=b[i]=0;
		for (int i=L;i<=mid;i++) a[i]=f[i],b[i]=g[i];
		NTT(a,0);NTT(b,0);
		for (int i=0;i<lim;i++) a[i]=(ll)a[i]*b[i]%MOD;
		NTT(a,1);
		for (int i=mid+1;i<=R;i++) f[i]=add(f[i],a[i]);
	}
	else
	{
		for (int i=0;i<lim;i++) a[i]=b[i]=0;
		for (int i=L;i<=mid;i++) a[i-L]=f[i];
		for (int i=0;i<=len;i++) b[i]=g[i];
		NTT(a,0);NTT(b,0);
		for (int i=0;i<lim;i++) a[i]=(ll)a[i]*b[i]%MOD;
		NTT(a,1);
		for (int i=mid+1;i<=R;i++) f[i]=add(f[i],a[i-L]);
		for (int i=0;i<lim;i++) a[i]=b[i]=0;
		for (int i=L;i<=mid;i++) a[i-L]=g[i];
		for (int i=0;i<=len;i++) b[i]=f[i];
		NTT(a,0);NTT(b,0);
		for (int i=0;i<lim;i++) a[i]=(ll)a[i]*b[i]%MOD;
		NTT(a,1);
		for (int i=mid+1;i<=R;i++) f[i]=add(f[i],a[i-L]);
	}
	
	solve(mid+1,R);
}
int main()
{
	rt[0][23]=qpow(3,119),rt[1][23]=inv(rt[0][23]);
	for (int i=22;i>=0;i--)
	{
		rt[0][i]=(ll)rt[0][i+1]*rt[0][i+1]%MOD;
		rt[1][i]=(ll)rt[1][i+1]*rt[1][i+1]%MOD;
	}
	scanf("%d",&n);
	f[1]=1;
	solve(1,n);
	int ans=f[n];
	printf("%d\n",ans);
	return 0;
}