图论计数
在组合数学中,图论计数(Graph Enumeration)是研究满足特定性质的图的计数问题的分支。生成函数、波利亚计数定理 与 符号化方法 和 OEIS 是解决这类问题时最重要的数学工具。图论计数可分为有标号和无标号两大类问题,大多数情况下1有标号版本的问题都比其对应的无标号问题更加简单,因此我们将先考察有标号问题的计数。
有标号树
即 Cayley 公式,参见 Prüfer 序列 一文,我们也可以使用 Kirchhoff 矩阵树定理 或 生成函数 和 拉格朗日定理 得到这一结果。
习题
有标号连通图
例题「POJ 1737」Connected Graph
例题 「POJ 1737」Connected Graph
题目大意:求有
这类问题最早出现于楼教主的男人八题系列中,我们设
移项得到
例题「集训队作业 2013」城市规划
例题 「集训队作业 2013」城市规划
题目大意:求有
对于数据范围更大的序列问题,往往我们需要构造这些序列的生成函数,以使用高效的多项式算法。
方法一:分治 FFT
上述的递推式可以看作一种自卷积形式,因而可以使用分治 FFT 进行计算,复杂度
方法二:多项式求逆
我们将上述递推式中的组合数展开,并进行变形:
构造多项式:
代换进上式得到
方法三:多项式 exp
另一种做法是使用 EGF 中多项式 exp 的组合意义,我们设有标号连通图和简单图序列的 EGF 分别为
使用 多项式 ln 解出
有标号欧拉图、二分图
例题「SPOJ KPGRAPHS」Counting Graphs
例题 「SPOJ KPGRAPHS」Counting Graphs
题目大意:求有
本题限制代码长度,因而无法直接使用多项式模板,但生成函数依然可以帮助我们进行分析。
连通图问题在之前的例题中已被解决,考虑欧拉图。注意到上述对连通图计数的几种方法,均可以在满足任意性质的有标号连通图进行推广。例如我们可以将连通图递推公式中的
我们将 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];
}
}
下面讨论有标号二分图计数,
我们设
接下来我们用两种不同的方法建立
方法一:算两次
我们设
比较两种
不难得到
方法二:连通化递推
方法二和方法三均使用连通二分图
注意到对于每个连通二分图,我们恰好有两种不同的染色方法,对应到两组不同的连通 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 理解上面的递推过程。
设
我们可以对等式两边分别进行求导并比较两边系数,以得到易于编码的递推公式,通过此题。 注意到做法二与做法三本质相同,且一般情况下做法三可以得到更优的时间复杂度。
参考代码
#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";
}
}
习题
- UOJ Goodbye Jihai D. 新年的追逐战
- BZOJ 3864. 大朋友和多叉树
- BZOJ 2863. 愤怒的元首
- Luogu P6295. 有标号 DAG 计数
- LOJ 6569. 仙人掌计数
- LOJ 6570. 毛毛虫计数
- Luogu P5434. 有标号荒漠计数
- Luogu P3343. [ZJOI2015] 地震后的幻想乡
- HDU 5279. YJC plays Minecraft
- Luogu P7364. 有标号二分图计数
- Luogu P5827. 点双连通图计数
- Luogu P5827. 边双连通图计数
- Luogu P6596. How Many of Them
- Luogu U152448. 有标号强连通图计数
- Project Euler 434. Rigid graphs
Riddell's Formula
上述关于 EGF 的 exp 的用法,有时又被称作 Riddell's formula for labeled graphs,生成函数的 欧拉变换,有时也被称为 Riddell's formula for unlabeled graphs,后者最早出现在欧拉对分拆数的研究中,除了解决图论计数问题之外,也在完全背包问题中出现。
对于给定序列
设
无标号树
例题「SPOJ PT07D」Let us count 1 2 3
例题 「SPOJ PT07D」Let us count 1 2 3
题目大意:求有 n 个结点的分别满足下列性质的树的方案数。
有根树
有标号情况以在前文中解决,下面考察无标号有根树,设其 OGF 为
取出系数即可。
无根树
考虑容斥,我们用有根树的方案中减去根不是重心的方案,并对
当
必然存在一棵子树大小
当
注意到当有两个重心的情况时,上面的过程只会减去一次,因此还需要减去
例题「Luogu P5900」无标号无根树计数
例题 「Luogu P5900」无标号无根树计数
题目大意:求有 n 个结点的无标号无根树的方案数(
对于数据范围更大的情况,做法同理,欧拉变换后使用多项式模板即可。
无标号简单图
例题「SGU 282. Isomorphism」Isomorphism
例题 「SGU 282. Isomorphism」Isomorphism
题目大意:求有 n 个结点的无标号完全图的边进行 m 染色的方案数。
注意到当 m = 2 时,所求对象就是无标号简单图 A000088,考察波利亚计数定理,
本题中置换群
考虑根据按照置换的循环结构进行分类,每种循环结构对应一种数的分拆,我们用 dfs() 生成分拆,那么问题即转化为求每一种分拆
考虑
这里
考虑
如果一条边关联的顶点处在同一个循环内,设该循环大小为
如果一条边关联的顶点处在两个不同的循环中,设分别为
参考代码
#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;
}
习题
- CodeForces 438 E. The Child and Binary Tree
- Luogu P5448. [THUPC2018] 好图计数
- Luogu P5818. [JSOI2011] 同分异构体计数
- Luogu P6597. 烯烃计数
- Luogu P6598. 烷烃计数
- Luogu P4128. [SHOI2006] 有色图
- Luogu P4727. [HNOI2009] 图的同构计数
- AtCoder Beginner Contest 222 H. Binary Tree
- AtCoder Beginner Contest 284 Ex. Count Unlabeled Graphs
- Luogu P4708. 画画
- Luogu P7592. 数树(2021 CoE-II E)
- Luogu P5206. [WC2019] 数树
参考资料与注释
- WC2015, 顾昱洲营员交流资料 Graphical Enumeration
- WC2019, 生成函数,多项式算法与图的计数
- Counting labeled graphs - Algorithms for Competitive Programming
- Graphical Enumeration Paperback, Frank Harary, Edgar M. Palmer
- The encyclopedia of integer sequences, N. J. A. Sloane, Simon Plouffe
- Combinatorial Problems and Exercises, László Lovász
- Graph Theory and Additive Combinatorics
-
也许无标号二叉树是一个反例,在结构简单的情况下,对应的置换群是恒等群(Identity Group),此时有标号版本可以直接通过乘以
得到。 ↩ -
粉兔的 blog 告诉我们,这个序列也可以使用 Chirp Z-Transform 优化。 ↩
创建日期: 2023年5月27日