跳至主要內容

[Hard] P5900 无标号无根树计数

xiezheyuan2025/1/11做题笔记生成函数图论计数约 2274 字大约 11 分钟

题目足够简要，无需简要题意，直接讲思路。

感觉符号化方法还没有普及，所以在这里介绍一下，大部分文段来自我自己在校内的《组合数学信息方法》讲稿，因此不构成侵权。

符号化方法

概述

在组合数学中，我们最常遇到的问题是，计算大小为 nn 的某种东西的数量，这种东西同时也要满足一些限制。我们称这种可以计数的对象称为组合对象。全体组合对象称为组合类。一个组合类的 OGF，通常为其计数序列（大小对应数量）的 OGF。

在实际应用中，这些对象往往具有一些特征，如：

  • 可以分割为若干个部分，每个部分对大小的贡献为 11。
  • 不同的部分通过经典的方法组合在一起，比如序列 SEQ\textsf{SEQ}、可重集合 MSET\textsf{MSET}、环 CYC\textsf{CYC} 等等。

符号化方法用来将组合对象建模，使得其容易用生成函数来表示，进而用多项式相关算法来求解。

基础无标号类

我们称 ϵ\epsilon 为中性对象，并认为 ∣ϵ∣=0|\epsilon|=0，表示其不会对大小产生贡献。中性对象的组合类 {ϵ}\{\epsilon\} 记作 E\mathcal{E}，其生成函数为 E(z)=1\mathbf{E}(z)=1。

我们称 ∙\bullet 为原子对象，并认为 ∣∙∣=1|\bullet|=1，表示其对大小产生贡献为 11。原子对象的组合类 {∙}\{\bullet\}，记作 Z\mathcal{Z}，其生成函数为 Z(z)=z\mathbf{Z}(z)=z。

对于两个组合对象的直积 C=A×B\mathcal{C}=\mathcal{A}\times \mathcal{B}，有 C(z)=A(z)B(z)\mathbf{C}(z)=\mathbf{A}(z)\mathbf{B}(z)，对于两个不交的组合对象的并 C=A∪B\mathcal{C}=\mathcal{A}\cup \mathcal{B}，有 C(z)=A(z)+B(z)\mathbf{C}(z)=\mathbf{A}(z)+\mathbf{B}(z)。

这两个结论是显然的，本质上是加法原理和乘法原理。

无标号序列构造

定义 SEQ(A)\textsf{SEQ}(\mathcal{A}) 表示 A\mathcal{A} 中所有元素组成的序列的组合类，其生成函数为：

Q(A)=∑i=0+∞Ai(z)=11−A(z)\mathbf{Q}(\mathcal{A})=\sum_{i=0}^{+\infty} \mathbf{A}^i(z)=\frac{1}{1-\mathbf{A}(z)}

例子：0101 组合类 B\mathcal{B} 的生成函数 B=2z\mathbf{B}=2z，则其序列的生成函数为 Q(B)=11−2z\mathbf{Q}(\mathcal{B})=\frac{1}{1-2z}。根据无穷级数的知识，可以得到 [zn]11−2z=2n[z^n]\frac{1}{1-2z}=2^n。

无标号可重集合构造

定义 MSET(A)\textsf{MSET}(\mathcal{A}) 表示 A\mathcal{A} 中所有元素组成的可重集合的组合类。

概念说清楚了，现在的问题是 MSET\textsf{MSET} 的生成函数是什么？

考虑一个显然正确，但是不够通用的生成函数：

M(A)=∏a∈AQ({a})=∏a∈A11−z∣a∣=∏i=1+∞(1−zi)−[zi]A(z)\mathbf{M}(\mathcal{A})=\prod_{a\in \mathcal{A}} \mathbf{Q}(\{a\})=\prod_{a\in \mathcal{A}}\frac{1}{1-z^{|a|}}=\prod_{i=1}^{+\infty} (1-z^i)^{-[z^i]\mathbf{A}(z)}

看到 ∏\prod 不好处理，对其取 ln⁡\ln 得：

ln⁡M(A)=∑i=1+∞(−[zi]A(z))ln⁡(1−zi)\ln \mathbf{M}(\mathcal{A}) = \sum_{i=1}^{+\infty} (-[z^i]\mathbf{A}(z))\ln(1-z^i)

根据 Taylor 展开的知识，我们有：

ln⁡(1−z)=−∑k=1+∞zkk\ln(1-z) = -\sum_{k=1}^{+\infty} \frac{z^k}{k}

代入，得：

ln⁡MSET(A)=∑i=1+∞([zi]A(z))∑k=1+∞zikk\ln \mathbf{MSET}(\mathcal{A}) = \sum_{i=1}^{+\infty} ([z^i]\mathbf{A}(z))\sum_{k=1}^{+\infty} \frac{z^{ik}}{k}

变形得：

ln⁡MSET(A)=∑k=1+∞1k∑i=1+∞([zi]A(z))zik=∑k=1+∞A(zk)k\ln \mathbf{MSET}(\mathcal{A}) = \sum_{k=1}^{+\infty}\frac{1}{k}\sum_{i=1}^{+\infty}([z^i]\mathbf{A}(z)) z^{ik} = \sum_{k=1}^{+\infty}\frac{\mathbf{A}(z^k)}{k}

取 exp⁡\exp 得：

MSET(A)=exp⁡(∑k=1+∞A(zk)k)\boxed{\mathbf{MSET}(\mathcal{A})=\exp\left(\sum_{k=1}^{+\infty}\frac{\mathbf{A}(z^k)}{k}\right)}

我们用 Exp\mathbf{Exp} 代替 MSET\mathbf{MSET} 表示无标号可重集合组合类的生成函数，并称 Exp\mathbf{Exp} 为 Euler 变换，如果你对 Pólya 定理熟悉，也许可以看出这是 Pólya 定理的扩展。

无标号有根树计数

第一部分

考虑定义 T\mathcal{T} 表示无标号有根树的组合类，则有：

T=Z×MSET(T)\mathcal{T}=\mathcal{Z}\times \textsf{MSET}(\mathcal{T})

其中原子对象组合类 Z\mathcal{Z} 表示树根，可重集合 MSET(T)\textsf{MSET}(\mathcal{T}) 表示子树。

构建 T\mathcal{T} 的生成函数 T\mathbf{T}，可以利用刚刚提到的 Euler 变换：

T(z)=z⋅Exp(T)=z⋅exp⁡(∑k=1+∞T(zk)k)\mathbf{T}(z)=z\cdot \mathbf{Exp}(\mathcal{T})=z\cdot \exp\left(\sum_{k=1}^{+\infty}\frac{\mathbf{T}(z^k)}{k}\right)

第二部分

两边取 ln⁡\ln：

ln⁡T(z)=ln⁡(z)+∑k=1+∞T(zk)k\ln \mathbf{T}(z)=\ln(z)+\sum_{k=1}^{+\infty}\frac{\mathbf{T}(z^k)}{k}

两边求导：

T′(z)T(z)=1z+∑k=1+∞T′(zk)zk−1\frac{\mathbf{T}'(z)}{\mathbf{T}(z)}=\frac{1}{z}+\sum_{k=1}^{+\infty} {\mathbf{T}'(z^k)}z^{k-1}

变形，得：

zT′(z)=T(z)+T(z)∑k=1+∞T′(zk)zkz\mathbf{T}'(z)=\mathbf{T}(z)+\mathbf{T}(z)\sum_{k=1}^{+\infty}{\mathbf{T}'(z^k)}z^{k}

令：

F(z)=∑k=1+∞T′(zk)zk\mathbf{F}(z)=\sum_{k=1}^{+\infty}{\mathbf{T}'(z^k)}z^{k}

得到：

[zn](zT′(z))=[zn]T(z)+[zn](T∗F)[z^n](z\mathbf{T'}(z))=[z^n]\mathbf{T}(z)+[z^n](\mathbf{T}\ast \mathbf{F})

由于 [zn−1]T′(z)=n⋅[zn]T(z)[z^{n-1}]\mathbf{T}'(z)=n\cdot [z^n]\mathbf{T}(z)，得到：

n⋅[zn]T(z)=[zn]T(z)+[zn](T∗F)n\cdot [z^n]\mathbf{T}(z)=[z^n]\mathbf{T}(z)+[z^n](\mathbf{T}\ast \mathbf{F})

变形得：

[zn]T(z)=1n−1[zn](T∗F)\boxed{[z^n]\mathbf{T}(z)=\frac{1}{n-1}[z^n](\mathbf{T}\ast \mathbf{F})}

接下来考虑 F\mathbf{F}，大力求导与化简，得：

F(z)=∑k=1+∞T′(zk)zk=∑k=1+∞zk∑i=1+∞([xi]T(z))i(zk)i−1=∑k=1+∞∑i=1+∞([xi]T(z))izki=∑p=1+∞zp∑i∣pi([xi]T(z))\begin{aligned} \mathbf{F}(z)&=\sum_{k=1}^{+\infty}{\mathbf{T}'(z^k)}z^k\\ &=\sum_{k=1}^{+\infty}z^k\sum_{i=1}^{+\infty} ([x^i]\mathbf{T}(z))i{(z^k)}^{i-1}\\ &=\sum_{k=1}^{+\infty}\sum_{i=1}^{+\infty}([x^i]\mathbf{T}(z))iz^{ki}\\ &=\sum_{p=1}^{+\infty}z^p\sum_{i\mid p} i([x^i]\mathbf{T}(z)) \end{aligned}

得到：

F(z)=∑p=1+∞zp∑i∣pi([xi]T(z))\boxed{\mathbf{F}(z)=\sum_{p=1}^{+\infty}z^p\sum_{i\mid p} i([x^i]\mathbf{T}(z))}

第三部分

你发现这玩意是一个半在线卷积的形式，我们可以考虑分治 FFT。

时间复杂度 O(nlog⁡2n)O(n\log^2 n)。

但是我悲伤地发现，我根本就不会半在线卷积！于是我滚去学了一下。

推荐 command_block 巨佬的博客。

无标号无根树计数

考虑钦定重心为根，则我们要去掉不以重心为根的情况。

先来考虑 n≡1(mod2)n\equiv 1\pmod{2}，此时重心是唯一的，由于有根树的重心，一定位于大小至少为 ⌊n2⌋+1\lfloor\frac{n}{2}\rfloor+1 的子树内，于是可以将这个重子树看成一组，其余子树与根一起看成一组，则根不为重心的方案数为：

∑k=⌊n2⌋+1n−1([zk]T(z))⋅([zn−k]T(z))\sum_{k=\lfloor\frac{n}{2}\rfloor+1}^{n-1} ([z^k]\mathbf{T}(z))\cdot ([z^{n-k}]\mathbf{T}(z))

最后是 2∣n2\mid n 的情况，此时重心可能为两个点，而根据重心的性质，这两个点必定相邻，且断掉中间的那条边后，形成的两个连通块大小一定是 n2\frac{n}{2}，每个连通块一共有 ([zn2])T(z)([z^{\frac{n}{2}}])\mathbf{T}(z) 种，所以要减去 (([zn2])T(z)2)\binom{([z^{\frac{n}{2}}])\mathbf{T}(z)}{2}。

最后写一下完整的式子：

([zn]T(z))−∑k=⌊n2⌋+1n−1([zk]T(z))⋅([zn−k]T(z))−[2∣n](([zn2])T(z)2)([z^n]\mathbf{T}(z)) - \sum_{k=\lfloor\frac{n}{2}\rfloor+1}^{n-1} ([z^k]\mathbf{T}(z))\cdot ([z^{n-k}]\mathbf{T}(z)) - [2\mid n]\binom{([z^{\frac{n}{2}}])\mathbf{T}(z)}{2}

总时间复杂度 O(nlog⁡2n)O(n\log^2 n) 可以通过本题。

代码

#include <bits/stdc++.h>
using namespace std;

const int mod = 998244353, g = 3, ginv = 332748118, inv2 = (mod + 1) >> 1;

const int N = 1e6 + 5;

int Add(int x, int y){ return (x + y) >= mod ? (x + y - mod) : (x + y); }
int Sub(int x, int y){ return (x - y) < 0 ? (x - y + mod) : (x - y); }
int Mul(int x, int y){ return 1ll * x * y % mod; }

int fastpow(int x, int y){
    int ret = 1;
    for(;y;y>>=1,x=Mul(x,x)){
        if(y & 1) ret = Mul(ret, x);
    }
    return ret;
}

int Inv(int x){ return fastpow(x, mod-2); }

struct Poly{
    vector<int> a;
    Poly(){}
    Poly(vector<int> _a):a(_a){}
    Poly(int siz){ a.resize(siz); }
    int& operator[](int x){ return a[x]; }
    int size(){ return a.size(); }
    Poly resize(int siz){ a.resize(siz, 0); return *this; }
    void ntt(bool flag = 1);
    static int init_ntt(int n, int m);
    static Poly conv(Poly a, Poly b, function<int(int, int)> f);
};

int r[N], wt[N];

void Poly::ntt(bool flag){
    int n = size();
    for(int i=0;i<n;i++){
        if(i < r[i]) swap(a[i], a[r[i]]);
    }
    for(int i=2;i<(n<<1);i<<=1){
        int w1 = fastpow(flag ? g : ginv, (mod - 1) / i);
        for(int j=0;j<n;j+=i){
            int w = 1;
            for(int k=0;k<(i>>1);k++,w=Mul(w,w1)){
                int x = a[j + k], y = Mul(w, a[(i >> 1) + j + k]);
                a[j + k] = Add(x, y), a[(i >> 1) + j + k] = Sub(x, y);
            }
        }
    }
    if(!flag){
        int inv = Inv(n);
        for(int i=0;i<n;i++) a[i] = Mul(a[i], inv);
    }
}

int Poly::init_ntt(int n, int m){
    int c = 1, l = 0;
    while(c < (n + m)) c <<= 1, l++;
    for(int i=0;i<=c;i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
    return c;
}

Poly Poly::conv(Poly a, Poly b, function<int(int, int)> f){
    int n = a.size(), m = b.size();
    int c = Poly::init_ntt(n, m);
    a.resize(c), b.resize(c);
    a.ntt(), b.ntt();
    for(int i=0;i<c;i++) a[i] = f(a[i], b[i]);
    a.ntt(0), a.resize(n + m + 1);
    return a;
}

Poly operator*(Poly a, Poly b){ return Poly::conv(a, b, Mul); }

int F[N], G[N], n;

void solve(int l, int r){
    if(l == r){
        F[l] = (l == 1) ? 1 : Mul(F[l], Inv(l - 1));
        for(int i=l;i<=n;i+=l) G[i] = Add(G[i], Mul(l, F[l]));
        return;
    }
    int mid = (l + r) >> 1;
    solve(l, mid);
    Poly a(mid - l + 1), b(r - l + 1);
    for(int i=l;i<=mid;i++) a[i - l] = F[i];
    for(int i=1;i<=(l==1?mid:r-l);i++) b[i] = G[i];
    a = a * b;
    for(int i=mid+1;i<=r;i++) F[i] = Add(F[i], a[i - l]);
    if(l != 1){
        a.resize(mid - l + 1), b.resize(r - l + 1);
        for(int i=1;i<=r-l;i++) b[i] = F[i];
        for(int i=l;i<=mid;i++) a[i - l] = G[i];
        a = a * b;
        for(int i=mid+1;i<=r;i++) F[i] = Add(F[i], a[i - l]);
    }
    solve(mid + 1, r);
}

signed main(){
    ios::sync_with_stdio(false);
    cin.tie(0);cout.tie(0);
    cin >> n;
    solve(1, n);
    int ans = F[n];
    for(int k=((n>>1)+1);k<n;k++) ans = Sub(ans, Mul(F[k], F[n - k]));
    if(!(n & 1)) ans = Sub(ans, Mul(Sub(F[n >> 1], 1), Mul(inv2, F[n >> 1])));
    cout << ans << '\n';
    return 0;
}

// Written by xiezheyuan