题目
分析
问题要求
由于矩阵树只能处理
,所以把式子变成
于是用矩阵树定理推广计算
,最后乘上
即为答案。
注意
会除零而出问题,给
减个
即可。
矩阵树定理推广:计算所有生成树边权乘积的和,即 其中 为一个有理数。
令 ,其中 ,那么 当边权 为整数的时候很好计算:把 之间连上 条边,然后计算生成树的个数就是原图所有生成树边权乘积和。
但是 很显然精度不够,所以我们用 行列式的可乘性 把 分配给行列式中每个元素,于是我们有了一个新的基尔霍夫矩阵:度为连出去所有边的边权和,边权为原边权的相反数。然后正常矩阵树定理得到的答案就是所求的了。
代码
#include <bits/stdc++.h>
const int MAXN = 50;
const double eps = 1e-5;
int N;
double Mat[MAXN + 5][MAXN + 5];
void AddEdge(int u, int v, double w) {
Mat[u][u] += w, Mat[v][v] += w;
Mat[u][v] -= w, Mat[v][u] -= w;
}
inline double Abs(const double &x) {
return x < 0 ? -x : x;
}
inline int Comp(const double &x, const double &y) {
if (Abs(x - y) < eps)
return 0;
return x > y ? 1 : -1;
}
double Det(int n) {
n--;
double ret = 1;
for (int i = 1; i <= n; i++) {
int p = i;
for (int j = i; j <= n; j++)
if (Comp(Mat[p][i], Mat[j][i]) < 0)
p = j;
if (p != i)
std::swap(Mat[i], Mat[p]);
if (!Comp(Mat[i][i], 0))
return 0;
ret *= Mat[i][i];
for (int j = i + 1; j <= n; j++) {
double tmp = Mat[j][i] / Mat[i][i];
for (int k = i; k <= n; k++)
Mat[j][k] -= Mat[i][k] * tmp;
}
}
return Abs(ret);
}
int main() {
scanf("%d", &N);
double All = 1.0;
for (int i = 1; i <= N; i++)
for (int j = 1; j <= N; j++) {
double p; scanf("%lf", &p);
if (i < j) {
if (Comp(p, 1) == 0)
p -= eps;
AddEdge(i, j, p / (1 - p));
All *= (1 - p);
}
}
printf("%.5f", Det(N) * All);
}