【UR #3】链式反应

Description

问一棵有 1 ~ n 个不同的点的树有多少种形态，点分黑白点，黑点一定为叶子节点，白点要么是叶子，要么下面可以接若干个黑点以及两个白点。

---

Solution

列列dp式：

[egin{aligned} f_n &=[n=1]+ frac{1}{2} sum_{t in G} sum_{i=1}^n C_{n-1}^{t} C_{n-t-1}^{i} f_{i} f_{n-t-i-1} \ &=[n=1]+frac 12(n-1)! sum_{t in G} frac{g_{n-t-1}}{t!}\ g_n&= sum_{i=1}^n frac{f_if_{n-i}}{i!(n-i)!} end{aligned} ]

那么就可以分治FFT解决啦！
我写的分治FFT是(O(frac {nlog^2n}{loglog n}))的，并且没有加上常数优化，但还是跑了uoj的rk1。

---

#include <bits/stdc++.h>

#define I inline
#define fi first
#define se second
#define LL long long
#define mp make_pair
#define reg register int
#define pii pair<int,int>
#define fo(i, a, b) for(int i = a; i <= b; i++)
#define fd(i, a, b) for(reg i = a; i >= b; i--)
#define ULL unsigned long long
#define cr const reg&
using namespace std;
const int inf = 2147483647;
const int mod = 998244353;
const int N = 1 << 18;
const int B = 16;

I int _max(cr x, cr y) {return x > y ? x : y;}
I int _min(cr x, cr y) {return x < y ? x : y;}
I LL read() {
	LL 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 << 3) + (x << 1) + (ch ^ 48), ch = getchar();
	return x * f;
}
I void ptt(LL x) {if(x >= 10) ptt(x / 10); putchar(x % 10 + '0');}
I void put(LL x) {x < 0 ? putchar('-'), ptt(-x) : ptt(x);}
I void pr1(LL x) {put(x), putchar(' ');}
I void pr2(LL x) {put(x), puts("");}

I int pow_mod(reg a, reg k) {reg ans = 1; for(; k; k >>= 1, a = (LL)a * a % mod) if(k & 1) ans = (LL)ans * a % mod; return ans;}

ULL p[N << 1], g[N], f[N];
int A[N << 1], R[N << 1], w[N << 1], jc[N], inv[N], iv[N];
int MEM1[N * 10], MEM2[N * 10], MEM3[N * 10], MEM4[N * 10];

I int Pre(cr n) {
	reg len = 1; for(; len <= n; len <<= 1);
	jc[0] = 1; fo(i, 1, n) jc[i] = (LL)jc[i - 1] * i % mod;
	inv[n] = pow_mod(jc[n], mod - 2); fd(i, n, 1) inv[i - 1] = (LL)inv[i] * i % mod;
	iv[0] = iv[1] = 1; fo(i, 2, n) iv[i] = (LL)iv[mod % i] * (mod - mod / i) % mod;
	for(reg i = 1; i < len; i <<= 1) {
		reg wn = pow_mod(3, (mod - 1) / (i << 1)), s = 1;
		fo(k, 0, i - 1) w[i + k] = s, s = (LL)s * wn % mod;
	} return len;
}
I int pre(cr n) {
	reg len = 1; for(; len <= n; len <<= 1);
	fo(i, 0, len - 1) R[i] = (R[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	return len;
}

I void DFT(int y[], cr len) {
	fo(i, 0, len - 1) p[R[i]] = y[i]; reg b;
	for(reg i = 1; i < len; i <<= 1) for(reg j = 0; j < len; j += i << 1)
		fo(k, 0, i - 1) b = p[i + j + k] * w[i + k] % mod, p[i + j + k] = p[j + k] - b + mod, p[j + k] = p[j + k] + b;
	fo(i, 0, len - 1) y[i] = p[i] % mod;
}
I void IDFT(int y[], cr len) {
	reverse(y + 1, y + len); DFT(y, len);
	reg hh = pow_mod(len, mod - 2); fo(i, 0, len - 1) y[i] = (LL)y[i] * hh % mod;
}

I void solve(cr l, cr r, int *MEMP1, int *MEMP2, int *MEMP3, int *MEMP4) {
	if(r - l + 1 <= 64) {
		fo(i, l, r) {
			if(i == 0) f[0] = g[0] = 0;
			else if(i == 1) f[i] = 1;
			else f[i] = (LL)f[i] * iv[i] % mod * iv[2] % mod;
			if(l == 0) fo(j, l, i - 1) g[i] = (g[i] + (LL)f[j] * f[i - j]) % mod;
			else fo(j, l, i - 1) g[i] = (g[i] + 2LL * f[j] * f[i - j]) % mod;
			fo(j, i + 1, r) f[j] = (f[j] + g[i] * A[j - i]) % mod;
		} return ;
	} int len = (r - l + 1) / B, ll = len << 1; pre(ll - 1);
	int *h1[B], *h2[B], *h3[B], *h4[B];
	fo(i, 0, B - 1) {
		h1[i] = MEMP1, MEMP1 += ll;
		h2[i] = MEMP2, MEMP2 += ll;
		h3[i] = MEMP3, MEMP3 += ll;
		h4[i] = MEMP4, MEMP4 += ll;
		fo(j, 0, ll - 1) h2[i][j] = h4[i][j] = 0;
	} if(l == 0) {
		fo(i, 0, B - 2) {
			fo(j, 0, ll - 1) h1[i][j] = A[j + i * len];
			DFT(h1[i], ll);
		}
	} fo(i, 0, B - 1) {
		IDFT(h2[i], ll);
		fo(j, 0, len - 1) f[l + i * len + j] = (f[l + i * len + j] + h2[i][j + len]) % mod;
		fo(j, 0, ll - 1) h1[B - 1][j] = 0;
		if(l == 0) fo(j, 0, i - 1) fo(k, 0, ll - 1) h1[B - 1][k] = (h1[B - 1][k] + (LL)h3[i - j - 1][k] * h4[j][k]) % mod;
		else fo(j, 0, i - 1) fo(k, 0, ll - 1) h1[B - 1][k] = (h1[B - 1][k] + 2LL * h3[i - j - 1][k] * h4[j][k]) % mod;
		IDFT(h1[B - 1], ll);
		fo(j, 0, len - 1) g[l + i * len + j] = (g[l + i * len + j] + h1[B - 1][j + len]) % mod;
		solve(l + i * len, l + (i + 1) * len - 1, MEMP1, MEMP2, MEMP3, MEMP4);
		if(i == B - 1) break;
		pre(ll - 1);
		if(l == 0) {
			fo(j, 0, ll - 1) h3[i][j] = h4[i][j] = 0;
			fo(j, 0, len - 1) h3[i][j] = h4[i][j] = f[l + i * len + j];
			if(i) {
				IDFT(h3[i - 1], ll);
				fo(j, 0, len - 1) h3[i - 1][j + len] = f[l + i * len + j];
				DFT(h3[i - 1], ll);
			} DFT(h3[i], ll), DFT(h4[i], ll);
		} else {
			fo(j, 0, ll - 1) h4[i][j] = 0;
			fo(j, 0, len - 1) h4[i][j] = f[l + i * len + j];
			DFT(h4[i], ll);
		} fo(j, 0, ll - 1) h1[B - 1][j] = 0;
		fo(j, 0, len - 1) h1[B - 1][j] = g[l + i * len + j];
		DFT(h1[B - 1], ll);
		fo(j, i + 1, B - 1) fo(k, 0, ll - 1) h2[j][k] = (h2[j][k] + (LL)h1[B - 1][k] * h1[j - i - 1][k]) % mod;
	} pre(ll - 1); if(l == 0) {
		IDFT(h3[B - 2], ll);
		fo(j, 0, len - 1) h3[B - 2][j + len] = f[l + (B - 1) * len + j];
		DFT(h3[B - 2], ll);
	}
}

char ss[N];

int main() {
	reg n = read(), len = Pre(n);
	scanf("%s", ss);
	fo(i, 0, n - 1) A[i] = ss[i] - '0';
	fd(i, n - 1, 0) A[i + 1] = A[i] * inv[i]; A[0] = 0;
	solve(0, len - 1, MEM1, MEM2, MEM3, MEM4);
	fo(i, 1, n) pr2((LL)f[i] * jc[i] % mod);
	return 0;
}