zjut 1423 （高斯消元 + 求期望）

题目链接：http://cpp.zjut.edu.cn/ShowProblem.aspx?ShowID=1423

此题为高斯消元求期望的入门题目， 做法就是对图中非X的点重新标号， 令Ei表示从该点出发到达终点的期望走过的步数，对于终点的E值显然为0， 最后的结果为起点的E值，方程也不难得到。 注意普通的高斯消元化成上三角阵用在这类题目可能会有问题， 因为会出现方程无解的情况， 这样在高斯消元最后一步回代的过程中会出现问题， lrj的白书中给出了一种化为对角阵的高斯消元法， 避免了这个问题， 细节详见代码。

#include <vector>
#include <list>
#include <map>
#include <set>
#include <deque>
#include <stack>
#include <bitset>
#include <algorithm>
#include <cstring>
#include <functional>
#include <numeric>
#include <utility>
#include <sstream>
#include <iostream>
#include <iomanip>
#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <ctime>

using namespace std;

const int N = 205;
const double eps = 1e-8;

#define fi first
#define se second


double A[N][N];
int idx[15][15];
bool inf[N];
int d[N];
char mat[15][15];
int n, m, c, st, ed;
pair<int, int> pts[N];
vector<int> adj[N];
int dx[] = {0, 0, -1, 1};
int dy[] = {-1, 1, 0, 0};

bool ok(int x, int y) {
	return x >= 0 && x < n && y >= 0 && y < m && idx[x][y] != -1;
}

void gauss() {
	int i, j, k, r;
	for (i = 0; i < c; i++) {
		r = i;
		for (j = i + 1; j < c; j++)
			if (fabs(A[j][i]) > fabs(A[r][i])) r = j;
		if (fabs(A[r][i]) < eps) continue;
		if (r != i) for (j = 0; j <= c; j++) swap(A[r][j], A[i][j]);

		for (k = 0; k < c; k++) if (k != i)
			for (j = c; j >= i; j--)
				A[k][j] -= A[k][i] / A[i][i] * A[i][j];
	}
}

void gao() {
	memset(idx, -1, sizeof(idx));
	for (int i = 0; i < n; i++)
		scanf("%s", mat[i]);
	c = 0;
	for (int i = 0; i < n; i++)
		for (int j = 0; j < m; j++) {
			if (mat[i][j] != 'X') {
				pts[c].fi = i, pts[c].se = j;
				idx[i][j] = c;
				if (mat[i][j] == 'D') st = c;
				else if (mat[i][j] == 'E') ed = c;
				c++;		
			} 
		}
	
	for (int i = 0; i < c; i++)
		fill(A[i], A[i] + c + 1, 0);

	for (int i = 0; i < c; i++) {
		adj[i].clear();
		if (i == ed) {
			A[i][i] = 1;
			continue;
		}
		int x = pts[i].fi, y = pts[i].se;
		int tmp = 0;
		for (int j = 0; j < 4; j++) {
			int tx = x + dx[j], ty = y + dy[j];
			if (ok(tx, ty)) {
				tmp++;
				adj[i].push_back(idx[tx][ty]);	
			}
		}
		if (tmp)
			A[i][i] = 1;
		for (int j = 0; j < adj[i].size(); j++) {
			int t = adj[i][j];
			A[i][t] = -1.0 / tmp;
		}
		A[i][c] = 1;
	}
	
	gauss();

	memset(inf, 0, sizeof(inf));

	for (int i = c - 1; i >= 0; i--) {
		if (fabs(A[i][i]) < eps && fabs(A[i][c]) > eps) inf[i] = 1;
		for (int j = i + 1; j < c; j++) 
			if (fabs(A[i][j]) > eps && inf[j]) inf[i] = 1;
	}

}

int main() {
	while (~scanf("%d%d", &n, &m)) {
		gao();
		if (inf[st] || A[st][c] / A[st][st] > 1000000) 
			puts("tragedy!");
		else if (fabs(A[st][st]) < eps)
			puts("tragedy!");
		else
			printf("%.2f\n", A[st][c] / A[st][st]);
	}	
	return 0;
}