传送门
首先,大家应该都能看出来这是矩阵树定理,然后大部分人应该就会把概率直接带进去算,然后就愉快地WA掉了(我当时就是这么想的,幸亏没交)
然后就来讲这个题的正解思路。
首先我们来看答案应该是怎样的:
ans=∑Tree∏(u,v)∈EP(u,v)∏(u,v)∉E(1−P(u,v))
然后我们来想一下怎么来构造这个答案:首先,我们直接矩阵树用高斯算出来的结果应该是这个:
now=∑Tree∏(u,v)∈EP(u,v)
那我们怎么让它变成答案那个样子呢?直观地,我们可以这样:
now=now∗∏(u,v)(1−P(u,v))
然后就变成了这样子:
now=∑Tree∏(u,v)∈EP(u,v)∗(1−P(u,v))∏(u,v)∉E(1−P(u,v))
现在我们发现now与答案有着中间那一部分的差距,那么如何消除那个差距呢?
我们可以注意到式子中的第一个 P(u,v) 其实就是矩阵的初值,那么也就是这个:
now=∑Tree∏(u,v)∈EAu,v∗(1−P(u,v))∏(u,v)∉E(1−P(u,v))
那么我们现在就是要让:
Au,v∗(1−P(u,v))→P(u,v)
想必我说到现在大家应该就懂了吧,我们只要这样:
Au,v=P(u,v)(1−P(u,v))
这样我们现在的now值就是答案了!
tmp=∏(u,v)(1−P(u,v))
now=tmp∗∑Tree∏(u,v)∈EP(u,v)(1−P(u,v))
还有就是几个技巧:
当矩阵中出现 |a|<eps 时就 a=eps
当矩阵中出现 |1−a|<eps 时就 a=1−eps
自己想一想为什么。
代码如下:
/**************************************************************
Problem: 3534
User: stone41123
Language: C++
Result: Accepted
Time:16 ms
Memory:1308 kb
****************************************************************/
#include<bits/stdc++.h>
#define ll long long
using namespace std;
int n;
double a[51][51];
double ans;
double eps=1e-8;
void gauss(){
for(int i=1;i<n;i++){
int mx=i;
for(int j=i+1;j<n;j++){
if(fabs(a[j][i])>fabs(a[mx][i]))mx=j;
}
if(mx!=i)for(int j=1;j<n;j++)swap(a[i][j],a[mx][j]);
for(int k=i+1;k<n;k++){
double mul=a[k][i]/a[i][i];
for(int j=i;j<n;j++){
a[k][j]-=a[i][j]*mul;
}
}
if(fabs(a[i][i])<eps){
ans=0;
return;
}
}
for(int i=1;i<n;i++){
ans*=a[i][i];
}
ans=fabs(ans);
}
int main(){
scanf("%d",&n);
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
scanf("%lf",&a[i][j]);
}
}
ans=1;
double tmp=1;
for(int i=1;i<=n;i++){
for(int j=1;j<=n;j++){
if(fabs(a[i][j])<eps)a[i][j]=eps;
if(fabs(1.0-a[i][j])<eps)a[i][j]=1-eps;
if(i<j)tmp*=1.0-a[i][j];
a[i][j]=a[i][j]/(1.0-a[i][j]);
}
}
for(int i=1;i<=n;i++){
a[i][i]=0;
for(int j=1;j<=n;j++){
if(i!=j)
a[i][i]-=a[i][j];
}
}
gauss();
ans*=tmp;
printf("%.10lf",ans);
return 0;
}