矩陣樹定理 Matrix Tree

　　
​　　矩陣樹定理主要用於圖的生成樹計數。
　　
　　看到給出圖求生成樹的這類問題就大概要往這方面想了。
　　
　　演算法會根據圖構造出一個特殊的基爾霍夫矩陣AA，接著根據矩陣樹定理，用AA計算出生成樹個數。
　　
　　
　　

1.無向圖的生成樹計數

　　
　　對於給定的可含重邊的連通無向圖GG，求其生成樹的個數。求法如下：
　　
　　定義度數矩陣DD：該矩陣僅在對角線上有值，Di,iDi,i表示ii號點的度數。對於圖中每一條無向邊(u,v)(u,v)，Du,uDu,u++，Dv,vDv,v++。
　　
　　定義鄰接矩陣CC：Ci,jCi,j表示ii到jj的邊數。對於圖中每一條無向邊(u,v)(u,v)，Cu,vCu,v++，Cv,uCv,u++。
　　
　　定義圖GG的基爾霍夫矩陣

A=D−CA=D−C。
　　
​　　矩陣樹定理：將AA去掉第ii行和第ii列（i∈[1,n]i∈[1,n]），將它當做一個行列式求解，則det(A)det(A)就是生成樹個數。
　　
　　　　
　　

2.有向圖的樹形圖計數

　　
　　對於有向圖，不存在“生成樹”的概念，但存在“樹形圖”的概念。有向圖中，若選定一個點作為樹根，能構造出一棵“樹”（包含n−1n−1條邊）使得根能到達任意節點，則這是一棵外向樹；若能構造出一棵“樹”使得任意節點能到達根，則這是一棵內向樹。
　　
　　定義度數矩陣DD：該矩陣僅在對角線上有值，Di,iDi,i表示ii號點的度數。
　　
　　　　　　　　　　　對於圖中每一條有向邊(u,v)(u,v)，若構造外向樹則Dv,vDv,v++；若構造內向樹則Du,uDu,u++。
　　
　　定義鄰接矩陣CC：Ci,jCi,j表示ii到jj的邊數。對於圖中每一條有向邊(u,v)(u,v)，Cu,vCu,v++。
　　
​　　定義圖GG的基爾霍夫矩陣

A=D−CA=D−C。
　　
　　矩陣樹定理：將AA去掉第ii行和第ii列（i∈[1,n]i∈[1,n]），將它當做一個行列式求解，則det(A)det(A)就是以ii為根的外向/內向樹形圖個數。很多時候我們會發現AA的對角線上某數為Ai,i=0Ai,i=0，刪去第ii行和第ii列可以幹掉0。只有這樣行列式才不等於0，其實也就是說只能從ii出發有解了。
　　
　　　
　　

3.細節

　　
　　求行列式的方法是：將行列式通過行列式初等變換消成上三角。此時對角線乘積即為行列式的值。
　　
　　注意，矩陣樹定理這一套演算法會考慮如何把所有的參與點構建成生成樹，所以編號不要跳躍。如果說有障礙之類的元素，千萬不要在矩陣中給它留一行一列，因為這一行一列都一定是0，演算法會嘗試將“障礙”構建進生成樹，最後只能得到無解。
　　
　　
　　
　　
　　

一些例題

BZOJ 4031

　　
　　傳送門在此
　　
　　這是一道無向圖生成樹計數的裸題，直接上基礎演算法即可。
　　
　　這裡就要注意障礙的處理了。我們應該對非障礙格子重新標號使得它們的編號連續，保證演算法正常進行。
　　
　　這題的模數真的噁心，不是質數，沒法用逆元。所以使用類輾轉相除法來消元，每行的消元從O(n)O(n)變成O(nlg)O(nlg)。
　　

#include <cstdio>
using namespace std;
const int N=10,MOD=1e9;
int n,m,id[N][N],idcnt,a[N*N][N*N];
char map[N][N];
inline void swap(int &x,int &y){x^=y^=x^=y;}
inline int plus(int x,int y){return (x+y)%MOD;}
inline int mul(int x,int y){return 1LL*x*y%MOD;}
inline bool ok(int x,int y){return 1<=x&&x<=n&&1<=y&&y<=m&&map[x][y]=='.';}
void addEdge(int u,int v){
    a[u][u]++; a[v][v]++;
    a[u][v]--; a[v][u]--;
}
int solve(){
    if(idcnt==1) return 1;
    int all=idcnt-1,res=1;
    for(int i=1;i<=all;i++)
        for(int j=i+1;j<=all;j++)
            while(a[j][i]){
                int t=a[i][i]/a[j][i];
                for(int k=i;k<=all;k++){
                    int q=plus(a[i][k],-mul(a[j][k],t));
                    a[i][k]=a[j][k];
                    a[j][k]=q;
                }
                res=-res;
            }
    for(int i=1;i<=all;i++) res=mul(res,a[i][i]);
    return (res+MOD)%MOD;
}
int main(){
    scanf("%d%d",&n,&m);
    for(int i=1;i<=n;i++) scanf("%s",map[i]+1);
    for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
            if(map[i][j]=='.') id[i][j]=++idcnt;
    for(int i=1;i<=n;i++)
        for(int j=1;j<=m;j++)
            if(ok(i,j)){
                if(ok(i+1,j))
                    addEdge(id[i][j],id[i+1][j]);
                if(ok(i,j+1))
                    addEdge(id[i][j],id[i][j+1]);
            }
    printf("%d\n",solve());
    return 0;
}

BZOJ 4894

　　
　　傳送門
　　
　　這是一道有向圖樹形圖計數。要求以1號點為根的外向樹形圖個數。
　　
　　按照上述做法直接寫即可。刪去AA的第1行第1列，因為1號點沒有入邊，若不刪第一行第一列行列式值為0，無法計算。
　　

#include <cstdio>
using namespace std;
const int N=305,MOD=1e9+7;
int n,a[N][N];
char str[N];
inline int mul(int x,int y){return 1LL*x*y%MOD;}
inline int plus(int x,int y){return (x+y)%MOD;}
inline void swap(int &x,int &y){x^=y^=x^=y;}
int ksm(int x,int y){
    int res=1;
    for(;y;x=mul(x,x),y>>=1)
        if(y&1) res=mul(res,x);
    return res;
}
int gaussian(){
    int res=1,size=n-1;
    for(int i=1;i<size;i++){
        if(!a[i][i]){
            int l;
            for(l=i+1;l<=size;l++)
                if(a[l][i]) break;
            if(l<=size&&a[l][i]){
                for(int j=i;j<=size;j++) swap(a[l][j],a[i][j]);
                res=-res;
            }
            else return 0;
        }
        for(int j=i+1;j<=size;j++){
            int t=mul(a[j][i],ksm(a[i][i],MOD-2));
            for(int k=i;k<=size;k++)         
                a[j][k]=plus(a[j][k],-mul(a[i][k],t));
        }
    }
    for(int i=1;i<=size;i++) res=mul(res,a[i][i]);
    return plus(res,MOD);
}
int main(){
    scanf("%d",&n);
    for(int i=1;i<=n;i++){
        scanf("%s",str+1);
        for(int j=1;j<=n;j++)
            if(str[j]=='1') 
                a[j][j]++,a[i][j]--;
    }
    for(int i=1;i<n;i++)
        for(int j=1;j<n;j++) a[i][j]=a[i+1][j+1];
    printf("%d\n",gaussian());
    return 0;
}

BZOJ 4596

　　
　　傳送門
　　
​　　我開始不會做了。
　　
​　　我們發現如果將所有公司提供的邊都加進圖中，然後求生成樹個數，是無法限制“每個公司至少要建一條”這個條件的。有的生成樹可能只有一部分公司參與，比如說某一種生成樹只含有SS集合的公司。
　　
　　如果僅加入SS集合所含公司的邊，我們發現這些答案也會被統計到。既然有重複統計，可以考慮消除嗎？
　　
​　　於是容斥的思想就體現出來了！
　　
　　我們列舉加入哪一些公司，分別求生成樹個數。記參與公司集合為SS時生成樹個數為f(S)f(S)，記有xx個公司參與時的生成樹總方案為中有個Ax=∑S中有x個f(S)中有個Ax=∑S中有x個f(S)，則

Ans=An−An−1+An−2−An−3......Ans=An−An−1+An−2−An−3......

​　　複雜度為O(2nn3)O(2nn3)，其實是可以跑的過的。
　　

#include <cstdio>
#include <vector>
#define mp make_pair
#define pb push_back
using namespace std;
typedef pair<int,int> pii;
const int N=18,MOD=1e9+7;
int n;
int a[N][N];
vector<pii> l[N];
vector<int> b[N];
inline bool in(int i,int j){return (i>>(j-1))&1;}
inline int plus(int x,int y){return (x+y)%MOD;}
inline int mul(int x,int y){return 1LL*x*y%MOD;}
inline void swap(int &x,int &y){x^=y^=x^=y;}
void clear_mat(){
    for(int i=1;i<=n;i++)
        for(int j=1;j<=n;j++) a[i][j]=0;
}
void add_company(int id){
    for(int i=0,sz=l[id].size();i<sz;i++){
        int u=l[id][i].first,v=l[id][i].second;
        a[u][u]++; a[v][v]++;
        a[u][v]--; a[v][u]--;
    }
}
int ksm(int x,int y){
    int res=1;
    for(;y;x=mul(x,x),y>>=1)
        if(y&1) res=mul(res,x);
    return res;
}
int gaussian(){
    int s=n-1,res=1;
    for(int i=1;i<=s;i++){
        if(!a[i][i]){
            int t;
            for(t=i+1;t<=s&&!a[t][i];t++);
            if(t>s) return 0;
            for(int j=i;j<=s;j++) swap(a[i][j],a[t][j]);
            res=-res;
        }
        int inv=ksm(a[i][i],MOD-2);
        for(int j=i+1;j<=s;j++){
            int t=mul(a[j][i],inv);
            for(int k=i;k<=s;k++)
                a[j][k]=plus(a[j][k],-mul(a[i][k],t));
        }
    }
    for(int i=1;i<=s;i++) res=mul(res,a[i][i]);
    return res;
}
int main(){
    scanf("%d",&n);
    for(int i=1,m;i<n;i++){
        scanf("%d",&m);
        for(int j=1,u,v;j<=m;j++){
            scanf("%d%d",&u,&v);
            l[i].pb(mp(u,v));   
        }
    }
    int all=1<<(n-1);
    for(int i=1;i<all;i++){
        int cnt=0;
        for(int j=1;j<=n;j++) 
            cnt+=in(i,j);
        b[cnt].pb(i);   
    }
    int ans=0;
    for(int i=n-1,r=1;i>=1;i--,r=-r)
        for(int j=0,sz=b[i].size();j<sz;j++){
            clear_mat();
            int st=b[i][j];
            for(int c=1;c<n;c++)
                if(in(st,c)) add_company(c);
            ans=plus(ans,gaussian()*r);             
        }
    ans=plus(ans,MOD);
    printf("%d\n",ans);
    return 0;
}

總結

　　
　　矩陣樹定理本身還是挺簡單的，但願自己不要忘得太快......
　　
　　但是要靈活運用（廢話）。
　　
​　　如果要深入透徹，我還是得研究一下矩陣樹定理的證明。不過就當一個大坑先留著吧。