题解 P4727 【[HNOI2009]图的同构记数】

· · 题解

一道Burnside题的流程大概是两步：

  1. 找出置换群<G,*>和目标集合S。
  2. 求出G中每一个置换在S中的不动点集合大小之和。

记一个置换t(t\in G)在S中的不动点集合是S_t=\{x\in S |x*t=x\}。

此题中，置换群G是所有的阶为n的置换，|G|=n!。目标集合是所有有标号的、点数为n的图的个数，因为有\frac{n(n-1)}{2}条可能的边，所以|S|=2^{\frac{n(n-1)}{2}}。

显然G太大了，不能拿出其中的每一个置换一个个算不动点集合大小。通过对小数据的手动计算可以发现，由于图的对称性，对于任意两个置换t_1,t_2\in G，|S_{t_1}|=|S_{t_2}|的充分条件是t_1,t_2的循环个数相同，而且把各个循环按照大小排序之后也分别相同。

所以，不妨把G中的元素以此为据分类，这样就可以使得每一类中的所有置换t的|S_t|都是相等的。这一类置换组成的集合用G_i表示；也就是说，G_i就是由各个循环分别相等的一些置换组成的一个G的子集。

实际上，每一个G_i中的置换们的循环们的大小，都是整数n的某一个划分——所谓划分就是指把n表示为任意多个整数之和，而且整数是无序的，也就是说如果n=3，那么2+1和1+2是同一种划分。n的划分的总方案数，就是放苹果问题中，苹果数为n而盘子个数从1到n的答案之和。n=60时，总划分数大约是百万级别（语出yyb）。

现在我们通过dfs得到了一个单调不下降的序列{a_1,a_2,...,a_m}，满足1\le m,a_i \le n,\sum_{i=1}^m a_i = n。这就是一个n的划分。这个序列对应了某一个G_i，表示G_i中置换的各个循环的大小。为求出这个G_i中的置换们的不动点集合大小之和，我们要求出两个值：G_i中元素的个数，以及那个相等的|S_t|是多少。

前一个值是：

\frac{n!}{(\prod a_i)(\prod b_x!)}

其中b_x = \sum_{i=1}^m [a_i = x]。

这个的证明可以用传统的排列组合知识解决。

后一个值较为复杂。后一个值要求的就是有多少种连边的方法，使得在图的点编号被置换一次之后，图的各个循环内部以及循环之间仍然是一模一样的。

我们需要再次把问题分成两部分：循环内和循环间。

对于m个循环中的一个大小为a_i的循环，不妨选取其中任意的一个点k作为考虑的核心。把a_i个点按照循环的顺序一字铺开，并把k放在第一个，如果k有一条边连到了它后面距离它为d(d\ge1)的点，那么每一个其他点都需要有一条边连到它后面距离它为d的点。所以从k出发，对于每一个跨度d，可以选择连一条跨度为d的边或者不连；而随即所有其他地方的连边情况也就都被确定了。

到这里，需要发现的一点是，因为边是无向的，所以如果k连出的边中，有某一条边跨度超过了\lfloor \frac{a_i}{2}\rfloor（记这个跨度为d），那么就相当于往前面连了一条跨度为a_i-d的边——故而k往后a_i-d条边的那个点也必定会有一条边连到k。这也就相当于是k往后连了一条跨度为a_i - d的边。综上所述，k往后连一条跨度为a_i-d的边和跨度为d的边是等价的，最后随之确定的整个循环中的边是相同的一些边，所以不能重复计算。

那么，每一个循环内部的方案数就是2^{\lfloor \frac{a_i}{2}\rfloor}。

（这也就是yyb题解中说的“一个n个点的循环只有\frac{n}{2}个边的循环”。）

现在考虑循环之间的连边方案数。从m个循环中挑出两个a_i,a_j来（i < j），考虑它们之间的连边。仍然从a_i循环中拿出一个点k，从a_j循环中拿出一个点l来作为考虑的核心，并把两个循环中的点都一字铺开。

如果k和l连了边，那么a_i循环中k的下一个也要和a_j循环中l的下一个连边，k往后两个的点也要和l往后两个的点连边......这样两边同时沿着各自的循环转，不知道转了多少圈之后，k点一共往a_j循环中连了多少条边？显然，k往a_j循环中每隔gcd(a_i,a_j)条边的距离就会连一条边过去，一共连了\frac{a_j}{gcd(a_i,a_j)}条。这样来说，a_j循环中的a_j个点，一共被分为了gcd(a_i,a_j)组，每一组包括了\frac{a_j}{gcd(a_i,a_j)}个点，这些点中，每相邻的两个点之间的距离都是gcd(a_i,a_j)。而且，每一组要么一起向k连边，要么都不向k连边。那么总连边方案自然就是2^{gcd(a_i,a_j)}。

综上所述，对于一个划分a_1,a_2,...a_m，这个划分对应的G_i中元素的共同的|S_t|值，就是：

\prod_{i=1}^m 2^{\lfloor \frac{a_i}{2}\rfloor} \prod _{i=1}^{m}\prod_{j=i+1}^m 2^{gcd(a_i,a_j)}

再把这个值乘上G_i的元素的个数：

\frac{n!}{(\prod a_i)(\prod b_x!)}\prod_{i=1}^m 2^{\lfloor \frac{a_i}{2}\rfloor} \prod _{i=1}^{m}\prod_{j=i+1}^m 2^{gcd(a_i,a_j)}

就是这个划分a_1,a_2,...a_m的答案了。

把所有划分的答案加起来之后，根据Burnside定理的公式，还要乘上一个\frac{1}{n!}才是答案。

而这个\frac{1}{n!}可以与上面式子中的n!约掉。所以如果你够懒，就直接算各个划分的这个式子的和：

\frac{1}{(\prod a_i)(\prod b_x!)}\prod_{i=1}^m 2^{\lfloor \frac{a_i}{2}\rfloor} \prod _{i=1}^{m}\prod_{j=i+1}^m 2^{gcd(a_i,a_j)}

把所有gcd(i,j)预处理出来会快一点。

/*
几个错误调试了好久。
*///using CRLF, UTF-8

#include <cstdio>
#include <algorithm>
#include <cstring>
#define pr printf
#define F(i, j, k) for(register int i = j, kkllkl = k; i <= kkllkl; ++i)
#define G(i, j, k) for(register int i = j, kkllkl = k; i >= kkllkl; --i)
#define clr(a, v) memset((a), v, sizeof(a))
using namespace std;
typedef long long ll;

#define isd(x) (('0' <= (x) && (x) <= '9') || (x) == '-')
int rd() {
    int ans = 0, sign = 1; char c = getchar();
    while(!isd(c)) c = getchar();
    if(c == '-') sign = -1, c = getchar();
    while(isd(c)) ans = (ans << 3) + (ans << 1) + c - '0', c = getchar();
    return sign == 1 ? ans : -ans;
}

#define chkmin(x, y) ((x) = min((x), (y)))
#define chkmax(x, y) ((x) = max((x), (y)))

/* --------------------------------------------------- */

const int N = 65, MOD = 997;
int Gcd(int x, int y) { return y ? Gcd(y, x % y) : x; }
int fsp(int a, int x) {
    int ans = 1;
    for(; x; x >>= 1, a = a * a % MOD)
        if(x & 1) ans = ans * a % MOD;
    return ans;
}

int n, frac[N], ifrac[N], gcd[N][N], p2[N], inv[N];

int a[N], m, b[N];
int solve() {
    // F(i, 1, m) pr("%d%c", a[i], " \n"[i == m]);
    int ret = 1;
    F(i, 1, m) ret = (ret * inv[a[i]]) % MOD;
    F(i, 1, n) ret = (ret * ifrac[b[i]]) % MOD;
    F(i, 1, m) ret = (ret * p2[a[i] / 2]) % MOD;
    F(i, 1, m) F(j, i + 1, m) ret = (ret * p2[gcd[a[i]][a[j]]]) % MOD;
    // pr("ret = %d\n", ret);
    return ret;
}

int ans = 0;
void dfs(int level, int sum) {
    if(sum == n) {
        m = level - 1;
        ans = (ans + solve()) % MOD;
        return ;
    }
    F(i, max(1, a[level - 1]), n - sum) {
        a[level] = i;
        ++b[i];
        dfs(level + 1, sum + i);
        --b[i];
    }
}

int main() {
    n = rd();
    F(i, 1, n) inv[i] = fsp(i, MOD - 2);
    frac[0] = 1;
    F(i, 1, n) frac[i] = frac[i - 1] * i % MOD;
    ifrac[n] = fsp(frac[n], MOD - 2);
    G(i, n - 1, 0) ifrac[i] = ifrac[i + 1] * (i + 1) % MOD;
    F(i, 1, n) F(j, 1, n) gcd[i][j] = Gcd(i, j);
    p2[0] = 1;
    F(i, 1, n) p2[i] = (p2[i - 1] * 2) % MOD;
    dfs(1, 0);
    pr("%d\n", ans);
    return 0;
}
/*
*/