题解：P12394 「RiOI-6」神曲
Register_int
·
·
题解
题解来自 @NaCly_fish。
设答案为 f_{n, m}，先考虑暴力如何写：先确定第 n 个区间，它将整个值域分为了三段 —— 剩下的 n - 1 个区间的值域都必须在这三段中，不能相交。那么枚举每一段的长度，再枚举对应段中有多少个区间，就有：
f_{n, m} = \sum_{l_1 + l_2 + l_3 = m} [l_2 \geq 1] \sum_{c_1 + c_2 + c_3 = n - 1} \binom{n - 1}{c_1, c_2, c_3} f_{c_1, l_1} f_{c_2, l_2} f_{c_3, l_3}
容易看到这是二元生成函数卷积的形式，设
F = \sum_{n = 0}^\infty \sum_{m = 0}^\infty \frac{x^ny^m}{n!}f_{n, m}
就可以得到方程（这里出现的 F - 1 是因为 l_2 \neq 0，要去除长度为 0 段的贡献）：
\frac{\partial F}{\partial x} = F (F - 1) F
虽然求的是偏导，但只与 x 有关，由此就能以反函数的形式求解。为了方便后续处理，设 G = F - 1，然后求出
x = \frac{1}{1 + G} + \ln \frac{-G}{1 + G} + C(y)
这里 \ln 中的东西常数项不为 1，难道我们算错了？当然不是，我们还没有确定 C(y) 是什么东西。可以令 C(y) = f(y) - \ln(-y)，其中 f(y) 是一个「正常的」函数，它可以在 y = 0 处做 Taylor 展开，得到对应的形式幂级数。
带回去，\ln 里面的东西常数项为零，可以正常计算了。那么只需要暴力算出 G 的低次项，就可以直接求出 f(y)，然后得到
x = \frac{1}{1 + G} + \ln \frac{G}{y(1 + G)} + y - 1
这个形式看起来可以 Lagrange 反演，但你别急。我们再稍微改写一下其形式：
y\text e^{x - y + 1} = \exp\left( \frac{1}{1 + G}\right) \frac{G}{1 + G}
设 H = G / (1 + G)，则：
y\text e^{x - y} = H\text e^{-H}
这让人想到什么？就是有标号有根树的 EGF T(x)！直接得到 H = T(y \text e^{x-y})。再带回 G 中就是
G = \frac{1}{1 - T(y \text e^{x - y})} - 1
根据一个经典结论（直接对 T(x) 的方程求导即可证明）：
x\frac{\text d T(x)}{\text dx} = \frac{1}{1 - T(x)} - 1
就能得到 G 的展开
G = \sum_{i\geq 1}^\infty \frac{i^i}{i!}(y \text e^{x - y})^i
这个东西就容易提取系数了：
f_{n,m} &= n![x^n y^m] G \\
&= \sum_{i = 1}^m n! \frac{i^i}{i!} [x^n y^{m-i}] \text e^{ix} \text e^{-iy} \\
&=\sum_{i = 1}^m \frac{i^{n + m}}{i!} \times \frac{(-1)^{m - i}}{(m - i)!} \\
&= \begin{Bmatrix} n + m \\ m \end{Bmatrix}
\end{aligned}
我们精心挑选了样例使得样例 oeis 不出来。
征集一个组合意义证法，找到了我 v 他 50。
感谢 @Watersphere 提供的组合意义证法！这部分请看他的题解。
问题转化为求一条对角线的第二类斯特林数。当然你可以贺 P8561 然后多项式多点求值，但是跑得太慢了。这里需要更快的方法，我们有：
\left\{\begin{matrix}n\\m\end{matrix}\right\}=\frac{n!}{m!}[x^n](e^x-1)^m
所以相当于对每个 m，要求
\frac{(n+m)!}{m!}[x^{n+m}](e^x-1)^m=\frac{(n+m)!}{m!}[x^n]\left(\frac{e^x-1}{x}\right)^m
当然你到这里就可以直接二元 bostan-mori 草上去了。不过虽然这玩意常数是真小，但是笔者测了一发这玩意跑了六秒。如果你硬卡过了的话那说明你是卡常大师，算你厉害。
怎么大家都是卡常大师？？？
我们当然希望可以拉反，但是 \frac{e^x-1}x 的常数项不为 0。怎么办呢？直接简单粗暴地把它变成 \frac{e^x-x-1}x+1 然后大力展开。设 F(x)=\frac{e^x-x-1}x，则原式为：
\begin{aligned}
&\frac{(n+m)!}{m!}[x^n](F+1)^m\\
=&\frac{(n+m)!}{m!}[x^n]\sum^m_i\binom miF^i\\
=&\frac{(n+m)!}{m!}\sum^m_i\binom mi[x^n]F^i
\end{aligned}
设 G=F^{<-1>}，施另类拉格朗日反演：
\begin{aligned}
=&\frac{(n+m)!}{m!}\sum^m_i\binom mi[x^{-i-1}]G'G^{-n-1}\\
=&\frac{(n+m)!}{m!}\sum^m_i\binom mi[x^{n-i}]G'\left(\frac x{G}\right)^{n+1}\\
=&(n+m)!\sum^m_i\frac1{(m-i)!}\left(\frac1{i!}[x^{n-i}]G'\left(\frac x{G}\right)^{n+1}\right)\\
\end{aligned}
明显是卷积形式，求出 G 即可。因为 F(G(x))=x，所以有：
\frac{e^G-G-1}G=x
可以牛顿迭代，总时间复杂度 O(n\log n)。
下面代码实际上是 O(n\log^2 n/\log\log n) 的。因为求 \exp 牛迭太慢了，用的是半在线卷积。其余部分是 O(n\log n)。
#include <bits/stdc++.h>
using ll = long long;
using namespace std;
/*
a big poly template by Aleph1022.
*/
poly poly::pow(int k) const {
int lead = a[0];
poly f = poly(vector<int>(a.begin(), a.begin() + deg() + 1)) * nt.inv(lead);
return (f.ln() * k).srcExp() * mpow(lead, (k % (mod - 1) + mod - 1) % (mod - 1));
}
poly iterate(int n) {
poly f = zeroes(1); f[1] = 2;
for (int len = 2; ; len <<= 1) {
f.resize(len + 2); poly g = f.srcExp(), h = g; h[1]--;
for (int i = 1; i <= len + 1; i++) sub(g[i], f[i]), sub(g[i], f[i - 1]);
g = g.shift(-1) * h.shift(-1).inv();
for (int i = 0; i <= len; i++) sub(f[i], g[i]);
if (len >= n) break;
}
return f.slice(n);
}
int n, m;
int main() {
scanf("%d%d", &n, &m);
poly g = iterate(n + 1), tg = g.shift(-1).pow(-n - 1);
g = (g.deriv() * tg).slice(n), reverse(g.a.begin(), g.a.end());
poly h = zeroes(m);
for (int i = 0; i <= m; i++) g[i] = (ll)g[i] * gifac[i] % mod, h[i] = gifac[i];
h = (h * g).slice(m);
for (int i = 1; i <= m; i++) printf("%lld ", (ll)h[i] * gfac[n + i] % mod);
}
火大了。赛后加强。等着我。