【BZOJ 1016】[JSOI2008]最小生成树计数 基尔霍夫矩阵||暴力

第一次写这道题的时候是直接暴力水过去的（貌似正解就是暴力），但是后来看到了关于基尔霍夫矩阵，回过来看感觉十分的有（e）趣（xin），就又写了一遍。

至于基尔霍夫矩阵，证明什么的都去见鬼吧。。。会用就好咯：

设基尔霍夫矩阵为C:1.C[i][i]为i这一个节点的度数

                                 2.C[i][j]表示邻接链表，数值是i和j之间有多少连接的边

基本性质：1.在矩阵前乘以一个系数，相当于对矩阵中每一个元素乘以一个系数

                   2.交换矩阵中任意两行，矩阵等于原来相反数（求解的时候会用到）

然后任意取出一行一列（一般是最后一行一列）得到一个n-1阶的主子式，他的值的绝对值就是原图生成树的个数。具体做法就是高斯消元将其转化为三角形行列式，他主对角线的乘积即为行列式的值。

最后为了避免高斯消元出现精度误差，我们把高斯线元中的除法变成欧几里得算法（就是gcd），这里是需要注意的（我在代码中标出）

至于这道题嘛，基本没什么了，除了一个基尔霍夫矩阵。。。

#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<vector>
#define maxn  50053
#define LL long long
#include<cstdlib>
using namespace std;
int n,m,head[maxn],tot,f[maxn],q[maxn],top,vis[maxn];
LL p=31011,ans=1,mat[105][105],a[105][105];
int find(int x){return x==f[x] ? x : f[x]=find(f[x]);}
int unio(int x,int y){
	int dx=find(x),dy=find(y);
	if(dx==dy)return false;
	f[dx]=f[dy];
	return true;
}
struct edge{
	int a,b,w;
	bool operator <(const edge& y)const{return w<y.w;}
}e[maxn],es[maxn];
vector<int>g[maxn];

void pre_act(){
	for(int i=1;i<=n;i++)f[i]=i;
	for(int i=1;i<=m;i++)scanf("%d%d%d",&e[i].a,&e[i].b,&e[i].w);
	sort(e+1,e+1+m);int cnt=0;
	for(int i=1;i<=m;i++){
		if(unio(e[i].a,e[i].b))cnt++;
	}
	if(cnt<n-1){puts("0");exit(0);}
}

void dfs(int u,int fa){
	q[++top]=u;vis[u]=1;
	for(int i=0;i<g[u].size();i++){
		if(!vis[g[u][i]])dfs(g[u][i],u);
	}
}

void calc(int x){
	top=0;
	dfs(x,x);
	for(int i=1;i<=top;i++){
		for(int j=1;j<=top;j++){
			a[i][j]=(p+mat[q[i]][q[j]])%p;
		}
	}
	
	for(int i=1;i<top&&a[i][i];i++){
		for(int j=i+1;j<top;j++){
			while(a[j][i]){//嗯，就是这里欧几里得 
				LL s=a[i][i]/a[j][i];
				for(int k=1;k<=n;k++){
					a[i][k]=(a[i][k]-a[j][k]*s%p+p)%p;swap(a[i][k],a[j][k]);
				}
				ans=(p-ans)%p;
			}
		}
		ans=ans*a[i][i]%p;
	}
	
}

void solve(){
	for(int i=1;i<=n;i++)f[i]=i;
	for(int j,i=1;i<=m;i=j+1){
		tot=0;
		for(j=i;e[j].w==e[i].w&&j<=m;j++);j--;
		for(int k=1;k<=n;k++)vis[k]=0;
		for(int k=i;k<=j;k++){
			int dx=find(e[k].a),dy=find(e[k].b);
			if(dx!=dy){
				es[++tot]=e[k];
				mat[dx][dx]++,mat[dy][dy]++,mat[dy][dx]--,mat[dx][dy]--;
				g[dx].push_back(dy),g[dy].push_back(dx);
			}
		}
		for(int k=1;k<=tot;k++)if(!vis[find(es[k].a)])calc(find(es[k].a));
		//清楚标记 
		for(int k=1;k<=tot;k++){
			int dx=find(es[k].a),dy=find(es[k].b);
			mat[dx][dx]--,mat[dy][dy]--,mat[dy][dx]++,mat[dx][dy]++;
			g[dx].clear(),g[dy].clear();
		}
		for(int k=1;k<=tot;k++)
			unio(es[k].a,es[k].b);
	}
}

int main(){
	scanf("%d%d",&n,&m);
	pre_act();
	solve();
	printf("%lld",ans);
	return 0;
}