题目
给出一个无向图,两个人初始在两个点上。当一个人在一个点\(i\)上的时候,每一次,他有\(p[i]\)的概率留在原位,有\(1-p[i]\)的概率等概率地选择直接连边的一个点走出去。当两个人在同一时刻走到同一个点,那么他们相遇,过程结束。现在求他们在每一个点相遇的概率。
分析
看似这是一个简单的概率dp,其实不然。发现图中如果出现环,那么dp的顺序无法保证。hjh说这和一个复杂电路中的电阻差不多。这题中概率无法直接dp,因为概率的性质需要容斥什么的,很麻烦,所以我们把它转化为简单的期望求解。
令\(f[x][y]\)表示两人一个在\(x\)一个在\(y\)这种状态的期望发生次数,那么就有可以通过四种情况的讨论列出方程:
- 从i走到x,y不动:\(f[i][y]*(1-p[i])/z[i]*p[y]\)
- 从j走到y,x不动:\(f[x][j]*(1-p[j])/z[j]*p[x]\)
- 从i走到x,从j走到y:\(f[i][j]*(1-p[i])*(1-p[j])/z[i]/z[j]\)
- 留在原地:\(f[x][y]*p[x]*p[y]\)
这里一定要注意,\(i\ne j\),因为相遇后不会继续走。
然后就开始列方程,最后一个问题是初始状态怎么处理。
我们计算的是期望经过的次数,而起点最开始就经过了,所以其实起点的\(f[x][y]\)等于上述四种情况再加上原有的1,所以移项到右边变成-1.
由于是矩阵,所以我们要把\(f[x][y]\)编个号,比如\((n-1)*x+y\)之类的,能用即可。
消元之后,我们就可以得到每个\(f[x][y]\)的值。由于到达终点的总期望是1次,所以这里的概率就等于期望,直接输出答案即可。
代码
#include<cstdio>
#include<algorithm>
#include<vector>
#include<cmath>
using namespace std;
const int maxn=21;
const int maxm=maxn*maxn;
typedef vector<int> vec;
typedef vec::iterator itt;
vec a[maxn];
int z[maxn],n,all;
double p[maxn];
double g[maxm][maxm];
void sw(double a[],double b[]) {
for (int i=1;i<=all+1;++i) swap(a[i],b[i]);
}
void elm() {
for (int i=1;i<all;++i) {
if (g[i][i]==0) {
for (int j=i+1;j<=all;++j) if (g[j][i]) {
sw(g[i],g[j]);
break;
}
}
for (int j=i+1;j<=all;++j) {
double tmp=g[j][i]/g[i][i];
for (int k=1;k<=all+1;++k) {
g[j][k]-=g[i][k]*tmp;
}
}
}
for (int i=all;i>1;--i) {
if (g[i][i]==0) {
for (int j=i-1;j>0;--j) if (g[j][i]) {
sw(g[i],g[j]);
break;
}
}
for (int j=i-1;j>0;--j) {
double tmp=g[j][i]/g[i][i];
for (int k=1;k<=all+1;++k) {
g[j][k]-=g[i][k]*tmp;
}
}
}
}
int get(int x,int y) {
return (x-1)*n+y;
}
int main() {
#ifndef ONLINE_JUDGE
freopen("test.in","r",stdin);
#endif
int m,s,t;
scanf("%d%d%d%d",&n,&m,&s,&t),all=n*n;
for (int i=1;i<=n;++i) a[i].push_back(i);
for (int i=1;i<=m;++i) {
int u,v;
scanf("%d%d",&u,&v);
a[u].push_back(v),a[v].push_back(u);
z[u]++,z[v]++;
}
for (int i=1;i<=n;++i) scanf("%lf",p+i);
for (int x=1;x<=n;++x) for (int y=1;y<=n;++y) {
int idn=get(x,y);
g[idn][idn]=-1;
for (itt it=a[x].begin();it!=a[x].end();++it) {
for (itt jt=a[y].begin();jt!=a[y].end();++jt) {
int i=*it,j=*jt,id=get(i,j);
if (i==j) continue;
if (i==x && j==y) g[idn][id]+=p[x]*p[y]; else
if (i==x) g[idn][id]+=p[x]*(1-p[j])/z[j]; else
if (j==y) g[idn][id]+=p[y]*(1-p[i])/z[i]; else
g[idn][id]+=(1-p[i])*(1-p[j])/z[i]/z[j];
}
}
}
g[(s-1)*n+t][all+1]=-1;
elm();
for (int i=1;i<=n;++i) printf("%.6lf ",g[(i-1)*n+i][all+1]/g[(i-1)*n+i][(i-1)*n+i]);
puts("");
}