图论计数

在组合数学中，图论计数（Graph Enumeration）是研究满足特定性质的图的计数问题的分支。生成函数、波利亚计数定理与符号化方法和 OEIS 是解决这类问题时最重要的数学工具。图论计数可分为有标号和无标号两大类问题，大多数情况下¹有标号版本的问题都比其对应的无标号问题更加简单，因此我们将先考察有标号问题的计数。

有标号树

即 Cayley 公式，参见 Prüfer 序列一文，我们也可以使用 Kirchhoff 矩阵树定理或生成函数和拉格朗日定理得到这一结果。

习题

Hihocoder 1047. Random Tree

有标号连通图

例题「POJ 1737」Connected Graph

题目大意：求有个结点的有标号连通图的方案数（ $n \leq 50$ ）。

这类问题最早出现于楼教主的男人八题系列中，我们设为个点有标号图的方案数，为待求序列。个点的图至多有 $\binom{n}{2}$ 条边，每条边根据其出现与否有两种状态，每种状态之间独立，因而有 $g_n = 2^{\binom{n}{2}}$ 。我们固定其中一个节点，枚举其所在连通块的大小，那么还需要从剩下的个节点中选择个节点组成一个连通块。连通块之外的节点可以任意连边，因而有如下递推关系：

\begin{align} \sum_{i=1}^{n} \binom{n-1}{i-1} c_i g_{n-i} &= g_n \\ c_n &= g_n - \sum_{i=1}^{n-1} \binom{n-1}{i-1} c_i g_{n-i} \end{align}

移项得到序列的递推公式，可以通过此题。

例题「集训队作业 2013」城市规划

题目大意：求有个结点的有标号连通图的方案数（ $n \leq 130000$ ）。

对于数据范围更大的序列问题，往往我们需要构造这些序列的生成函数，以使用高效的多项式算法。

方法一：分治 FFT

上述的递推式可以看作一种自卷积形式，因而可以使用分治 FFT 进行计算，复杂度 $O(n\log^2n)$ 。

方法二：多项式求逆

我们将上述递推式中的组合数展开，并进行变形：

\begin{align} \sum_{i=1}^{n} \binom{n-1}{i-1} c_i g_{n-i} &= g_n \\ \sum_{i=1}^{n} \frac{c_i}{(i-1)!} \frac{g_{n-i}}{(n-i)!} &= \frac{g_n}{(n-1)!} \end{align}

构造多项式：

\begin{align} C(x) &= \sum_{n=1} \frac{c_n}{(n-1)!} x^n \\ G(x) &= \sum_{n=0} \frac{g_n}{n!} x^n \\ H(x) &= \sum_{n=1} \frac{g_n}{(n-1)!} x^n \end{align}

代换进上式得到，使用多项式求逆后再卷积解出即可。

方法三：多项式 exp

另一种做法是使用 EGF 中多项式 exp 的组合意义，我们设有标号连通图和简单图序列的 EGF 分别为和，那么它们将有下列关系：

\begin{align} \exp(C(x)) &= G(x) \\ C(x) &= \ln(G(x)) \end{align}

使用多项式 ln 解出即可。

有标号欧拉图、二分图

例题「SPOJ KPGRAPHS」Counting Graphs

题目大意：求有个结点的分别满足下列性质的有标号图的方案数（ $n \leq 1000$ ）。

连通图 A001187。
欧拉图 A033678。
二分图 A047864。

本题限制代码长度，因而无法直接使用多项式模板，但生成函数依然可以帮助我们进行分析。

连通图问题在之前的例题中已被解决，考虑欧拉图。注意到上述对连通图计数的几种方法，均可以在满足任意性质的有标号连通图进行推广。例如我们可以将连通图递推公式中的，从任意图替换成满足顶点度数均为偶数的图，此时得到的即为欧拉图。

我们将 POJ 1737 的递推过程封装成连通化函数，

void ln(Int C[], Int G[]) {
  for (int i = 1; i <= n; ++i) {
    C[i] = G[i];
    for (int j = 1; j <= i - 1; ++j)
      C[i] -= binom[i - 1][j - 1] * C[j] * G[i - j];
  }
}

前两问即可轻松解决：

for (int i = 1; i <= n; ++i) G[i] = pow(2, binom[i][2]);
ln(C, G);
for (int i = 1; i <= n; ++i) G[i] = pow(2, binom[i - 1][2]);
ln(E, G);

注意到这里的连通化递推过程其实等价于对其 EGF 求多项式 ln，同理我们也可以写出逆连通化函数，它等价于对其 EGF 求多项式 exp。

void exp(Int G[], Int C[]) {
  for (int i = 1; i <= n; ++i) {
    G[i] = C[i];
    for (int j = 1; j <= i - 1; ++j)
      G[i] += binom[i - 1][j - 1] * C[j] * G[i - j];
  }
}

下面讨论有标号二分图计数，

我们设表示 n 个结点的二分图方案数，表示个结点对结点进行 2 染色，满足相同颜色的结点之间不存在边的图的方案数。枚举其中一种颜色节点的数量，有²：

g_n = \sum_{i=0}^{n} \binom{n}{i}2^{i(n-i)}

接下来我们用两种不同的方法建立与之间的关系。

方法一：算两次

我们设 $c_{n, k}$ 表示有 k 个连通分量的二分图方案数，那么不难得到如下关系：

\begin{align} b_n &= \sum_{i=1}^{n} c_{n, i} \\ g_n &= \sum_{i=1}^{n} c_{n, i} 2^i \end{align}

比较两种的表达式，展开得：

\begin{align} \sum_{i=0}^{n} \binom{n}{i}2^{i(n-i)} &= \sum_{i=1}^{n} c_{n, i} 2^i \\ c_{n, i} &= \sum_{i=0}{n-1} \binom{n-1}{i-1} c_{n, 1}c_{n-i,k-1} \end{align}

不难得到的递推关系，复杂度，进一步使用容斥原理，可以优化到通过本题。

方法二：连通化递推

方法二和方法三均使用连通二分图 A001832 来建立与之间的桥梁。

注意到对于每个连通二分图，我们恰好有两种不同的染色方法，对应到两组不同的连通 2 染色图，因而对进行连通化，得到的序列恰好是的两倍，而则由进行逆连通化得到。

因此：

for (int i = 1; i <= n; ++i) {
  G[i] = 0;
  for (int j = 0; j < i + 1; ++j) G[i] += binom[i][j] * pow(2, j * (i - j));
}
ln(B1, G);
for (int i = 1; i <= n; ++i) B1[i] /= 2;
exp(B, B1);

两种递推的过程复杂度均为，可以通过本题。

方法三：多项式 exp

我们注意到也可以使用 EGF 理解上面的递推过程。

设为的 EGF，为的 EGF，为的 EGF，应用做法二的方法，我们有：

\begin{align} G(x) &= \exp(2B1(x)) \\ B(x) &= \exp(B1(x)) \\ &= \exp(\frac{\ln{G(x)}}{2}) \\ &= \sqrt{G} \end{align}

我们可以对等式两边分别进行求导并比较两边系数，以得到易于编码的递推公式，通过此题。注意到做法二与做法三本质相同，且一般情况下做法三可以得到更优的时间复杂度。

\begin{align} B_n^2 &= G \\ 2B_nB_n' &= G' \end{align}

参考代码

#include <iostream>
using namespace std;

using LL = long long;
constexpr int MOD = int(1e9) + 7;

// <<= '2. Number Theory .,//{
namespace NT {
void INC(int& a, int b) {
  a += b;
  if (a >= MOD) a -= MOD;
}

int sum(int a, int b) {
  a += b;
  if (a >= MOD) a -= MOD;
  return a;
}

void DEC(int& a, int b) {
  a -= b;
  if (a < 0) a += MOD;
}

int dff(int a, int b) {
  a -= b;
  if (a < 0) a += MOD;
  return a;
}

void MUL(int& a, int b) { a = (LL)a * b % MOD; }

int pdt(int a, int b) { return (LL)a * b % MOD; }

int _I(int b) {
  int a = MOD, x1 = 0, x2 = 1, q;
  while (1) {
    q = a / b, a %= b;
    if (!a) return x2;
    DEC(x1, pdt(q, x2));

    q = b / a, b %= a;
    if (!b) return x1;
    DEC(x2, pdt(q, x1));
  }
}

void DIV(int& a, int b) { MUL(a, _I(b)); }

int qtt(int a, int b) { return pdt(a, _I(b)); }

int pow(int a, LL b) {
  int c(1);
  while (b) {
    if (b & 1) MUL(c, a);
    MUL(a, a), b >>= 1;
  }
  return c;
}

template <class T>
T pow(T a, LL b) {
  T c(1);
  while (b) {
    if (b & 1) c *= a;
    a *= a, b >>= 1;
  }
  return c;
}

template <class T>
T pow(T a, int b) {
  return pow(a, (LL)b);
}

struct Int {
  int val;

  operator int() const { return val; }

  Int(int _val = 0) : val(_val) {
    val %= MOD;
    if (val < 0) val += MOD;
  }

  Int(LL _val) : val(_val) {
    _val %= MOD;
    if (_val < 0) _val += MOD;
    val = _val;
  }

  Int& operator+=(const int& rhs) {
    INC(val, rhs);
    return *this;
  }

  Int operator+(const int& rhs) const { return sum(val, rhs); }

  Int& operator-=(const int& rhs) {
    DEC(val, rhs);
    return *this;
  }

  Int operator-(const int& rhs) const { return dff(val, rhs); }

  Int& operator*=(const int& rhs) {
    MUL(val, rhs);
    return *this;
  }

  Int operator*(const int& rhs) const { return pdt(val, rhs); }

  Int& operator/=(const int& rhs) {
    DIV(val, rhs);
    return *this;
  }

  Int operator/(const int& rhs) const { return qtt(val, rhs); }

  Int operator-() const { return MOD - *this; }
};

}  // namespace NT

using namespace NT;

constexpr int N = int(1e3) + 9;
Int binom[N][N], C[N], E[N], B[N], B1[N], G[N];
int n;

void ln(Int C[], Int G[]) {
  for (int i = 1; i <= n; ++i) {
    C[i] = G[i];
    for (int j = 1; j <= i - 1; ++j)
      C[i] -= binom[i - 1][j - 1] * C[j] * G[i - j];
  }
}

void exp(Int G[], Int C[]) {
  for (int i = 1; i <= n; ++i) {
    G[i] = C[i];
    for (int j = 1; j <= i - 1; ++j)
      G[i] += binom[i - 1][j - 1] * C[j] * G[i - j];
  }
}

int main() {
  cin.tie(nullptr)->sync_with_stdio(false);

  n = 1000;
  for (int i = 0; i < n + 1; ++i) {
    binom[i][0] = 1;
    for (int j = 0; j < i; ++j)
      binom[i][j + 1] = binom[i - 1][j] + binom[i - 1][j + 1];
  }

  for (int i = 1; i <= n; ++i) G[i] = pow(2, binom[i][2]);
  ln(C, G);
  for (int i = 1; i <= n; ++i) G[i] = pow(2, binom[i - 1][2]);
  ln(E, G);
  for (int i = 1; i <= n; ++i) {
    G[i] = 0;
    for (int j = 0; j < i + 1; ++j) G[i] += binom[i][j] * pow(2, j * (i - j));
  }
  ln(B1, G);
  for (int i = 1; i <= n; ++i) B1[i] /= 2;
  exp(B, B1);

  int T;
  cin >> T;
  while (T--) {
    cin >> n;
    cout << "Connected: " << C[n] << '\n'
         << "Eulerian: " << E[n] << '\n'
         << "Bipartite: " << B[n] << "\n\n";
  }
}

习题

Riddell's Formula

上述关于 EGF 的 exp 的用法，有时又被称作 Riddell's formula for labeled graphs，生成函数的欧拉变换，有时也被称为 Riddell's formula for unlabeled graphs，后者最早出现在欧拉对分拆数的研究中，除了解决图论计数问题之外，也在完全背包问题中出现。

对于给定序列，和对应的 OGF ，定义的欧拉变换为：

\begin{align} \mathcal{E}(A(x)) &= \prod_{i} (1-x^i)^{-a_i} \\ &= \exp (\sum_{i} \frac{A(x^i)}{i}) \end{align}

设 $\mathcal{E}(A(x))$ 的各项系数为，定义辅助数组 $c_i = \sum_{d|n} d a_d$ ，则有递推公式

n b_n = c_n + \sum_{i=1}^{n-1} c_i b_{n-i}

无标号树

例题「SPOJ PT07D」Let us count 1 2 3

题目大意：求有 n 个结点的分别满足下列性质的树的方案数。

有标号有根树 A000169。
有标号无根树 A000272。
无标号有根树 A000081。
无标号无根树 A000055。

有根树

有标号情况以在前文中解决，下面考察无标号有根树，设其 OGF 为，应用欧拉变换，可得：

F(x) = x\mathcal{E}(F(x))

取出系数即可。

无根树

考虑容斥，我们用有根树的方案中减去根不是重心的方案，并对的奇偶性进行讨论。

当是奇数时：

必然存在一棵子树大小 $\geq \left\lceil \frac{n}{2}\right\rceil$ ，枚举这棵子树的大小有。

g_n = f_n - \sum_{i=\left\lceil\frac{n}{2}\right\rceil}^{n-1} f_i f_{n-i}

当是偶数时：

注意到当有两个重心的情况时，上面的过程只会减去一次，因此还需要减去

g_n = f_n - \sum_{i=\left\lceil\frac{n}{2}\right\rceil}^{n-1} f_i f_{n-i} - \binom{f_{\frac{n}{2}}}{2}

例题「Luogu P5900」无标号无根树计数

题目大意：求有 n 个结点的无标号无根树的方案数（ $n \leq 200000$ ）。

对于数据范围更大的情况，做法同理，欧拉变换后使用多项式模板即可。

无标号简单图

例题「SGU 282. Isomorphism」Isomorphism

题目大意：求有 n 个结点的无标号完全图的边进行 m 染色的方案数。

注意到当 m = 2 时，所求对象就是无标号简单图 A000088，考察波利亚计数定理，

\frac{1}{|G|}\sum_{g\in G} m^{c(g)}

本题中置换群为顶点的阶对称群生成的边集置换群，但暴力做法的枚举量为，无法通过此题。

考虑根据按照置换的循环结构进行分类，每种循环结构对应一种数的分拆，我们用 dfs() 生成分拆，那么问题即转化为求每一种分拆所对应的置换数目和每一类置换中的循环个数，答案为

\frac{1}{|G|} \sum_{p \in P} w(p) m^{c(p)}

考虑，每一个分拆对应一个循环排列，同时同一种大小的分拆之间的顺序无关，因而我们有：

w(p) = \frac{n!}{\prod_{i}(p_i)\prod_{i}(q_i!)}

这里表示大小为的分拆在中出现的次数。

考虑，所影响的点集的循环即为，但题目考察的是边染色，所以还需要考察点置换所生成的边置换，

如果一条边关联的顶点处在同一个循环内，设该循环大小为，那么边所生成的循环数恰好为 $\left\lfloor \frac{p_i}{2} \right\rfloor$ 。

如果一条边关联的顶点处在两个不同的循环中，设分别为 ,，每个循环节的长度均为 $\operatorname{lcm}(p_i,p_j)$ ，因而边所生成的循环数恰好为 $\frac{p_i p_j}{\operatorname{lcm}(p_i,p_j)} = \gcd(p_i, p_j)$ 。

参考代码

#include <iostream>
#include <vector>
using namespace std;

using LL = long long;
int MOD = int(1e9) + 7;

namespace NT {
void INC(int& a, int b) {
  a += b;
  if (a >= MOD) a -= MOD;
}

int sum(int a, int b) {
  a += b;
  if (a >= MOD) a -= MOD;
  return a;
}

void DEC(int& a, int b) {
  a -= b;
  if (a < 0) a += MOD;
}

int dff(int a, int b) {
  a -= b;
  if (a < 0) a += MOD;
  return a;
}

void MUL(int& a, int b) { a = (LL)a * b % MOD; }

int pdt(int a, int b) { return (LL)a * b % MOD; }

int _I(int b) {
  int a = MOD, x1 = 0, x2 = 1, q;
  while (1) {
    q = a / b, a %= b;
    if (!a) return x2;
    DEC(x1, pdt(q, x2));

    q = b / a, b %= a;
    if (!b) return x1;
    DEC(x2, pdt(q, x1));
  }
}

void DIV(int& a, int b) { MUL(a, _I(b)); }

int qtt(int a, int b) { return pdt(a, _I(b)); }

int pow(int a, LL b) {
  int c(1);
  while (b) {
    if (b & 1) MUL(c, a);
    MUL(a, a), b >>= 1;
  }
  return c;
}

template <class T>
T pow(T a, LL b) {
  T c(1);
  while (b) {
    if (b & 1) c *= a;
    a *= a, b >>= 1;
  }
  return c;
}

template <class T>
T pow(T a, int b) {
  return pow(a, (LL)b);
}

struct Int {
  int val;

  operator int() const { return val; }

  Int(int _val = 0) : val(_val) {
    val %= MOD;
    if (val < 0) val += MOD;
  }

  Int(LL _val) : val(_val) {
    _val %= MOD;
    if (_val < 0) _val += MOD;
    val = _val;
  }

  Int& operator+=(const int& rhs) {
    INC(val, rhs);
    return *this;
  }

  Int operator+(const int& rhs) const { return sum(val, rhs); }

  Int& operator-=(const int& rhs) {
    DEC(val, rhs);
    return *this;
  }

  Int operator-(const int& rhs) const { return dff(val, rhs); }

  Int& operator*=(const int& rhs) {
    MUL(val, rhs);
    return *this;
  }

  Int operator*(const int& rhs) const { return pdt(val, rhs); }

  Int& operator/=(const int& rhs) {
    DIV(val, rhs);
    return *this;
  }

  Int operator/(const int& rhs) const { return qtt(val, rhs); }

  Int operator-() const { return MOD - *this; }
};

}  // namespace NT

using namespace NT;

constexpr int N = int(5e1) + 9;
Int Fact[N];
vector<vector<int>> Partition;
vector<int> cur;
int n, m;

void gen(int n = ::n, int s = 1) {
  if (!n) {
    Partition.push_back(cur);
  } else if (n >= s) {
    cur.push_back(s);
    gen(n - s, s);
    cur.pop_back();
    gen(n, s + 1);
  }
}

Int w(const vector<int> P) {
  Int z = Fact[n];
  int c = 0, l = P.front();

  for (auto p : P) {
    z /= p;
    if (p != l) {
      z /= Fact[c];
      l = p;
      c = 1;
    } else {
      ++c;
    }
  }

  z /= Fact[c];
  return z;
}

int gcd(int x, int y) { return y ? gcd(y, x % y) : x; }

int c(const vector<int> P) {
  int z = 0;
  for (int i = 0; i < P.size(); ++i) {
    z += P[i] / 2;
    for (int j = 0; j < i; ++j) z += gcd(P[i], P[j]);
  }
  return z;
}

int main() {
  cin >> n >> m >> MOD;
  Fact[0] = 1;
  for (int i = 1; i <= n; ++i) Fact[i] = Fact[i - 1] * i;

  gen();

  Int res = 0;
  for (auto P : Partition) {
    res += w(P) * pow(m, c(P));
  }
  res /= Fact[n];
  cout << res << endl;
}

习题

参考资料与注释

也许无标号二叉树是一个反例，在结构简单的情况下，对应的置换群是恒等群（Identity Group），此时有标号版本可以直接通过乘以得到。 ↩
粉兔的 blog 告诉我们，这个序列也可以使用 Chirp Z-Transform 优化。 ↩

本页面最近更新：2024/10/5 17:45:55，更新历史
发现错误？想一起完善？在 GitHub 上编辑此页！
本页面贡献者：ComeIntoCalm, GoodCoder666, lychees, megakite, Tiphereth-A
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用