用高斯消元法(gauss)求解一类期望问题
http://acm.zjut.edu.cn/ShowProblem.aspx?ShowID=1423
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可达这个条件,不然就会导致无解的情况。
1 #include<iostream> 2 #include<string> 3 #include<cmath> 4 #include<algorithm> 5 using namespace std; 6 #define MAXN 800 7 8 double mat[MAXN][MAXN]; 9 char map[40][40]; 10 int num[40][40]; 11 12 struct node 13 { 14 int x; 15 int y; 16 }; 17 node Star,end; 18 int n,m,len; 19 int di[4][2]={1,0,0,1,-1,0,0,-1}; 20 int cut; 21 22 int isok(int x,int y) 23 { 24 return x >= 0 && x < n && y >= 0 && y < m && map[x][y]!='X'; 25 } 26 27 node queue[2*MAXN]; 28 int head,tail; 29 30 void fillnum(int x,int y) //对地图进行重新标号,标号从0开始 31 { 32 num[x][y]=++cut; 33 tail=head=0; 34 node e,q; 35 e.x=x;e.y=y; 36 queue[tail++]=e; 37 while(head<tail) 38 { 39 e=queue[head++]; 40 for(int i=0;i<4;i++) 41 { 42 q.x=e.x+di[i][0]; 43 q.y=e.y+di[i][1]; 44 if(isok(q.x,q.y) && num[q.x][q.y]==-1) 45 { 46 num[q.x][q.y]=++cut; 47 queue[tail++]=q; 48 } 49 } 50 } 51 } 52 53 bool gauss(int n) 54 { 55 int i, j, row, idx; 56 double buf, maxx; 57 for(row = 0; row < n; row ++) 58 { 59 for(maxx = 0, i = row; i < n; i ++) 60 { 61 if(maxx < fabs(mat[i][row])) 62 { 63 maxx = fabs(mat[i][row]); 64 idx = i; 65 } 66 } 67 if(maxx == 0) 68 continue; 69 if(idx != row) 70 { 71 for(i = row; i <= n; i ++) 72 swap(mat[row][i], mat[idx][i]); 73 } 74 for(i = row + 1; i < n; i ++) 75 { 76 if(fabs(mat[i][row])<1e-8) 77 continue; 78 buf = mat[i][row] / mat[row][row]; 79 for(j = row; j <= n; j ++) 80 mat[i][j] -= buf * mat[row][j]; 81 } 82 } 83 for(i = n - 1; i >= 0; i --) 84 { 85 for(j = i + 1; j < n; j ++) 86 mat[i][n] -= mat[i][j] * mat[j][j]; 87 mat[i][i] = mat[i][n] / mat[i][i]; 88 } 89 return true; 90 } 91 92 int main() 93 { 94 int i,j,k,x,y; 95 freopen("D:\\in.txt","r",stdin); 96 while(scanf("%d%d",&n,&m)!=EOF) 97 { 98 len=0; 99 for(i=0;i<n;i++) 100 { 101 scanf("%s",map[i]); 102 for(j=0;j<m;j++) 103 { 104 if(map[i][j]=='D') 105 { 106 Star.x=i; 107 Star.y=j; 108 } 109 else if(map[i][j]=='E') 110 { 111 end.x=i; 112 end.y=j; 113 } 114 } 115 } 116 memset(num,255,sizeof(num)); 117 cut=-1; 118 fillnum(Star.x,Star.y); 119 if(num[end.x][end.y]==-1) 120 { 121 printf("tragedy!\n"); 122 continue; 123 } 124 memset(mat,0,sizeof(mat)); 125 for(i=0;i<n;i++) //建立N元一次方程组 126 { 127 for(j=0;j<m;j++) 128 { 129 if(num[i][j]!=-1) 130 { 131 int now=num[i][j]; 132 int count=0; 133 for(k=0;k<4;k++) 134 { 135 x=i+di[k][0];y=j+di[k][1]; 136 if(isok(x,y)) 137 { 138 mat[now][num[x][y]]=-1; 139 count++; 140 } 141 mat[now][now]=count; 142 mat[now][cut+1]=count; 143 } 144 } 145 } 146 } 147 x=num[end.x][end.y]; 148 memset(mat[x],0,sizeof(mat[x])); 149 mat[x][x]=1; 150 if(gauss(cut+1)) 151 { 152 if(mat[num[Star.x][Star.y]][num[Star.x][Star.y]]<=1000000) 153 printf("%0.2lf\n",mat[num[Star.x][Star.y]][num[Star.x][Star.y]]); 154 else 155 printf("tragedy!\n"); 156 } 157 else 158 { 159 printf("tragedy!\n"); 160 } 161 } 162 return 0; 163 }
http://acm.zjut.edu.cn/ShowProblem.aspx?ShowID=1317
Description:
m个人位于正m边形的顶点上,彼此抛掷飞盘。他们共有两个飞盘,且开始时这两个飞盘位于相距为n的两个人的手中(相邻两个人相距为1,依此类推)。在每次抛掷时两个飞盘被同时抛出,飞盘都以1/2的概率被抛到掷飞盘的人左边相邻的人,1/2的概率被抛到右边相邻的人。此过程一直进行,直到两个飞盘被掷到同一个人手中,求此抛掷飞盘的游戏平均情况下(期望)会在抛掷几次后结束。
Input:
每行有两个整数m (2<m<=100),n (0 < n < m)。
这题我们以两个飞盘的距离为状态进行转移。
那么E[n]=E[n+2]/4+E[n-2]/4+E[n]/2+1,化简成:2E[n]-E[n+2]-E[n-2]=4。
首先对于两个飞盘给定的起始距离n,我们可以先搜索一下可否到达状态0,如果不行,则直接输出INF。在搜索的过程中,顺便把重新标号也进行了。为什么这题也要重新标号呢?
因为该题中,假设m是偶数,那么对于任意的n,n+1和n-1都是不可达的状态。请一定记得,如果方程组中有不可达的点的话,就会导致无解的情况!
接下来对每一个可达的状态建立一个如上的方程,最后用高斯消元法解除E[No[n]]即可!
1 #include<iostream> 2 #include<string> 3 #include<cmath> 4 #include<stdio.h> 5 #include<memory.h> 6 #include<algorithm> 7 using namespace std; 8 #define maxn 200 9 10 double mat[maxn][maxn]; 11 int num[maxn]; 12 int n,m; 13 int cut; 14 15 int getnum(int x) //获取x的状态,对于一个环来说,x和m-x是一样的,我们取小的 16 { 17 if(x<0) 18 x=-x; 19 if(x>n/2) 20 x=n-x; 21 return x; 22 } 23 24 void dfs(int x) //重标号 25 { 26 num[x]=cut++; 27 int y=getnum(x+2); 28 if(num[y]==-1) 29 dfs(y); 30 y=getnum(x-2); 31 if(num[y]==-1) 32 dfs(y); 33 } 34 35 bool gauss(int n) 36 { 37 int i, j, row, idx; 38 double buf, maxx; 39 for(row = 0; row < n; row ++) 40 { 41 for(maxx = 0, i = row; i < n; i ++) 42 { 43 if(maxx < fabs(mat[i][row])) 44 { 45 maxx = fabs(mat[i][row]); 46 idx = i; 47 } 48 } 49 if(maxx == 0) return false; 50 if(idx != row) 51 { 52 for(i = row; i <= n; i ++) 53 swap(mat[row][i], mat[idx][i]); 54 } 55 for(i = row + 1; i < n; i ++) 56 { 57 buf = mat[i][row] / mat[row][row]; 58 for(j = row; j <= n; j ++) 59 mat[i][j] -= buf * mat[row][j]; 60 } 61 } 62 for(i = n - 1; i >= 0; i --) 63 { 64 for(j = i + 1; j < n; j ++) 65 mat[i][n] -= mat[i][j] * mat[j][j]; 66 mat[i][i] = mat[i][n] / mat[i][i]; 67 } 68 return true; 69 } 70 71 int main() 72 { 73 int i,len; 74 freopen("D:\\in.txt","r",stdin); 75 while(scanf("%d%d",&n,&m)==2) 76 { 77 memset(num,255,sizeof(num)); 78 if(m>n/2) 79 m=n-m; 80 cut=0; 81 dfs(m); 82 // for(i=0;i<=n/2;i++) 83 // cout<<num[i]<<" "; 84 // cout<<endl; 85 if(num[0]==-1) 86 { 87 printf("INF\n"); 88 continue; 89 } 90 memset(mat,0,sizeof(mat)); 91 len=n/2; 92 int now,j; 93 for(i=0;i<=len;i++) 94 { 95 if(num[i]!=-1) 96 { 97 now=num[i]; 98 mat[now][now]=2; 99 mat[now][cut]=4; 100 j=getnum(i-2); 101 mat[now][num[j]]-=1; 102 j=getnum(i+2); 103 mat[now][num[j]]-=1; 104 } 105 } 106 j=num[0]; 107 memset(mat[j],0,sizeof(mat[j])); 108 mat[j][j]=1; 109 // for(i=0;i<cut;i++) 110 // { 111 // for(j=0;j<=cut;j++) 112 // cout<<mat[i][j]<<" "; 113 // cout<<endl; 114 // } 115 if(gauss(cut)) 116 { 117 printf("%0.2lf\n",mat[num[m]][num[m]]); 118 } 119 else 120 { 121 printf("INF\n"); 122 } 123 } 124 return 0; 125 }