容斥和反演 2
复杂容斥系数
容斥原理的核心是，通过设置一系列容斥系数，使得每种合法方案被恰好统计 \(1\) 次，不合法方案被恰好统计 \(0\) 次。
如果我们计算一个集合 \(T\)，就会把所有 \(T_0\supseteq T\)或 \(T_0\subseteq T\)）的 \(T_0\) 统计一遍，那么我们可以列一个表，手动递推容斥系数。
给定 \(n\) 个平面上的等腰直角三角形，三点分别为 \((x_i,y_i),(x_i+r_i,y_i),(x_i,y_i+r_i),~(r_i>0)\)。询问恰好被奇数个三角形覆盖的部分的面积。
\(n\le 10\)
暴力容斥不难做到 \(O(3^n)\)。这里讲如何手动递推容斥系数并做到 \(O(n2^n)\)。
设 \(S\) 为一个三角形的集合，\(f_S\) 为 \(S\) 内三角形交集的面积。设 \(w_i\) 为 \(f_S,~|S|=i\) 的容斥系数。那么一块恰好被集合 \(S_0\) 覆盖的三角形会被统计：
\[
\sum_{\varnothing\subset S\subseteq S_0}w_{|S|}=\sum_{i=1}^{|S_0|}\binom{|S|}{i}w_i
\]
而我们希望上式满足
\[
\sum_{\varnothing\subset S\subseteq S_0}w_{|S|}=\sum_{i=1}^{|S_0|}\binom{|S|}{i}w_i=tar_{|S|}=|S|\bmod 2
\]
注意到：
- \(w_1\) 满足 \(w_1=tar_1=1\)；
- \(w_2\) 满足 \(2w_1+w_2=tar_2=0\)，解得 \(w_2=-2\)；
- \(w_3\) 满足 \(3w_1+3w_2+w_3=tar_3=1\)，解得 \(w_3=4\)；
- \(w_4\) 满足 \(4w_1+6w_2+4w_3+w_4=tar_4=0\)，解得 \(w_4=-8\)；
- ...
对于这种容斥，手动计算前 \(4\) 项大部分可以发现规律，然后再尝试证明即可。我们可以列出一个表格，然后手动递推即可。这种方法兼顾严谨、快速、通用。
| \(\vert S\vert\) |
\(1\) |
\(2\) |
\(3\) |
\(4\) |
\(5\) |
... |
| \(tar\) |
\(1\) |
\(0\) |
\(1\) |
\(0\) |
\(1\) |
... |
| \(w\) |
\(1\) |
\(-2\) |
\(4\) |
\(-8\) |
\(16\) |
... |
时间复杂度 \(O(n2^n)\)。
代码
| #include<iostream>
#include<iomanip>
#define int long long
using namespace std;
const int N = 12;
int n, ans;
int px[N], py[N], pr[N];
int cnt[1 << N], w[N];
int calc(int s) {
int x = 0, y = 0, r = 3e6 + 10;
for(int i = 0; i < n; i++) {
if(s >> i & 1) {
int nx = max(x, px[i]), ny = max(y, py[i]);
int nr = min(x + y + r - nx - ny, px[i] + py[i] + pr[i] - nx - ny);
x = nx, y = ny, r = nr;
if(r < 0) return 0;
}
}
return r * r;
}
signed main() {
cin >> n;
for(int i = 0; i < n; i++) cin >> px[i] >> py[i] >> pr[i];
w[1] = 1;
for(int i = 2; i <= n; i++) w[i] = w[i - 1] * (-2);
for(int i = 1; i < (1 << n); i++) {
cnt[i] = cnt[i >> 1] + (i & 1);
ans += w[cnt[i]] * calc(i);
}
cout << fixed << setprecision(1) << (double)ans / 2 << '\n';
return 0;
}
|
反射容斥
在网格图上，你需要从 \((0,0)\) 走到 \((2n,0)\)，每次你可以向右上或者右下走一步，不能走到 \(y\) 轴下面，也不能走到 \(y=m\) 上面（两者都可以接触）。问方案数。
设两条直线，一条叫 \(a\)，一条叫 \(b\)。维护一个序列，当你越过某条直线的时候，在序列末尾添加 a 或者 b。那么每种不合法的方案都对应这样一个字符串，例如 aabab，a 等。将多个连续的相同字符缩成一个，这样就只有形如 ababa.. 这样以 \(a\) 开头的和 babab.. 这样以 \(b\) 开头的字符串。
考虑如何对一个字符串计数对应的路径数。可以使用经典 trick 在碰到直线的时候把整张地图都沿着直线对称，这样最终走下来会得到一个显然在合法范围之外的目的地镜像，每个走到这个位置的路径都对应一个非法序列包含给定字符串的路径。
注意到会算重，例如减去 a 和 b 时，会把 ab 和 ba 分别减去两遍，所以需要加回来。那么可以递推得到一个容斥系数：
* 1
a -1
b -1
ab 1
ba 1
aba -1
bab -1
abab 1
baba 1
...
递推到最后，因为目的地镜像太远，走到它的路径数为 \(0\)，就可以停止了。时间复杂度应该为 \(O(\frac{n}{d})\)，\(n\) 是到目的地的距离，\(d\) 是两条直线之间的距离。
代码
| // 260224 T1
#include<iostream>
using namespace std;
typedef long long ll;
const int N = 2e6 + 10;
const int MOD = 1e9 + 7;
inline void add(int &a, int b) { a += b; (a >= MOD) && (a -= MOD); }
int T, n, m, ans;
int inv[N], fac[N], ifac[N];
inline int C(int a, int b) { return (ll)fac[a] * ifac[a - b] % MOD * ifac[b] % MOD; }
int main() {
#ifndef LOCAL
freopen("counting.in", "r", stdin);
freopen("counting.out", "w", stdout);
#endif
inv[1] = fac[0] = ifac[0] = 1;
for(int i = 2; i < N; i++) inv[i] = MOD - (ll)inv[MOD % i] * (MOD / i) % MOD;
for(int i = 1; i < N; i++) fac[i] = (ll)fac[i - 1] * i % MOD;
for(int i = 1; i < N; i++) ifac[i] = (ll)ifac[i - 1] * inv[i] % MOD;
cin >> T >> T;
while(T--) {
cin >> n >> m;
ans = C(2 * n, n);
for(int y = 0; abs(y) <= 2 * n; ) {
y = -2 - y;
if(abs(y) > 2 * n) break;
add(ans, MOD - C(2 * n, n + y / 2));
y = 2 * m + 2 - y;
if(abs(y) > 2 * n) break;
add(ans, C(2 * n, n + y / 2));
}
for(int y = 0; abs(y) <= 2 * n; ) {
y = 2 * m + 2 - y;
if(abs(y) > 2 * n) break;
add(ans, MOD - C(2 * n, n + y / 2));
y = -2 - y;
if(abs(y) > 2 * n) break;
add(ans, C(2 * n, n + y / 2));
}
cout << ans << '\n';
}
return 0;
}
|
给定一棵 \(n\) 个节点的无根树。称两条边是相邻的，当且仅当它们有一个公共端点。
考虑一个 dfs 的过程：初始选择一条边 \(s\) 作为起点，每次找到与当前边相邻且未被访问过的边，递归进入访问。如果不存在，则回溯。每次递归儿子时，都在两条边之间连一条“树边”，最后得到一棵以边为节点的搜索树。搜索树是无根且不区分儿子顺序的。
给定 \(k\) 条关键边，问：可以从某一条关键边出发得到的搜索树 有多少个。答案对 \(10^9+7\) 取模。
\(n\le 10^5\)
首先发现一个点的所有邻边在搜索树上形成一条链，并且顺序（几乎）任意排列。先考虑从一条关键边出发会得到什么样的搜索树。不难发现每个节点的邻边中，必须有一个待在链头，因此答案为 \(\prod (deg_i-1)!\)。
再考虑多条关键边的情况。发现从不同的关键边出发，可能得到同一棵搜索树。而给定一棵搜索树，合法的起始边会形成一条链（也就是一个固定的集合 \(S_0\)，任选一个子集都能统计到它）。因此如果我们直接枚举关键边的集合 \(S\)，再计数能同时被它们得到的搜索树数量，则会算重。因此考虑容斥，对于关键边集合 \(S\)，设 \(g_S\) 表示能被关键边集合 \(S\) 中的每条边同时得到的搜索树数量。给 \(g_S\) 赋予 \((-1)^{|S|+1}\) 的容斥系数，即可使每棵合法的搜索树恰好被统计 \(1\) 次：
| \(\vert S\vert\) |
\(0\) |
\(1\) |
\(2\) |
\(3\) |
\(4\) |
\(5\) |
... |
| \(tar\) |
\(0\) |
\(1\) |
\(1\) |
\(1\) |
\(1\) |
\(1\) |
... |
| \(w\) |
\(0\) |
\(1\) |
\(-1\) |
\(1\) |
\(-1\) |
\(1\) |
... |
至于如何求 \(g_S\)，我们发现链上的点的贡献变为了 \((deg_i-2)!\)。因此对于链上的点都给答案乘以 \(\frac{1}{deg_i-1}\) 即可。这里可以用树上 dp 维护。
代码
| #include<iostream>
#define int long long
using namespace std;
const int N = 1e5 + 10;
const int MOD = 1e9 + 7;
struct Edge {
int v, next;
} pool[2 * N];
int ne = 1, head[N];
void addEdge(int u, int v) {
pool[++ne] = {v, head[u]};
head[u] = ne;
}
int inv[N], fact[N];
int T;
int n, k, ans;
int deg[N], dep[N], key[N];
int f[N];
void dfs1(int u, int fa) {
dep[u] = dep[fa] + 1;
for(int e = head[u]; e; e = pool[e].next) {
int v = pool[e].v;
if(v == fa) continue;
dfs1(v, u);
}
}
void dfs2(int u, int fa) {
f[u] = 0;
for(int e = head[u]; e; e = pool[e].next) {
int v = pool[e].v;
if(v == fa) continue;
dfs2(v, u);
if(!key[v]) {
ans = (ans + f[u] * f[v] % MOD) % MOD;
f[u] = (f[u] + f[v] * inv[deg[u] - 1] % MOD) % MOD;
} else {
ans = (ans + MOD - f[v]) % MOD;
ans = (ans + MOD - 1) % MOD;
ans = (ans + MOD - f[u]) % MOD;
f[u] = (f[u] + MOD - inv[deg[u] - 1]) % MOD;
}
}
}
void clear() {
for(int i = 0; i <= n + 5; i++) head[i] = key[i] = deg[i] = 0;
ne = 1; ans = 0;
}
signed main() {
inv[1] = fact[0] = 1;
for(int i = 2; i < N; i++) inv[i] = MOD - inv[MOD % i] * (MOD / i) % MOD;
for(int i = 1; i < N; i++) fact[i] = fact[i - 1] * i % MOD;
cin >> T >> T;
while(T--) {
cin >> n >> k;
clear();
for(int i = 1; i <= n - 1; i++) {
int u, v;
cin >> u >> v;
++deg[u], ++deg[v];
addEdge(u, v);
addEdge(v, u);
}
dfs1(1, 0);
for(int i = 1; i <= k; i++) {
int e, u, v;
cin >> e;
u = pool[2 * e].v, v = pool[2 * e + 1].v;
if(dep[v] < dep[u]) swap(u, v);
key[v] = 1;
}
dfs2(1, 0);
ans = MOD - ans;
for(int i = 1; i <= n; i++) ans = ans * fact[deg[i] - 1] % MOD;
cout << ans << '\n';
}
return 0;
}
|
DAG 子图计数
给定一个有向图，请问有多少个边集是无环的。
\(n\le 17\)
考虑状压，设 \(f_S\) 表示点集 \(S\) 的答案。对于 DAG，我们可以枚举它的顶层节点（入度为 \(0\)）组成的集合 \(T\)，然后从 \(f_{S-T}\) 转移。然而这样做会将一个 DAG 算重多次：考虑拆贡献，对于一个 DAG，设它的顶层节点集合为 \(T_0\)；那么只要我们枚举的 \(T\) 是真实的 \(T_0\) 的子集，就会将这个 DAG 计算一遍。
可以和上一题一样，构造出容斥系数 \((-1)^{|T|+1}\)，使得每个 \(|T_0|\) 恰好被统计一次：
| \(\vert T\vert\) |
\(1\) |
\(2\) |
\(3\) |
\(4\) |
\(5\) |
... |
| \(tar\) |
\(1\) |
\(1\) |
\(1\) |
\(1\) |
\(1\) |
... |
| \(w\) |
\(1\) |
\(-1\) |
\(1\) |
\(-1\) |
\(1\) |
... |
或者形式化的：
\[
\sum_{\varnothing\subset T\subseteq S}(-1)^{|T|+1}=\sum_{i=1}^{|S|}(-1)^{i+1}\binom{|S|}{i}=1
\]
因此转移为：
\[
f_S=\sum_{\varnothing\subset T\subseteq S}(-1)^{|T|+1}f_{S-T}2^{E(T,S-T)}
\]
至于如何在 \(O(3^n)\) 内计算出两个点集之间的边数，考虑刷表法，从当前集合 \(S\) 转移到 \(S+T\)。这样我们只需在 \(S\) 固定的情况下，动态维护 \(E(S,T)\) 即可。
给定一个 \(n\) 个节点的有向图，询问它有多少个边集是强联通的（所有节点之间都相互可达）。
\(n\le 15\)
设 \(ans_S\) 表示点集 \(S\) 构成恰好一个 SCC 的方案数。考虑正难则反，统计图缩点之后的 DAG 包含大于 \(1\) 个节点的情况。这里可以使用经典的 DAG 容斥方法，枚举 DAG 的顶层节点集合 \(T\)：
\[
ans_S=2^{E(S,S)}-\sum_{\varnothing\subset T\subseteq S}{(-1)^{?}f(T)2^{E(T,S-T)+E(S-T,S-T)}}
\]
其中 \(f_{T}\) 是将点集 \(T\) 划分为若干 SCC 且它们之间没有边的方案数。但是此时我们没法确定上式中的容斥系数，因为 \(-1\) 的指数应该等于 \(T\) 中 SCC 的数量加一。因此我们将容斥系数放到 \(f\) 中计算。考虑如何计算 \(f\)，枚举它的划分方案：
\[
f_{S}=-\sum_{\{x\}\subseteq T\subseteq S}ans_Tf_{S-T}
\]
其中钦定一个元素 \(x\) 在 \(T\) 中是为了防止算重。那么 \(ans\) 的转移可以写成：
\[
ans_S=2^{E(S,S)}+\sum_{\varnothing\subset T\subseteq S}{f(T)2^{E(T,S-T)+E(S-T,S-T)}}
\]
然而还会出现一个问题：\(f_S\) 和 \(ans_S\) 之间似乎存在循环依赖的问题。但马上我们就会注意到，\(ans_S\) 其实不应该统计 \(T=S\) 且 \(T\) 恰好形成一个 SCC 的情况。因此计算 \(f_S\) 时先忽略 \(T=S\) 的项，然后再计算 \(ans_S\)，最后再让 \(f_S\) 加上 \(ans_S\) 的贡献，去转移别人。通过刷表法可以动态维护 \(E(T,S-T)\)，时间复杂度 \(O(3^n)\)。
代码
| #include<iostream>
#define int long long
using namespace std;
const int N = 15;
const int MOD = 1e9 + 7;
int n, m;
int adj[1 << N];
int pw[300], cnt[1 << N];
int f[1 << N], ans[1 << N];
int ess[1 << N], est[1 << N];
signed main() {
pw[0] = 1;
for(int i = 1; i < 300; i++) pw[i] = pw[i - 1] * 2 % MOD;
for(int i = 1; i < (1 << N); i++) cnt[i] = cnt[i >> 1] + (i & 1);
cin >> n >> m;
for(int i = 1; i <= m; i++) {
int u, v;
cin >> u >> v;
--u, --v;
adj[1 << u] |= 1 << v;
}
for(int s = 1; s < (1 << n); s++) {
for(int i = 0; i < n; i++) {
if(!(s >> i & 1)) continue;
for(int j = 0; j < n; j++) {
if((s >> j & 1) && (adj[1 << i] >> j & 1)) ++ess[s];
}
}
}
f[0] = ans[0] = 1;
for(int s = 1; s < (1 << n); s++) {
int x = s & -s, is = (~s) & ((1 << n) - 1);
for(int t = (s - 1) & s; t > 0; t = (t - 1) & s) {
if(~t & x) continue;
f[s] = (f[s] + MOD - ans[t] * f[s - t] % MOD) % MOD;
}
ans[s] = (ans[s] + f[s]) % MOD;
ans[s] = (pw[ess[s]] + ans[s]) % MOD;
f[s] = (f[s] + MOD - ans[s]) % MOD;
for(int tt = (is - 1) & is; tt > 0; tt = (tt - 1) & is) {
int t = is ^ tt;
est[t] = est[t - (t & -t)] + cnt[adj[t & -t] & s];
}
est[is] = est[is - (is & -is)] + cnt[adj[is & -is] & s];
for(int t = is; t > 0; t = (t - 1) & is) {
ans[s + t] = (ans[s + t] + f[s] * pw[est[t]] % MOD * pw[ess[t]] % MOD);
}
}
cout << ans[(1 << n) - 1] << '\n';
return 0;
}
|
发现切成恰好 \(k\) 段的方案数是固定的，可以一次二项式反演求出。考虑如何判定一个划分方案是否合法，我们把每一段的子集缩成一个点，要求形成的图还是一个 DAG，并且段内的顺序是子集的一个拓扑序。集合拓扑序可以直接状压排列 \(O(n2^n)\) 求出。然后直接跑 DAG 容斥即可。
具体的，设 \(f_S\) 表示集合 \(S\) 的拓扑序数量，\(g_S\) 表示将 \(S\) 划分为若干不能相互到达的集合 \(T_{1\sim k}\)，集合内部任意按照拓扑序排列，乘以容斥系数 \((-1)^{k+1}\)，\(h_S\) 表示集合 \(S\) 划分为若干集合，集合缩完点之后形成 DAG 的方案数。转移是容易的。时间复杂度 \(O(3^nn^2)\)，小常数可以通过。
代码
| #include<iostream>
using namespace std;
typedef long long ll;
const int N = 15;
const int N2 = 40;
const int MOD = 1e9 + 7;
struct Edge {
int u, v;
} e[N * N];
inline void add(int &a, int b) { a += b; (a >= MOD) && (a -= MOD); }
int n, m, k, ans;
int to[N], rec[1 << N];
int f[1 << N], g[N + 1][1 << N], dp[N + 1][1 << N];
int h[N2];
int inv[N2], fac[N2], ifac[N2];
inline int A(int a, int b) { return (ll)fac[a] * ifac[a - b] % MOD; }
inline int C(int a, int b) { return (ll)fac[a] * ifac[a - b] % MOD * ifac[b] % MOD; }
int main() {
inv[1] = fac[0] = ifac[0] = 1;
for(int i = 2; i < N2; i++) inv[i] = MOD - (ll)inv[MOD % i] * (MOD / i) % MOD;
for(int i = 1; i < N2; i++) fac[i] = (ll)fac[i - 1] * i % MOD;
for(int i = 1; i < N2; i++) ifac[i] = (ll)ifac[i - 1] * inv[i] % MOD;
cin >> n >> m >> k;
for(int i = 1; i <= m; i++) {
int u, v; cin >> u >> v; --u, --v;
to[u] |= (1 << v);
e[i] = (Edge){u, v};
}
for(int i = 0; i < n; i++) rec[1 << i] = to[i];
for(int cnt = 1; cnt <= n; cnt++) {
for(int i = 1; i <= m; i++) rec[1 << e[i].u] |= rec[1 << e[i].v];
}
for(int s = 1; s < (1 << n); s++) {
if(s != (s & -s)) rec[s] = rec[s ^ (s & -s)] | rec[s & -s];
}
for(int i = 0; i < n; i++) f[1 << i] = 1;
for(int s = 1; s < (1 << n); s++) {
if(s == (s & -s)) continue;
for(int i = 0; i < n; i++) {
if(~s >> i & 1) continue;
if(rec[1 << i] & s) continue;
add(f[s], f[s ^ (1 << i)]);
}
}
for(int s = 1; s < (1 << n); s++) {
g[1][s] = f[s];
for(int i = 2; i <= __builtin_popcount(s); i++) {
int x = s & -s;
for(int t = (s - 1) & s; t; t = (t - 1) & s) {
if(~t & x) continue;
if(!(rec[s - t] & t) && !(rec[t] & s - t)) add(g[i][s], MOD - (ll)f[t] * g[i - 1][s - t] % MOD);
}
}
}
dp[0][0] = 1;
for(int s = 1; s < (1 << n); s++) {
dp[1][s] = f[s];
for(int i = 2; i <= __builtin_popcount(s); i++) {
for(int t = s; t; t = (t - 1) & s) {
if(!(rec[s - t] & t))
for(int j = 1; j <= i; j++) {
add(dp[i][s], (ll)dp[i - j][s - t] * g[j][t] % MOD);
}
}
}
}
for(int i = 1; i <= n; i++) h[i] = A(i + k, k);
for(int i = 2; i <= n; i++) {
for(int j = 1; j < i; j++) {
add(h[i], MOD - (ll)h[j] * C(i - 1, j - 1) % MOD);
}
}
for(int i = 1; i <= n; i++) add(ans, (ll)dp[i][(1 << n) - 1] * h[i] % MOD * fac[i] % MOD);
ans = (ll)ans * ifac[n] % MOD;
ans = (ll)ans * ((ll)ifac[n + k] * fac[n] % MOD) % MOD;
cout << ans << '\n';
return 0;
}
|
给定一个 \(n\) 个点 \(m\) 条边的无向图，将每条无向边拆成两条有向边，得到一个有向图。这个有向图的每条有向边都有 \(\frac 12\) 的概率保留，问剩下的图的最小外向生成树大小和原图的 mst 相同的概率。
\(n\le 15\)，保证原图无重边自环
特殊性质 C：所有边边权相同
可以先把概率计数转化成方案数计数（也可以不转化）。
考虑性质 C，此时问题等价于剩下的图存在一棵外向生成树。考虑对最终的图进行缩点，那么图存在一棵外向生成树当且仅当恰好有一个入度为 \(0\) 的 SCC。后文称这个 SCC 为原图的根集。先跑一遍主旋律，得到 \(f_S\) 表示集合 \(S\) 恰好形成一个 SCC 的方案数，然后钦定一些入度为 \(0\) 的 SCC，稍微推一下容斥系数即可。
考虑原问题，记原无向图为 \(G\)，最终得到的有向图为 \(G'\)。考察 Kruskal 的过程，对于一值域内的数字 \(w\)，我们保留 \(G\) 和 \(G'\) 中边权小于等于 \(w-1\) 的边，分别记为 \(G_{w-1},G'_{w-1}\)，那么 \(G'\) 合法的一必要条件是 \(G_{w-1},G'_{w-1}\) 两个图的弱连通性相同。由此我们可以对每个 \(w\) 知道加入了边权为 \(w\) 的边之后哪些连通块发生了合并。考虑还差了什么，我们将 \(G'_{w-1}\) 中的连通块视为点，全集 \(U\) 设为 \(G'_{w}\) 中的某一个连通块，那么加入的这些边权为 \(w\) 的边必须要存在一棵外向生成树子集连接这些点，同时这棵外向生成树每条边的末端必须都在原先连通块的根集内。
由此我们可以得到一种大概的思路：考察 \(G'_{w-1}\to G'_{w}\) 合并的过程，大概也是一个 C 性质。但是我们需要处理“生成树每条边末端都必须在原先连通块的根集内”的限制，因此考虑将根集写进 dp 状态里。设 \(ext_S\) 表示点集 \(S\) 原本所处的连通块集合的并，这样就可以用一维状态 \(S\) 表示点集 \(S\) 为 \(ext_S\) 内所有连通块根集的并，同时记录了根集和连通块集合（连通块集合对应 C 性质子问题的点集）的状态。
考虑如何跑 C 性质。还是先考虑 \(G'_{w}\) 中的一个连通块，并将 \(G'_{w-1}\) 中的连通块视为点。
- 设 \(f_S\) 表示 \(S\) 是 \(ext_S\) 所有连通块根集的并，且 \(S\) 所有点处于同一个 SCC 中，\(S\) 导出子图的方案数。
- 设 \(g_S\) 表示 \(S\) 是 \(ext_S\) 所有连通块根集的并，将 \(S\) 划分为若干相互不连边的 SCC，带上 \((-1)^{k}\) 容斥系数的加权方案数之和（\(k\) 表示 SCC 数量）。
- 设 \(dp_S\) 表示 \(S\) 恰好是当前连通块的根集的方案数。
在推转移式之前，有必要说明几点：
- \(G'_w\) 中一个连通块的根集，一定是 \(G'_{w-1}\) 中那些连通块根集并集的子集。根集并不是一个极大的强联通子图，但如果将小连通块看成点，那么根集是一个 SCC。即使加入的一些边权为 \(w\) 的边将小连通块中不在根集内的点和根集连在了同一个 SCC 内，这些点也不属于根集（否则会破坏最小外向生成树）。
- 对于末端不属于小连通块根集的边，可以随便选；
- 对于两端属于同一个小连通块的边，可以随便选；
考虑转移，大致和主旋律一样。然而发现 \(f_S\) 正难则反的过程中，需要保证小连通块内部是合法的，因此总方案数并不是 4^(导出子图边数)，而是 \(ext_S\) 中所有连通块的 \(dp\) 值的乘积，再乘以互相随意连边的方案数。
- 记 \(h_S\) 表示 \(S\) 是 \(ext_S\) 所有连通块根集的并，小连通块内部合法，再随意连边的方案数。
记 \(x=lowbit(S)\)，\(T=ext_x\cap S\)（随意一个连通块的根集），有转移：
\[
h_S=h_{S-T}dp_T2^{\small E(ext(S-T),ext(T))}
\]
由此可以写出 \(f,g\) 的转移：
\[
f_S=h_S-\sum_{T\subseteq S}[ext(T)\cap ext(S-T)=\varnothing]g_Th_{S-T}2^{E(ext(T),ext(S-T))+E(ext(S-T),ext(T)-T)}\\
g_S=-f_S-\sum_{\{x\}\subseteq T\subseteq S}[ext(T)\cap ext(S-T)=\varnothing]f_Tg_{S-T}2^{E(ext(S-T),ext(T)-T)+E(ext(T),ext(S-T)-(S-T))}
\]
初始时，对于 \(S\subseteq ext(lowbit(S))\)，\(f_S=dp_S\)。在尝试写出 \(dp\) 的转移时，还需要用到 \(sh_S\) 表示不钦定根集，集合 \(S\) 导出子图的方案数。其转移也是容易的（第二行 \(T\) 也是随意一个连通块的根集）：
\[
\begin{cases}
sh_S=\sum_{T\subseteq S}h_T&,S=ext(lowbit(S))\\
sh_S=sh_T+sh_{S-T}&,S=ext(S)
\end{cases}
\]
那么 \(dp_S\) 的转移如下（钦定集合 \(T\) 是额外入度为 \(0\) 的 SCC 的点集，容斥系数为 \((-1)^{k}\)，剩下的点随便连边即可）：
\[
dp_S=\sum_{T\subseteq U-S}[ext(S)\cap ext(T)=\varnothing]g_Tf_Ssh_{U-ext(S+T)}2^{E(ext(S+T),U-ext(S+T))+E(U-ext(S+T),ext(S+T)-(S+T))}2^{E(ext(S),ext(T)-T)+E(ext(T),ext(S)-S)}
\]
最终答案为 \(\sum_S dp_S\)。至于时间复杂度，显然不超过 \(O(n3^n)\)。可以加一点小优化：只有连通块发生变化时，才更新内部的 dp 值。这样的话，每下降一层时间就至少除以 \(3\)，因此时间复杂度 \(O(3^n)\)。
代码
| #include<iostream>
#include<vector>
#define ll long long
#define subseteq(s, U) for(int ss = (U - 1) & U, s = U ^ ss; ss != U; ss = (ss - 1) & U, s = U ^ ss)
#define subseteq_empty(s, U) for(int ss = U | (1 << N), s = 0; ss != U; ss = (ss - 1) & U, s = U ^ ss)
using namespace std;
const int N = 15;
const int M = N * N + 10;
const int MOD = 1e9 + 7;
struct Edge {
int u, v;
};
int T, n, m; ll ans1, ans2, ans;
int ext[1 << N], net[1 << N];
int ei[1 << N];
ll f[1 << N], g[1 << N], h[1 << N], sh[1 << N], dp[1 << N];
int vis[1 << N];
ll pw2[M], pw4[M];
vector<Edge> edg[M];
inline int E(int s, int t) {
return ei[s + t] - ei[s] - ei[t];
}
void clear() {
ans1 = 1;
ans2 = 0;
for(int i = 1; i <= m; i++) edg[i].clear();
for(int i = 0; i < (1 << n); i++) ext[i] = net[i] = ei[i] = f[i] = g[i] = h[i] = sh[i] = dp[i] = vis[i] = 0;
f[0] = g[0] = h[0] = sh[0] = 1;
}
void work() {
for(int i = 1; i <= m; i++) {
int u, v, w;
cin >> u >> v >> w;
--u, --v;
edg[w].push_back({u, v});
}
for(int i = 0; i < N; i++) {
net[1 << i] = ext[1 << i] = 1 << i;
dp[1 << i] = 1;
}
for(int w = 1; w <= m; w++) {
for(Edge e : edg[w]) {
int u = e.u, v = e.v;
if(ext[1 << u] == ext[1 << v]) { ans1 = ans1 * 4 % MOD; continue; }
if(net[1 << u] == net[1 << v]) continue;
int tmp = net[1 << u], res = net[1 << u] + net[1 << v];
while(tmp) {
net[tmp & -tmp] = res;
tmp -= tmp & -tmp;
}
tmp = net[1 << v];
while(tmp) {
net[tmp & -tmp] = res;
tmp -= tmp & -tmp;
}
}
for(int i = 0; i < n; i++) {
if(net[1 << i] != ext[1 << i] && !vis[net[1 << i]]) {
int U = net[1 << i];
vis[U] = 1;
subseteq(s, U) sh[s] = 0;
subseteq(s, U) {
ext[s] = ext[s - (s & -s)] | ext[s & -s];
h[s] = dp[s];
sh[ext[s]] = (sh[ext[s]] + h[s]) % MOD;
ei[s] = 0;
for(Edge e : edg[w]) {
if(((s >> e.u) & 1) && ((s >> e.v) & 1)) ei[s]++;
}
}
subseteq(s, U) {
if(!h[s]) {
int t = ext[s & -s] & s;
h[s] = h[s - t] * dp[t] % MOD * pw4[E(ext[s - t], ext[t])] % MOD;
}
}
subseteq(s, U) {
if(s != ext[s]) continue;
int t = ext[s & -s] & s;
sh[s] = sh[s - t] * sh[t] % MOD * pw4[E(s - t, t)] % MOD;
}
subseteq(s, U) {
f[s] = h[s]; g[s] = 0;
int x = s & -s;
subseteq(t, s) {
if(ext[t] & ext[s - t]) continue;
if(!(t & x)) continue;
if(t == s) continue;
g[s] = (g[s] - f[t] * g[s - t] % MOD * pw2[E(ext[s - t], ext[t] - t) + E(ext[t], ext[s - t] - (s - t))] % MOD + MOD) % MOD;
}
subseteq(t, s) {
if(ext[t] & ext[s - t]) continue;
f[s] = (f[s] + g[t] * h[s - t] % MOD * pw2[E(ext[t], ext[s - t]) + E(ext[s - t], ext[t] - t)] % MOD + MOD) % MOD;
}
g[s] = (g[s] - f[s] + MOD) % MOD;
}
subseteq(s, U) {
dp[s] = sh[U - ext[s]] % MOD * pw2[E(ext[s], U - ext[s]) + E(U - ext[s], ext[s] - s)] % MOD;
subseteq(t, U - s) {
if(ext[s] & ext[t]) continue;
dp[s] = (dp[s] + g[t] * sh[U - ext[s + t]] % MOD * pw2[E(ext[s + t], U - ext[s + t]) + E(U - ext[s + t], ext[s + t] - (s + t))] % MOD * pw2[E(ext[s], ext[t] - t) + E(ext[t], ext[s] - s)] % MOD) % MOD;
}
dp[s] = dp[s] * f[s] % MOD;
}
subseteq(s, U) {
ext[s] = U;
}
}
}
}
subseteq(s, (1 << n) - 1) {
ans2 = (ans2 + dp[s]) % MOD;
}
ans = ans1 * ans2 % MOD;
for(int i = 1; i <= m; i++) ans = ans * 250000002 % MOD;
cout << ans << '\n';
}
int main() {
pw2[0] = pw4[0] = 1;
for(int i = 1; i < M; i++) pw2[i] = pw2[i - 1] * 2 % MOD, pw4[i] = pw4[i - 1] * 4 % MOD;
cin >> T >> T;
while(T--) {
cin >> n >> m;
clear();
work();
}
return 0;
}
|