【UOJ424】count(笛卡尔树,DP,生成函数,矩阵快速幂)

首先可以发现两个序列 A , B A,B A,B 同构当且仅当它们的笛卡尔树同构。

那么可以考虑枚举笛卡尔树,然后判断它能否构成满足题目条件的序列。

发现一棵笛卡尔树满足条件当且仅当它有 n n n 个节点~~(废话)~~,而且它的最长左链长度不超过 m − 1 m-1 m1

定义一条链的左链长度为这条链上向左的边数,一棵树的最长左链长度为这棵树所有从根到叶子的链的左链长度的最大值。

那么就可以 DP 了:设 f n , m f_{n,m} fn,m 表示 n n n 个点、最长左链长度+1不超过 m m m 的笛卡尔树个数。

这里定义状态的时候把最长左链长度+1是为了方便初始状态定义,因为理论上只有 1 1 1 个点的笛卡尔树的最长左链长度为 0 0 0,只有 0 0 0 个点的笛卡尔树的最长左链长度为 − 1 -1 1

枚举根节点的左子树大小,容易得到转移:
f n , 0 = [ n = 0 ] f n , 1 = 1 f n , m = ∑ i = 0 n − 1 f i , m − 1 f n − 1 − i , m \begin{aligned} f_{n,0}&=[n=0]\\ f_{n,1}&=1\\ f_{n,m}&=\sum_{i=0}^{n-1}f_{i,m-1}f_{n-1-i,m} \end{aligned} fn,0fn,1fn,m=[n=0]=1=i=0n1fi,m1fn1i,m
暴力转移是 O ( n 3 ) O(n^3) O(n3) 的。

考虑优化,设 F m F_m Fm f ∗ , m f_{*,m} f,m 的生成函数,有:
F m = x F m − 1 F m F m = 1 1 − x F m − 1 \begin{aligned} F_m&=xF_{m-1}F_m\\ F_m&=\dfrac{1}{1-xF_{m-1}} \end{aligned} FmFm=xFm1Fm=1xFm11
初始状态:
F 0 = 1 F 1 = 1 1 − x \begin{aligned} F_0&=1\\ F_1&=\dfrac{1}{1-x} \end{aligned} F0F1=1=1x1
于是得到了 O ( n m log ⁡ n ) O(nm\log n) O(nmlogn) 的做法。

观察递推式,发现可以写成矩阵快速幂的形式:设 F m − 1 = A B F_{m-1}=\dfrac{A}{B} Fm1=BA,那么 F m = B B − x A F_{m}=\dfrac{B}{B-xA} Fm=BxAB

于是有了 O ( n log ⁡ n log ⁡ m ) O(n\log n\log m) O(nlognlogm) 的做法。

进一步地,设 F m = A m B m F_m=\dfrac{A_m}{B_m} Fm=BmAm,可以归纳证明 A m A_m Am m − 1 m-1 m1 次多项式, B m B_m Bm m m m 次多项式。

这意味着我们可以先 NTT 一遍,把 O ( n ) O(n) O(n) 个点值带进去求矩阵快速幂,再把求出来的点值 NTT 回来,而这过程中不会发生循环卷积的 BUG。

时间复杂度 O ( n log ⁡ n + n log ⁡ m ) O(n\log n+n\log m) O(nlogn+nlogm)

#include<bits/stdc++.h>

#define LN 20
#define N 100010

using namespace std;

namespace modular
{
	const int mod=998244353;
	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;}
}using namespace modular;

inline int poww(int a,int b)
{
	int ans=1;
	while(b)
	{
		if(b&1) ans=mul(ans,a);
		a=mul(a,a);
		b>>=1;
	}
	return ans;
}

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;
}

int n,m;
int A[N<<2],B[N<<2],C[N<<2],D[N<<2],invD[N<<2];
int rev[N<<2],w[LN][N<<2][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++,g=mul(g,gn),ig=mul(ig,ign))
			w[bit][j][0]=g,w[bit][j][1]=ig;
	}
}

void NTT(int *a,int limit,int opt)
{
	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][j][opt],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++)
			a[i]=mul(a[i],tmp);
	}
}

void getinv(int *f,int *g,int n)
{
    static int ff[N<<2];
    assert(f[0]);
    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(g[i],dec(2,mul(ff[i],g[i])));
        NTT(g,limit,-1);
        for(int i=now;i<limit;i++) g[i]=0;
    }
    for(int i=n;i<now;i++) g[i]=0;
    for(int i=0;i<now;i++) ff[i]=0;
}

struct Matrix
{
	int a[2][2];
	Matrix(){a[0][0]=a[0][1]=a[1][0]=a[1][1]=0;};
	void init(){a[0][0]=a[1][1]=1;}
}trans[N<<2],st[N<<2];

Matrix mul(Matrix a,Matrix b)
{
	Matrix c;
	for(int i=0;i<2;i++)
		for(int j=0;j<2;j++)
			for(int k=0;k<2;k++)
				Add(c.a[i][j],mul(a.a[i][k],b.a[k][j]));
	return c;
}

Matrix poww(Matrix a,int b)
{
	Matrix ans;
	ans.init();
	while(b)
	{
		if(b&1) ans=mul(ans,a);
		a=mul(a,a);
		b>>=1;
	}
	return ans;
}

int main()
{
	n=read(),m=read();
	if(n<m)
	{
		puts("0");
		return 0;
	}
	int limit=1;
	while(limit<=(n<<1)) limit<<=1;
	init(limit);
	A[1]=1;
	NTT(A,limit,1);
	B[0]=1;
	NTT(B,limit,1);
	for(int i=0;i<limit;i++)		
	{
		trans[i].a[0][0]=0,trans[i].a[0][1]=1,trans[i].a[1][0]=dec(0,A[i]),trans[i].a[1][1]=1;
		st[i].a[0][0]=st[i].a[1][0]=B[i];
		Matrix res=mul(poww(trans[i],m),st[i]);
		C[i]=res.a[0][0],D[i]=res.a[1][0];
	}
	NTT(D,limit,-1);
	getinv(D,invD,n+1);
	NTT(invD,limit,1);
	for(int i=0;i<limit;i++) C[i]=mul(C[i],invD[i]);
	NTT(C,limit,-1);
	printf("%d\n",C[n]);
	return 0;
}
  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值