分析
要求的就是每个节点的儿子数都不超过 $3$ 的无标号有根树数。
设其 OGF 为 $F(x)$,考虑列出一个关于 $F(x)$ 的方程。
不妨规定 $[x^0]F(x) = 1$,那么可以看成每个节点的儿子数都恰好为 $3$。
直接得出 $F(x) = 1 + F(x) ^ 3$ 显然是错误的,因为我们同构的情况被重复计算了。
因此可以考虑 Burnside 引理。对于 $(1, 2, 3)$,所有情况都是不动点;对于 $(1, 3, 2)$、$(3, 2, 1)$、$(2, 1, 3)$,在置换中成环的两棵子树同构的情况是不动点;对于 $(3, 1, 2)$、$(2, 3, 1)$,三棵子树都同构的情况是不动点。
这样子我们就可以列出方程
$$
F(x) = 1 + x \frac{F(x) ^ 3 + 3 F(x ^ 2) F(x) + 2 F(x ^ 3)}{6}
$$
考虑使用牛顿迭代解出 $F(x)$。设 $G(F(x)) = 1 + x \frac{F(x) ^ 3 + 3 F(x ^ 2) F(x) + 2 F(x ^ 3)}{6} - F(x) = 0$。
假设当前已知 $F_0(x) = F(x) \bmod x ^ {\frac{n}{2}}$,要求 $F(x) \bmod x ^ n$。那么事实上 $F(x ^ 2) \bmod{x ^ n}$ 和 $F(x ^ 3) \bmod x ^ n$ 已知,因此可以将其看作常量 $A(x)$ 和 $B(x)$。这样子有
$$
G'(F(x)) = x \frac{F(x) ^ 2 + 3 A(x)}{6} - 1
$$
套牛顿迭代的式子
$$
F(x) = F_0(x) - \frac{G(F_0(x))}{G'(F_0(x))}
$$
这样子就可以求出 $F(x)$ 了。时间复杂度 $\mathcal{O}(n\log n)$。
代码
// ====================================
// author: M_sea
// website: https://m-sea-blog.com/
// ====================================
#include <bits/stdc++.h>
#define file(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
#define debug(...) fprintf(stderr, __VA_ARGS__)
using namespace std;
typedef long long ll;
int read() {
int X = 0, w = 1; char c = getchar();
while (c < '0' || c > '9') { if (c == '-') w = -1; c = getchar(); }
while (c >= '0' && c <= '9') X = X * 10 + c - '0', c = getchar();
return X * w;
}
const int N = 262144 + 10;
const int mod = 998244353, inv2 = 499122177;
int qpow(int a, int b) {
int c = 1;
for (; b; b >>= 1, a = 1ll * a * a % mod)
if (b & 1) c = 1ll * c * a % mod;
return c;
}
namespace Poly {
int r[N];
void NTT(int *A, int n, int op) {
for (int i = 0; i < n; ++i)
if (i < r[i]) swap(A[i], A[r[i]]);
for (int i = 1; i < n; i <<= 1) {
int rot = qpow(op == 1 ? 3 : 332748118, (mod - 1) / (i << 1));
for (int j = 0; j < n; j += i << 1)
for (int k = 0, w = 1; k < i; ++k, w = 1ll * w * rot % mod) {
int x = A[j + k], y = 1ll * w * A[j + k + i] % mod;
A[j + k] = (x + y) % mod, A[j + k + i] = (x - y + mod) % mod;
}
}
if (op == -1) {
int inv = qpow(n, mod - 2);
for (int i = 0; i < n; ++i) A[i] = 1ll * A[i] * inv % mod;
}
}
int NTT_init(int n) {
int lim = 1, l = 0;
for (; lim < n; lim <<= 1, ++l);
for (int i = 0; i < lim; ++i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
return lim;
}
void PolyInv(int *F, int *G, int n) {
static int A[N], B[N];
if (n == 1) { G[0] = qpow(F[0], mod - 2); return; }
PolyInv(F, G, n >> 1);
for (int i = 0; i < n; ++i) A[i] = F[i], B[i] = G[i];
int lim = NTT_init(n << 1);
NTT(A, lim, 1), NTT(B, lim, 1);
for (int i = 0; i < lim; ++i) A[i] = 1ll * A[i] * B[i] % mod * B[i] % mod;
NTT(A, lim, -1);
for (int i = 0; i < n; ++i) G[i] = (2ll * G[i] - A[i] + mod) % mod;
for (int i = 0; i < lim; ++i) A[i] = B[i] = 0;
}
}
void Newton(int *F, int n) {
static int F0[N], A[N], B[N], X[N], Y[N], iY[N];
F[0] = 1;
for (int lim = 2; lim <= n; lim <<= 1) {
for (int i = 0; i < lim; ++i) F0[i] = F[i];
for (int i = 0; i < lim; i += 2) A[i] = F[i / 2];
for (int i = 0; i < lim; i += 3) B[i] = F[i / 3];
Poly::NTT_init(lim << 1);
Poly::NTT(F0, lim << 1, 1), Poly::NTT(A, lim << 1, 1);
for (int i = 0; i < lim << 1; ++i)
X[i] = (1ll * F0[i] * F0[i] % mod * F0[i] + 3ll * A[i] * F0[i]) % mod;
Poly::NTT(X, lim << 1, -1), Poly::NTT(F0, lim << 1, -1), Poly::NTT(A, lim << 1, -1);
for (int i = 0; i < lim; ++i) X[i] = (X[i] + 2ll * B[i]) % mod;
for (int i = lim - 1; i; --i) X[i] = X[i - 1];
X[0] = 6;
for (int i = 0; i < lim; ++i) X[i] = (X[i] + 1ll * (mod - 6) * F0[i]) % mod;
Poly::NTT(F0, lim << 1, 1);
for (int i = 0; i < lim << 1; ++i) Y[i] = 3ll * F0[i] * F0[i] % mod;
Poly::NTT(Y, lim << 1, -1), Poly::NTT(F0, lim << 1, -1);
for (int i = 0; i < lim; ++i) Y[i] = (Y[i] + 3ll * A[i]) % mod;
for (int i = lim - 1; i ; --i) Y[i] = Y[i - 1];
Y[0] = mod - 6;
for (int i = lim; i < lim << 1; ++i) X[i] = Y[i] = 0;
Poly::PolyInv(Y, iY, lim);
Poly::NTT_init(lim << 1);
Poly::NTT(X, lim << 1, 1), Poly::NTT(iY, lim << 1, 1);
for (int i = 0; i < lim << 1; ++i) X[i] = 1ll * X[i] * iY[i] % mod;
Poly::NTT(X, lim << 1, -1);
for (int i = 0; i < lim; ++i) F[i] = (F0[i] - X[i] + mod) % mod;
for (int i = 0; i < lim << 1; ++i) F0[i] = A[i] = B[i] = X[i] = Y[i] = iY[i] = 0;
}
}
int F[N], G[N], H[N];
int main() {
int n = read(), lim = 1;
for (; lim <= n; lim <<= 1);
Newton(F, lim);
Poly::NTT_init(lim << 1);
Poly::NTT(F, lim << 1, 1);
for (int i = 0; i < lim << 1; ++i) G[i] = 1ll * F[i] * F[i] % mod;
Poly::NTT(G, lim << 1, -1), Poly::NTT(F, lim << 1, -1);
for (int i = 0; i < lim; i += 2) G[i] = (G[i] + F[i / 2]) % mod;
for (int i = 0; i < lim; ++i) G[i] = 1ll * G[i] * inv2 % mod;
for (int i = lim - 1; i; --i) G[i] = G[i - 1];
G[0] = 0;
for (int i = lim; i < lim << 1; ++i) G[i] = 0;
Poly::NTT(G, lim << 1, 1);
for (int i = 0; i < lim << 1; ++i) H[i] = 1ll * G[i] * G[i] % mod;
Poly::NTT(H, lim << 1, -1), Poly::NTT(G, lim << 1, -1);
for (int i = 0; i < lim; i += 2) H[i] = (H[i] + G[i / 2]) % mod;
for (int i = 0; i < lim; ++i) H[i] = 1ll * H[i] * inv2 % mod;
for (int i = 2; i <= n; ++i) printf("%d\n", H[i]);
return 0;
}