题目链接:戳我
变元矩阵树定理,具体可以去看蒟蒻学习笔记中[树上问题]这一篇的整理QAQ
原式子为\(\sum_T\prod_{(u,v)\in T}p_{u,v}\prod_{(u,v)\notin T}(1-p_{u,v})\)
\(=\prod_{(u,v)\in G}(1-p_{u,v})\sum_T\prod_{(u,v)\in T}\frac{p_{u,v}}{1-p_{u,v}}\)
所以直接上变元矩阵树定理就可以了.
但是要注意我们中间涉及到除法操作,所以记得把概率为1的设置成为1-eps,概率为0的设置成eps
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define MAXN 110
#define eps 1e-8
using namespace std;
int n,p;
double all=1.0;
double a[MAXN][MAXN];
inline double solve()
{
double cur_ans=1.0;
for(int i=2;i<=n;i++)
{
int pos=i;
for(int j=i+1;j<=n;j++)
if(fabs(a[j][i])>fabs(a[pos][i]))
pos=j;
if(pos!=i) swap(a[i],a[pos]),cur_ans*=-1.0;
for(int j=i+1;j<=n;j++)
{
double t=a[j][i]/a[i][i];
for(int k=i;k<=n;k++) a[j][k]-=t*a[i][k];
}
}
for(int i=2;i<=n;i++)
cur_ans*=a[i][i];
return cur_ans;
}
int main()
{
#ifndef ONLINE_JUDGE
freopen("ce.in","r",stdin);
#endif
scanf("%d",&n);
for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
{
scanf("%lf",&a[i][j]);
if(a[i][j]<eps) a[i][j]=eps;
if(a[i][j]>1.0-eps) a[i][j]=1.0-eps;
if(i<j) all*=(1.0-a[i][j]);
a[i][j]=a[i][j]/(1.0-a[i][j]);
}
}
/*
for(int i=1;i<=n;i++)
{
for(int j=1;j<=n;j++)
printf("%.3lf ",a[i][j]);
printf("\n");
}
*/
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
if(i!=j)
a[i][i]-=a[i][j];
for(int i=1;i<=n;i++)
for(int j=1;j<=n;j++)
a[i][j]=-a[i][j];
//rintf("all=%.3lf\n",all);
printf("%.6lf\n", all*solve( ));
return 0;
}