【学习笔记】矩阵树定理

imily's notebook

2020-07-17

upd on 2021.3.26: 感性证明一下！？：

设点数为 $n$。

我们都知道 kirchhoff 矩阵 $[i, j]$ 位置上的值表示 $i$ 走到 $j$ 的走法个数。减去度数矩阵是减去自环的走法。

然后开始容斥，容斥掉所有轮换个数 $> 1$ 的排列。而「$n$ - 轮换个数」与「排列逆序对数」的奇偶性相同，证明你就考虑给排列排序时每个轮换依次排，交换一对时逆序对的奇偶性也会变，总共交换「$n$ - 轮换个数」次，最后排列逆序对为 $0$。所以直接算行列式就好啦。

对容斥系数的另外一种理解是，$(-1)^{n - 环数} = (-1)^{n + 环数} = (-1)^{偶环数}$，所以图的邻接矩阵的行列式 $= \sum\limits_{G 的所有分解方式，分解为若干环} (-1)^{偶环数}$

OIWiki上的矩阵树相关

行列式

$n \times n$ 的矩阵 $A$ 的行列式可以理解为行或列向量所构成的超平行多面体的有向面积或有向体积，是一个标量。

矩阵 $A$ 的行列式用 $det(A)$ 表示。

$det(A) = \sum\limits_{j_1 j_2 \cdots j_n} (-1)^{\sigma(j_1 j_2 \cdots j_n)} a_{1, j_1} a_{1, j_2} \cdots a_{n, j_n}$

其中 $\sigma(j_1 j_2 \cdots j_n)$ 表示 $j_1 j_2 \cdots j_n$ 的逆序对数。

行列式初等变换(最下面)

行列互换，行列式值不变
行列式一行或一列的公因子可以提出去
行列式中若有某一行是两组数的和，则该行列式等于两个行列式的和
两行或两列交换，行列式变号（区别于第一条）
两行或两列相同，行列式为 0
两行或两列成比例，行列式为 0

高斯消元化为上三角矩阵：

运用了行列式的初等变换。

实数：直接处理
模意义：1. 模为质可用逆元 2. 模不为质用辗转相除法，复杂度会多一个 log（详情见例题代码）

基尔霍夫矩阵

A 为邻接矩阵，D 为度数矩阵，Kirchhoff矩阵为 K = D - A

行列式 a[i, i] 记录点 i 度数，a[i, j] 表示 i, j 之间边数的相反数。

具体实现的话，设 kirchhoff 矩阵为 a，若存在边 (u, v) 则 a[u, u]++, a[v, v]++, a[u, v]—, a[v, u]—

Kirchhoff矩阵每行内数的和和每列内数的和都为 0，所以行列式为 0

矩阵树定理

用于求解图上的生成树个数。

无向图生成树个数 = Kirchhoff 矩阵任何一个 N - 1 阶主子式的行列式的绝对值.

有向图相关见这里（由于本人过于菜连有向图生成树是啥都不会先咕咕后续做到题了再来填坑！）

可能是，去掉根所在的那行那列，度数矩阵外向树就是入度，内向树就是出度。

无向图，谁是根并不重要所以随便去掉哪一行&列都可以；

有向图，删去指定的根所在的行和列，求剩下的矩阵行列式即可。

如果图不连通，那么任意 N - 1 阶主子式为 0。证明：如果图不连通，那么每个连通块内的点构成的矩阵仍然是 Kirchhoff 矩阵，而连通块不止一个，所以去掉第 i 行第 i 列之后，一定有一个连通块仍然是 Kirchhoff 矩阵，也就是行列式为 0 。

变元矩阵树定理

求所有生成树总边积的和。行列式 a[i, i] 记录点 i 边权和，a[i, j] 表示 i, j 之间边权的相反数。

例题1. [HEOI2015]-小Z的房间

模板题，注意模数非质

code

#include <bits/stdc++.h>
#define rep(i, x, y) for (int i = x; i <= y; i++)
using namespace std;

const int mod = 1e9;
const int dir[][2] = {{0, 1}, {1, 0}};
typedef long long ll;
ll n, m, a[100][100], id[10][10], idx;
char s[10][10];

int main() {
    cin >> n >> m;
    rep(i, 1, n) {
        scanf("%s", s[i] + 1);
        rep(j, 1, m)
            if (s[i][j] == '.') id[i][j] = ++idx;
    }
    rep(i, 1, n) {
        rep(j, 1, m) {
            if (id[i][j]) {
                rep(k, 0, 1) {
                    int tx = i + dir[k][0], ty = j + dir[k][1];
                    if (tx > n || ty > m || !id[tx][ty]) continue;
                    int x = id[i][j], y = id[tx][ty];
                    a[x][x]++, a[y][y]++, a[x][y]--, a[y][x]--;
                }
            }
        }
    }
    ll ans = 1;
    rep(i, 1, idx - 1) {  // n - 1 阶主子式
        rep(j, i + 1, idx - 1) {
            while (a[j][i]) {  // 在模意义下用辗转相除法，类比 gcd(a, b) = gcd(b, a % b) 直到 b = 0，用第 i 行消第 j 行
                ll d = a[i][i] / a[j][i];
                rep(k, i, idx - 1)
                    a[i][k] = (a[i][k] - d * a[j][k] % mod + mod) % mod;
                swap(a[i], a[j]), ans = -ans;
            }
        }
        ans = (ans * a[i][i] % mod + mod) % mod;
    }
    printf("%lld\n", ans);
    return 0;
}

例题2. [SHOI2016]-黑暗前的幻想乡

容斥 + 矩阵树定理

code

#include <bits/stdc++.h>
#define rep(i, x, y) for (int i = x; i <= y; i++)
using namespace std;

const int mod = 1e9 + 7;
typedef long long ll;
typedef pair<int, int> pii;
ll n, a[20][20], res;
vector<pii> g[20];

int main() {
    cin >> n;
    rep(i, 1, n - 1) {
        int m, x, y; scanf("%d", &m);
        while (m--) {
            scanf("%d%d", &x, &y);
            g[i].push_back(make_pair(x, y));
        }
    }
    rep(s, 0, (1 << (n - 1)) - 1) {
        int cnt = 0;
        memset(a, 0, sizeof(a));
        rep(i, 1, n) {
            if (s & (1 << (i - 1))) {
                ++cnt;
                for (int j = 0; j < g[i].size(); j++) {
                    int x = g[i][j].first, y = g[i][j].second;
                    a[x][x]++, a[y][y]++, a[x][y]--, a[y][x]--;
                }
            }
        }
        ll ans = 1;
        rep(i, 1, n - 1) {
            rep(j, i + 1, n - 1) {
                while (a[j][i]) {
                    ll d = a[i][i] / a[j][i];
                    rep(k, i, n - 1)
                        a[i][k] = (a[i][k] - d * a[j][k] % mod + mod) % mod;
                    swap(a[i], a[j]), ans = -ans;
                }
            }
            ans = (ans * a[i][i] % mod + mod) % mod;
        }
        res = (res + ((n - 1 - cnt) & 1 ? -1ll : 1ll) * ans + mod) % mod;  // !!!
    }
    printf("%lld\n", res);
    return 0;
}

例题3. [SDOI2014]-重建

推柿子 + 变元矩阵树定理，显然边权跟概率有关

$\sum\limits_{T}(\prod\limits_{e \in T} P_e \prod\limits_{e \not\in T}(1-P_e))$

$\sum\limits_T( \prod\limits_{e \in T} P_e \frac{ \prod\limits_e (1 - P_e) }{ \prod\limits_{e \in T} (1 - P_e) } )$

$\prod\limits_e (1 - P_e) \sum \prod\limits_{e \in T} \frac{P_e}{1 - P_e}$

code

#include <bits/stdc++.h>
#define rep(i, x, y) for (int i = x; i <= y; i++)
using namespace std;

const int N = 55;
const double eps = 1e-9;
int n;
double D[N][N], ans = 1.0;

double solve() {
    bool ff = 0;
    double ret = 1.0;
    int N = n - 1;
    rep(i, 1, N) {
        int mx = i;
        rep(j, i + 1, N)
            if (D[mx][i] < D[j][i]) mx = j;
        if (mx != i) ff ^= 1, swap(D[mx], D[i]);
        if (D[i][i] > -eps && D[i][i] < eps) return 0.0;
        rep(j, i + 1, N) {
            double t = D[j][i] / D[i][i];
            rep(k, i, N) D[j][k] -= t * D[i][k];
        }
        ret *= D[i][i];
    }
    if (ff) ret *= -1.0;
    return ret;
}

int main() {
    cin >> n;
    rep(i, 1, n) rep(j, 1, n) cin >> D[i][j];
    rep(i, 1, n)
        rep(j, 1, n) {
            if (fabs(D[i][j]) < eps) D[i][j] = eps;
            else if (fabs(1.0 - D[i][j]) < eps) D[i][j] = 1.0 - eps;
            if (j > i) ans *= (1 - D[i][j]);
            D[i][j] /= (1.0 - D[i][j]);
        }
    rep(i, 1, n)
        rep(j, 1, n)
            if (i != j)
                D[i][i] += D[i][j], D[i][j] *= -1;
    printf("%.5lf\n", ans * solve());
    return 0;
}

例题4. [JSOI2008]-最小生成树计数

最小生成树性质：对于所有最小生成树，每种边权的边数相同；且对于所有生成树，某种权值的边连完后图的联通性相同

所以可以分别处理每种权值，乘起来。

算同种边权的边的贡献，由于具有相同权值的边不超过 10 条，暴搜也可以过。。（2^10 很稳的！

正解是矩阵树。注意，同种边权的边联通性相同并不等于连通，因此要连一些桥，并不影响矩阵树算答案。

code

#include <bits/stdc++.h>
#define rep(i, x, y) for (int i = x; i <= y; i++)
using namespace std;

const int mod = 31011, N = 105, M = 1e3 + 10;
typedef long long ll;
ll n, m, ans = 1, num;
ll fa[N], fat[N], id[M << 1], a[N][N];
struct node {
    ll u, v, w;
    bool operator < (const node &a) const { return w < a.w; }
}e[M], t[M];

int getfa(int x) {
    return fa[x] == x ? x : fa[x] = getfa(fa[x]);
}
int getfa2(int x) {
    return fat[x] == x ? x : fat[x] = getfa2(fat[x]);
}

ll Gauss(ll n) {
    ll ret = 1;
    rep(i, 1, n - 1) {
        rep(j, i + 1, n - 1) {
            while (a[j][i]) {
                ll d = a[i][i] / a[j][i];
                rep(k, i, n - 1)
                    a[i][k] = (a[i][k] - d * a[j][k] % mod + mod) % mod;
                swap(a[i], a[j]), ret = -ret;
            }
        }
        ret = (ret * a[i][i] % mod + mod) % mod;
    }
    return ret;
}

void calc(int l, int r) {
    int cnt = 0;
    memset(a, 0, sizeof(a));
    rep(i, l, r) {
        t[i] = e[i];
        int u = getfa(t[i].u), v = getfa(t[i].v);
        t[i].u = u, t[i].v = v;  // 类似于缩点
        if (u == v) continue;
        id[++cnt] = u, id[++cnt] = v;
    }
    sort(id + 1, id + cnt + 1);
    cnt = unique(id + 1, id + cnt + 1) - id - 1;
    rep(i, 1, cnt) fat[i] = i;
    rep(i, l, r) {
        if (t[i].u == t[i].v) continue;

        int u = getfa(t[i].u), v = getfa(t[i].v);
        if (u != v) --num, fa[u] = v;

        u = lower_bound(id + 1, id + cnt + 1, t[i].u) - id;
        v = lower_bound(id + 1, id + cnt + 1, t[i].v) - id;
        a[u][u]++, a[v][v]++, a[u][v]--, a[v][u]--;

        u = getfa2(u), v = getfa2(v);
        if (u != v) fat[u] = v;
    }
    rep(i, 2, cnt) {  // 有可能不连通，那么连一些桥，并不影响矩阵树算答案
        int u = getfa2(i), v = getfa2(i - 1);
        if (u == v) continue;
        a[u][u]++, a[v][v]++, a[u][v]--, a[v][u]--;
        fat[u] = v;
    }
    ans = ans * Gauss(cnt) % mod;
}

int main() {
    cin >> n >> m;
    rep(i, 1, m)
        scanf("%lld%lld%lld", &e[i].u, &e[i].v, &e[i].w);
    sort(e + 1, e + m + 1);
    rep(i, 1, n) fa[i] = i;
    num = n;
    for (int i = 1, j = 1; i <= m; i = j) {
        while (j <= m && e[i].w == e[j].w) ++j;
        for (j = i; j <= m; ++j) if (e[i].w != e[j].w) break;
        if (i + 1 < j) {
            calc(i, j - 1);
        } else {
            int u = getfa(e[i].u), v = getfa(e[i].v);
            if (u != v) fa[u] = v;
            --num;
        }
    }
    if (num > 1) puts("0");  // 判定最小生成树的存在
    else printf("%lld\n", ans);
    return 0;
}