BZOJ 3534 [Sdoi2014]重建

Posted 2020-10-24 ws_zzy

题解：矩阵树定理

邻接矩阵-度数矩阵（期望下）

求出来的行列式为所有(生成树边权乘积)的和

每条边边权化为 c/(1-c),最后乘上π(1-c),对1边权特殊处理一下

问题:矩阵树定理不熟，不会证明

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
using namespace std;
const int maxn=60;
const double eps=1e-9;

double tmp=1;

inline int dcmp(double x){
	if(fabs(x)<eps)return 0;
	if(x>0)return 1;
	else return -1;
}
int n;

double A[maxn][maxn];
double G[maxn][maxn];

double Gauss(){
	--n;
	for(int j=1;j<=n;++j){
		int maxline=j;
		for(int i=j+1;i<=n;++i){
			if(fabs(A[i][j])>fabs(A[maxline][j]))maxline=i;
		}
		if(dcmp(A[maxline][j])==0)return 0.0;
		if(maxline!=j)for(int i=j;i<=n;++i)swap(A[j][i],A[maxline][i]);
		for(int i=j+1;i<=n;++i){
			for(int k=j+1;k<=n;++k){
				A[i][k]-=A[i][j]*A[j][k]/A[j][j];
			}
			A[i][j]=0;
		}
	}
	double ret=1;
	for(int i=1;i<=n;++i)ret*=A[i][i];
	return ret*tmp;
}

int main(){
	scanf("%d",&n);
	for(int i=1;i<=n;++i){
		for(int j=1;j<=n;++j){
			scanf("%lf",&G[i][j]);
			if(dcmp(G[i][j]-1)==0)G[i][j]-=eps;
			if(i<j)tmp*=(1-G[i][j]);
			G[i][j]/=(1-G[i][j]);
		}
	}
	for(int i=1;i<=n;++i){
		double sum=0;
		for(int j=1;j<=n;++j){
			A[i][j]=G[i][j];sum+=A[i][j];
		}
//		A[i][i]=-sum;
		A[i][i]=-sum;
	}
	printf("%.9f\n",fabs(Gauss()));
	return 0;
}