【BZOJ3534】【SDOI2014】重建

【题目链接】

• 点击打开链接
  

【思路要点】

• 无向图的生成树，考虑矩阵树定理。
• 若一条边出现的概率为\(P_i\)，则令其权值为\(\frac{P_i}{1-P_i}\)。
• 构造基尔霍夫矩阵，令其去掉一行一列的行列式为\(Det\)。
• 则答案为\(Det*\prod_{i\in E}(1-P_i)\)。
• 注意到可能有\(P_i=1\)，此时令\(P_i=1-10^{-8}\)即可。
• 时间复杂度\(O(N^3)\)。

【代码】

#include<bits/stdc++.h>
using namespace std;
const int MAXN = 105;
const long double eps = 1e-15;
template <typename T> void chkmax(T &x, T y) {x = max(x, y); }
template <typename T> void chkmin(T &x, T y) {x = min(x, y); } 
template <typename T> void read(T &x) {
	x = 0; int f = 1;
	char c = getchar();
	for (; !isdigit(c); c = getchar()) if (c == '-') f = -f;
	for (; isdigit(c); c = getchar()) x = x * 10 + c - '0';
	x *= f;
}
template <typename T> void write(T x) {
	if (x < 0) x = -x, putchar('-');
	if (x > 9) write(x / 10);
	putchar(x % 10 + '0');
}
template <typename T> void writeln(T x) {
	write(x);
	puts("");
}
int n; long double a[MAXN][MAXN];
long double solve(int n) {
	long double f = 1;
	for (int i = 1; i <= n; i++) {
		for (int j = i; j <= n; j++)
			if (fabs(a[j][i]) > eps) {
				if (i != j) f = -f;
				swap(a[j], a[i]);
				break;
			}
		for (int j = i + 1; j <= n; j++) {
			long double tmp = a[j][i] / a[i][i];
			for (int k = 1; k <= n; k++)
				a[j][k] -= a[i][k] * tmp;
		}
	}
	for (int i = 1; i <= n; i++)
		f *= a[i][i];
	return f;
}
int main() {
	read(n);
	for (int i = 1; i <= n; i++)
	for (int j = 1; j <= n; j++) {
		scanf("%Lf", &a[i][j]);
		if (a[i][j] == 1) a[i][j] -= 1e-8;
	}
	long double sum = 1;
	for (int i = 2; i <= n; i++)
	for (int j = 1; j <= i - 1; j++) {
		sum *= (1 - a[i][j]);
		a[i][j] = a[j][i] = a[i][j] / (1 - a[i][j]);
	}
	for (int i = 1; i <= n; i++) {
		long double ans = 0;
		for (int j = 1; j <= n; j++)
			ans += a[i][j], a[i][j] = -a[i][j];
		a[i][i] = ans;
	}
	printf("%.10Lf\n", solve(n - 1) * sum);
	return 0;
}