CF 24 D. Broken robot

D. Broken robot

链接、

题意：

　　一个方格，从(x,y)出发，等价的概率向下，向左，向右，不动。如果在左右边缘上，那么等价的概率不动，向右/左，向下。走到最后一行即结束。求期望结束的步数。

分析：

　　因为不能往上走，所以行与行之间存在转移，即上一行转移到下一行。

　　同一行内的位置可以互相转移，所以可以对每一行内进行高斯消元，那么复杂度是$O(n^4)$，但是发现高斯消元的矩阵中每行只有三个位置有数，这个矩阵叫三对角矩阵，观察这个矩阵，发现可以O(n)消元。复杂度$O(n^2)$

代码：

#include<cstdio>
#include<algorithm>
#include<cstring>
#include<iostream>
#include<cmath>
#include<cctype>
#include<set>
#include<queue>
#include<vector>
#include<map>
using namespace std;
typedef long long LL;

inline int read() {
    int x=0,f=1;char ch=getchar();for(;!isdigit(ch);ch=getchar())if(ch=='-')f=-1;
    for(;isdigit(ch);ch=getchar())x=x*10+ch-'0';return x*f;
}

const int mod = 998244353, N = 1005;
const double f31 = 1.0 / 3.0, f32 = 2.0 / 3.0, f41 = 1.0 / 4.0, f43 = 3.0 / 4.0, f38 = 8.0 / 3.0, f23 = 3.0 / 2.0;

double dp[N][N], a[N][N];
int n, m, x, y; 

int ksm(int a,int b) {
    int res = 1;
    while (b) {
        if (b & 1) res = 1ll * res * a % mod;
        a = 1ll * a * a % mod;
        b >>= 1;
    }
    return res;
}

void Make(int i) {
    a[1][1] = f32;
    a[1][2] = -f31;
    a[1][m + 1] = f31 * dp[i + 1][1] + 1;
    a[m][m - 1] = -f31;
    a[m][m] = f32;
    a[m][m + 1] = f31 * dp[i + 1][m] + 1;
    for (int j = 2; j < m; ++j) {
        a[j][j - 1] = -f41;
        a[j][j + 1] = -f41;
        a[j][j] = f43;
        a[j][m + 1] = f41 * dp[i + 1][j] + 1;
    }
}
void pr(int l,int r) {
    for (int i = l; i <= r; ++i) {
        for (int j = 1; j <= m + 1; ++j)
            printf("% .2lf ", a[i][j]);
        puts("");
    }
    puts("");
}
void solve(int i) {
    a[1][1] = 1; a[1][2] *= f23; a[1][m + 1] *= f23;
    for (int j = 2; j < m; ++j) {
        a[j][j - 1] = 0;
        a[j][j + 1] /= (a[j][j] + f41 * a[j - 1][j]);
        a[j][m + 1] += f41 * a[j - 1][m + 1];
        a[j][m + 1] /= (a[j][j] + f41 * a[j - 1][j]);
        a[j][j] = 1;
    }
    a[m][m - 1] = 0;
    a[m][m + 1] += (f31 * a[m - 1][m + 1]);
    a[m][m + 1] /= (a[m][m] + f31 * a[m - 1][m]);
    a[m][m] = 1;

    dp[i][m] = a[m][m + 1];
    for (int j = m - 1; j >= 1; --j) 
        dp[i][j] = a[j][m + 1] - a[j][j + 1] * dp[i][j + 1];
}
int main() {
    n = read(), m = read(), x = read(), y = read();
    if (m == 1) {
        printf("%.10lf", 2.0 * (n - x)); return 0;
    }
    for (int i = n - 1; i >= x; --i) {
        Make(i);
        solve(i);
    }
    printf("%.10lf", dp[x][y]);
    return 0;
}