@算法 - [email protected] matrix - tree 定理(矩阵树定理)

目录

  • @0 - 参考资料@
  • @0.5 - 你所需要了解的线性代数知识@
  • @1 - 定理主体@
    • @证明 part - [email protected]
    • @证明 part - [email protected]
    • @证明 part - [email protected]
    • @证明 part - 4@
  • @2 - 一些简单的推广@
  • @3 - 例题与应用@

@0 - 参考资料@

MoebiusMeow 的讲解(超喜欢这个博主的!)

网上找的另外一篇讲解

@0.5 - 你所需要了解的线性代数知识@

什么是矩阵?

什么是高斯消元?这个虽然与主题无关,但是求解行列式需要用。

什么是行列式?行列式采用排列的方式定义。但是为了简洁(我懒得写),我们直接使用拉普拉斯展开。

设一个 n*n 矩阵 A:

\[A = \begin{bmatrix}
a_{11} & a_{12} & \cdots & a_{1n} \\
a_{21} & a_{22} & \cdots & a_{2n} \\
\vdots & \vdots & \ddots & \vdots \\
a_{n1} & a_{n2} & \cdots & a_{nn}
\end{bmatrix}\]

对第 i 行(i 可以是任意的)进行展开,得到 A 的行列式 det(A)为:

\[\det(A) = \sum_{j=1}^{n}(-1)^{i+j}a_{ij}M_{ij}\]

其中 \(M_{ij}\) 表示去掉第 i 行第 j 列后剩下的矩阵的行列式,称为 \(a_{ij}\) 的余子式。

行列式的一些性质(我真的懒得写证明)

(1)矩阵 A 的行列式等于其转置的行列式(可以理解为行和列是等价的)。

(2)把矩阵 A 的某一行变为 k 倍,行列式变为 k 倍。

(3)把矩阵 A 的第 i 行加第 j 行的若干倍,行列式不变。

(4)交换两行,行列式变号(正变负,负变正)。

(5)上/下三角矩阵的行列式等于主对角线元素乘积。

(1)表示行列式中行和列是等价的。

(2)、(3)和(4)都是关于初等行变换的。

(5)提供了一种计算特殊矩阵行列式的方法。

因此,我们可以通过高斯消元将一个矩阵消成三角矩阵 O(n^3) 求行列式。

@1 - 定理主体@

先考虑最简单的:无重边无自环的无向图 G。(简称三无图?)

我们定义这个无向图 G 的度数矩阵 \(D\)。只有 \(D_{ii}\) 为点 i 的度数,其他项为 0。

我们定义这个无向图 G 的邻接矩阵 \(A\)。\(A_{ij}\) 为 1 当且仅当点 i 与点 j 有连边。

我们定义这个无向图 G 的基尔霍夫矩阵 \(L\)。\(L = D - A\)。

matrix - tree 定理:这个无向图的生成树个数等于任意一个主对角线上的元素的余子式。

@证明 part - [email protected]

我们先来证明一个简单的引理:任意一个基尔霍夫矩阵的行列式总等于零。

观察到基尔霍夫矩阵,其每一行每一列的和等于零。

因此我们可以将其他行加到第一行,第一行就变成全零的。

再把第一行展开就可以证明了。

@证明 part - [email protected]

再证明这样一个引理:如果 G 不连通,主对角线上的余子式总等于零。

我们将基尔霍夫矩阵 L 作这样的变换:将每一个连通块的点都放在一起。

每次交换 i, j 两行,再交换 i, j 两列,就可以将点 i 与点 j 的位置对换,且行列式两次变号后还是保持不变。

这样最终基尔霍夫矩阵长这样:

\[L = \begin{bmatrix}A_1 & 0 & \cdots & 0\\ 0 & A_2 & \cdots & 0 \\ \vdots & \vdots& \ddots & \vdots \\ 0 & 0 & \cdots & A_k \end{bmatrix}\]

这个时候,\(\det(L)=\det(A_1)\det(A_2)\dots\det(A_k)\)。至于为什么直观感受一下,对角线的乘积嘛。

其实证明也比较简单,将所有 A 高斯消元消成三角矩阵就可以了。

因为不连通,所以 \(2 \leq k\)。

我们的余子式中至少包含一个完整的 A 矩阵。

而每一个矩阵 A 都是基尔霍夫矩阵。所以最终的余子式一定等于零。

@证明 part - [email protected]

接下来我们证明最后一个引理:如果 G 是一棵树,主对角线上的余子式总等于一。

我们先把要展开的那一行的点作为根结点,然后以根为起点进行 bfs。

按照 bfs 遍历的逆序将结点排列(最后遍历的放最前面)。

一样还是每次交换 i, j 两行再交换 i, j 两列交换 i, j 两点的位置。

然后,我们从左往右扫描每一列。每次扫描到某一列,我们将这一列加到这个结点的父亲那一列上(如果它的父亲是根结点就什么也不做)。

通过归纳法可以证明,扫描到第 i 列的时候,第 i 列的第 1~i-1 行全是 0,第 i 行为 1,第 i+1~n 行上除了它父亲那一行为 -1 其他都为 0。

直观理解也可以理解。

于是最后剩下的矩阵是三角矩阵且主对角线上都是 1。

行列式自然就等于 1 了。

@证明 part - 4@

要证明矩阵树定理,我们要先知道 Cauchy - Binet定理

对于 n*m 矩阵 A 与 m*n 矩阵 B (m < n):

\[\det(AB)=\sum_{|s|=m}\det(A_{s*})\det(B_{s*})=\sum_{|s|=m}\det(A_{s*}*B_{s*})\]

表示我们选择 A 中的 m 行构成的 m*m 的方阵乘上选择 B 中相应的 m 列构成的 m*m 的方阵之和等于 AB 的行列式。

该定理的证明超出了我的理解范围。

然后我们开始我们正式的证明:

我们构造 m*n 矩阵 \(B\),其中 m 是边数,n 是点数。

如果第 i 条边为 \((u_i, v_i)\),矩阵中就有 \(b_{iu_i}=1, b_{iv_i}=-1\)(因为这条边是无向边,所以哪一个是 -1 哪一个是 1 没有什么实际影响,只要保证一个为 1 另一个为 -1 就可以了。)

然后,有 \(L = B^TB\)。为什么呢?可以把矩阵乘法的式子列出来看。

再根据上面那个定理,就有:

\[M_{ii}=\det(B^T_iB_i)=\sum_{|s|=n-1}\det((B_i^T)_{s*}(B_i)_{*s})\]

右边那个其实就是从 m 条边里面选择 n-1 条边构成的图,展开的余子式。根据我们的几个引理,只有当这个图是树的时候,会对答案产生贡献且贡献为 1。所以余子式等于生成树个数。

@2 - 一些简单的推广@

首先最关注的应该是重边自环

根据上面的推导过程,重边是没有影响。

自环的话,根据我们构造的 B 矩阵,你把 Bij++ 又把 Bij-- 过后,实际上这个自环就不会造成任何影响了。但是实现上可能要特殊处理一下。

然后就是有向图的推广。有向图的生成树,我们这里只考虑以结点 i 为根的外向生成树。其他情况可以类比得出。

考虑我们证明树的余子式总为 1 时,我们是用一个点的度数去消去它父亲的邻接矩阵。因此,我们将邻接矩阵改成只记出边的矩阵,将度数矩阵改为只记入度的矩阵(当然因为矩阵的行列是等价的所以你也不需要特别去弄清楚谁是出边谁是入边,多试试就可以了)。

注意我们如果是有向图,不能随便展开行,展开的那一行就是根。

再然后,就是其他类型的生成树统计而不仅仅是个数统计。最基础的应该是每条边带边权,求不同情况生成树的边权乘积的总和。可以发现边权为 1 时就是我们的正常个数统计,因此我们可以直接将该加 1 减 1 的地方换成边权即可。

还有一个高斯消元上的细节问题。因为我们要求解的是行列式,所以用不着使用将某一行变 k 倍的初等行变换,让每一行的第一个非零元变成 1。

因此就有一个问题,是否我们可以在不一定存在逆元的模数(即非质数)意义下作行列式的求解。

我们知道使用辗转相除法,在 log 次就可以使得某一个数为 0。

因此我们对两行使用辗转相除法,在 O(n^3log n) 的时间复杂度实现非质数模数的行列式求解。

@3 - 例题与应用@

一开始竟然把 gauss 拼写成了 guess,果然太久没有写过这东西了……。

首先是无向图的生成树计数模板题(可能建图不是那么模板?)。

bzoj - 4031 : [HEOI2015]小Z的房间

#include<cstdio>
#include<algorithm>
using namespace std;
const int MAXN = 9;
const int MAXV = 81;
const int MOD = int(1E9);
const int dx[] = {0, 0, 1, -1};
const int dy[] = {1, -1, 0, 0};
int mat[MAXV][MAXV];
int gauss(int r, int c) {
    int nwr = 0, nwc = 0, f = 1;
    for(int i=0;i<r;i++)
        for(int j=0;j<c;j++)
            mat[i][j] = (mat[i][j] + MOD)%MOD;
    while( nwr < r && nwc < c ) {
/*
        for(int i=0;i<r;i++)
            for(int j=0;j<c;j++)
                printf("%9d%c", mat[i][j], (j + 1 == c) ? ‘\n‘ : ‘ ‘);
        puts("");
*/
        int mxr = nwr;
        for(int i=nwr+1;i<r;i++)
            if( mat[i][nwc] > mat[mxr][nwc] ) mxr = i;
        if( mat[mxr][nwc] ) {
            if( mxr != nwr ) {
                for(int i=nwc;i<c;i++)
                    swap(mat[nwr][i], mat[mxr][i]);
                f *= -1;
            }
            for(int i=nwr+1;i<r;i++)
                while( mat[i][nwc] ) {
                    if( mat[nwr][nwc] < mat[i][nwc] ) {
                        for(int j=nwc;j<c;j++)
                            swap(mat[nwr][j], mat[i][j]);
                        f *= -1; continue;
                    }
                    else {
                        int x = mat[nwr][nwc]/mat[i][nwc];
                        for(int j=nwc;j<c;j++)
                            mat[nwr][j] = (mat[nwr][j] + MOD - 1LL*x*mat[i][j]%MOD)%MOD;
                    }
                }
            nwr++;
        }
        nwc++;
    }
    int ans = 1;
    for(int i=0;i<r;i++)
        ans = 1LL*ans*mat[i][i]%MOD;
    return (ans*f + MOD)%MOD;
}
int num[MAXN][MAXN];
char str[MAXN];
int main() {
    int n, m, cnt = 0;
    scanf("%d%d", &n, &m);
    for(int i=0;i<n;i++) {
        scanf("%s", str);
        for(int j=0;j<m;j++) {
            if( str[j] == ‘.‘ )
                num[i][j] = (cnt++);
            else num[i][j] = -1;
        }
    }
    for(int i=0;i<n;i++)
        for(int j=0;j<m;j++) {
            if( num[i][j] == -1 ) continue;
            for(int k=0;k<4;k++) {
                int x = i + dx[k], y = j + dy[k];
                if( x < 0 || y < 0 || x >= n || y >= m ) continue;
                if( num[x][y] == -1 ) continue;
                mat[num[i][j]][num[i][j]]++, mat[num[i][j]][num[x][y]]--;
            }
        }
    printf("%d\n", gauss(cnt - 1, cnt - 1));
}

然后是有向图生成树计数模板题。

bzoj - 4894:天赋/bzoj - 5297: [Cqoi2018]社交网络

两者代码细节上有些不同,这里给出 5297 的代码。

#include<cstdio>
#include<algorithm>
using namespace std;
const int MAXN = 250;
const int MOD = 10007;
int pow_mod(int b, int p) {
    int ret = 1;
    while( p ) {
        if( p & 1 ) ret = ret*b%MOD;
        b = b*b%MOD;
        p >>= 1;
    }
    return ret;
}
int mat[MAXN][MAXN];
int gauss(int r, int c) {
    int nwr = 1, nwc = 1, f = 1;
    for(int i=1;i<r;i++)
        for(int j=1;j<c;j++)
            mat[i][j] = (mat[i][j] + MOD)%MOD;
    while( nwr < r && nwc < c ) {
        int mxr = nwr;
        for(int i=nwr+1;i<r;i++)
            if( mat[i][nwc] > mat[mxr][nwc] )
                mxr = i;
        if( mat[mxr][nwc] ) {
            if( nwr != mxr ) {
                for(int i=nwc;i<c;i++)
                    swap(mat[nwr][i], mat[mxr][i]);
                f *= -1;
            }
            int inv = pow_mod(mat[nwr][nwc], MOD-2);
            for(int i=nwr+1;i<r;i++) {
                int x = 1LL*inv*mat[i][nwc]%MOD;
                for(int j=nwc;j<c;j++)
                    mat[i][j] = (mat[i][j] + MOD - 1LL*x*mat[nwr][j]%MOD)%MOD;
            }
            nwr++;
        }
        nwc++;
    }
    int ans = 1;
    for(int i=1;i<r;i++)
        ans = 1LL*ans*mat[i][i]%MOD;
    return (ans*f + MOD)%MOD;
}
int main() {
    int n, m; scanf("%d%d", &n, &m);
    for(int i=1;i<=m;i++) {
        int a, b; scanf("%d%d", &a, &b); a--, b--;
        mat[a][a]++, mat[a][b]--;
    }
    printf("%d\n", gauss(n, n));
}

然后是一道不是特别模板的题。

bzoj - 3534: [Sdoi2014]重建

给出每条边的出现概率,最后出现的图是树的概率。

可以发现这个概率其实等于所有树边出现的概率*所有非树边不出现的概率。

但是因为我们求解生成树只能统计树边。所以我们稍微转换一下,变成:所有树边出现的概率*所有边不出现的概率/所有树边不出现的概率。

然后就可以跑带权的矩阵树定理了。

注意涉及到除法,除数可能为 0(某条边一定为 0)。我们给它微小扰动一下(给一定出现的边的概率减去 eps)。

最后……如果你输出 5 位小数都是过不了的,反正我直接输出 10 位小数。

#include<cmath>
#include<cstdio>
#include<algorithm>
using namespace std;
const int MAXN = 50;
const double EPS = 1E-7;
double mat[MAXN][MAXN];
double gauss(int r, int c) {
    int nr = 0, nc = 0; double f = 1;
    while( nr < r && nc < c ) {
        int mxr = nr;
        for(int i=nr+1;i<r;i++)
            if( fabs(mat[i][nc]) > fabs(mat[mxr][nc]) ) mxr = i;
        if( fabs(mat[mxr][nc]) > EPS ) {
            if( nr != mxr ) {
                for(int j=nc;j<c;j++)
                    swap(mat[nr][j], mat[mxr][j]);
                f *= -1;
            }
            for(int i=nr+1;i<r;i++) {
                double x = mat[i][nc]/mat[nr][nc];
                for(int j=nc;j<c;j++)
                    mat[i][j] = (mat[i][j] - x*mat[nr][j]);
            }
            nr++;
        }
        nc++;
    }
    double ans = 1;
    for(int i=0;i<r;i++)
        ans = ans*mat[i][i];
    return ans*f;
}
int main() {
    double f = 1;
    int n; scanf("%d", &n);
    for(int i=0;i<n;i++)
        for(int j=0;j<n;j++) {
            double x; scanf("%lf", &x);
            if( x == 1 ) x = 1 - EPS;
            mat[i][i] += x/(1-x), mat[i][j] -= x/(1-x);
            if( i < j ) f *= (1-x);
        }
    printf("%.10lf\n", gauss(n-1, n-1)*f);
}

然后可能不是模板题的一道题。

bzoj - 4596: [Shoi2016]黑暗前的幻想乡

给出 n-1 个公司以及每个公司能够修出来的公路,现在让你统计有多少种修公路的方案,使得每个公司都能负责其中的恰好一条公路且 n 个点连通。

注意如果两种方案即使长得一样,修建某一条公路的公司如果不一样我们也认为是不同的方案。n<= 17。

初看毫无思路。但是如果我们消去每个公司恰好负责一条路这个限制就是一个很简单的生成树计数了。

因此,也许我们可以容斥?

一个方案不合法,则至少有一个公司没有参与到修建公路的工程中。

我们用总方案 - 至少有一个公司没参与 + 至少有两个公司没参与 - ...。

二进制枚举,暴力重构矩阵跑高斯消元。

#include<vector>
#include<cstdio>
#include<iostream>
#include<algorithm>
using namespace std;
typedef pair<int, int> pii;
const int MOD = int(1E9) + 7;
const int MAXN = 17;
int pow_mod(int b, int p) {
    int ret = 1;
    while( p ) {
        if( p & 1 ) ret = 1LL*ret*b%MOD;
        b = 1LL*b*b%MOD;
        p >>= 1;
    }
    return ret;
}
int A[MAXN][MAXN];
int gauss(int r, int c) {
    int nwr = 0, nwc = 0, f = 1;
    for(int i=0;i<r;i++)
        for(int j=0;j<c;j++)
            A[i][j] = (A[i][j] + MOD)%MOD;
    while( nwr < r && nwc < c ) {
        int mxr = nwr;
        for(int i=nwr+1;i<r;i++)
            if( A[i][nwc] > A[mxr][nwc] ) mxr = i;
        if( A[mxr][nwc] ) {
            if( mxr != nwr )
                swap(A[mxr], A[nwr]), f *= -1;
            for(int i=nwr+1;i<r;i++)
                while( A[i][nwc] ) {
                    if( A[nwr][nwc] ) {
                        int x = A[i][nwc]/A[nwr][nwc];
                        for(int j=nwc;j<c;j++)
                            A[i][j] = (A[i][j] + MOD - 1LL*x*A[nwr][j]%MOD)%MOD;
                    }
                    swap(A[i], A[nwr]), f *= -1;
                }
            nwr++;
        }
        nwc++;
    }
    int ans = 1;
    for(int i=0;i<r;i++)
        ans = 1LL*ans*A[i][i]%MOD;
    return ans*f;
}
int bits[1<<MAXN];
vector<pii>edge[MAXN];
int main() {
    int N; scanf("%d", &N);
    for(int i=0;i<N-1;i++) {
        int m, x, y; scanf("%d", &m);
        for(int j=1;j<=m;j++) {
            scanf("%d%d", &x, &y);
            edge[i].push_back(make_pair(x-1, y-1));
        }
    }
    int ans = 0, t = (1<<(N-1));
    for(int s=0;s<t;s++) {
        for(int i=0;i<N;i++)
            for(int j=0;j<N;j++)
                A[i][j] = 0;
        for(int i=0;i<N-1;i++)
            if( (1<<i) & s ) {
                for(int j=0;j<edge[i].size();j++) {
                    int u = edge[i][j].first, v = edge[i][j].second;
                    A[u][u]++, A[v][v]++, A[u][v]--, A[v][u]--;
                }
            }
        bits[s] = bits[s>>1] + (s & 1);
        int f = ( (N - bits[s] - 1) & 1 ) ? -1 : 1;
        ans = (ans + f*gauss(N-1, N-1))%MOD;
    }
    printf("%d\n", (ans + MOD)%MOD);
}

原文地址:https://www.cnblogs.com/Tiw-Air-OAO/p/10280720.html

时间: 01-16

@算法 - [email protected] matrix - tree 定理(矩阵树定理)的相关文章

【算法】Matrix - Tree 矩阵树定理 &amp; 题目总结

最近集中学习了一下矩阵树定理,自己其实还是没有太明白原理(证明)类的东西,但想在这里总结一下应用中的一些细节,矩阵树定理的一些引申等等. 首先,矩阵树定理用于求解一个图上的生成树个数.实现方式是:\(A\)为邻接矩阵,\(D\)为度数矩阵,则基尔霍夫(Kirchhoff)矩阵即为:\(K = D - A\).具体实现中,记 \(a\) 为Kirchhoff矩阵,则若存在 \(E(u, v)\) ,则\(a[u][u] ++, a[v][v] ++, a[u][v] --, a[v][u] --\

luoguP3317 [SDOI2014]重建 变元矩阵树定理 + 概率

首先,我们需要求的是 $$\sum\limits_{Tree} \prod\limits_{E \in Tree} E(u, v) \prod\limits_{E \notin Tree} (1 - E(u, v))$$ 我们知道变元矩阵树定理 ---> 不知道请见此 我们自然希望要求和的事物只跟生成树的边有关 因此考虑把$\prod\limits_{E \notin Tree} (1 - E(u, v))$转化为$\prod\limits_{E} (1 - E(u, v)) * \frac{1

@雅礼集训01/10 - [email&#160;protected] matrix

目录 @[email protected] @[email protected] @[email protected] @[email protected] @[email protected] 给定一个矩阵.求它的所有子矩阵中本质不同的行的个数之和. input 第一行,两个正整数 n, m. 第二行,n * m 个正整数,第 i 个数表示 A[i/m][i mod m]. 保证 n * m <= 10^5, 1 <= A[i][j] <= 10^9 output 输出一个非负整数表示

矩阵树定理速证

凯莱公式: spanning_trees_num( G ) = spanning_trees_num( G - e ) + spanning_trees_num( G · e ) 矩阵树定理: G 相应的拉普拉斯矩阵(度矩阵 - 邻接矩阵)L( G )   删除随意一行一列得到的行列式的值det( L*( G ) ) 即生成树的个数,即spanning_trees_num( G ) = det( L*( G ) ) 证: 归纳如果 spanning_trees_num( G - e ) = de

[spoj104][Highways] (生成树计数+矩阵树定理+高斯消元)

In some countries building highways takes a lot of time... Maybe that's because there are many possiblities to construct a network of highways and engineers can't make up their minds which one to choose. Suppose we have a list of cities that can be c

CSU 1805 Three Capitals(矩阵树定理+Best定理)

http://acm.csu.edu.cn/csuoj/problemset/problem?pid=1805 题意: A和B之间有a条边,A和G之间有b条边,B和G之间有c条边.现在从A点出发走遍所有的边,然后再回到A点,问一共有多少种方法. 思路: 16年湖南省赛题目,这道题目是求欧拉回路的个数,和生成树的计数有一定的联系. 首先给出神奇的Best定理,这是什么鬼定理,反正查不到什么有关该定理的文章... $ec(G)=t_s(G)\cdot deg(s)! \cdot \prod_{v\i

【BZOJ4031】[HEOI2015]小Z的房间 矩阵树定理

[BZOJ4031][HEOI2015]小Z的房间 Description 你突然有了一个大房子,房子里面有一些房间.事实上,你的房子可以看做是一个包含n*m个格子的格状矩形,每个格子是一个房间或者是一个柱子.在一开始的时候,相邻的格子之间都有墙隔着. 你想要打通一些相邻房间的墙,使得所有房间能够互相到达.在此过程中,你不能把房子给打穿,或者打通柱子(以及柱子旁边的墙).同时,你不希望在房子中有小偷的时候会很难抓,所以你希望任意两个房间之间都只有一条通路.现在,你希望统计一共有多少种可行的方案.

【bzoj2467】[中山市选2010]生成树 矩阵树定理

题目描述 有一种图形叫做五角形圈.一个五角形圈的中心有1个由n个顶点和n条边组成的圈.在中心的这个n边圈的每一条边同时也是某一个五角形的一条边,一共有n个不同的五角形.这些五角形只在五角形圈的中心的圈上有公共的顶点.如图0所示是一个4-五角形圈. 现在给定一个n五角形圈,你的任务就是求出n五角形圈的不同生成树的数目.还记得什么是图的生成树吗?一个图的生成树是保留原图的所有顶点以及顶点的数目减去一这么多条边,从而生成的一棵树. 注意:在给定的n五角形圈中所有顶点均视为不同的顶点. 输入 输入包含多

【bzoj4894】天赋 矩阵树定理

题目描述 小明有许多潜在的天赋,他希望学习这些天赋来变得更强.正如许多游戏中一样,小明也有n种潜在的天赋,但有一些天赋必须是要有前置天赋才能够学习得到的.也就是说,有一些天赋必须是要在学习了另一个天赋的条件下才能学习的.比如,要想学会"开炮",必须先学会"开枪".一项天赋可能有多个前置天赋,但只需习得其中一个就可以学习这一项天赋.上帝不想为难小明,于是小明天生就已经习得了1号天赋-----"打架".于是小明想知道学习完这n种天赋的方案数,答案对1