BZOJ 3270: 博物馆 [概率DP 高斯消元]

时间：2022-09-16 06:24:34

http://www.lydsy.com/JudgeOnline/problem.php?id=3270

题意：一张无向图，一开始两人分别在$x$和$y$，每一分钟在点$i$不走的概率为$p[i]$,走的话等概率走到相邻的点，求两人在每个点相遇的概率
对于100%的数据有 n <= 20，n-1 <= m <= n(n-1)/2

因为两个人，所以状态肯定要二元组呀
$f(i,j)$表示一人在$i$另一人在$j$的概率，转移方程：
$f(i,j)=f(i,j)p_ip_j-\sum\limits_{(x,i)\ ,\ (y,j) \in E}{f(x,i)p_jt_x+f(i,y)p_it_y+f(x,y)t_xt_y}$
$t_i=\frac{1-p_i}{d_i}$
然后有环没法$DP$.....
高斯消元解方程...$n^2$个方程$n^2$个变量

注意：
$1.$ $f(i,i)$不能走到其他啦(也不能不走啦)
$2.$ 一开始$f(x,y)$的概率除了走出来的左面还要加上$1$，因为一开始一定在$(x,y)$

#include <iostream>

#include <cstdio>

#include <cstring>

#include <algorithm>

#include <cmath>

using namespace std;

typedef long long ll;

const int N=,M=;

inline int read(){

    char c=getchar();int x=;

    while(c<''||c>''){c=getchar();}

    while(c>=''&&c<=''){x=x*+c-'';c=getchar();}

    return x;

}

int n,m,x,y,u,v,nn;

int d[N];

double p[N],a[N][N];

inline int id(int i,int j){return (i-)*n+j;}

struct edge{

    int v,ne;

}e[M<<];

int h[N],cnt=;

inline void ins(int u,int v){

    cnt++;

    e[cnt].v=v;e[cnt].ne=h[u];h[u]=cnt;

    cnt++;

    e[cnt].v=u;e[cnt].ne=h[v];h[v]=cnt;

}

double t[N];

void buildEquation(){

    for(int i=;i<=n;i++) t[i]=(-p[i])/d[i];

    for(int i=;i<=n;i++)

        for(int j=;j<=n;j++){

            int num=id(i,j);double *g=a[num];

            if(i==j) g[num]=;

            else g[num]=-p[i]*p[j];

            int x,y;

            for(int k=h[i];k;k=e[k].ne){

                x=e[k].v;

                if(x!=j) g[id(x,j)]=-p[j]*t[x];

            }

            for(int k=h[j];k;k=e[k].ne){

                y=e[k].v;

                if(y!=i) g[id(i,y)]=-p[i]*t[y];

            }

            for(int k=h[i];k;k=e[k].ne)

                for(int l=h[j];l;l=e[l].ne){

                    x=e[k].v,y=e[l].v;

                    if(x!=y) g[id(x,y)]=-t[x]*t[y];

                }

        }

    a[id(x,y)][nn+]=;

}

void GaussElimination(int n){

    for(int i=;i<=n;i++){

        int r=i;

        for(int j=i+;j<=n;j++)

            if(abs(a[j][i])>abs(a[r][i])) r=j;

        if(r!=i) for(int k=;k<=n+;k++) swap(a[i][k],a[r][k]);

        for(int j=i+;j<=n;j++){

            double t=a[j][i]/a[i][i];

            for(int k=i;k<=n+;k++) a[j][k]-=a[i][k]*t;

        }

    }

    for(int i=n;i>=;i--){

        for(int j=n;j>i;j--) a[i][n+]-=a[j][n+]*a[i][j];

        a[i][n+]/=a[i][i];

    }

}

int main(){

    freopen("in","r",stdin);

    n=read();m=read();nn=n*n;

    x=read();y=read();

    for(int i=;i<=m;i++)

        u=read(),v=read(),ins(u,v),d[u]++,d[v]++;

    for(int i=;i<=n;i++) scanf("%lf",&p[i]);

    buildEquation();

    GaussElimination(nn);

    for(int i=;i<=n;i++) printf("%.6lf ",a[id(i,i)][nn+]);

}

相关文章