题面在这里
题意
给定一个(n imes m)的格网,从((1,1))出发,每一格((i,j))往上下左右移动的概率已经给出,询问到达((n,m))的期望步数
数据范围
[n,mle40
]
sol
列方程高斯消元求解是显然的
[E_{i,j}=1+p1_{i,j} imes E_{i+1,j}+p2_{i,j} imes E_{i,j+1}+p3_{i,j} imes E_{i-1,j}+p4_{i,j} imes E_{i,j-1}
]
但是直接套高斯消元复杂度为(O(n^3m^3)),会(TLE)
于是考虑优化
因为每项方程的系数只有(4)个,可以知道最后得到的矩阵满足一定的规律
比如当(n=4,m=3)时，系数矩阵如下
[left[
egin{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}
ight] imes
left[
egin{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}
ight]=
left[
egin{matrix}
1\
1\
1\
1\
1\
1\
1\
1\
1\
1\
1\
0\
end{matrix}
ight]
]
(建议缩小屏幕后食用)
在这个矩阵中,每行的系数都占据了((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
",g[1][n*m+1]);
}
return 0;
}