【BZOJ3534】【Luogu P3317】 [SDOI2014]重建 变元矩阵树，高斯消元

题解看这里，主要想说一下以前没见过的变元矩阵树还有前几个题见到的几个小细节。

邻接矩阵是可以带权值的。求所有生成树边权和的时候我们有一个基尔霍夫矩阵，是度数矩阵减去邻接矩阵。而所谓变元矩阵树实际上就是把度数矩阵和邻接矩阵带权化，也就是度数矩阵变成该点连接的所有边的权值和，邻接矩阵变成边权矩阵，剩下的依然是求一个行列式。变元矩阵树求的是所有可能生成树的边权之积。

值得注意的点：

• 交换两行，行列式取反。在\(double\)存矩阵的时候可以最后取对角线乘积的绝对值，但如果答案要取膜就需要套上一个辗转相除来解这个矩阵，这时就要在交换行时更新答案，对答案取反处理。

• 求行列式的时候要随便去掉一行和一列，比如去掉最后一行和最后一列就可以。可以传一个\(n-1\)进去避免写错。

• 推式子也很重要。矩阵树定理维护的东西是可以转化为一个式子的，有时候要把它抽象出来。

#include <bits/stdc++.h>
using namespace std;

const int N = 50;
const double eps = 1e-8;

int n; double k = 1, p[N][N], mat[N][N];

double gauss (int n) {
	double ret = 1;
	for (int i = 1; i <= n; ++i) {
		int besti = i;
		for (int j = i; j <= n; ++j) {
			if (fabs (mat[besti][i]) < fabs (mat[j][i])) {
				besti = j;
			}
		}
		if (i != besti) {
			ret = -ret;
			swap (mat[i], mat[besti]);
		}
		for (int j = i + 1; j <= n; ++j) {
			if (fabs (mat[j][i]) > eps) {
				double d = mat[j][i] / mat[i][i];
				for (int k = i; k <= n; ++k) {
					mat[j][k] -= mat[i][k] * d;
				}
			}
		}
		ret *= mat[i][i];
	}
	return ret;
}

int main () {
	cin >> n;
	for (int i = 1; i <= n; ++i) {
		for (int j = 1; j <= n; ++j) {
			cin >> p[i][j];
			if (i != j) {
				p[i][j] = max (p[i][j], 0.0 + eps);
				p[i][j] = min (p[i][j], 1.0 - eps);
				if (i < j) k *= (1 - p[i][j]);	
				mat[i][j] -= p[i][j] / (1.0 - p[i][j]);
			}
		}
	}
	for (int i = 1; i <= n; ++i) {
		double res = 0.0;
		for (int j = 1; j <= n; ++j) if (i != j) {
			res -= mat[i][j];
		}
		mat[i][i] = res;
	}
//	cout << k << endl;;
	cout << k * gauss (n - 1) << endl;
}