【loj6538】烷基计数 加强版 加强版 Burnside引理+多项式牛顿迭代

别问我为啥突然刷了道OI题，也别问我为啥花括号不换行了...

题目描述

求含 (n) 个碳原子的本质不同的烷基数目模 (998244353) 的结果。(1le nle 10^5) 。

---

题解

Burnside引理+多项式牛顿迭代

不考虑同构的话，很容易想到dp方程 (egin{cases}f_0=1\f_i=sumlimits_{j+k+l+1=i}f_jf_kf_lend{cases}) 。

考虑同构，可以通过容斥原理，大力讨论一下容斥系数。一个更简单的方法是考虑Burnside引理，即：等价类的数目等于每个置换下不动点数目的平均值。

• 对于置换 ((1,2,3)) ，所有组合都是不动点；
• 对于置换 ((1,3,2)) 、((2,1,3)) 和 ((3,2,1)) ，“(2+1)” 的组合是不动点；
• 对于置换 ((2,3,1)) 和 ((3,1,2)) ，只有 “(3)” 的组合是不动点。

于是新的dp方程为 (egin{cases}f_0=1\f_i=frac{sumlimits_{j+k+l+1=i}f_jf_kf_l+sumlimits_{2j+k+1=i}3f_jf_k+sumlimits_{3j+1=i}2f_j}6end{cases}) 。

(n) 这么大肯定不能直接dp，考虑多项式解法，则dp方程的多项式形式为 (F(x)=xcdotfrac{F^3(x)+3F(x)F(x^2)+2F(x^3)}6+1) 。

由于出现了三次方和 (F(x^2)) 、(F(x^3)) 项，因此这个方程难以直接解出。

考虑牛顿迭代，则当我们已知 (F(x)mod x^n) 时，(F(x^2)mod x^{2n}) 和 (F(x^3)mod x^{2n}) 就已经是已知量，在迭代时可以当作常量处理。

记 (S(x)=F(x^2)) ，(C(x)=F(x^3)) ，则我们要迭代的方程就是 (G(T(x))=xcdotfrac{T^3(x)+3S(x)T(x)+2C(x)}6-T(x)+1) 的零点 (G(F(x))=0) 。

又因为 (G'(T(x))=xcdotfrac{3T^2(x)+3S(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^3(x)+3S(x)F_0(x)+2C(x))-6F_0(x)+6}{3F_0^2(x)+3S(x)-6}, ]

其中 (F_0(x)) 是上次迭代所得的多项式。

最后的答案就是 (F(x)[n]) 。

时间复杂度 (O(nlog n)) 。

#include <cstdio>
#include <algorithm>
#define N 262155
#define mod 998244353
using namespace std;
typedef long long ll;
ll F[N];
inline ll qpow(ll x , ll y) {
	ll ans = 1;
	while(y) {
		if(y & 1) ans = ans * x % mod;
		x = x * x % mod , y >>= 1;
	}
	return ans;
}
inline void ntt(ll *A , int n , ll flag) {
	int i , j , k;
	for(k = i = 0 ; i < n ; i ++ ) {
		if(i < k) swap(A[i] , A[k]);
		for(j = (n >> 1) ; (k ^= j) < j ; j >>= 1);
	}
	for(k = 2 ; k <= n ; k <<= 1) {
		ll wn = qpow(3 , (mod - 1) / k * flag);
		for(i = 0 ; i < n ; i += k) {
			ll w = 1 , t;
			for(j = i ; j < i + (k >> 1) ; j ++ , w = w * wn % mod)
				t = w * A[j + (k >> 1)] % mod , A[j + (k >> 1)] = (A[j] - t + mod) % mod , A[j] = (A[j] + t) % mod;
		}
	}
	if(flag == mod - 2) {
		ll t = qpow(n , flag);
		for(i = 0 ; i < n ; i ++ ) A[i] = A[i] * t % mod;
	}
}
inline void inv(ll *A , ll *B , int n) {
	static ll T[N];
	int i , j;
	for(i = 0 ; i < (n << 1) ; i ++ ) B[i] = 0;
	B[0] = qpow(A[0] , mod - 2);
	for(i = 2 ; i <= n ; i <<= 1) {
		for(j = 0 ; j < i ; j ++ ) T[j] = A[j] , T[j + i] = 0;
		ntt(T , i << 1 , 1) , ntt(B , i << 1 , 1);
		for(j = 0 ; j < (i << 1) ; j ++ ) B[j] = B[j] * (2 - T[j] * B[j] % mod + mod) % mod;
		ntt(B , i << 1 , mod - 2);
		for(j = i ; j < (i << 1) ; j ++ ) B[j] = 0;
	}
}
inline void solve(int n) {
	static ll G[N] , H[N] , S[N] , C[N] , T[N];
	int i , j;
	F[0] = 1;
	for(i = 2 ; i <= n ; i <<= 1) {
		for(j = 0 ; j < i ; j ++ ) T[j] = F[j] , S[j] = C[j] = T[j + i] = S[j + i] = C[j + i] = 0;
		for(j = 0 ; j < i ; j += 2) S[j] = F[j / 2];
		for(j = 0 ; j < i ; j += 3) C[j] = F[j / 3];
		ntt(T , i << 1 , 1) , ntt(S , i << 1 , 1);
		for(j = 0 ; j < (i << 1) ; j ++ ) G[j] = T[j] * (T[j] * T[j] % mod + 3 * S[j]) % mod , H[j] = 3 * (T[j] * T[j] + S[j]) % mod;
		ntt(G , i << 1 , mod - 2) , ntt(H , i << 1 , mod - 2);
		for(j = i ; j < (i << 1) ; j ++ ) G[j] = H[j] = 0;
		for(j = i - 1 ; j ; j -- ) G[j] = ((G[j - 1] + 2 * C[j - 1] - 6 * F[j]) % mod + mod) % mod , H[j] = H[j - 1];
		G[0] = 0 , H[0] = mod - 6;
		ntt(G , i << 1 , 1) , inv(H , T , i) , ntt(T , i << 1 , 1);
		for(j = 0 ; j < (i << 1) ; j ++ ) G[j] = G[j] * T[j] % mod;
		ntt(G , i << 1 , mod - 2);
		for(j = 0 ; j < i ; j ++ ) F[j] = (F[j] - G[j] + mod) % mod;
	}
}
int main() {
	int n , len = 1;
	scanf("%d" , &n);
	while(len <= n) len <<= 1;
	solve(len);
	printf("%lld
" , F[n]);
	return 0;
}