0%

生成树计数

(本帖排版有问题,回头需要重新整理)

生成树计数就是统计一张图中一共有多少种构造生成树的方案。

大概要用到组合数学等等的数学知识。

问题的引入

以下内容均来自NOI2007国家集训队论文 周冬 《生成树的计数及其应用》:


Matrix-Tree定理(Kirchhoff矩阵-树定理)。Matrix-Tree定理是解决生成树计数问题最有力的武器之一。它首先于1847年被Kirchhoff证明。在介绍定理之前,我们首先明确几个概念:

1、G的度数矩阵D[G]是一个n*n的矩阵,并且满足:当i≠j时,dij=0;当i=j时,dij等于vi的度数。

2、G的邻接矩阵A[G]也是一个n*n的矩阵, 并且满足:如果vi、vj之间有边直接相连,则aij=1,否则为0。

我们定义G的Kirchhoff矩阵(也称为拉普拉斯算子)C[G]为C[G]=D[G]-A[G],则Matrix-Tree定理可以描述为:

G的所有不同的生成树的个数等于其Kirchhoff矩阵C[G]任何一个n-1阶主子式的行列式的绝对值。

所谓n-1阶主子式,就是对于r(1≤r≤n),将C[G]的第r行、第r列同时去掉后得到的新矩阵,用Cr[G]表示。

在证明Matrix-Tree定理的过程中,我们使用了无向图G的关联矩阵B。B是一个n行m列的矩阵,行对应点而列对应边。B满足,如果存在一条边e={vi,vj},那在e所对应的列中,vi和vj所对应的那两行,一个为1、另一个为-1,其他的均为0。至于谁是1谁是-1并不重要。

接下来,我们考察BBT。容易证明,(BBT)ij等于B的第i行与第j列的点积。所以,当i=j时,(BBT)ij等于与vi相连的边的个数,即vi的度数;而当i≠j时,如果有一条连接vi、vj的边,则(BBT)ij等于-1,否则等于0。这与Kirchhoff矩阵的定义完全相同!因此,我们得出:C=BBT!也就是说,我们可以将C的问题转化到BBT上来。


这里没办法帖公式,只好省去详细的证明过程了,详见论文原文吧:http://files.cnblogs.com/jcf94/Counting_Spanning_Tree.rar


例题有:

SPOJ P104 Highways

题目就是裸的生成树计数,用Matrix_Tree定理可以直接解出。

然后贴一个Kuang神的模板:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
/* ***********************************************
MYID : Chen Fan
LANG : G++
PROG : Count_Spaning_Tree_From_Kuangbin
************************************************ */

#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <iostream>
#include <math.h>

using namespace std;
const double eps = 1e-8;
const int MAXN = 110;

int sgn(double x)
{
if(fabs(x) < eps)return 0;
if(x < 0)return -1;
else return 1;
}

double b[MAXN][MAXN];
double det(double a[][MAXN],int n)
{
int i, j, k, sign = 0;
double ret = 1;
for(i = 0;i < n;i++)
for(j = 0;j < n;j++) b[i][j] = a[i][j];
for(i = 0;i < n;i++)
{
if(sgn(b[i][i]) == 0)
{
for(j = i + 1; j < n;j++)
if(sgn(b[j][i]) != 0) break;
if(j == n)return 0;
for(k = i;k < n;k++) swap(b[i][k],b[j][k]);
sign++;
}
ret *= b[i][i];
for(k = i + 1;k < n;k++) b[i][k]/=b[i][i];

for(j = i+1;j < n;j++)
for(k = i+1;k < n;k++) b[j][k] -= b[j][i]*b[i][k];
}
if(sign & 1)ret = -ret;
return ret;
}

double a[MAXN][MAXN];
int g[MAXN][MAXN];

int main()
{
int T;
int n,m;
int u,v;
scanf("%d",&T);
while(T--)
{
scanf("%d%d",&n,&m);
memset(g,0,sizeof(g));
while(m--)
{
scanf("%d%d",&u,&v);
u--;v--;
g[u][v] = g[v][u] = 1;
}
memset(a,0,sizeof(a));
for(int i = 0;i < n;i++)
for(int j = 0;j < n;j++)
if(i != j && g[i][j])
{
a[i][i]++;
a[i][j] = -1;
}
double ans = det(a,n-1);
printf("%.0lf\n",ans);
}
return 0;
}

最小生成树计数:

在ACM2012金华赛区的网赛中曾经出现过这么一道题,现在是HDU 4408,统计最小生成树的个数。

方法是Kruskal+Matrix_Tree定理,也就是把上面的算法加上Kruskal:

来自:http://blog.csdn.net/jarily/article/details/8902402 的算法解释:

  • 算法引入:
  • 给定一个含有N个结点M条边的无向图,求它最小生成树的个数t(G);
  • 算法思想:
  • 抛开“最小”的限制不看,如果只要求求出所有生成树的个数,是可以利用Matrix-Tree定理解决的;
  • Matrix-Tree定理此定理利用图的Kirchhoff矩阵,可以在O(N3)时间内求出生成树的个数;
  • kruskal算法:
  • 将图G={V,E}中的所有边按照长度由小到大进行排序,等长的边可以按照任意顺序;
  • 初始化图G’为{V,Ø},从前向后扫描排序后的边,如果扫描到的边e在G’中连接了两个相异的连通块,则将它插入G’中;
  • 最后得到的图G’就是图G的最小生成树;
  • 由于kruskal按照任意顺序对等长的边进行排序,则应该将所有长度为L0的边的处理当作一个阶段来整体看待;
  • 令kruskal处理完这一个阶段后得到的图为G0,如果按照不同的顺序对等长的边进行排序,得到的G0也是不同;
  • 虽然G0可以随排序方式的不同而不同,但它们的连通性都是一样的,都和F0的连通性相同(F0表示插入所有长度为L0的边后形成的图);
  • 在kruskal算法中的任意时刻,并不需要关注G’的具体形态,而只要关注各个点的连通性如何(一般是用并查集表示);
  • 所以只要在扫描进行完第一阶段后点的连通性和F0相同,且是通过最小代价到达这一状态的,接下去都能找到最小生成树;
  • 经过上面的分析,可以看出第一个阶段和后面的工作是完全独立的;
  • 第一阶段需要完成的任务是使G0的连通性和F0一样,且只能使用最小的代价;
  • 计算出这一阶段的方案数,再乘上完成后续事情的方案数,就是最终答案;
  • 由于在第一个阶段中,选出的边数是一定的,所有边的长又都为L0;
  • 所以无论第一个阶段如何进行代价都是一样的,那么只需要计算方案数就行了;
  • 此时Matrix-Tree定理就可以派上用场了,只需对F0中的每一个连通块求生成树个数再相乘即可;

同样,贴一个bin神模板:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
/* ***********************************************
MYID : Chen Fan
LANG : G++
PROG : Counting_MST_From_Kuangbin_HDU4408
************************************************ */

#include <map>
#include <stack>
#include <queue>
#include <math.h>
#include <vector>
#include <string>
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <iostream>
#include <algorithm>
#define N 405
#define M 4005
#define E
#define inf 0x3f3f3f3f
#define dinf 1e10
#define linf (LL)1<<60
#define LL long long
#define clr(a,b) memset(a,b,sizeof(a))
using namespace std;

LL mod;
struct Edge
{
int a,b,c;
bool operator<(const Edge & t)const
{
return c<t.c;
}
}edge[M];
int n,m;
LL ans;
int fa[N],ka[N],vis[N];
LL gk[N][N],tmp[N][N];
vector<int>gra[N];

int findfa(int a,int b[]){return a==b[a]?a:b[a]=findfa(b[a],b);}

LL det(LL a[][N],int n)
{
for(int i=0;i<n;i++)for(int j=0;j<n;j++)a[i][j]%=mod;
long long ret=1;
for(int i=1;i<n;i++)
{
for(int j=i+1;j<n;j++)
while(a[j][i])
{
LL t=a[i][i]/a[j][i];
for(int k=i;k<n;k++)
a[i][k]=(a[i][k]-a[j][k]*t)%mod;
for(int k=i;k<n;k++)
swap(a[i][k],a[j][k]);
ret=-ret;
}
if(a[i][i]==0)return 0;
ret=ret*a[i][i]%mod;
//ret%=mod;
}
return (ret+mod)%mod;
}

int main()
{
while(scanf("%d%d%I64d",&n,&m,&mod)==3)
{

if(n==0 && m==0 && mod==0)break;

memset(gk,0,sizeof(gk));
memset(tmp,0,sizeof(tmp));
memset(fa,0,sizeof(fa));
memset(ka,0,sizeof(ka));
memset(tmp,0,sizeof(tmp));

for(int i=0;i<N;i++)gra[i].clear();
for(int i=0;i<m;i++)
scanf("%d%d%d",&edge[i].a,&edge[i].b,&edge[i].c);
sort(edge,edge+m);
for(int i=1;i<=n;i++)fa[i]=i,vis[i]=0;
int pre=-1;
ans=1;
for(int h=0;h<=m;h++)
{
if(edge[h].c!=pre||h==m)
{
for(int i=1;i<=n;i++)
if(vis[i])
{
int u=findfa(i,ka);
gra[u].push_back(i);
vis[i]=0;
}
for(int i=1;i<=n;i++)
if(gra[i].size()>1)
{
for(int a=1;a<=n;a++)
for(int b=1;b<=n;b++)
tmp[a][b]=0;
int len=gra[i].size();
for(int a=0;a<len;a++)
for(int b=a+1;b<len;b++)
{
int la=gra[i][a],lb=gra[i][b];
tmp[a][b]=(tmp[b][a]-=gk[la][lb]);
tmp[a][a]+=gk[la][lb];tmp[b][b]+=gk[la][lb];
}
long long ret=(long long)det(tmp,len);
ret%=mod;
ans=(ans*ret%mod)%mod;
for(int a=0;a<len;a++)fa[gra[i][a]]=i;
}
for(int i=1;i<=n;i++)
{
ka[i]=fa[i]=findfa(i,fa);
gra[i].clear();
}
if(h==m)break;
pre=edge[h].c;
}
int a=edge[h].a,b=edge[h].b;
int pa=findfa(a,fa),pb=findfa(b,fa);
if(pa==pb)continue;
vis[pa]=vis[pb]=1;
ka[findfa(pa,ka)]=findfa(pb,ka);
gk[pa][pb]++;gk[pb][pa]++;
}
int flag=0;
for(int i=2;i<=n&&!flag;i++)if(ka[i]!=ka[i-1])flag=1;
ans%=mod;
printf("%I64d\n",flag?0:ans);
}
return 0;
}

然后,最后这一题:

HDU4305

综合性相当强的一题:

首先需要用计算几何的知识对平面上的点构图,形成图之后用Matrix-Tree定理求解生成树的个数。

但是题目数据比较大,对行列式求值的时候,需要用到高斯消元求上三角阵,行列式的值等于对角线元素的积。
由于是整数然后再mod。消元时需要求最小公倍数。还需要拓展欧几里得算法求逆元。

日后等有人能集我们队三者之大成之后,应该就可以挑战这道题了。