码迷,mamicode.com
首页 > 其他好文 > 详细

@hdu - 6426@ Problem A.Alkane

时间:2020-04-01 23:28:08      阅读:62      评论:0      收藏:0      [点我收藏+]

标签:hdu   cpp   gcc   namespace   ons   大小   rip   ade   区间   


@description@

求包含 n 个碳的烷烃与烷基的同分异构体个数 mod 998244353。

如果你没学过有机化学,你可以认为烷烃是 n 个点且每个点度数 <= 4 的无根树;烷基是 n 个点且每个点儿子个数 <= 3 的有根树。

原题传送门。

@solution@

先考虑有根树的情况,设 \(A(x)\) 为烷基个数对应的生成函数(包括空树)。

显然直接写 \(A(x) = A^3(x) + 1\) 是不对的,因为儿子之间是无序的。

考虑使用 burnside 定理将有序对个数转成无序对个数,其中置换群就是所有大小为 3 的置换。

那么可以得到 \(A(x)\) 的真正转移式:

\[A(x) = \frac{A^3(x) + 3A(x)A(x^2) + 2A(x^3)}{6}x + 1 \]

至于系数为什么是那样的,考虑所有置换,有 1 个拥有 3 个循环的置换;有 3 个拥有 2 个循环的置换;有 2 个拥有 1 个循环的置换。

然后虽然这个式子可以牛顿迭代(把 \(A(x^2), A(x^3)\) 当作常量推导),不过分治 fft 不用动脑子。

虽然可能大家觉得分治 fft 很显然但是我自闭了好久才找到了合理的解释。
具体来说,如果分治区间 [L, R] 中的 L = 0,则直接算 [L, mid] 的 3 次方;否则 [L, mid] 恰好在 3 次方中出现 1 次(分治的过程可以解构成线段树),把贡献乘 3 即可。

可以先提交loj6538看一下自己写的算法正确性


然后是烷烃个数。不难想到以重心为根转成有根树,但是如果按重心的定义做(最大子树 <= n/2),这样只能解决单次询问。

考虑一个经典等式:点数 - 边数 = 1。进一步地有 ∑点数 - ∑边数 = 树个数。

记 p 表示一棵树内点等价类个数(两个点等价当且仅当以两个点为根的有根树同构);同理记 q 表示一棵树内边等价类个数。最后特别地,记 s 表示是否有对称边(即是否存在双重心且双重心等价)。

则:\(p - q + s = 1\)

至于为什么,一样以重心为根。如果 s = 0,则重心是独一无二,其他点等价类对应该点的父边等价类;如果 s = 1,少一个点等价类贡献,填上 s 就行了。

那么自然有 \(\sum p - \sum q + \sum s = ans\)。我们只需要分别求出 \(\sum p, \sum q, \sum s\) 的生成函数 \(P(x), Q(x), S(x)\) 即可。

显然 \(S(x) = A(x^2)\),不多解释。

可以得到 \(P(x)\) 其实就是每个点度数 <= 4 的有根树。它和 \(A(x)\) 就差在树根的儿子个数 <= 4。
于是熟练地运用 burnside 引理可以得到:

\[P(x) = \frac{A^4(x) + 6A^2(x)A(x^2) + 3A^2(x^2) + 8A(x)A(x^3) + 6A(x^4)}{24}x + 1 \]

不过这一次就没有必要进行分治 fft 了。

\(Q(x)\),可以发现它是由两个非空子树拼起来得到。继续简单地运用 burnside 引理得到:

\[Q(x) = \frac{(A(x^2) - 1) + (A(x) - 1)^2}{2}x \]

这样本题就可以得到解决了。

@accepted code@

#pragma GCC optimize (2)
#pragma G++ optimize (2)

#include <cstdio>
#include <algorithm>
using namespace std;

const int MAXN = 100000;
const int MOD = 998244353;
const int INV2 = (MOD + 1) / 2;
const int INV3 = (MOD + 1) / 3;
const int INV4 = 1LL*INV2*INV2%MOD;
const int INV6 = (MOD + 1) / 6;
const int INV8 = 1LL*INV4*INV2%MOD;
const int INV24 = 1LL*INV8*INV3%MOD;

inline int add(int x, int y) {return (x + y >= MOD ? x + y - MOD : x + y);}
inline int sub(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;}

int pow_mod(int b, int p) {
	int ret = 1;
	for(int i=p;i;i>>=1,b=mul(b,b))
		if( i & 1 ) ret = mul(ret, b);
	return ret;
}

int w[20], iw[20], iv[1<<20];
void ntt(int *A, int n, int type) {
	for(int i=0,j=0;i<n;i++) {
		if( i < j ) swap(A[i], A[j]);
		for(int k=(n>>1);(j^=k)<k;k>>=1);
	}
	for(int i=1,s=2,t=1;s<=n;i++,s<<=1,t<<=1) {
		int u = (type == 1 ? w[i] : iw[i]);
		for(int j=0;j<n;j+=s) {
			for(int k=0,p=1;k<t;k++,p=mul(p,u)) {
				int x = A[j + k], y = mul(A[j + k + t], p);
				A[j + k] = add(x, y), A[j + k + t] = sub(x, y);
			}
		}
	}
	if( type == -1 ) {
		for(int i=0;i<n;i++)
			A[i] = mul(A[i], iv[n]);
	}
}
int length(int n) {
	int len; for(len = 1; len < n; len <<= 1);
	return len;
}

int t1[8*MAXN + 5], t2[8*MAXN + 5], t3[8*MAXN + 5], t4[8*MAXN + 5];
int a[MAXN + 5], a2[MAXN + 5], a3[MAXN + 5], a4[MAXN + 5];
void solve(int l, int r) {
	if( l == r ) {
		a[l] = (l == 0 ? 1 : add(a[l], mul(a3[l - 1], INV3)));
		if( 2*l <= MAXN ) a2[2*l] = a[l];
		if( 3*l <= MAXN ) a3[3*l] = a[l];
		if( 4*l <= MAXN ) a4[4*l] = a[l];
		return ;
	}
	int m = (l + r) >> 1;
	solve(l, m);
	
	if( l == 0 ) {
		int n1 = m - l + 1, n2 = r - l + 1, len = length((n2 - 1) + (n1 - 1) + 1);
		for(int i=0;i<len;i++) t2[i] = t3[i] = 0;
		for(int i=l;i<=m;i++) t2[i-l] = a[i];
		for(int i=0;i<n2;i++) t3[i] = a2[i];
		ntt(t2, len, 1), ntt(t3, len, 1);
		for(int i=0;i<len;i++)
			t4[i] = add(mul(mul(t2[i], t2[i]), mul(t2[i], INV6)), mul(mul(t2[i], t3[i]), INV2));
		ntt(t4, len, -1);
		for(int i=m+1;i<=r;i++) a[i] = add(a[i], t4[i - l - 1]);
	}
	else {
		int n1 = m - l + 1, n2 = r - l + 1, len = length(2*(n2 - 1) + 1);
		for(int i=0;i<len;i++) t1[i] = t2[i] = t3[i] = 0;
		for(int i=l;i<=m;i++) t2[i-l] = a[i];
		for(int i=0;i<n2;i++) t1[i] = a[i], t3[i] = a2[i];
		ntt(t1, len, 1), ntt(t2, len, 1), ntt(t3, len, 1);
		for(int i=0;i<len;i++) t4[i] = mul(mul(t2[i], add(mul(t1[i], t1[i]), t3[i])), INV2);
		ntt(t4, len, -1);
		for(int i=m+1;i<=r;i++) a[i] = add(a[i], t4[i - l - 1]);
	}
	
	solve(m + 1, r);
}
int f[MAXN + 5];
void init() {
	for(int i=0;i<20;i++) {
		w[i] = pow_mod(3, (MOD - 1) / (1 << i));
		iw[i] = pow_mod(w[i], MOD - 2);
		iv[1 << i] = pow_mod(1 << i, MOD - 2);
	}
	
	solve(0, MAXN);
	
	int len = length(4*(MAXN - 1) + 1);
	for(int i=0;i<len;i++) t1[i] = t2[i] = t3[i] = t4[i] = 0;
	for(int i=0;i<MAXN;i++) t1[i] = a[i], t2[i] = a2[i], t3[i] = a3[i];
	ntt(t1, len, 1), ntt(t2, len, 1), ntt(t3, len, 1);
	for(int i=0;i<len;i++) {
		int t = mul(t1[i], t1[i]);
		t4[i] = add(mul(mul(t, t), INV24), mul(mul(t, t2[i]), INV4));
		t4[i] = add(t4[i], mul(mul(t2[i], t2[i]), INV8));
		t4[i] = add(t4[i], mul(mul(t1[i], t3[i]), INV3));
	}
	ntt(t4, len, -1);
	for(int i=1;i<=MAXN;i++)
		f[i] = add(mul(a4[i-1], INV4), t4[i-1]);
	
	for(int i=0;i<len;i++) t1[i] = sub(t1[i], 1), t1[i] = mul(t1[i], t1[i]);
	ntt(t1, len, -1);
	for(int i=1;i<=MAXN;i++)
		f[i] = sub(f[i], mul(add(a2[i], t1[i]), INV2));
	
	for(int i=1;i<=MAXN;i++) f[i] = add(f[i], a2[i]);
	f[0] = 1;
}

int read() {
	int x = 0, ch = getchar();
	while( ch < ‘0‘ || ch > ‘9‘ ) ch = getchar();
	while( ‘0‘ <= ch && ch <= ‘9‘ ) x = 10*x + ch - ‘0‘, ch = getchar();
	return x;
}

void write(int x) {
	if( !x ) return ;
	write(x / 10), putchar(x % 10 + ‘0‘);
}

int main() {
	init(); int T = read();
	while( T-- ) {
		int n = read();
		if( f[n] ) write(f[n]); else putchar(‘0‘);
		putchar(‘ ‘);
		if( a[n] ) write(a[n]); else putchar(‘0‘);
		puts("");
	}
}

@details@

为什么 hdu 上只是把编译器从 C++ 换成 G++ 就可以避免超时了啊。为什么啊。

话说我找到最早讨论这个问题的网站竟然是 2011 年的贴子,真是惊了.jpg。

@hdu - 6426@ Problem A.Alkane

标签:hdu   cpp   gcc   namespace   ons   大小   rip   ade   区间   

原文地址:https://www.cnblogs.com/Tiw-Air-OAO/p/12616707.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!