zujt ZJUT 1423 地下迷宫 高斯消元法 概率DP

参考：

http://www.cnblogs.com/kuangbin/archive/2012/10/02/2710606.html

第二个高斯消元解概率dp，1ac

/*

地下迷宫

Description:
由于山体滑坡，DK被困在了地下蜘蛛王国迷宫。为了抢在DH之前来
到TFT，DK必须尽快走出此迷宫。此迷宫仅有一个出口，而由于大BOSS
的力量减弱影响到了DK，使DK的记忆力严重下降，他甚至无法记得他
上一步做了什么。所以他只能每次等概率随机的选取一个方向走。
当然他不会选取周围有障碍的地方走。如DK周围只有两处空地，则每
个都有1/2的概率。现在要求他平均要走多少步可以走出此迷宫。

Input:
先是一行两个整数N, M(1<=N, M<=10)表示迷宫为N*M大小，
然后是N行，每行M个字符，'.'表示是空地，'E’表示出口，'D’表示DK，'X’表示障碍。
Output:
如果DK无法走出或要超过1000000步才能走出，输出tragedy!，
否则输出一个实数表示平均情况下DK要走几步可以走出迷宫，四舍五入到小数点后两位。
Sample Input:
1 2
ED
3 3
D.X
.X.
X.E
Sample Output:
1.00
tragedy!


首先对地图节点重新标号。假设E[i]表示DK从i点开始走出迷宫的期望值。
那么E[i]=(E[a1]+E[a2]+E[a3]+...+E[an])/n+1，其中a1...an是i的相邻节点。
那么对于每一个DK可达的节点来说，都可以为它建立这样的一个方程。现
在假设DK可达的点有N个，那么我们最终将会得到N元一次方程组。最后
利用高斯消元解出E[No[S]]。其中S是DK的起点，No[S]是重标号后的起点
这里要重点注意的是，我们联立方程的时候，一定要注意DK可达这个条件，
不然就会导致无解的情况。

*/

//#pragma warning (disable: 4786)
//#pragma comment (linker, "/STACK:16777216")
//HEAD
#include <cstdio>
#include <ctime>
#include <cstdlib>
#include <cstring>
#include <queue>
#include <string>
#include <set>
#include <stack>
#include <map>
#include <cmath>
#include <vector>
#include <iostream>
#include <algorithm>
using namespace std;
//LOOP
#define FE(i, a, b) for(int i = (a); i <= (b); ++i)
#define FED(i, b, a) for(int i = (b); i>= (a); --i)
#define REP(i, N) for(int i = 0; i < (N); ++i)
#define CLR(A,value) memset(A,value,sizeof(A))
//INPUT
#define RI(n) scanf("%d", &n)
#define RII(n, m) scanf("%d%d", &n, &m)
#define RIII(n, m, k) scanf("%d%d%d", &n, &m, &k)
#define RS(s) scanf("%s", s)
typedef long long LL;
const int INF = 1000000007;


const double eps = 1e-9;
const int MAXN=110;
double a[MAXN][MAXN],x[MAXN];///方程的左边的矩阵和等式右边的值，求解之后x存的就是结果
int equ,var;///方程数和未知数个数


int Gauss()
{
    int i,j,k,col,max_r;
    for(k=0,col=0;k<equ&&col<var;k++,col++)
    {
        max_r=k;
        for(i=k+1;i<equ;i++)
          if(fabs(a[i][col])>fabs(a[max_r][col]))
            max_r=i;
        if(fabs(a[max_r][col])<eps)return 0;
        if(k!=max_r)
        {
            for(j=col;j<var;j++)
              swap(a[k][j],a[max_r][j]);
            swap(x[k],x[max_r]);
        }
        x[k]/=a[k][col];
        for(j=col+1;j<var;j++)a[k][j]/=a[k][col];
        a[k][col]=1;
        for(i=0;i<equ;i++)
          if(i!=k)
          {
              x[i]-=x[k]*a[i][k];
              for(j=col+1;j<var;j++)a[i][j]-=a[k][j]*a[i][col];
              a[i][col]=0;
          }
    }
    return 1;
}


char ch[12][12];
int num[12][12];
int ind[12][12];
int n, m;
int cnt;
int di[] = {0, 1, 0, -1};
int dj[] = {1, 0, -1, 0};
int si, sj;
int ei, ej;


bool check(int x, int y)
{
    if (x >= 0 && x < n && y >= 0 && y < m) return true;
    return false;
}
void bfs(int si, int sj)
{
    queue<int>q;
    CLR(num, 0);
    CLR(ind, -1);
    cnt = 0;
    q.push(si * 10 + sj);
    ind[si][sj] = cnt++;
    while (!q.empty())
    {
        int u = q.front(); q.pop();
        int ui = u / 10;
        int uj = u % 10;
        REP(r, 4)
        {
            int vi = ui + di[r];
            int vj = uj + dj[r];
            if (check(vi, vj) && ch[vi][vj] != 'X')
            {
                num[ui][uj]++;///
                if (ind[vi][vj] == -1)
                {
                    ind[vi][vj] = cnt++;
                    q.push(vi * 10 + vj);
                }
            }
        }
    }
}


void build()
{
    equ = var = cnt;
    CLR(a, 0);
    CLR(x, 0);
    REP(i, n) REP(j, m)
    {
        if (ind[i][j] != -1)
        {
            int nn = ind[i][j];
            if (i == ei && j == ej)
            {
                a[nn][nn] = 1;
                x[nn] = 0;
                continue;
            }

            a[nn][nn] = 1;
            REP(r, 4)
            {
                int vi = i + di[r];
                int vj = j + dj[r];
                if (check(vi, vj) && ch[vi][vj] != 'X')
                {
                    a[nn][ind[vi][vj]] -= 1.0 / num[i][j];
                }
            }
            x[nn] = 1;
        }
    }


}
int main ()
{
    while (~RII(n, m))
    {
        REP(i, n)
        {
            RS(ch[i]);
            int l = strlen(ch[i]);
            REP(j, m)
            {
                if (ch[i][j] == 'D')
                {
                    si = i;
                    sj = j;
                }
                else if (ch[i][j] == 'E')
                {
                    ei = i;
                    ej = j;
                }
            }
        }
        bfs(si, sj);
        if (ind[ei][ej] == -1)
        {
            puts("tragedy!");
            continue;
        }
        build();
        if (Gauss())
        {
            if (x[ind[si][sj]] <= 1000000)
                printf("%.2lf\n", x[ind[si][sj]]);
            else puts("tragedy!");
        }
        else puts("tragedy!");
    }
    return 0;
}