P3317 [SDOI2014]重建（Matrix-tree+期望）

P3317 [SDOI2014]重建

详情看这位神犇的blog

剩下的注释在code里吧.......

#include<iostream>
#include<cstdio>
#include<cstring>
#include<cmath>
using namespace std;
typedef double db;
const db eps=1e-7;
db d[55][55],ans=1.0;
int n;
db det(){//矩阵树定理板子，计算行列式
    db re=1.0,div;
    for(int i=1;i<n;++i){
        int x=i;
        for(int j=i+1;j<n;++j)
            if(d[j][i]>d[x][i]) x=j;
        if(x!=i) re=-re,swap(d[x],d[i]);
        for(int j=i+1;j<n;++j){
            div=d[j][i]/d[i][i];
            for(int u=i;u<n;++u)
                d[j][u]-=d[i][u]*div;
        }re*=d[i][i];
        if(fabs(re)<eps) break;
    }return re;
}
int main(){
    scanf("%d",&n);
    for(int i=1;i<=n;++i)
        for(int j=1;j<=n;++j){
            scanf("%lf",&d[i][j]);
            if(i==j) continue;
            if(fabs(d[i][j])<eps) d[i][j]=eps;
            if(fabs(1.0-d[i][j])<eps) d[i][j]=1.0-eps;//防止分母出现inf的情况
            if(i<j) ans*=(1.0-d[i][j]);//不要重复算
            d[i][j]=d[i][j]/(1.0-d[i][j]);
        }
    for(int i=1;i<=n;++i){//构造基尔霍夫矩阵
        d[i][i]=0.0;
        for(int j=1;j<=n;++j)
            if(i!=j) d[i][i]+=d[i][j],d[i][j]=-d[i][j];//度数矩阵+=边权，邻接矩阵为负边权
    }ans*=det();
    printf("%.10lf",ans);
    return 0;
}