Loj 6538 烷基计数 加强版 加强版

Posted on | Views:

Burnside 引理 + 牛顿迭代.

就是求每个节点的儿子个数 $\le 3$ 的有根树数目.

记答案的生成函数为 $F(x)$ ,直接 $F(x)=1+F(x)^3$ 显然是错的,因为可能出现同构的情况.

考虑 Burnside 引理,对 $3$ 个儿子的置换有 $6$ 种.

对于 $(1,2,3)​$ 这 $1​$ 种置换,每种方案都是不动点.

对于 $(1,3,2),(3,2,1),(2,1,3)$ 这 $3$ 种置换,只有在置换中成环的两棵子树相同的方案才是不动点.

对于 $(2,3,1),(3,1,2)$ 这 $2​$ 种置换,只有三棵子树都相同的方案才是不动点.

根据 Burnside 引理,本质不同的方案数是各个置换的不动点数目平均数,得到
$$
F(x)=1+x\frac{F(x)^3+3F(x^2)F(x)+2F(x^3)}{6}
$$
考虑牛顿迭代解出这个方程,设 $G(F(x))=x\frac{F(x)^3+3F(x^2)F(x)+2F(x^3)}{6}-F(x)+1$ .

若当前已知了 $F(x)\bmod x^n$ ,要求出 $F(x)\bmod x^{2n}$ ,可以发现 $F(x^2)\bmod x^{2n},F(x^3)\bmod x^{2n}$ 也已知了.

于是可以直接把它们看做常量,记 $A(x)=F(x^2)\bmod x^{2n},B(x)=F(x^3)\bmod x^{2n}$ .

注意是对 $F(x)​$ 求导,所以 $x​$ 也可以看做常量.

得到 $G’(F(x))=x\frac{3F(x)^2+3A(x)}{6}-1$ ,代入牛顿迭代公式 $F(x)=F_0(x)-\frac{G(F_0(x))}{G’(F_0(x))}$ ,
$$
F(x)=F_0(x)-\frac{x(F_0(x)^3+3A(x)F_0(x)+2B(x))-6F_0(x)+6}{x(3F_0(x)^2+3A(x))-6}
$$
时间复杂度 $O(n\log n)​$ .

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
//%std
#include<bits/stdc++.h>
#define endl '\n'
using namespace std;
typedef long long ll;
inline int read()
{
	int out=0,fh=1;
	char jp=getchar();
	while ((jp>'9'||jp<'0')&&jp!='-')
		jp=getchar();
	if (jp=='-')
		fh=-1,jp=getchar();
	while (jp>='0'&&jp<='9')
		out=out*10+jp-'0',jp=getchar();
	return out*fh;
}
const int P=998244353,G=3;
int add(int a,int b)
{
	return (a+b>=P)?(a+b-P):(a+b);
}
void inc(int &a,int b)
{
	a=add(a,b);
}
int mul(int a,int b)
{
	return 1LL * a * b % P;
}
int fpow(int a,int b)
{
	int res=1;
	while(b)
	{
		if(b&1)
			res=mul(res,a);
		a=mul(a,a);
		b>>=1;
	}
	return res;
}
const int MAXN=4e5+10;
int rev[MAXN],omega[MAXN],inv[MAXN],curn=0;
void init(int n)
{
	if(curn==n)
		return;
	for(int i=0;i<n;++i)
		rev[i]=(rev[i>>1]>>1)|((i&1)*(n>>1));
	for(int l=2;l<=n;l<<=1)
	{
		omega[l]=fpow(G,(P-1)/l);
		inv[l]=fpow(omega[l],P-2);
	}
	curn=n;
}
void DFT(int *a,int n,bool invflag)
{
	init(n);
	for(int i=0;i<n;++i)
		if(i<rev[i])
			swap(a[i],a[rev[i]]);
	for(int l=2;l<=n;l<<=1)
	{
		int m=l>>1;
		int gi=omega[l];
		if(invflag)
			gi=inv[l];
		for(int *p=a;p!=a+n;p+=l)
		{
			int g=1;
			for(int i=0;i<m;++i)
			{
				int t=mul(g,p[i+m]);
				p[i+m]=add(p[i],P-t);
				p[i]=add(p[i],t);
				g=mul(g,gi);
			}
		}
	}
	if(invflag)
	{
		int invn=fpow(n,P-2);
		for(int i=0;i<n;++i)
			a[i]=mul(a[i],invn);
	}
}
void NTT(int *A,int *B,int *C,int lenA,int lenB)
{
	int lenC=lenA+lenB-1,n=1;
	while(n<lenC)
		n<<=1;
	static int a[MAXN],b[MAXN];
	copy(A,A+lenA,a);
	fill(a+lenA,a+n,0);
	copy(B,B+lenB,b);
	fill(b+lenB,b+n,0);
	DFT(a,n,false);
	DFT(b,n,false);
	for(int i=0;i<n;++i)
		C[i]=mul(a[i],b[i]);
	DFT(C,n,true);
}
void PolyInverse(int *A,int *B,int N) // B=A^(-1)
{
	int n=1;
	while(n<N)
		n<<=1;
	static int res[MAXN],tmp[MAXN];
	fill(A+N,A+n,0);
	res[0]=fpow(A[0],P-2);
	for(int i=2;i<=n;i<<=1)
	{
		NTT(A,res,tmp,i,i);
		NTT(tmp,res,tmp,i,i);
		for(int j=0;j<i;++j)
			res[j]=add(mul(2,res[j]),P-tmp[j]);
	}
	copy(res,res+N,B);
}
void solve(int *F,int N)
{
	int n=1;
	while(n<N)
		n<<=1;
	static int A[MAXN],B[MAXN],F0[MAXN],up[MAXN],down[MAXN];
	F[0]=1;
	for(int i=2;i<=n;i<<=1)
	{
		for(int j=0;j<i;++j)
			F0[j]=F[j];
		for(int j=0;j<i;j+=2)
			A[j]=F0[j/2];
		for(int j=0;j<i;j+=3)
			B[j]=F0[j/3];
			
		NTT(F0,F0,up,i,i);
		for(int j=0;j<i;++j)
			inc(up[j],mul(3,A[j]));
		NTT(F0,up,up,i,i);
		for(int j=0;j<i;++j)
			inc(up[j],mul(2,B[j]));
		for(int j=i-1;j>=1;--j)
			up[j]=up[j-1];
		up[0]=0;
		for(int j=0;j<i;++j)
			inc(up[j],P-mul(6,F0[j]));
		inc(up[0],6);
			
		NTT(F0,F0,down,i,i);
		for(int j=0;j<i;++j)
			down[j]=add(mul(3,down[j]),mul(3,A[j]));
		for(int j=i-1;j>=1;--j)
			down[j]=down[j-1];			
		down[0]=P-6;
		PolyInverse(down,down,i);
		
		NTT(up,down,up,i,i);
		for(int j=0;j<i;++j)
			F[j]=add(F0[j],P-up[j]);
	}
}
int F[MAXN];
int main()
{
	int n=read();
	solve(F,n+1);
	cout<<F[n]<<endl;
	return 0;
}