只有我是暴力缩点搞的吗。。。
令$E$为生成树边集,显然我们要求的就是$\sum\limits_E\left(\prod\limits_{e\in E}p_e\prod\limits_{e\notin E}(1-p_e)\right)$
令$P=\prod(1-p_e)$,即原图所有边的$1-p_e$之积。上面的式子就成了$P\sum\limits_E\left(\prod\limits_{e\in E}\dfrac{p_e}{1-p_e}\right)$
直接变元矩阵树定理搞它就完事了。。。吗?
若$p_e=1$,会出现分母为$0$的情况。其他题解里都是改为$1-eps$,精度误差可以接受。
没想到于是暴力缩点。$p_e=1$说明生成树中一定有这条边,用并查集缩成一个点。
对于两端点在同一并查集的$e$,若$p_e\neq 1$,当成自环不用管它;若$p_e=1$,说明一定会出来一个环,输出0
。然而我一开始没判过了,大概是数据水了。
代码:
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <cstring>
#define maxn 55
#define inf 0x3f3f3f3f
using namespace std;
inline int read(){
int x=0,y=0;
char ch=getchar();
while(ch<'0'||ch>'9'){if(ch=='-')y=1;ch=getchar();}
while(ch>='0'&&ch<='9')x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
return y?-x:x;
}
double a[maxn][maxn],e[maxn][maxn];
int cnt,fa[maxn],id[maxn];
inline void add(int x,int y,double l){a[x][x]+=l,a[y][y]+=l,a[x][y]-=l,a[y][x]-=l;}
int find(int x){return fa[x]==x?x:fa[x]=find(fa[x]);}
double Gauss(){
double ans=1.0;
for(register int i=2;i<=cnt;++i){
int t=i;
for(register int j=i+1;j<=cnt;++j)if(fabs(a[j][i])>fabs(a[t][i]))t=j;
if(t!=i)swap(a[t],a[i]),ans=-ans;
for(register int j=i+1;j<=cnt;++j)
for(register int k=cnt;k>=i;--k)
a[j][k]-=a[i][k]*a[j][i]/a[i][i];
ans*=a[i][i];
}
return ans;
}
int main(){
int n=read();
double P=1.0;
for(register int i=1;i<=n;++i)fa[i]=i;
for(register int i=1;i<=n;++i)
for(register int j=1;j<=n;++j){
scanf("%lf",&e[i][j]);
if(i<j&&e[i][j]==1.0){
int u=find(i),v=find(j);
if(u==v){puts("0.000000");return 0;}
fa[u]=v;
}
}
for(register int i=1;i<=n;++i){
if(!id[find(i)])id[find(i)]=++cnt;
id[i]=id[find(i)];
}
for(register int i=1;i<=n;++i)
for(register int j=1;j<=n;++j)
if(i<j&&e[i][j]!=1.0)add(id[i],id[j],e[i][j]/(1.0-e[i][j])),P*=1.0-e[i][j];
printf("%.8lf\n",Gauss()*P);
}