bzoj1778 驱逐猪猡 [高斯消元+概率DP]

in 题目 with 0 comment

Description

奶牛们建立了一个随机化的臭气炸弹来驱逐猪猡。猪猡的文明包含1到N一共N个猪城。这些城市由M条由两个不同端点$A_j$和$B_j$ $(1 \le A_j \le N; 1 \le B_j \le N)$ 表示的双向道路连接。保证城市1至少连接一个其它的城市。一开始臭气弹会被放在城市1。每个小时(包括第一个小时),它有P/Q $(1 \le P \le 1,000,000 , 1 \le Q \le 1,000,000)$ 的概率污染它所在的城市。如果这个小时内它没有污染它所在的城市,那麽它随机地选择一条道路,在这个小时内沿着这条道路走到一个新的城市。可以离开这个城市的所有道路被选择的概率均等。因为这个臭气弹的随机的性质,奶牛们很困惑哪个城市最有可能被污染。给定一个猪猡文明的地图和臭气弹在每个小时内爆炸的概率。计算每个城市最终被污染的概率。如下例,假设这个猪猡文明有两个连接在一起的城市。臭气炸弹从城市1出发,每到一个城市,它都有1/2的概率爆炸。 1--2 可知下面这些路径是炸弹可能经过的路径(最后一个城市是臭气弹爆炸的城市):
1: 1
2: 1-2
3: 1-2-1
4: 1-2-1-2
5: 1-2-1-2-1 ...
要得到炸弹在城市1终止的概率,我们可以把上面的第1,第3,第5……条路径的概率加起来,(也就是上表奇数编号的路径)。上表中第k条路径的概率正好是$(1/2)^k$,也就是必须在前k-1个回合离开所在城市(每次的概率为1 - 1/2 = 1/2)并且留在最后一个城市(概率为1/2)。所以在城市1结束的概率可以表示为$1/2 + (1/2)^3 + (1/2)^5 + \dots$ 当我们无限地计算把这些项一个个加起来,我们最后会恰好得到2/3,也就是我们要求的概率,大约是0.666666667。这意味着最终停留在城市2的概率为1/3,大约为0.333333333。


Input


Output


Sample Input

2 1 1 2
1 2


Sample Output

0.666666667
0.333333333


Hint

$ 2 \le N \le 300 $
$ 1 \le M \le 44,850 $


Solution

我们构造一个nn的转移矩阵E[ ][ ],其中E[i][j]表示$i\rightarrow j$转移的概率,易得$E[i][j]=\frac{(1-P/Q)}{deg[i]}$
然后我们构造一个1
n的单位向量T[],表示当前走到每个点的概率,初始时$T[


Code

#include<bits/stdc++.h>
#define maxn 501
#define maxm 50010
#define eps 1e-7
// #define DEBUG
using namespace std;
int n,m,p,q;
double rate;
int edge[maxn][maxn];
double a[maxn][maxn];

inline int dcmp(double x)
{
  if(abs(x)<=eps) return 0;
  else return x>0?1:-1;
}

inline int read()
{
  char ch;
  int read=0,sign=1;
  do
    ch=getchar();
  while((ch<'0'||ch>'9')&&ch!='-');
  if(ch=='-') ch=getchar(),sign=-1;
  while(ch>='0'&&ch<='9')
  {
    read=read*10+ch-'0';
    ch=getchar();
  }
  return sign*read;
}

void change()
{
  for(int i=1;i<=n;++i)
    for(int j=1;j<=n;++j)
    {
      if(edge[i][j]) a[j][i]=(rate-1)/edge[i][0];
      if(i==j) a[i][j]+=1;
    }
  a[1][n+1]=rate;
}

void prework()
{
  n=read(),m=read(),p=read(),q=read();
  rate=double(p)/double(q);
  for(int i=1,u,v;i<=m;++i)
  {
    int u=read();
    int v=read();
    edge[u][v]=1;
    edge[v][u]=1;
    edge[u][0]++;
    edge[v][0]++;
  }
  change();
#ifdef DEBUG
  for(int i=1;i<=n;++i){
    for(int j=1;j<=n+1;++j)
      printf("%.3lf ",a[i][j]);
    printf("\n");}
#endif
}

void gauss()
{
  for(int i=1;i<=n;++i)
  {
    if(!dcmp(a[i][i]))
    {
      int k=i;
      while(!dcmp(a[k][k])) k++;
      swap(a[k],a[i]);
    }
    for(int j=i+1;j<=n;++j)
    {
      for(int k=i+1;k<=n+1;++k)
        a[j][k]-=a[i][k]*(a[j][i]/a[i][i]);
      a[j][i]=0;
    }
  }
  for(int i=n;i>1;--i)
  {
    a[i][n+1]/=a[i][i];
    for(int j=1;j<i;++j)
    {
      a[j][n+1]-=a[i][n+1]*a[j][i];
      a[j][i]=0;
    }
  }
  for(int i=1;i<=n;++i) printf("%.9lf\n",a[i][n+1]);
}

int main()
{
  prework();
  gauss();
  return 0;
}

扫描二维码,在手机上阅读!
Responses