bzoj 1778: [Usaco2010 Hol]Dotp 驱逐猪猡

Description

奶牛们建立了一个随机化的臭气炸弹来驱逐猪猡。猪猡的文明包含1到N (2 <= N <= 300)一共N个猪城。这些城市由M (1 <= M <= 44,850)条由两个不同端点A_j和B_j (1 <= A_j<= N; 1 <= B_j <= N)表示的双向道路连接。保证城市1至少连接一个其它的城市。一开始臭气弹会被放在城市1。每个小时(包括第一个小时),它有P/Q (1 <= P <=1,000,000; 1 <= Q <= 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 + ...。当我们无限地计算把这些项一个个加起来,我们最后会恰好得到2/3,也就是我们要求的概率,大约是0.666666667。这意味着最终停留在城市2的概率为1/3,大约为0.333333333。

Input

* 第1行: 四个由空格隔开的整数: N, M, P, 和 Q * 第2到第M+1行: 第i+1行用两个由空格隔开的整数A_j和B_j表示一条道路。

Output

* 第1到第N行: 在第i行,用一个浮点数输出城市i被摧毁的概率。误差不超过10^-6的答桉会 被接受(注意这就是说你需要至少输出6位有效数字使得答桉有效)。
solution
再不写题解可能又忘了....
设f[i]为在i爆炸的概率   与i相连的点为j
f[i]=   f[j]/(p/q)      *     (1-p/q)          *   1/outdegree[j]        *        p/q
         在j点的概率      j点不爆炸的概率         不爆炸走到j点                   在i点爆炸
    =  f[j] * (1-p/q) * (1/outdegree[j])   (化简之后的式子)
高斯解方程即可
  1 #include<cstdio>
  2 #include<cstring>
  3 #include<iostream>
  4 #define mem(a,b) memset(a,b,sizeof(a))
  5 using namespace std;
  6 const double qqq=1e-8;
  7  
  8 double abss(double x)
  9 {
 10     if(x<0)
 11       return -x;
 12     return x;
 13 }
 14 
 15 void swap(double &x,double &y)
 16 {
 17     double temp;
 18     temp=x;
 19     x=y;
 20     y=temp;
 21 }
 22  
 23 struct son
 24 {
 25     int v,next;
 26 };
 27 son a1[100001];
 28 int first[100001];
 29 int e;
 30  
 31 void addbian(int u,int v)
 32 {
 33     a1[e].v=v;
 34     a1[e].next=first[u];
 35     first[u]=e++;
 36 }
 37  
 38 int n,m,p,q;
 39 int u,o;
 40 double degree[301];
 41 double a[301][301];
 42 double ans[301];
 43 double b[301];
 44  
 45 void out11()
 46 {
 47     for(int i=1;i<=n;i++)
 48       printf("%lf ",degree[i]);
 49     printf("\n");
 50     for(int i=1;i<=n;i++)
 51     {
 52         for(int j=1;j<=n;j++)
 53           printf("%.3lf ",a[i][j]);
 54         printf("%.3lf\n",b[i]);
 55     }
 56     printf("\n");
 57 }
 58  
 59 void gaos()
 60 {
 61     int now=1;
 62     for(int k=1;k<n;++k,++now)
 63     {
 64         int order=k;
 65         for(int i=k+1;i<=n;i++)
 66           if(abss(a[order][k])<abss(a[i][k]))
 67             order=i;
 68         if(order!=k)
 69         {
 70             for(int i=k;i<=n;i++)
 71               swap(a[now][i],a[order][i]);
 72             swap(b[now],b[order]);
 73         }
 74         if(!a[now][k])
 75         {
 76             --now;
 77             continue;
 78         }
 79         double temp;
 80         for(int i=now+1;i<=n;i++)
 81         {
 82             temp=a[i][k]/a[now][k];
 83             for(int j=k;j<=n+1;j++)
 84               a[i][j]-=temp*a[k][j];
 85             b[i]-=temp*b[k];
 86         }
 87     }
 88     
 89     for(int i=n;i>=1;--i)
 90     {
 91         for(int j=n;j>=i+1;--j)
 92           b[i]-=a[i][j]*ans[j];
 93         ans[i]=b[i]/a[i][i];
 94     }
 95     for(int k=1;k<=n;++k)
 96         printf("%.9lf\n",ans[k]);
 97 }
 98  
 99 int main(){
100     freopen("dotp.in","r",stdin);
101     freopen("dotp.out","w",stdout);
102     mem(first,-1);
103     scanf("%d%d%d%d",&n,&m,&p,&q);
104     for(int i=1;i<=m;i++)
105     {
106         scanf("%d%d",&u,&o);
107         addbian(u,o);
108         addbian(o,u);
109         ++degree[o];
110         ++degree[u];
111     }
112     for(int i=1;i<=n;i++)
113       a[1][i]=1;
114     b[1]=1;
115     int temp;
116     double qq=1-(double)p/q;
117     for(int i=2;i<=n;i++)
118     {
119         a[i][i]=1;
120         for(int j=first[i];j!=-1;j=a1[j].next)
121         {
122             temp=a1[j].v;
123             a[i][temp]=-(1.0/degree[temp])*qq;
124         }
125         b[i]=0;
126     }
127     gaos();
128     
129     //while(1);
130     return 0;
131 }
code
posted @ 2017-08-14 19:53  A_LEAF  阅读(146)  评论(0编辑  收藏  举报