【P2109 [NOI2007]】 生成树计数 题解（状压 dp + 矩阵快速幂）
状压 dp + 矩阵快速幂优化。
感觉题解区几篇题解说得云里雾里的……也没有一定的证明……
Update
- 2023年6月10日：贴上代码，修改错别字。
Solution
Part1 dp 状态设计
发现 \(k\) 的范围很小，具有很强的指向性——状压 dp。
结合题面，发现一个点 \(i\) 只可能和 \(j\in (i-k,i)\) 的点连边，即 \(i-k+1\) 之前的点我们根本不用考虑它和 \(i\) 的连接情况，我们只需要保证 \([1,i-k]\) 这些点连通即可。
故可以设计状态 \(f(i,sta)\) 表示当前考虑到点 \(i\)，且 \((i-k,i]\) 这 \(k\) 个点的连通状态为 \(sta\) 时，生成树的方案数。
因为最后要求的是生成树，而生成树的必要条件是 \([1,n]\) 都是连通的。所以我们每次从 \(f(i-1,sta')\) 转移到 \(f(i, sta)\) 时，都要保证 \(i-k\) 和 \(j\in (i-k,i]\) 中的任何一个点是连通的。
为什么这样就可以保证 \(f(n,sta_{final})\) 的方案数（其中 \(sta_{final}\) 代表 \((n-k,n]\) 都连通的状态）都是合法的（即 \([1,n]\) 都是连通的）？
此处口胡一个反证法。易知，若 \([1,n]\) 不连通，即有若干个连通块，那我们一定可以找到某个连通块中，编号最大的节点 \(i\)，满足它和 \((i,i+k]\) 一定都不连通。故我们如果保证了转移时每一个 \(i-k\) 和 \(j\in (i-k,i]\) 中的任何一个点是连通的，那这个图最终一定就是一个连通图。
Part2 状态压缩
然后考虑状态中 \(sta\) 的表述。
其余众多题解都提及了，我们可以使用最小表示法，将每种 \(k\) 点连接状态都表示出来。具体地，给同一个连通块中的点都赋上同一数值即可，不同连通块数值不同即可。
对于每种连接状态，预处理 dfs() 出来即可，状态总数不超 \(60\)。
Part3 dp 状态转移
考虑状态转移。
记 \(g[sta'][sta]\) 表示从状态 \(sta'\) 转移到 \(sta\) 的方案数。即，对于一个待转移的点 \(i+1\)，它可以通过 \(g[sta'][sta]\) 种方案，和 \(j\in [i - k+2,i]\) 的点连边，使得原先 \([i-k+1,i]\) 这些点构成的状态 \(sta'\)，能转变为现在 \((i-k+1,i+1]\) 这些点构成的状态 \(sta'\)。
如何保证最终答案的合法性？易知生成树一定是一个无环的连通图。连通性的保证在 Part1 中已经讲解了，现在保证图中一定不出现环。所以我们在两两状态转移的时候一定要注意不能出现环，否则当前状态不合法。
那么存在一个十分显然的转移：
然后仔细观察，会发现这个形式实际上是矩阵乘法的形式。故我们在求出了 \(g[][]\) 之后，可以使用矩阵快速幂求出最终方案数。
进一步地，不难得出可以省去 \(f\) 的第一维。因矩阵快速幂的使用，\(f\) 可简化为 \(f(1,sta)\)。
Part4 \(f\) 数组初始化
注意，\(f\) 不能直接采用 dfs() 中暴力推出来的状态种类进行初始化。
原因是不同的连边方式可能导致同样的连通状态。
故，我们这样初始化 \(f\)：对于前 \(k\) 个点，一共有 \(\dbinom{k}{2}\) 条边，而对于每条边，我们可以选择连或不连，故一共有 \(2^{\frac{k\times (k-1)}{2}}\) 种连边方式。对于每种连边方式，我们直接暴力连边，使用并查集维护连通性，即可得到最终连通状态。
总复杂度可过。
Code
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
#define rep(i, a, b) for(int i = a; i <= b; ++i)
#define Rep(i, a, b) for(int i = a; i < b; ++i)
const int mod = 65521;
ll n;
int id[600005], st[53], p[25], v[25];
int fa[25], tot, k;
struct mtr{
ll a[53][53];
inline mtr(){
memset(a, 0, sizeof a);
}
inline ll* operator [](int x){ return a[x];}
inline mtr operator *(mtr b) const{
mtr c;
rep(i, 1, tot) rep(j, 1, tot) rep(k, 1, tot)
c[i][k] += a[i][j] * b[j][k];
rep(i, 1, tot) rep(j, 1, tot) c[i][j] %= mod;
return c;
}
}f, tr;
inline void dfs(int x){
if(x > k){
int cnt, sta; cnt = sta = 0;
rep(i, 1, k) v[i] = 0;
rep(i, 1, k) if(!v[p[i]]) v[p[i]] = ++cnt;
rep(i, 1, k) sta = sta * 10 + v[p[i]];
if(!id[sta]) id[sta] = ++tot, st[tot] = sta;
return;
} rep(i, 1, k) p[x] = i, dfs(x + 1);
}
inline int fnd(int x){
return fa[x] == x ? x : fa[x] = fnd(fa[x]);
}
inline void init(int s){
rep(i, 1, k) fa[i] = i;
for (int i = 1, c = 0; i <= k; ++i)
for (int j = 1; j < i; ++j, ++c)
if(s & (1 << c)){
if(fnd(i) == fnd(j)) return;
fa[fnd(i)] = fnd(j);
}
int cnt = 0, sta = 0;
rep(i, 1, k) v[i] = 0, p[i] = fnd(i);
rep(i, 1, k) if(!v[p[i]]) v[p[i]] = ++cnt;
rep(i, 1, k) sta = sta * 10 + v[p[i]];
++f[1][id[sta]];
}
inline void add(int fr, int s){
rep(i, 1, k) fa[i] = i;
fa[k + 1] = k + 1, v[k + 1] = 0;
for(int st = fr, i = k; i; --i) p[i] = st % 10 + k, st /= 10;
rep(i, 1, k) rep(j, i + 1, k)
if(p[i] == p[j] and fnd(i) != fnd(j))
fa[fnd(j)] = fnd(i);
rep(i, 1, k) if(s & (1 << (i - 1))){
if(fnd(i) == fnd(k + 1)) return;
fa[fnd(k + 1)] = fnd(i);
} int flg = 0;
rep(i, 1, k) if(fnd(1) == fnd(i + 1)) flg = 1;
if(!flg) return;
int cnt = 0, sta = 0;
rep(i, 1, k) v[i] = 0, p[i] = fnd(i + 1);
rep(i, 1, k) if(!v[p[i]]) v[p[i]] = ++cnt;
rep(i, 1, k) sta = sta * 10 + v[p[i]];
++tr[id[fr]][id[sta]];
}
int main(){
scanf("%d%lld", &k, &n), n -= k, dfs(1);
int s = 1 << (k * (k - 1) >> 1);
Rep(i, 0, s) init(i);
s = 1 << k;
rep(i, 1, tot) Rep(j, 0, s) add(st[i], j);
for(; n; n /= 2, tr = tr * tr) if(n & 1)
f = f * tr;
printf("%lld\n", f[1][1]);
return 0;
}
Thanks for reading.
浙公网安备 33010602011771号