闲话 22.12.9

闲话

我说我今天主要刷的是网络流题你信嘛（
主要是在上网课 然后老师讲课：“一说荒漠就会想到仙人掌”
就在群里发这题
浅看了一下发现不会 然后写了

是……辣个男人！他又来了（

有标号荒漠计数

P5434

我们定义一张图是荒漠，当且仅当其每个极大联通分量都是仙人掌。
给定 \(n\)，求有多少 \(n\) 个点的不同的荒漠。结果对 \(998244353\) 取模。

\(n \le 10^5\)。

首先要计数有标号荒漠，可以先计数有标号仙人掌。最后作 \(\exp\) 就能得到答案。
注意答案是无根的，但是仙人掌显然是有根好统计一些。因此我们记 \(G\) 为有标号有根仙人掌计数的生成函数。考虑将一系列根连成环接在指定的一个根上，并由环没有方向，有

\[\begin{aligned} G & = x\exp\left(G + \frac 12 \sum_{i \ge 2} G^i\right) \\ & = x\exp\left(G + \frac 12 \left(\frac{1}{1 - G} - (1 + G) \right) \right) \\ & = x\exp\left(G + \frac{G^2}{2(1 - G)} \right) \\ & = x\exp\left( \frac{2G - G^2}{2(1 - G)} \right) \end{aligned}\]

无根考虑答案即为 \([x^n]\exp\left(xG'(x)\right)\)。
由于 \(\exp\) 里面东西太多，我们考虑记

\[h(x) = \frac{2x - x^2}{2(1-x)} \]

设 \(F\) 为答案计数的生成函数，有 \(F = x\left(\exp(h(x))\right)'\)。
到这一步的话自然想到 Lagrange 反演。设 \(F = g(f(x))\)，有 \(G = x\left(g(f(x))\right)' = x\ f'(x)\ g'(f(x))\)。现在只需要把另一部分化成同一形式就行了。

\[\begin{aligned} G' & = \exp(h(G)) + x\left(\exp(h(G))\right)' \\ & = \exp(h(G)) + x\ G'\ h'(G)\ \exp(h(G)))' \\ & = \exp(h(G))(1 + x\ G'\ h'(G)) \\ xG'& = x\exp(h(G))(1 + x\ G'\ h'(G)) \\ & = G(1 + x\ G'\ h'(G)) \\ & = G + x\ G\ G'\ h'(G) \\ G & = xG' - x\ G\ G'\ h'(G) \\ & = xG'(1 - G \ h'(G) ) \end{aligned}\]

自然导出了 \(f = G, g'(x) = 1 - xh'(x)\)。记 \(H(x) = \int h(x) dx\)，有

\[\begin{aligned} g(x) & = \int1 - xh'(x) \text dx \\ & = x - \int xh'(x) \text dx \\ & = x - \left(x h(x) - \int h(x) \text dx\right) \\ & = x - x h(x) + H(x) \end{aligned}\]

这样我们需要的所有拼图都齐全了。拉格朗日反演得到

\[[x^n] e^{F(x)} = [x^n] e^{g(f(x))} = \frac 1n [x^{n-1}] g'(x) e^{g(x)} \left(\frac {x}{f^{⟨-1⟩}(x)}\right) \]

于是可以分别讨论三个（两个？）部分如何做。

\(g'(x)\)

显然

下面需要 \(h'(x)\)。根据 Wolfram|Alpha，

\[h'(x) = \frac{x^2 - 2x - 2}{2(1 - x)^2} \]

于是

\[\begin{aligned} g'(x) = 1 - xh'(x) = 1 - \frac{x^3 - 2x^2 - 2x}{2(1 - x)^2} = \frac{-\frac 12 x^3 + 2x^2 - 3x + 1}{(1 - x)^2} \end{aligned}\]

\(x / f^{⟨-1⟩}(x)\)

观察 \(f^{⟨-1⟩}\) 是个啥。其一是 \(f^{⟨-1⟩} (f(x)) = x\)，其二是 \(f = x\exp\left( \frac{2f- f^2}{2(1 - f)} \right)\)。我们发现，假设 \(f^{⟨-1⟩}\) 中是连乘形式，且其中有 \(x\) 项，将 \(f(x)\) 带入后，只需要将后面的 \(\exp\) 消除即可。
不妨考虑 \(f^{⟨-1⟩}(x) = x\exp\left(\frac{x^2 - 2x}{2(1 - x)}\right)\)，容易发现此函数满足需要。

因此我们能得到 \(x / f^{⟨-1⟩}(x) = \exp\left(\frac{2x - x^2}{2(1 - x)}\right)\)。
因为上面还有一个 \(n\) 次方，该部分的贡献是 \(\exp\left(\frac{n(2x - x^2)}{2(1 - x)}\right)\) 。

\(e^{g(x)}\)

下面需要 \(H(x)\)。根据 Wolfram|Alpha，

\[H(x) = \frac 14\left((x-1)^2 - 2\ln(1- x)\right) \]

于是

\[\begin{aligned} g(x) & = x - x h(x) + H(x) \\ & = x + \left(\frac{x - 1}2\right)^2 - \frac 12 \ln (1 - x) - \frac{x(2x - x^2)}{2(1 - x)} \\ & = \frac{2x + x^2} 4 - \frac 12 \ln (1 - x) - \frac{x(2x - x^2)}{2(1 - x)} \\ e^{g(x)} & = \exp\left( \frac{2x + x^2} 4 - \frac 12 \ln (1 - x) - \frac{x(2x - x^2)}{2(1 - x)} \right) \\ [后两项] & = \exp\left( \frac{2x + x^2} 4 - \frac 12 \ln (1 - x) - \frac{x(2x - x^2)}{2(1 - x)} + \frac{n(2x - x^2)}{2(1 - x)} \right) \\ [后两项] & = (1 - x)^{- 1/2} \exp\left( \frac{2x + x^2} 4 - \frac{x(2x - x^2)}{2(1 - x)} + \frac{n(2x - x^2)}{2(1 - x)} \right) \\ [后两项] & = (1 - x)^{- 1/2} \exp\left( \frac{2x + x^2} 4 + (n - x)\frac{(2x - x^2)}{2(1 - x)} \right) \\ [整个的] & = \frac{-\frac 12 x^3 + 2x^2 - 3x + 1}{(1 - x)^{2/5}} \exp\left( \frac{2x + x^2} 4 + (n - x)\frac{(2x - x^2)}{2(1 - x)} \right) \end{aligned}\]

根据化简过程中不很简的部分，这玩意微分有限，存在线性递推，因此可以 \(O(n) / O(\sqrt n \log n)\) 地算。

具体地，我们首先考察 \(a(x) = \exp\left( \frac{2x + x^2} 4 + (n - x)\frac{(2x - x^2)}{2(1 - x)} \right)\)。求导后化简可以得到

\[2(1 - x)^2 a'(x) = \left((2n + 1) - (2n + 5) x + (n + 4)x^2 - x^3\right) a(x) \]

过程 dirty，不细说了（

随后考察 \(b(x) = (1 - x)^{-5/2} a(x)\)。同理有

\[2(1 - x)^2 b'(x) = \left((2n + 6) - (2n + 10) x + (n + 4)x^2 - x^3\right) b(x) \]

这自然导出了递推方案。答案即为

\[[x^n / n!]\ b(x) \left(-\frac 12 x^3 + 2x^2 - 3x + 1\right) \]

初值抄鰰的（
代码全屏好看（

code

#include <bits/stdc++.h>
using namespace std; using pii = pair<int,int>; using vi = vector<int>; using ll = long long; 
using ull = unsigned long long; using db = double; using ld = long double; using lll = __int128_t;
mt19937 rnd(chrono::steady_clock::now().time_since_epoch().count());
template <typename T> T rand(T l, T r) { return uniform_int_distribution<T>(l, r)(rnd); }
template <typename T1, typename T2> T1 min(T1 a, T2 b) { return a < b ? a : b; }
template <typename T1, typename T2> T1 max(T1 a, T2 b) { return a > b ? a : b; }
#define file(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
#define ufile(x) 
#define rep(i,s,t) for (register int i = (s), i##_ = (t) + 1; i < i##_; ++ i)
#define pre(i,s,t) for (register int i = (s), i##_ = (t) - 1; i > i##_; -- i)
#define Aster(i, s) for (int i = head[s], v; i; i = e[i].next)
#define all(s) s.begin(), s.end()
#define eb emplace_back
#define pb pop_back
#define em emplace
const int N = 1e5 + 10, mod = 998244353, ans[] = {1, 1, 2, 8, 57};
const int inf = 0x3f3f3f3f;
const ll infll = 0x3f3f3f3f3f3f3f3fll;
int n, inv[N], fac;

signed main() {
    ios::sync_with_stdio(false), cin.tie(0), cout.tie(0);
    cin >> n;
    if (n <= 4) cout << ans[n], exit(0);
    inv[0] = inv[1] = fac = 1;
    rep(i,2,n) inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod, fac = 1ll * fac * i % mod;
    vector<int> coef{ 0, 2 * n + 6, mod - (2 * n + 10), n + 4, mod - 1 },
    g{ 1, n + 3, ((1ll * (n + 7) * n >> 1) + 5) % mod, ((1ll * (n + 12) * n % mod + 42) * n % mod + 38) * inv[6] % mod };
    g.resize(n + 1);
    rep(i,4,n) {
        g[i] = (1ll * coef[1] * g[i - 1] % mod + 1ll * coef[2] * g[i - 2] % mod +
                1ll * coef[3] * g[i - 3] % mod + 1ll * coef[4] * g[i - 4] % mod) * inv[2] % mod;
        g[i] = (g[i] + 1ll * ((i - 1) << 1) * g[i - 1] % mod - 1ll * (i - 2) * g[i - 2] % mod + mod) % mod * inv[i] % mod;
    }
    fac = 1ll * fac * (g[n] - 3ll * g[n - 1] % mod + 2ll * g[n - 2] - 1ll * inv[2] * g[n - 3] % mod + 2 * mod) % mod;
    cout << (fac + mod) % mod << endl;
}