以下 “二叉树” 均默认为有根无标号但区分左右儿子的二叉树。
设 \(h_{n,k}\) 表示 \(n,k\) 的答案,有:
\[h_{n,k}=\sum_{i=0}^{n-1}\left(h_{i,k}\cdot f_{n-i-1}+f_{i}\cdot h_{n-i-1,k}+\sum_{j=0}^{k}\binom{k}{j}g_{n-i-1,k-j}\sum_{j'=0}^j\binom{j}{j'}g_{i,j'}\right) \]其中 \(f_n\) 表示 \(n\) 个点的二叉树数量。\(g_{n,k}\) 表示 \(n\) 个点的所有二叉树中,所有叶子节点到根的路径所经点数的 \(k\) 次方之和。
然后后面那一堆组合数是 \((a+b+1)^k\) 展开成 \(\sum_{j=0}^k\binom{k}{j}b^{k-j}\sum_{j'=0}^j\binom{j}{j'}a^{j'}\) 形式的结果。
考虑先把 \(f\) 和 \(g\) 推出来:
\[f_0=1,\ f_n=\sum_{i=0}^{n-1}f_if_{n-1-i} \]设 \(F\) 为 \(f\) 的 OGF,有:
\[\begin{aligned} F&=1+xF^2\\ xF^2-F+1&=0\\ &=\frac{1\pm\sqrt{1-4x}}{2x}\\ &=\frac{1-\sqrt{1-4x}}{2x} \end{aligned} \]然后对于 \(g\):
\[g_{1,k}=1,\ g_{n,k}=2\sum_{i=0}^{n-1}f_{n-1-i}\sum_{j=0}^k\binom{k}{j}g_{i,j} \]设 \(G\) 为 \(g\) 的 EGF(\([x^n][t^k]G=\frac{g_{n,k}}{k!}\)),有:
\[\begin{aligned} G&=xe^t+2xFGe^t\\ (1-2xFe^t)G&=xe^t\\ G&=\frac{xe^t}{1-2xFe^t}\\ &=\frac{xe^t}{1-(1-\sqrt{1-4x})e^t} \end{aligned} \]回到最上面那条式子。设 \(H\) 为 \(h\) 的 EGF,有:
\[\begin{aligned} H&=2xHF+xG^2e^t+xe^t\\ (1-2xF)H&=xG^2e^t+xe^t\\ H&=\frac{xG^2e^t+xe^t}{1-2xF}\\ H&=\frac{xe^t}{\sqrt{1-4x}}+\frac{x^3e^{3t}}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})e^t)^2}\\ \end{aligned} \]我们要算的即为:
\[R(t)=[x^n]\frac{xe^t}{\sqrt{1-4x}}+\frac{x^3e^{3t}}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})e^t)^2} \]模 \(t^{m+1}\) 的结果。
比较神奇的做法是进行变量代换,我们设 \(y=e^t-1\),并把 \(R(t)\) 表示成关于 \(y\) 的多项式 \(R'(y)\),然后求出 \(R'(y)\) 再把 \(y=e^t-1\) 代入求出 \(R(t)\)。
注意这里为什么取 \(y=e^t-1\) 而不是 \(y=e^t\),因为我们在将 \(y=e^t\) 代入 \(R'(y)\) 时,\(y\) 的很高次项也有可能影响到 \(R(t)\) 的低次项。而若我们用 \(y=e^t-1\) 代入,\(e^t-1\) 常数项为 \(0\),代入后 \(R'(y)\) 中 \(y\) 的大于 \(m\) 次项都不会影响到 \(R(t)\) 中 \(t\) 的小于等于 \(m\) 次项,这样我们就可以只需要用 \(R'(y)\) 的前 \(m\) 次项的系数即可。
由于 \(e^t=y+1\),可得:
\[\begin{aligned} R'(y)=R(t)&=[x^n]\frac{xe^t}{\sqrt{1-4x}}+\frac{x^3e^{3t}}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})e^t)^2}\\ &=[x^n]\frac{x(y+1)}{\sqrt{1-4x}}+\frac{x^3(y+1)^3}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})(y+1))^2}\\ \end{aligned} \]设 \(C(n,m)=[x^n](\sqrt{1-4x})^m=\dbinom{\frac{m}{2}}{n}(-4)^n\),代进去可得:
\[\begin{aligned} R'(y)&=[x^n]\frac{x(y+1)}{\sqrt{1-4x}}+\frac{x^3(y+1)^3}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})(y+1))^2}\\ &=(y+1)C(n-1,-1)+(y+1)^3[x^{n-3}]\frac{1}{\sqrt{1-4x}(1-(1-\sqrt{1-4x})(y+1))^2}\\ \end{aligned} \]设 \(z=\sqrt{1-4x}\),然后:
\[\begin{aligned} R'(y)&=(y+1)C(n-1,-1)+(y+1)^3[x^{n-3}]\frac{1}{z(1-(1-z)(y+1))^2}\\ &=(y+1)C(n-1,-1)+(y+1)^3[x^{n-3}]\frac{1}{z(zy+z-y)^2}\\ &=(y+1)C(n-1,-1)+(y+1)^3[x^{n-3}]z^{-1}\left(\sum_{i\geq 0}\frac{(1-z)^i}{z^{i+1}}y^i\right)^2\\ &=(y+1)C(n-1,-1)+(y+1)^3[x^{n-3}]\sum_{i\geq 0}(i+1)\frac{(1-z)^i}{z^{i+3}}y^i\\ &=(y+1)C(n-1,-1)+(y+1)^3\sum_{i\geq 0}(i+1)y^i\sum_{j=0}^{i}\binom{i}{j}[x^{n-3}]\frac{(-z)^j}{z^{i+3}}\\ &=(y+1)C(n-1,-1)+(y+1)^3\sum_{i\geq 0}(i+1)y^i\sum_{j=0}^{i}\binom{i}{j}(-1)^jC(n-3,j-i-3)\\ \end{aligned} \]后面那玩意就能卷积求了。现在我们求出了 \(R'(y)=R(t)\)。接下来我们要求多项式复合 \(R(t)=R'(e^t-1)\)。
我们可以先求出 \(A(y)=R'(y-1)\),一次减法卷积即可求出。
然后再把 \(e^t\) 代入,即要求 \(A(e^t)=\sum_i a_ie^{it}\)。分治 NTT 即可。
总时间复杂度 \(O(n+m\log ^2m)\)。
#include<bits/stdc++.h>
#define N 100010
#define V 20000010
using namespace std;
namespace modular
{
const int mod=998244353,inv2=(mod+1)>>1;
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
inline int dec(int x,int y){return x-y<0?x-y+mod:x-y;}
inline int mul(int x,int y){return 1ll*x*y%mod;}
inline void Add(int &x,int y){x=x+y>=mod?x+y-mod:x+y;}
inline void Dec(int &x,int y){x=x-y<0?x-y+mod:x-y;}
inline void Mul(int &x,int y){x=1ll*x*y%mod;}
inline int poww(int a,int b){int ans=1;for(;b;Mul(a,a),b>>=1)if(b&1)Mul(ans,a);return ans;}
}using namespace modular;
inline int read()
{
int x=0,f=1;
char ch=getchar();
while(ch<'0'||ch>'9')
{
if(ch=='-') f=-1;
ch=getchar();
}
while(ch>='0'&&ch<='9')
{
x=(x<<1)+(x<<3)+(ch^'0');
ch=getchar();
}
return x*f;
}
namespace Poly
{
using modular::add;using modular::mul;using modular::dec;
const int LN=20,NN=N<<2;
typedef vector<int> poly;
vector<int> w[LN][2];
void init(int limit)
{
for(int bit=0,mid=1;mid<limit;bit++,mid<<=1)
{
int len=mid<<1;
int gn=poww(3,(mod-1)/len);
int ign=poww(gn,mod-2);
int g=1,ig=1;
for(int j=0;j<mid;j++,Mul(g,gn),Mul(ig,ign))
w[bit][0].push_back(g),w[bit][1].push_back(ig);
}
}
void NTT(int *a,int limit,int opt)
{
static int rev[NN];
opt=(opt<0);
for(int i=0;i<limit;i++)
rev[i]=(rev[i>>1]>>1)|((i&1)*(limit>>1));
for(int i=0;i<limit;i++)
if(i<rev[i]) swap(a[i],a[rev[i]]);
for(int bit=0,mid=1;mid<limit;bit++,mid<<=1)
{
for(int i=0,len=mid<<1;i<limit;i+=len)
{
for(int j=0;j<mid;j++)
{
int x=a[i+j],y=mul(w[bit][opt][j],a[i+mid+j]);
a[i+j]=add(x,y),a[i+mid+j]=dec(x,y);
}
}
}
if(opt)
{
int tmp=poww(limit,mod-2);
for(int i=0;i<limit;i++) Mul(a[i],tmp);
}
}
poly bmul(const poly &a,const poly &b)
{
const int sa=a.size(),sb=b.size();
poly c(sa+sb-1);
for(int i=0;i<sa;i++)
for(int j=0;j<sb;j++)
Add(c[i+j],mul(a[i],b[j]));
return c;
}
poly mul(const poly &a,const poly &b)
{
static int A[NN],B[NN];
const int sa=a.size(),sb=b.size();
for(int i=0;i<sa;i++) A[i]=a[i];
for(int i=0;i<sb;i++) B[i]=b[i];
int limit=1;
while(limit<(sa+sb-1)) limit<<=1;
NTT(A,limit,1),NTT(B,limit,1);
for(int i=0;i<limit;i++) Mul(A[i],B[i]);
NTT(A,limit,-1);
poly c(sa+sb-1);
for(int i=0;i<sa+sb-1;i++) c[i]=A[i];
for(int i=0;i<limit;i++) A[i]=B[i]=0;
return c;
}
poly dmul(const poly &a,poly b)//a[i]b[j]->c[i-j]
{
const int sa=a.size(),sb=b.size();
reverse(b.begin(),b.end());
poly c=mul(a,b),d(sa);
for(int i=0;i<sa;i++) d[i]=c[sb+i-1];
return d;
}
poly add(const poly &a,const poly &b)
{
const int sa=a.size(),sb=b.size();
poly c(max(sa,sb));
for(int i=0;i<sa;i++) c[i]=a[i];
for(int i=0;i<sb;i++) Add(c[i],b[i]);
return c;
}
poly dec(const poly &a,const poly &b)
{
const int sa=a.size(),sb=b.size();
poly c(max(sa,sb));
for(int i=0;i<sa;i++) c[i]=a[i];
for(int i=0;i<sb;i++) Dec(c[i],b[i]);
return c;
}
poly getinv(const poly &f,int n)
{
static int ff[NN],g[NN];
g[0]=poww(f[0],mod-2);
int now=2;
for(;now<(n<<1);now<<=1)
{
int limit=now<<1;
for(int i=0;i<now;i++) ff[i]=f[i];
NTT(ff,limit,1),NTT(g,limit,1);
for(int i=0;i<limit;i++)
g[i]=mul(dec(2,mul(ff[i],g[i])),g[i]);
NTT(g,limit,-1);
for(int i=now;i<limit;i++) g[i]=0;
}
poly res(g,g+n);
for(int i=0;i<now;i++) ff[i]=g[i]=0;
return res;
}
}using Poly::poly;
int n,m;
int fac[V],ifac[V];
int C1[N];//C1[m] = C[n-3,-m] = Binomial[(-m)/2,n-3]*(-4)^(n-3)
poly C;
int Binom(int n,int m)
{
if(n<m||n<0||m<0) return 0;
return mul(mul(fac[n],ifac[m]),ifac[n-m]);
}
void comp1()
{
poly a(m+1),b(m+1);
for(int i=0;i<=m;i++) a[i]=mul(C[i],fac[i]);
for(int i=0;i<=m;i++) b[i]=((i&1)?dec(0,ifac[i]):ifac[i]);
C=Poly::dmul(a,b);
for(int i=0;i<=m;i++) Mul(C[i],ifac[i]);
}
poly fz[N<<2],fm[N<<2];
void solve(int k,int l,int r)
{
if(l==r)
{
fz[k]=poly{C[l]};
fm[k]=poly{1,dec(0,l)};
return;
}
int mid=(l+r)>>1,lc=k<<1,rc=lc|1;
solve(lc,l,mid),solve(rc,mid+1,r);
fz[k]=Poly::add(Poly::mul(fz[k<<1],fm[k<<1|1]),Poly::mul(fm[k<<1],fz[k<<1|1]));
fm[k]=Poly::mul(fm[k<<1],fm[k<<1|1]);
}
void comp2()
{
solve(1,0,m);
C=Poly::mul(fz[1],Poly::getinv(fm[1],m+1));
for(int i=0;i<=m;i++) Mul(C[i],ifac[i]);
}
int main()
{
n=read(),m=read();
int limit=1;
while(limit<=(m<<1|1)) limit<<=1;
Poly::init(limit);
const int t=max(n<<1,m); fac[0]=1;
for(int i=1;i<=t;i++) fac[i]=mul(fac[i-1],i);
ifac[t]=poww(fac[t],mod-2);
for(int i=t;i>=1;i--) ifac[i-1]=mul(ifac[i],i);
C1[3]=mul(ifac[n-3],poww(mod-4,n-3));
for(int i=1,v=dec(inv2,2);i<=n-3;i++,Dec(v,1)) Mul(C1[3],v);
for(int i=5,vl=dec(inv2,2),vr=dec(inv2,n-1);i<=m+3;i+=2,Dec(vl,1),Dec(vr,1)) C1[i]=mul(C1[i-2],mul(vr,poww(vl,mod-2)));
C1[4]=mul(ifac[n-3],poww(mod-4,n-3));
for(int i=1,v=dec(0,2);i<=n-3;i++,Dec(v,1)) Mul(C1[4],v);
for(int i=6,vl=dec(0,2),vr=dec(0,n-1);i<=m+3;i+=2,Dec(vl,1),Dec(vr,1)) C1[i]=mul(C1[i-2],mul(vr,poww(vl,mod-2)));
poly A(m+1),B(m+1);
for(int i=0;i<=m;i++) A[i]=((i&1)?dec(0,ifac[i]):ifac[i]);
for(int i=0;i<=m;i++) B[i]=mul(ifac[i],C1[-(-i-3)]);
C=Poly::mul(A,B); C.resize(m+1);
for(int i=0;i<=m;i++) Mul(C[i],fac[i+1]);
C=Poly::bmul(C,poly{1,3,3,1}); C.resize(m+1);
Add(C[0],Binom(2*n-2,n-1)); Add(C[1],Binom(2*n-2,n-1));
comp1();
comp2();
for(int i=0;i<=m;i++) printf("%d ",mul(C[i],fac[i]));
return 0;
}
/*
3 5
*/
标签:frac,4x,int,永无,sqrt,GF,XSY4375,aligned,sum
From: https://www.cnblogs.com/ez-lcw/p/16842986.html