Description
今天DZY 想要玩一个古老的游戏。他在一个有n 个房间并有m 个走廊互相连接的大迷宫里(每个走廊都允许双向通行)。你可以认为所有房间都被走廊直接或间接连接。
DZY 在迷宫里迷路了。现在他在第一房间并且有k 条命。他将会按如下所述行动:
•首先,他会随机抽取一条从他现在所处房间出发的走廊。每个抽取范围内的走廊选中的机率相等。
•然后他会沿着走廊走到走廊的另一端,并且回到第一步重复这个过程。
迷宫中的一些房间里面埋着陷阱。第一房间明显没有陷阱,第n 号房间明确地有一个陷阱。每次DZY 进入这些有陷阱的房间,他都会失去一条命。现在,DZY 知道如果他恰好有两条命时进入了第n 号房间,那么首先他会失去一条命,但是然后他会开启一个福利关卡。他想要知道他开启福利关卡的机率到底为多少。请帮助他。
Solution
只要把从点1到每个点的概率求出,再求出两两黑点之间的概率,最后再用矩阵乘法即可。
于是这题就是高斯消元;
先列出原方程,
设#d_i#表示点i的度数,
bi
为点i被经过的期望次数,因为黑点到了就停止,所以对于黑点来说,期望经过次数和到达概率是一样的;
首先每个点都可以从每个与它相邻的点走来,概率是走过来的点的
1/d
;
如果从点1出发,那么与它相邻的点都要加上(
1/d1
)这个常数,有因为点1不是黑点,所以可以经过多次,那么原来1的方程就表示走出去以后有走回来的期望次数;
如果从黑点出发,那么于上面一样,与它相邻的点都要加上一个常数,当然可能从这个黑点出去走了一堆白点又走回来;
这样复杂度是 O(n4)
观察方程,发现把未知数移到左边后,每次左边都不会变了,变得是右边的常数,
所以我们可以一开始就求出每个b是由哪几个方程右边的常数的几倍加起来,每次变就直接代入即可,
复杂度: O(n3+log(m)∗n3)
Code
#include <iostream>
#include <cstdio>
#include <cstdlib>
#define fo(i,a,b) for(int i=a;i<=b;i++)
#define fod(i,a,b) for(int i=a;i>=b;i--)
#define efo(i,q) for(int i=A[q];i;i=B[i][0])
#define abs(q) ((q)>0?(q):-(q))
using namespace std;
typedef double db;
const int N=505,M=100500;
int read(int &n)
{
char ch=' ';int q=0,w=1;
for(;(ch!='-')&&((ch<'0')||(ch>'9'));ch=getchar());
if(ch=='-')w=-1,ch=getchar();
for(;ch>='0' && ch<='9';ch=getchar())q=q*10+ch-48;n=q*w;return n;
}
int m,n,m1;
int B[M*2][2],A[N],B0,B1[N];
int b1[N],zx[N];
db s[N][N],t[N][N],di[N][N];
void link(int q,int w)
{
B1[q]++,B1[w]++;
B[++B0][0]=A[q];A[q]=B0,B[B0][1]=w;
B[++B0][0]=A[w];A[w]=B0,B[B0][1]=q;
}
void chenans()
{
fo(i,1,m)fo(j,1,m)t[i][j]=0;
fo(k,1,m)fo(i,1,m)fo(j,1,m)t[i][j]+=s[i][k]*di[k][j];
fo(i,1,m)fo(j,1,m)s[i][j]=t[i][j];
}
void chen()
{
fo(i,1,m)fo(j,1,m)t[i][j]=0;
fo(k,1,m)fo(i,1,m)fo(j,1,m)t[i][j]+=di[i][k]*di[k][j];
fo(i,1,m)fo(j,1,m)di[i][j]=t[i][j];
}
db JC()
{
while(m1)
{
if(m1&1)chenans();
chen();m1>>=1;
}
db ans=0;
fo(i,1,m)ans+=s[i][m];
return ans;
}
db b[N][N],c[N][N];
int z[N];
void GS()
{
fo(i,1,n)
{
b[i][i]=c[i][i]=1;
efo(j,i)if(!zx[B[j][1]])b[i][B[j][1]]-=1.0/B1[B[j][1]];
}
fo(i,1,n)
{
int q=i;
fo(j,i+1,n)if(abs(b[j][i])>abs(b[q][i]))q=j;
if(q!=i)
{
fo(j,1,n)swap(b[i][j],b[q][j]),swap(c[i][j],c[q][j]);
}
fo(j,1,n)if(i!=j)
{
db t=b[j][i]/b[i][i];
fo(k,1,n)b[j][k]-=t*b[i][k],c[j][k]-=t*c[i][k];
}
}
fo(i,1,n)fo(j,1,n)c[i][j]/=b[i][i];
efo(i,1)z[B[i][1]]++;
fo(i,1,n)if(zx[i])
{
fo(j,1,n)
s[zx[i]][zx[i]]+=c[i][j]*((db)z[j]/B1[1]);
}
fo(I,1,n)if(zx[I]&&B1[I])
{
fo(i,1,n)z[i]=0;
efo(i,I)z[B[i][1]]++;
fo(i,1,n)if(zx[i])
{
fo(j,1,n)di[zx[I]][zx[i]]+=c[i][j]*((db)z[j]/B1[I]);
}
}
}
int main()
{
freopen("games.in","r",stdin);
freopen("games.out","w",stdout);
int q,w;
read(n),read(m),read(m1);
m1-=2;if(m1<0){printf("0\n");return 0;}
fo(i,1,n)if(read(q))b1[++b1[0]]=i,zx[i]=b1[0];
fo(i,1,m)read(q),read(w),link(q,w);
GS();
m=b1[0];
printf("%.12lf\n",JC());
return 0;
}