bzoj3270-博物馆

题目

给出一个无向图，两个人初始在两个点上。当一个人在一个点\(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("");
}