HDU--2262（高斯消元，概率DP）
2015-04-13 22:49:30
思路：参考了爱神的博客A的... 也改善了下高斯的模板。
首先给每个格子编号，假设某个格子（设为第 i 个）到终点的期望为Ei，考虑第 i 格的所有可行后继节点（可行节点的意思是这些后继节点可以到达终点）
记为En1,En2....Enk，总共k个，可列方程：Ei = (En1+En2+En3+ ... + Enk) / k + 1
转化为：En1+En2+En3+ ... + k*Eni = -k
如此这样可以列出 n*m 条方程，高斯消元求出主元即可。终点的方程稍特殊：1*Eed = 0
关于求可行后继节点，可以从每个终点（注意，可能有多个出口）开始bfs求出。
1 #include <cstdio> 2 #include <cstring> 3 #include <cstdlib> 4 #include <cmath> 5 #include <vector> 6 #include <map> 7 #include <set> 8 #include <stack> 9 #include <queue> 10 #include <string> 11 #include <iostream> 12 #include <algorithm> 13 using namespace std; 14 15 #define MEM(a,b) memset(a,b,sizeof(a)) 16 #define REP(i,n) for(int i=0;i<(n);++i) 17 #define FOR(i,a,b) for(int i=(a);i<=(b);++i) 18 #define getmid(l,r) ((l) + ((r) - (l)) / 2) 19 #define MP(a,b) make_pair(a,b) 20 21 typedef long long ll; 22 typedef pair<int,int> pii; 23 const int INF = (1 << 30) - 1; 24 const int MAXN = 20; 25 const int MAXM = 250; 26 const double eps = 1e-9; 27 28 int n,m,top; 29 char g[MAXN][MAXN]; 30 double A[MAXM][MAXM]; 31 bool can[MAXN][MAXN]; 32 int dir[4][2] = {{-1,0},{1,0},{0,-1},{0,1}}; 33 34 void Bfs(){ 35 MEM(can,false); 36 queue<int> x,y; 37 while(!x.empty()) x.pop(); 38 while(!y.empty()) y.pop(); 39 REP(i,n) REP(j,m) if(g[i][j] == '$'){ 40 x.push(i); 41 y.push(j); 42 can[i][j] = true; 43 } 44 while(!x.empty()){ 45 int curx = x.front(); x.pop(); 46 int cury = y.front(); y.pop(); 47 for(int k = 0; k < 4; ++k){ 48 int tx = curx + dir[k][0]; 49 int ty = cury + dir[k][1]; 50 if(tx >= 0 && tx < n && ty >= 0 && ty < m && g[tx][ty] != '#'){ 51 if(can[tx][ty] == false){ 52 can[tx][ty] = true; 53 x.push(tx); 54 y.push(ty); 55 } 56 } 57 } 58 } 59 } 60 61 bool Gauss(){ 62 int i,col; 63 for(i = 0,col = 0; col < top; ++i,++col){ 64 int r = i; 65 for(int j = i + 1; j < top; ++j) 66 if(fabs(A[j][col]) > fabs(A[r][col])) r = j; 67 if(r != i) for(int j = col; j <= top; ++j) swap(A[r][j],A[i][j]); 68 if(fabs(A[i][col]) < eps){ 69 --i; 70 continue; 71 } 72 for(int j = top; j >= col; --j) A[i][j] /= A[i][col]; 73 for(int k = 0; k < top; ++k) if(k != i) 74 for(int j = top; j >= col; --j) 75 A[k][j] -= A[k][col] * A[i][j]; 76 } 77 for(; i < top; ++i) if(!(fabs(A[i][top]) < eps)) return false; 78 return true; 79 } 80 81 int main(){ 82 while(scanf("%d%d",&n,&m) != EOF){ 83 top = n * m; 84 REP(i,n) scanf("%s",g[i]); 85 REP(i,top) REP(j,top + 1) A[i][j] = 0.0; 86 Bfs(); 87 int sx,sy; 88 REP(i,n) REP(j,m) if(g[i][j] == '@'){ 89 sx = i; 90 sy = j; 91 break; 92 } 93 REP(i,n) REP(j,m) if(g[i][j] != '#'){ 94 int a = i * m + j; 95 if(g[i][j] == '$'){ 96 A[a][top] = 0; 97 A[a][a] = 1; 98 continue; 99 } 100 int cnt = 0; 101 REP(k,4){ 102 int tx = i + dir[k][0]; 103 int ty = j + dir[k][1]; 104 if(tx >= 0 && tx < n && ty >= 0 && ty < m 105 && g[tx][ty] != '#' && can[tx][ty]){ 106 A[a][tx * m + ty] = 1; 107 cnt++; 108 } 109 } 110 A[a][a] = A[a][top] = -1.0 * cnt; 111 } 112 if(can[sx][sy] && Gauss()){ 113 //找主元 114 for(int i = 0; i < top; ++i){ 115 if(fabs(A[i][sx * m + sy] - 1.0) < eps){ 116 printf("%.6f\n",A[i][top]); 117 break; 118 } 119 } 120 } 121 else printf("-1\n"); 122 } 123 return 0; 124 }
