Bzoj3270:博物館:概率與期望,高斯消元
阿新 • • 發佈:2019-02-15
題目連結:博物館
我們用id[i][j]代表一個人到了i另一個人在j的狀態
假設id[i][j]代表的狀態可以一步走到id[x][y]代表的狀態,那麼id[x][y]一步也可以走到id[i][j]
所以狀態之間的轉移形成了一個環,這時候要用高斯消元來解決一下了
設a[id[x][y]][id[i][j]]為從狀態id[i][j]轉移到id[x][y]的概率來列個方程組
假設現在在id[i][j]這個狀態上,接下來有這麼幾種情況
1:兩個人同時停留在原點,a[id[i][j]][id[i][j]]=p[i]*p[j];
2:第一個人走到了x(前提是i->x有邊存在),a[id[x][j]][id[i][j]]+=(1-p[i])/du[i]*p[j];
3:第二個人走到了y,a[id[i][y]][id[i][j]]+=p[i]*(1-p[j])/du[j];
4:第一個人走到了x,第二個人走到了y,a[id[x][y]][id[i][j]]+=(1-p[i])/du[i]*(1-p[j])/du[j];
這樣我們列出了n*n個方程,形如a[id[x][y]][id[i][j]]=ax+by+cz+...
然後移項即可,只有一個特例即初始點,他的概率要+1,所以移項後常數項為-1,其餘的都為0
消元即可
#include<cmath> #include<queue> #include<vector> #include<cstdio> #include<cstdlib> #include<iostream> #include<algorithm> #define pb push_back using namespace std; const int maxn=21; const int maxm=410; int n,m,id[maxn][maxn],cnt=0,x,y,du[maxn]; double a[maxm][maxm],equ[maxm],p[maxn]; vector<int>nxt[maxn]; void guass(){ for (int i=1;i<=cnt;++i){ int tmp=i; for (int j=i+1;j<=cnt;++j) if (abs(a[j][i])>abs(a[tmp][i])) tmp=j; if (tmp!=i) for (int j=1;j<=cnt+1;++j) swap(a[i][j],a[tmp][j]); for (int j=i+1;j<=cnt;++j){ double N=a[j][i]/a[i][i]; for (int k=i;k<=cnt+1;++k) a[j][k]-=a[i][k]*N; } } equ[cnt]=a[cnt][cnt+1]/a[cnt][cnt]; for (int i=cnt-1;i>=1;--i){ double N=a[i][cnt+1]; for (int j=i+1;j<=cnt;++j) N-=a[i][j]*equ[j]; equ[i]=N/a[i][i]; } for (int i=1;i<=n;++i) printf("%.6lf ",equ[id[i][i]]); } int main(){ scanf("%d%d%d%d",&n,&m,&x,&y); for (int i=1;i<=m;++i){ int u,v; scanf("%d%d",&u,&v); nxt[u].pb(v); nxt[v].pb(u); du[u]++; du[v]++; } 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]=++cnt;if(i!=j)a[cnt][cnt]=p[i]*p[j];} for (int i=1;i<=n;++i) for (int j=1;j<=n;++j) if (i!=j){ for (int k=0;k<nxt[i].size();++k){ int tx=nxt[i][k]; a[id[tx][j]][id[i][j]]+=(1-p[i])/du[i]*p[j]; } for (int k=0;k<nxt[j].size();++k){ int ty=nxt[j][k]; a[id[i][ty]][id[i][j]]+=p[i]*(1-p[j])/du[j]; } for (int k=0;k<nxt[i].size();++k) for (int q=0;q<nxt[j].size();++q){ int tx=nxt[i][k],ty=nxt[j][q]; a[id[tx][ty]][id[i][j]]+=(1-p[i])*(1-p[j])/du[i]/du[j]; } } for (int i=1;i<=n;++i) for (int j=1;j<=n;++j) a[id[i][j]][id[i][j]]-=1.0; a[id[x][y]][cnt+1]=-1; guass(); }