题面在这里
题意
给定一个\(n\times m\)的格网,从\((1,1)\)出发,每一格\((i,j)\)往上下左右移动的概率已经给出,询问到达\((n,m)\)的期望步数
数据范围
\[n,m\le40
\]
\]
sol
列方程高斯消元求解是显然的
\[E_{i,j}=1+p1_{i,j}\times E_{i+1,j}+p2_{i,j}\times E_{i,j+1}+p3_{i,j}\times E_{i-1,j}+p4_{i,j}\times E_{i,j-1}
\]
\]
但是直接套高斯消元复杂度为\(O(n^3m^3)\),会\(TLE\)
于是考虑优化
因为每项方程的系数只有\(4\)个,可以知道最后得到的矩阵满足一定的规律
比如当\(n=4,m=3\)时，系数矩阵如下
\[\left[
\begin{matrix}
1&-p2_{1,1}&0&0&-p1_{1,1}&0&0&0&0&0&0&0\\
-p4_{1,2}&1&-p2_{1,2}&0&0&-p1_{1,2}&0&0&0&0&0&0\\
0&-p4_{1,3}&1&-p2_{1,3}&0&0&-p1_{1,3}&0&0&0&0&0\\
0&0&-p4_{1,4}&1&0&0&0&-p1_{1,4}&0&0&0&0\\
-p3_{2,1}&0&0&0&1&-p2_{2,1}&0&0&-p1_{2,1}&0&0&0\\
0&-p3_{2,2}&0&0&-p4_{2,2}&1&-p2_{2,2}&0&-p1_{2,2}&0&0&0\\
0&0&-p3_{2,3}&0&0&-p4_{2,3}&1&-p2_{2,3}&0&-p1_{2,3}&0\\
0&0&0&-p3_{2,4}&0&0&-p4_{2,4}&1&-p2_{2,4}&0&-p1_{2,4}\\
0&0&0&0&-p3_{3,1}&0&0&0&1&-p2_{3,1}&0&0\\
0&0&0&0&0&-p3_{3,2}&0&0&-p4_{3,2}&1&-p2_{3,2}&0\\
0&0&0&0&0&0&-p3_{3,3}&0&0&-p4_{3,3}&1&-p2_{3,3}\\
0&0&0&0&0&0&0&0&0&0&0&1\\
\end{matrix}
\right]\times
\left[
\begin{matrix}
E_{1,1}\\
E_{1,2}\\
E_{1,3}\\
E_{1,4}\\
E_{2,1}\\
E_{2,2}\\
E_{2,3}\\
E_{2,4}\\
E_{3,1}\\
E_{3,2}\\
E_{3,3}\\
E_{3,4}\\
\end{matrix}
\right]=
\left[
\begin{matrix}
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
0\\
\end{matrix}
\right]
\]
\begin{matrix}
1&-p2_{1,1}&0&0&-p1_{1,1}&0&0&0&0&0&0&0\\
-p4_{1,2}&1&-p2_{1,2}&0&0&-p1_{1,2}&0&0&0&0&0&0\\
0&-p4_{1,3}&1&-p2_{1,3}&0&0&-p1_{1,3}&0&0&0&0&0\\
0&0&-p4_{1,4}&1&0&0&0&-p1_{1,4}&0&0&0&0\\
-p3_{2,1}&0&0&0&1&-p2_{2,1}&0&0&-p1_{2,1}&0&0&0\\
0&-p3_{2,2}&0&0&-p4_{2,2}&1&-p2_{2,2}&0&-p1_{2,2}&0&0&0\\
0&0&-p3_{2,3}&0&0&-p4_{2,3}&1&-p2_{2,3}&0&-p1_{2,3}&0\\
0&0&0&-p3_{2,4}&0&0&-p4_{2,4}&1&-p2_{2,4}&0&-p1_{2,4}\\
0&0&0&0&-p3_{3,1}&0&0&0&1&-p2_{3,1}&0&0\\
0&0&0&0&0&-p3_{3,2}&0&0&-p4_{3,2}&1&-p2_{3,2}&0\\
0&0&0&0&0&0&-p3_{3,3}&0&0&-p4_{3,3}&1&-p2_{3,3}\\
0&0&0&0&0&0&0&0&0&0&0&1\\
\end{matrix}
\right]\times
\left[
\begin{matrix}
E_{1,1}\\
E_{1,2}\\
E_{1,3}\\
E_{1,4}\\
E_{2,1}\\
E_{2,2}\\
E_{2,3}\\
E_{2,4}\\
E_{3,1}\\
E_{3,2}\\
E_{3,3}\\
E_{3,4}\\
\end{matrix}
\right]=
\left[
\begin{matrix}
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
1\\
0\\
\end{matrix}
\right]
\]
(建议缩小屏幕后食用)
在这个矩阵中,每行的系数都占据了\((2m+1)\)的长度,且以\(f_{i,i}\)为中心
因此我们在高斯消元的时候,只需要消除后\(m\)行中的\(m\)个系数
枚举每一行是\(O(nm)\)的,因此总复杂度即为\(O(nm^3)\)
于是就可以过了......
代码
#include<bits/stdc++.h>
#include<algorithm>
#include<iostream>
#include<cstdlib>
#include<iomanip>
#include<cstring>
#include<complex>
#include<vector>
#include<cstdio>
#include<string>
#include<bitset>
#include<cmath>
#include<queue>
#include<stack>
#include<map>
#include<set>
#define mp make_pair
#define pb push_back
#define RG register
#define il inline
using namespace std;
typedef unsigned long long ull;
typedef vector<int>VI;
typedef long long ll;
typedef double dd;
const dd eps=1e-10;
const int mod=1e9+7;
const int N=40;
const int M=50010*2;
il ll read(){
RG ll data=0,w=1;RG char ch=getchar();
while(ch!='-'&&(ch<'0'||ch>'9'))ch=getchar();
if(ch=='-')w=-1,ch=getchar();
while(ch<='9'&&ch>='0')data=data*10+ch-48,ch=getchar();
return data*w;
}
il void file(){
freopen("a.in","r",stdin);
freopen("a.out","w",stdout);
}
int n,m;
dd p[5][N+10][N+10],g[N*N+10][N*N+10];
il void gauss(int nm){
for(RG int i=1;i<=nm;i++){
for(RG int j=i;j<=min(nm,i+m);j++)
if(g[j][i]>eps){swap(g[i],g[j]);break;}
for(RG int j=i+1;j<=min(nm,i+m);j++){
g[j][nm+1]-=g[i][nm+1]*g[j][i]/g[i][i];
for(RG int k=min(nm,i+m);k>=i;k--)
g[j][k]-=g[i][k]*g[j][i]/g[i][i];
}
}
for(RG int i=nm;i;i--){
for(RG int j=i+1;j<=min(nm,i+m);j++)
g[i][nm+1]-=g[j][nm+1]*g[i][j];
g[i][nm+1]/=g[i][i];
}
}
int main()
{
while(1){
n=read();m=read();if(!n)return 0;
for(RG int k=1;k<=4;k++)
for(RG int i=1;i<=n;i++)
for(RG int j=1;j<=m;j++)
scanf("%lf",&p[k][i][j]);
for(RG int i=1;i<=n*m;i++)
for(RG int j=1;j<=n*m+1;j++)
g[i][j]=0;
for(RG int i=1;i<=n;i++)
for(RG int j=1;j<=m;j++){
RG int x=(i-1)*m+j;g[x][x]=g[x][n*m+1]=1;
if(i<n)g[x][i*m+j]=-p[1][i][j];
if(j<m)g[x][(i-1)*m+j+1]=-p[2][i][j];
if(i>1)g[x][(i-2)*m+j]=-p[3][i][j];
if(j>1)g[x][(i-1)*m+j-1]=-p[4][i][j];
g[x][n*m+1]=1;
}
g[n*m][(n-1)*m]=g[n*m][n*m-1]=0;g[n*m][n*m+1]=0;
gauss(n*m);
printf("%.6lf\n",g[1][n*m+1]);
}
return 0;
}