Matrix - Tree 定理学习笔记

前置芝士

Matrix - Tree 定理，中文名：矩阵树定理。

度数矩阵

对于无向图，度数矩阵 $D$ 的定义是：

D_{i,j}= \begin{cases} deg_i &\text{if } i=j \\ 0 &\text{if } i\not =j \\ \end{cases}

对于有向图，分入度矩阵 $D^{in}$ 和出度矩阵 $D^{out}$，定义如下：

D^{in}_{i,j}= \begin{cases} in_i &\text{if } i=j \\ 0 &\text{if } i\not =j \\ \end{cases}

D^{out}_{i,j}= \begin{cases} out_i &\text{if } i=j \\ 0 &\text{if } i\not =j \\ \end{cases}

邻接矩阵

对于无向图，邻接矩阵 $A$ 的定义为：

$$A_{i,j}=A_{j,i}=\mathrm{\#}e(i,j)$$

其中 $\mathrm{\#}e(i,j)$ 表示 $i,j$ 两点之间所连边数。

对于有向图，邻接矩阵 $A$ 的定义为：

$$A_{i,j}=\mathrm{\#}e(i,j)$$

其中 $\mathrm{\#}e(i,j)$ 表示从点 $i$ 指向 $j$ 的边数。

拉普拉斯矩阵

说到拉普拉斯就不得不提到一位恶魔系拉普拉斯妖学妹

对于无向图，拉普拉斯矩阵 $L$ 的定义为：

$$L=D-A$$

对于有向图，入度拉普拉斯矩阵 $L^{in}$ 和出度拉普拉斯矩阵 $L^{out}$ 的定义为：

$$L^{in}=D^{in}-A$$

$$L^{out}=D^{out}-A$$

Matrix - Tree 定理

对于无向图，其生成树的个数为（$k\in [1,n]$)：

$$\text{det}{L}_{\left[n\right]\setminus \left\{k\right\},\left[n\right]\setminus \left\{k\right\}}$$

其中 $\text{det}$ 表示行列式求值，$L_{\left[n\right]\setminus \left\{k\right\},\left[n\right]\setminus \left\{k\right\}}$ 表示将矩阵 $L$ 删除第 $k$ 行与第 $k$ 列后得到的矩阵。

对于有向图，其以 $k$ 为根的内向生成树（根向生成树）的个数为：

$$\text{det}{L}^{{out}_{\left[n\right]\setminus \left\{k\right\},\left[n\right]\setminus \left\{k\right\}}}$$

其以 $k$ 为根的外向生成树（叶向生成树）的个数为：

$$\text{det}{L}^{{in}_{\left[n\right]\setminus \left\{k\right\},\left[n\right]\setminus \left\{k\right\}}}$$

带权拓展

问题：对于一个带权图，定义一个生成树的权值为其所有边权之积，求其所有生成树权值之和。

这里我们修改一下度数矩阵的定义即可，不在存储边数，而是所有边边权之和，按照这种方式得到拉普拉斯矩阵后其余和上面的式子是一模一样的。

例题

P6178 【模板】Matrix-Tree 定理

带权拓展的模板。注意，这里有向图生成树的定义是以 $1$ 为根的外向树，所以我们要看入度，并删除拉普拉斯矩阵的第一行、第一列。

代码：

C++

#include
using namespace std;

#define int long long
const int N=305;
const int mod=1e9+7;
int gh[N][N];

int det(int n){
    int res=1,f=1;
    for(int i=1;i>n>>m>>t;
    for(int i=1,u,v,w;i>u>>v>>w;
        u--,v--;
        if(t) (gh[v][v]+=w)%=mod,(((gh[u][v]-=w)%=mod)+=mod)%=mod;
        else (gh[v][v]+=w)%=mod,(((gh[u][v]-=w)%=mod)+=mod)%=mod,
             (gh[u][u]+=w)%=mod,(((gh[v][u]-=w)%=mod)+=mod)%=mod;
    }
    cout

`#include
using namespace std;

#define int long long
const int N=305;
const int mod=1e9+7;
int gh[N][N];

int det(int n){
    int res=1,f=1;
    for(int i=1;in;i++){
        for(int j=i+1;jn;j++){
            while(gh[j][i]){
                int t=gh[i][i]/gh[j][i];
                for(int k=i;kn;k++) (((gh[i][k]-=gh[j][k]*t)%=mod)+=mod)%=mod;
                swap(gh[i],gh[j]);
                f=-f;
            }
        }
        if(!gh[i][i]) return 0;
        (res*=gh[i][i])%=mod;
    }
    return ((res*f)%mod+mod)%mod;
}

signed main(){
    int n,m,t;
    cin>>n>>m>>t;
    for(int i=1,u,v,w;im;i++){
        cin>>u>>v>>w;
        u--,v--;
        if(t) (gh[v][v]+=w)%=mod,(((gh[u][v]-=w)%=mod)+=mod)%=mod;
        else (gh[v][v]+=w)%=mod,(((gh[u][v]-=w)%=mod)+=mod)%=mod,
             (gh[u][u]+=w)%=mod,(((gh[v][u]-=w)%=mod)+=mod)%=mod;
    }
    coutdet(n-1);
    return 0;
}`

P3317 [SDOI2014] 重建

推一下式子：

\begin{aligned} P&=\sum_{T \text{ 是生成树}} \left( \prod_{e \in T} p_e \prod_{e \notin T}(1 - p_e)\right) \\ &=\left( \prod_{e \in E}(1 - p_e) \right) \cdot \sum_{T \text{ 是生成树}} \left( \prod_{e \in T}\frac{p_e}{1 - p_e}\right) \end{aligned}

后面就是矩阵树定理的板子了，代码：

C++

#include
using namespace std;

#define double long double
const int N=105;
const double eps=1e-9;
double gh[N][N],res=1,f=1;

void det(int n){
    for(int i=1;ifabs(gh[now][i])) now=j;
        swap(gh[i],gh[now]);
        if(i!=now) f=-f;
        if(gh[i][i]>n;
    for(int u=1;u>w;
        if(w>1-eps) w=1-eps;
        gh[v][v]+=w/(1-w),gh[u][v]-=w/(1-w);
        if(u>v) res*=(1-w);
    }
    det(n-1);
    cout

`#include
using namespace std;

#define double long double
const int N=105;
const double eps=1e-9;
double gh[N][N],res=1,f=1;

void det(int n){
    for(int i=1;in;i++){
        int now=i;
        for(int j=i+1;jn;j++) if(fabs(gh[j][i])>fabs(gh[now][i])) now=j;
        swap(gh[i],gh[now]);
        if(i!=now) f=-f;
        if(gh[i][i]eps){
            res=0;
            return;
        }
        for(int j=i+1;jn;j++){
            if(gh[j][i]==0) continue;
            double t=gh[j][i]/gh[i][i];
            for(int k=i;kn;k++) gh[j][k]-=t*gh[i][k];
        }
        res*=gh[i][i];
    }
}

signed main(){
    int n;
    cin>>n;
    for(int u=1;un;u++) for(int v=1;vn;v++){
        double w;
        cin>>w;
        if(w>1-eps) w=1-eps;
        gh[v][v]+=w/(1-w),gh[u][v]-=w/(1-w);
        if(u>v) res*=(1-w);
    }
    det(n-1);
    coutfixedsetprecision(5)res*f;
    return 0;
}`

P4455 [CQOI2018] 社交网络

外向生成树计数板子题。

代码：

C++

#include
using namespace std;

#define int long long
const int N=305;
const int mod=1e4+7;
int gh[N][N];

int det(int n){
    int res=1,f=1;
    for(int i=1;i>n>>m;
    for(int i=1,u,v,w;i>v>>u;
        u--,v--;
        (gh[v][v]+=1)%=mod,(((gh[u][v]-=1)%=mod)+=mod)%=mod;
    }
    cout

`#include
using namespace std;

#define int long long
const int N=305;
const int mod=1e4+7;
int gh[N][N];

int det(int n){
    int res=1,f=1;
    for(int i=1;in;i++){
        for(int j=i+1;jn;j++){
            while(gh[j][i]){
                int t=gh[i][i]/gh[j][i];
                for(int k=i;kn;k++) (((gh[i][k]-=gh[j][k]*t)%=mod)+=mod)%=mod;
                swap(gh[i],gh[j]);
                f=-f;
            }
        }
        if(!gh[i][i]) return 0;
        (res*=gh[i][i])%=mod;
    }
    return ((res*f)%mod+mod)%mod;
}

signed main(){
    int n,m;
    cin>>n>>m;
    for(int i=1,u,v,w;im;i++){
        cin>>v>>u;
        u--,v--;
        (gh[v][v]+=1)%=mod,(((gh[u][v]-=1)%=mod)+=mod)%=mod;
    }
    coutdet(n-1);
    return 0;
}`

P4336 [SHOI2016] 黑暗前的幻想乡

考虑容斥原理：

$$\mathrm{所}\mathrm{有}\mathrm{公}\mathrm{司}\mathrm{都}\mathrm{没}\mathrm{有}\mathrm{被}\mathrm{漏}\mathrm{掉}\mathrm{的}\mathrm{生}\mathrm{成}\mathrm{树}\mathrm{方}\mathrm{案}\mathrm{数}=\sum _{S\subseteq U}(-1{)}^{\left|S\right|}\times (S\mathrm{内}\mathrm{的}\mathrm{公}\mathrm{司}\mathrm{被}\mathrm{漏}\mathrm{掉}\mathrm{的}\mathrm{生}\mathrm{成}\mathrm{树}\mathrm{方}\mathrm{案}\mathrm{数})$$

之后就是矩阵树定理了。

代码：

C++

#include
using namespace std;

#define int long long
const int N=20;
const int mod=1e9+7;
vector> ve[N];
int gh[N][N],ans;

int det(int n){
    int res=1,f=1;
    for(int i=1;i>n;
    for(int i=1;i>m;
        for(int j=1;j>u>>v;
            ve[i].push_back({u,v});
        }
    }
    for(int s=0;s>i)&1)) continue;
            cnt++;
            for(auto[u,v]:ve[i+1]) (gh[u][u]+=1)%=mod,(gh[v][v]+=1)%=mod,(gh[u][v]+=mod-1)%=mod,(gh[v][u]+=mod-1)%=mod;
        }
        ans=((ans+((n-cnt-1)&1?-1:1)*det(n-1))%mod+mod)%mod;
    }
    cout

`#include
using namespace std;

#define int long long
const int N=20;
const int mod=1e9+7;
vectorpairint,int>> ve[N];
int gh[N][N],ans;

int det(int n){
    int res=1,f=1;
    for(int i=1;in;i++){
        for(int j=i+1;jn;j++){
            while(gh[j][i]){
                int t=gh[i][i]/gh[j][i];
                for(int k=i;kn;k++) (((gh[i][k]-=gh[j][k]*t)%=mod)+=mod)%=mod;
                swap(gh[i],gh[j]);
                f=-f;
            }
        }
        if(!gh[i][i]) return 0;
        (res*=gh[i][i])%=mod;
    }
    return ((res*f)%mod+mod)%mod;
}

signed main(){
    int n;
    cin>>n;
    for(int i=1;in;i++){
        int m;
        cin>>m;
        for(int j=1;jm;j++){
            int u,v;
            cin>>u>>v;
            ve[i].push_back({u,v});
        }
    }
    for(int s=0;s(1(n-1));s++){
        memset(gh,0,sizeof gh);
        int cnt=0;
        for(int i=0;in-1;i++){
            if(!((s>>i)&1)) continue;
            cnt++;
            for(auto[u,v]:ve[i+1]) (gh[u][u]+=1)%=mod,(gh[v][v]+=1)%=mod,(gh[u][v]+=mod-1)%=mod,(gh[v][u]+=mod-1)%=mod;
        }
        ans=((ans+((n-cnt-1)&1?-1:1)*det(n-1))%mod+mod)%mod;
    }
    coutans;
    return 0;
}`

P2143 [JSOI2010] 巨额奖金 & P4208 [JSOI2008] 最小生成树计数

考虑所有的最小生成树的边权序列排序后都是一致的，所以我们先跑一遍最小生成树，对于这个最小生成树的每一种边权 $w$，我们在新图中先把原最小生成树中边权不为 $w$ 的边上的点进行缩点，再按照原图中所有边权为 $w$ 的边建图，跑矩阵树求生成树个数即可。

代码：

C++

#include
using namespace std;

#define int long long
const int N=105,M=1005;
const int mod=31011;
int n,m,ans=1;
int gh[N][N],nd[N];
map mp;

struct edge{
    int u,v,w;
    bool operator mst;

struct DSU{
    int fa[N],siz[N];
    void init(){
        for(int i=1;isiz[fy]) swap(fx,fy);
        return fa[fy]=fx,1;
    }
}dsu;

int det(int n){
    int res=1,f=1;
    for(int i=1;i>n>>m;
    for(int i=1;i>e[i].u>>e[i].v>>e[i].w;
    if(!kk()){
        cout

`#include
using namespace std;

#define int long long
const int N=105,M=1005;
const int mod=31011;
int n,m,ans=1;
int gh[N][N],nd[N];
mapint,bool> mp;

struct edge{
    int u,v,w;
    bool operatorconst edge b)const{
        return wb.w;
    }
}e[M];

vectoredge> mst;

struct DSU{
    int fa[N],siz[N];
    void init(){
        for(int i=1;in;i++) fa[i]=i,siz[i]=1;
    }
    int find(int x){
        while(x!=fa[x]) x=fa[x]=fa[fa[x]];
        return x;
    }
    bool merge(int x,int y){
        int fx=find(x),fy=find(y);
        if(fx==fy) return 0;
        if(siz[fx]>siz[fy]) swap(fx,fy);
        return fa[fy]=fx,1;
    }
}dsu;

int det(int n){
    int res=1,f=1;
    for(int i=1;in;i++){
        for(int j=i+1;jn;j++){
            while(gh[j][i]){
                int t=gh[i][i]/gh[j][i];
                for(int k=i;kn;k++) (((gh[i][k]-=gh[j][k]*t)%=mod)+=mod)%=mod;
                swap(gh[i],gh[j]);
                f=-f;
            }
        }
        if(!gh[i][i]) return 0;
        (res*=gh[i][i])%=mod;
    }
    return ((res*f)%mod+mod)%mod;
}

bool kk(){
    int res=0;
    dsu.init();
    sort(e+1,e+1+m);
    for(int i=1;im;i++){
        if(!dsu.merge(e[i].u,e[i].v)) continue;
        mp[e[i].w]=1,res++;
        mst.push_back(e[i]);
        if(res==n-1) return 1;
    }
    return 0;
}

signed main(){
    cin>>n>>m;
    for(int i=1;im;i++) cin>>e[i].u>>e[i].v>>e[i].w;
    if(!kk()){
        cout0;
        return 0;
    }
    for(auto[val,tag]:mp){
        int idx=0;
        dsu.init();
        memset(gh,0,sizeof gh);
        for(auto[u,v,w]:mst){
            if(w==val) continue;
            dsu.merge(u,v);
        }
        for(int i=1;in;i++) if(dsu.find(i)==i) nd[i]=++idx;
        for(int i=1;in;i++) if(dsu.find(i)!=i) nd[i]=nd[dsu.find(i)];
        for(int i=1;im;i++){
            if(e[i].w!=val) continue;
            int x=nd[e[i].u],y=nd[e[i].v];
            (gh[x][x]+=1)%=mod,(gh[y][y]+=1)%=mod,
            (gh[x][y]+=mod-1)%=mod,(gh[y][x]+=mod-1)%=mod;
        }
        (((ans*=det(idx-1))%=mod)+=mod)%=mod;
    }
    coutans;
    return 0;
}`