zoukankan      html  css  js  c++  java
  • 【bzoj3270】博物馆

    同样是高斯消元,我写的版本就受到了歧视

    我怎么又犯把 $j$ 打成 $i$ 这种 $sb$ 错误


    题意

    一张无向图,两个人分别从 $s_1$ 号点和 $s2$ 号点开始,每轮两人都会同时进行一次以下操作:有 $p_i$ 的概率留在原地,剩下的 $1-p_i$ 的概率等概率随机选择一条出边走到连向的点。问两人在每个点相遇的概率(边上不相遇)。

    $nle 20$

    题解

    首先得知道,这种题的概率就是期望,因为在点 $i$ 相遇($1le ile n$)的期望和是 $1$,总概率也是 $1$。

    期望和概率的一个重要差别就是期望可以不在 $[0,1]$ 范围内,之后转化的时候记住了。

    这种题都有一种特性,就是感觉可以拓扑排序 $dp$。

    但拓扑排序只适用于 $DAG$(有向无环图),对于有环图(包括无向图),没法确定 $dp$ 转移顺序。

    好了,现在我们降低自己的智商,回到小学状态,想想在小学时你怎么做这种难题

    如果你学过小学奥数(没学过也无所谓吧),应该能想到列方程。

    那方程是什么?就是我们的 $dp$ 式子,因为我们只是不确定转移顺序,但式子肯定是对的。

    众所周知,方程是互相之间都有约束性的,解必须同时满足所有方程,所以可以解决没法确定转移顺序的情况。

    所以这是高斯消元入门题。

     

    由于 $nle 20$,时间空间都随便用,先不用管会不会爆。

    考虑正常 $dp$,设 $f(i,j)$ 表示两人分别在点 $i,j$ 的概率。

    那 $f(i,j)$ 可以从它自己所有直接相连的点转移过来。

    我们再枚举哪个点转 $x$ 移到点 $i$,哪个点 $y$ 转移到点 $j$。

    分类讨论,

    当 $i=x,space j=y$ 时,从 $f(x,y)$ 转移到 $f(i,j)$ 的概率是 $p_x imes p_y$;

    当 $i=x,space j≠y$ 时,从 $f(x,y)$ 转移到 $f(i,j)$ 的概率是 $p_x imes (1-p_y)/du_y$;

    当 $i≠x,space j=y$ 时,从 $f(x,y)$ 转移到 $f(i,j)$ 的概率是 $(1-p_x) imes p_y/du_x$;

    当 $i≠x,space j≠y$ 时,从 $f(x,y)$ 转移到 $f(i,j)$ 的概率是 $(1-p_x) imes (1-p_y)/du_x/du_y$;

    其中 $du_x$ 表示点 $x$ 的度(就是有多少个点与它直接相连),注意点 $x$ 自己不算。

    我们把所有的 $f(i,j)space |space 1le ile n,space 1le jle n$ 都看作一个未知数,这样方程组共有 $n^2$ 个未知数,并且我们知道一个未知数表示的二元组(就是两人的位置)到另一个未知数表示的二元组的概率,也就是转移系数(其实就是系数,比如状态 $A$ 有 $k$ 的概率转移到 $B$,方程中就是 $B=kA$)是多少,这样就可以把所有转移系数列成矩阵,解方程组了。

    一行方程大概长这样(区分了一下 $i,x$ 和 $j,y$) $$f(i,j) space =space f(i,j) imes ... + f(i,y) imes ... + f(x,j) imes ... + f(x,y) imes ...$$

    把等号左边的 $f(i,j)$ 移到右边得到 $$0 space =space f(i,j) imes ... + f(i,y) imes ... + f(x,j) imes ... + f(x,y) imes ... - f(i,j)$$

    也就是对于转移到的状态 $f(i,j)$,它对应的那行方程 等号右边的 $f(i,j)$ 那项的系数先减 $1$。

    这个移项处理是高斯消元的终点,把未知数统一移到一边,常数统一移到另一边,省得再判等号两边的未知数。

    起点作为被转移项时,其方程的等号左边是 $-1$,因为两人一开始就都在起点,所以已经期望有 $1$ 次到达这种状态,等号右边(所有未知数那边)就加了 $1$,移项到另一边就变成减 $1$。

    再注意一点,两人终止的状态不再往其它状态转移,即两人的转移来源点 $x$ 和 $y$ 不能相等。

      1 #include<bits/stdc++.h>
      2 #define rep(i,x,y) for(int i=x;i<=y;++i)
      3 #define dwn(i,x,y) for(int i=x;i>=y;--i)
      4 #define rep_e(i,u) for(int i=hd[u];i;i=e[i].nxt)
      5 #define ll long long
      6 #define N 22
      7 #define eps 1e-7
      8 using namespace std;
      9 inline int read(){
     10     int x=0; bool f=1; char c=getchar();
     11     for(;!isdigit(c);c=getchar()) if(c=='-') f=0;
     12     for(; isdigit(c);c=getchar()) x=(x<<3)+(x<<1)+(c^'0');
     13     if(f) return x;
     14     return 0-x;
     15 }
     16 int n,m,a,b;
     17 double p[N],dp[N*N][N*N],ans[N*N];
     18 struct edge{int v,nxt;}e[(N*N)<<1];
     19 int hd[N],cnt,du[N];
     20 inline void add(int u,int v){e[++cnt]=(edge){v,hd[u]}, hd[u]=cnt;}
     21 inline int getCnt(int x,int y){return (x-1)*n+y;}
     22 bool Gauss(int n){
     23     rep(i,1,n){
     24         int r;
     25         //cout<<"wow:"<<i<<endl; 
     26         for(r=i; r<=n && fabs(dp[i][r])<eps; ++r);
     27         
     28         if(fabs(dp[i][r])<eps) return 0;
     29         /*
     30             rep(i,1,n){
     31                 rep(j,1,n+1) cout<<dp[i][j]<<' ';
     32                 cout<<endl;
     33             }*/
     34         if(i^r) swap(dp[i],dp[r]);
     35         /*
     36             rep(i,1,n){
     37                 rep(j,1,n+1) cout<<dp[i][j]<<' ';
     38                 cout<<endl;
     39             }*/
     40         double div=dp[i][i];
     41         //cout<<div<<endl;
     42         rep(j,i,n+1) dp[i][j]/=div;
     43         rep(j,i+1,n){
     44             div=dp[j][i];
     45             rep(k,i,n+1)
     46                 dp[j][k]-=dp[i][k]*div;
     47         }
     48     }
     49     /*
     50             rep(i,1,n){
     51                 rep(j,1,n+1) cout<<dp[i][j]<<' ';
     52                 cout<<endl;
     53             }*/
     54     ans[n]=dp[n][n+1];
     55     dwn(i,n-1,1){
     56         ans[i]=dp[i][n+1];
     57         rep(j,i+1,n) ans[i]-=ans[j]*dp[i][j];
     58     }
     59     /*
     60     rep(i,1,n){
     61         ans[i]=dp[fd][n+1]/dp[fd][fd];
     62     }*/ 
     63     return 1;
     64 }
     65 int main(){
     66     n=read(),m=read(),a=read(),b=read();
     67     int u,v,all=n*n;
     68     rep(i,1,n) add(i,i);
     69     rep(i,1,m)
     70         u=read(),v=read(),
     71         ++du[u],++du[v],
     72         add(u,v),add(v,u);
     73     rep(i,1,n)
     74         cin>>p[i];
     75     dp[getCnt(a,b)][all+1]=-1;
     76     rep(i,1,n)
     77         rep(j,1,n){
     78             int cnt1=getCnt(i,j);
     79             --dp[cnt1][cnt1]; //i=j时自己不能转移到自己,所以自己给自己的转移的期望次数-1 
     80             rep_e(ii,i)
     81                 rep_e(jj,j){
     82                     int x=e[ii].v, y=e[jj].v;
     83                     //cout<<"try:"<<i<<' '<<j<<' '<<x<<' '<<y<<endl;
     84                     if(x^y){
     85                         int cnt2=getCnt(x,y);
     86                         if(i==x && j==y) dp[cnt1][cnt2]+=p[i]*p[j];
     87                         else if(i==x) dp[cnt1][cnt2]+=p[i]*(1-p[y])/du[y];
     88                         else if(j==y) dp[cnt1][cnt2]+=(1-p[x])*p[j]/du[x];
     89                         else dp[cnt1][cnt2]+=(1-p[x])*(1-p[y])/du[x]/du[y];
     90                         //cout<<i<<' '<<j<<' '<<x<<' '<<y<<' '<<cnt1<<' '<<cnt2<<' '<<p[x]<<' '<<p[y]<<' '<<du[x]<<' '<<du[y]<<' '<<dp[cnt1][cnt2]<<endl;
     91                     }
     92                 }
     93         }
     94     
     95     if(!Gauss(all)){printf("I AK IOI!
    "); return -1;}
     96     rep(i,1,n) printf("%.6lf ",ans[getCnt(i,i)]);
     97     return 0;
     98 }
     99 /*
    100 3 2 1 3
    101 1 2
    102 2 3
    103 0.5
    104 0.5
    105 0.5
    106 显然,ans[1]和ans[3]应该相等,否则就是程序写错了 
    107 */
    View Code

    怎么网上那么多人写的代码都一个样(包括高斯消元部分),这么不求上进的吗?

    我写的高斯消元跟他们不一样,一开始结果不对,我还以为我的高斯消元是错的。

    后来调对了,才意识到只是被鄙视了而已,因为我跟其他人写的都不一样。

    然后呢?

  • 相关阅读:
    MySQL之数据库结构优化
    MySQL之索引
    Spring之单元测试
    Spring之IOC容器加载初始化的方式
    LeetCode之Sort List
    [译]Java 垃圾回收的监控和分析
    [译]Java垃圾回收器的类型
    [译]Java垃圾回收是如何工作的
    [译]Java 垃圾回收介绍
    JSP之项目路径问题(${pageContext.request.contextPath},<%=request.getContextPath()%>以及绝对路径获取)
  • 原文地址:https://www.cnblogs.com/scx2015noip-as-php/p/bzoj3270.html
Copyright © 2011-2022 走看看