bzoj3270 博物馆(期望+高斯消元)
有一天Petya和他的朋友Vasya在进行他们众多旅行中的一次旅行,他们决定去参观一座城堡博物馆。这座博物馆有着特别的样式。它包含由m条走廊连接的n间房间,并且满足可以从任何一间房间到任何一间别的房间。
两个人在博物馆里逛了一会儿后两人决定分头行动,去看各自感兴趣的艺术品。他们约定在下午六点到一间房间会合。然而他们忘记了一件重要的事:他们并没有选好在哪儿碰面。等时间到六点,他们开始在博物馆里到处乱跑来找到对方(他们没法给对方打电话因为电话漫游费是很贵的)
不过,尽管他们到处乱跑,但他们还没有看完足够的艺术品,因此他们每个人采取如下的行动方法:每一分钟做决定往哪里走,有Pi 的概率在这分钟内不去其他地方(即呆在房间不动),有1-Pi 的概率他会在相邻的房间中等可能的选择一间并沿着走廊过去。这里的i指的是当期所在房间的序号。在古代建造是一件花费非常大的事,因此每条走廊会连接两个不同的房间,并且任意两个房间至多被一条走廊连接。
两个男孩同时行动。由于走廊很暗,两人不可能在走廊碰面,不过他们可以从走廊的两个方向通行。(此外,两个男孩可以同时地穿过同一条走廊却不会相遇)两个男孩按照上述方法行动直到他们碰面为止。更进一步地说,当两个人在某个时刻选择前往同一间房间,那么他们就会在那个房间相遇。
两个男孩现在分别处在a,b两个房间,求两人在每间房间相遇的概率。
第一行包含四个整数,n表示房间的个数;m表示走廊的数目;a,b (1 ≤ a, b ≤ n),表示两个男孩的初始位置。
之后m行每行包含两个整数,表示走廊所连接的两个房间。
之后n行每行一个至多精确到小数点后四位的实数 表示待在每间房间的概率。
题目保证每个房间都可以由其他任何房间通过走廊走到。
输出一行包含n个由空格分隔的数字,注意最后一个数字后也有空格,第i个数字代表两个人在第i间房间碰面的概率(输出保留6位小数)
注意最后一个数字后面也有一个空格
Sample Input
2 1 1 2
1 2
0.5
0.5
1 2
0.5
0.5
Sample Output
0.500000 0.500000
对于100%的数据有 n <= 20,n-1 <= m <= n(n-1)/2
(可以列个概率方程并暴力矩乘,然鹅并不知道什么时候收敛,精度会出现误差)
设两个人在$(x,y)$的期望经过次数为$f(x,y)$,$(x,y)$转移到$(r,w)$的概率为$a[(r,w)][(x,y)]$,
点$i$连$in[i]$条边,在点$i$停留的概率为$p[i]$
$a[(x,y)][(x,y)]=p[x]*p[y]$
$a[(x,y+1)][(x,y)]=p[x]*(1-p[y])/in[y]$
$a[(x+1,y)][(x,y)]=p[y]*(1-p[x])/in[x]$
$a[(x+1,y+1)][(x,y)]=(1-p[y])/in[y]*(1-p[x])/in[x]$
那么对于每个$f(x,y)$,我们都能得出一个方程:
$f(x,y)= \sum f(r,w)*a[(x,y)][(r,w)]$($r=x$或存在边$(x,r)$,$w$同理)
酱紫我们就得出了$n^2$个方程。
发现方程组存在环的关系,而对于$n$元1次方程组,常规做法是套个高斯消元上去
把上面的方程移项一下
$ \sum f(r,w)*a[(x,y)][(r,w)]-f(x,y)=0$
蓝后跑高斯消元就好辣
最后对于$f(a,b)$,因为是初始点,所以要期望次数+1
终点只走一次,所以期望经过次数就等于概率
#include<iostream> #include<cstdio> #include<cstring> #include<vector> using namespace std; typedef double db; const db eps=1e-9; inline db Fabs(db x){return x<0?-x:x;} #define N 420 int n,m,w,A,B,id[N][N]; db p[N],a[N][N],sol[N]; vector <int> g[N]; void draw(int x,int y){ int e=id[x][y],lx=g[x].size(),ly=g[y].size();//lx,ly:x,y的连边数 for(int i=0;i<lx;++i){ int r=g[x][i]; a[id[r][y]][e]+=(1-p[x])/(1.0*lx)*p[y]; } for(int i=0;i<ly;++i){ int r=g[y][i]; a[id[x][r]][e]+=(1-p[y])/(1.0*ly)*p[x]; } for(int i=0;i<lx;++i) for(int j=0;j<ly;++j){ int r1=g[x][i],r2=g[y][j]; a[id[r1][r2]][e]+=(1-p[x])/(1.0*lx)*(1-p[y])/(1.0*ly); } } void gauss(){ for(int i=1,x=1;i<=w;x=++i){ for(int j=i+1;j<=w;++j) if(Fabs(a[j][i])>Fabs(a[x][i])) x=j; if(x!=i) swap(a[x],a[i]); if(Fabs(a[i][i])<eps) continue; for(int j=1;j<=w;++j){ if(i==j) continue; db div=a[j][i]/a[i][i]; for(int u=i;u<=w+1;++u) a[j][u]-=a[i][u]*div; } } for(int i=w;i;--i){ sol[i]=a[i][w+1]; for(int j=w;j>i;--j) sol[i]-=sol[j]*a[i][j]; sol[i]/=a[i][i]; } for(int i=1;i<=n;++i) printf("%.6f ",sol[id[i][i]]+eps);//防止输出-0.00 } int main(){ scanf("%d%d%d%d",&n,&m,&A,&B); for(int i=1,u,v;i<=m;++i){ scanf("%d%d",&u,&v); g[u].push_back(v); g[v].push_back(u); } for(int i=1;i<=n;++i) scanf("%lf",&p[i]); for(int i=1;i<=n;++i) for(int j=1;j<=n;++j){ id[i][j]=++w; if(i!=j) a[w][w]=p[i]*p[j]; } for(int i=1;i<=n;++i) for(int j=1;j<=n;++j) if(i!=j) draw(i,j); for(int i=1;i<=w;++i) a[i][i]-=1.0;//移项系数为-1 a[id[A][B]][w+1]=-1.0;//起点+1 gauss(); return 0; }