Codeforces 113D Museum (期望DP + 高斯消元)

题目链接

题目大意：
有$N$个房间，其中两个房间 $A$,$B$ 里各有一个人。有$M$对房间通过双向边连通。每一个时刻，每个人有一定概率$p_i$留在 $i$ 号房间，或者有$1-p_i$的概率离开，等概率前往相邻房间。两人一旦处于同一房间则立即停下来。求对于所有的房间 $i$，两人在 $i$ 号房间相遇的概率是多少
PS:
1.在过道内相遇不算
2.所有房间保证连通

解法：
考虑DP，设$f[i][j]$ 为一个人在房间 $i$ ，同时另一个人在 房间 $j$ 的状态出现期望，转移情况显然为：

$$f\left[i\right]\left[j\right]=P_i\cdot P_j\cdot f\left[i\right]\left[j\right]+P_i\cdot \sum \frac{1-P_k}{\mathrm{deg}\left[k\right]}\cdot f\left[i\right]\left[k\right]+P_j\cdot \sum \frac{1-P_k}{\mathrm{deg}\left[k\right]}\cdot f\left[k\right]\left[j\right]+\sum \sum \frac{1-P_k}{\mathrm{deg}\left[k\right]}\cdot \frac{1-P_h}{\mathrm{deg}\left[h\right]}\cdot f\left[k\right]\left[h\right]$$

一共有$n^2$ 个变量，所以复杂度为 $O(n^6)$

又因为终止状态$(i,i)$的出现次数只能为$1/0$ 所以终止状态的概率即为期望

注意：
1.初始化$f[A][B]=1$
2. 终止状态$(i,i)$不能转移到其他状态，但可以从其他状态转移过来

Code：

#include <bits/stdc++.h>
using namespace std;
const long double eps = 1e-8;
const int MX = 25 + 7;
int n,m,A,B,cnt;
double p[MX];

inline int F(int x,int y){return (x - 1) * n + y;}
struct Edge{
	int v,next;
}e[2 * MX * MX];
int head[MX],ecnt,deg[MX];
void add(int u,int v){
	deg[u]++;
	e[++ecnt].v = v;
	e[ecnt].next = head[u];
	head[u] = ecnt;
}

double f[MX*MX][MX*MX];

void build(){
	cnt = n * n;
	f[F(A,B)][cnt+1] = 1;
	for(int i = 1;i <= n;++i){
		for(int j = 1;j <= n;++j){
			int u = F(i,j);
			f[u][u] = 1;
			if(i != j) f[u][u] -= (p[i] * p[j]);
			for(int k = head[j];k;k = e[k].next){
				int v = e[k].v;
				if(v != i) f[u][F(i,v)] -= p[i] * (1 - p[v]) / deg[v];
			}
			for(int k = head[i];k;k = e[k].next){
				int v = e[k].v;
				if(v != j) f[u][F(v,j)] -= p[j] * (1 - p[v]) / deg[v];
			}
			for(int k = head[i];k;k = e[k].next){
				for(int h = head[j];h;h = e[h].next){
					int v1 = e[k].v, v2 = e[h].v;
					if(v1 != v2) f[u][F(v1,v2)] -= (1 - p[v1]) * (1 - p[v2]) / (1.0 * deg[v1] * deg[v2]);
				}
			}
		}
	}
}
double ans[MX*MX];
void Gauss(){
	for(int i = 1;i <= cnt;++i){
		int r = i;
		for(int j = i + 1;j <= cnt;++j) if(fabs(f[r][i]) < fabs(f[j][i])) r = j;
		if(fabs(f[r][i]) < eps) continue;
		if(i != r) swap(f[i],f[r]);
		double div = f[i][i];
		for(int j = i;j <= cnt + 1;++j) f[i][j] /= div;
		for(int j = i + 1;j <= cnt;++j){
			div = f[j][i];
			for(int k = i;k <= cnt + 1;++k)
				f[j][k] -= f[i][k] * div;
		}
	}
	ans[cnt] = f[cnt][cnt+1];
	for(int i = cnt - 1;i;--i){
		ans[i] = f[i][cnt+1];
		for(int j = i + 1;j <= cnt;++j)
			ans[i] -= (f[i][j] * ans[j]);
	}
}
int main(){
	scanf("%d %d %d %d",&n,&m,&A,&B);
	for(int i = 1;i <= m;++i){
		int u,v;scanf("%d %d",&u,&v);
		add(u,v);add(v,u);
	}
	for(int i = 1;i <= n;++i) scanf("%lf",&p[i]);
	build();Gauss();
	for(int i = 1;i <= n;++i) printf("%.6f%c", ans[F(i,i)],(i == n) ? '\n' : ' ');
	return 0;
}