sss

<更新提示>

---

<正文>

Broken Robot

Description

你作为礼物收到一个非常聪明的机器人走在矩形板上。不幸的是，你明白它已经破碎并且行为相当奇怪（随机）。该板由N行和M列单元组成。机器人最初位于第i行和第j列的某个单元格中。然后在每一步，机器人都可以去另一个细胞。目的是走到最底层（N.排。机器人可以停留在当前单元格中，向左移动，向右移动或移动到当前单元格下方的单元格。如果机器人位于最左侧的列中，则它不能向左移动，如果它位于最右侧的列中，则它不能向右移动。在每一步中，所有可能的动作都是同样可能的。返回预期的步数以到达最下面一行。

Input Format

在第一行中你将被提供两个空间隔开的整数N和M(1≤N,M≤?1000)。在第二行中你将得到另外两个空间隔开的整数i和j(1≤i≤N,1≤j≤M) ——初始行的数目和初始列的数量。注意，（1,1）是板的左上角，（N， M）是右下角。

Output Format

在自身的一行上输出预期的步数，小数点后至少有4位数。

Sample Input

10 14
5 14

Sample Output

18.0038068653

解析

显然，这道题看起来很像一道期望(dp)，那么我们就用期望(dp)的套路设一个状态试试：(f[i][j])代表机器人从((i,j))走到最后一行的步数期望值。

我们可以根据移动规则很容易列出(dp)方程：

(1.) 当机器人处于第一列时：(f[i][1]=frac{1}{3}(f[i][1]+f[i][2]+f[i+1][1])+1)

(2.) 当机器人处于最后一列时：(f[i][m]=frac{1}{3}(f[i][m]+f[i][m-1]+f[i+1][m])+1)

(3.) 当机器人处于中间列时：(f[i][j]=frac{1}{4}(f[i][j]+f[i][j-1]+f[i][j+1]+f[i+1][j])+1)

我们发现，在行维度上，这个状态转移方程时没有问题的，可以倒序枚举每一行作为阶段，来进行转移。

但是，在同一行中，这个状态转移方程并不满足无后效性这一动态规划基本原则，于是我们决定使用高斯消元算法来解方程。

总体上，我们还是以行号为阶段，倒序进行转移。在第(i)行行内，我们将(i+1)行的状态看为常数，剩下的状态看做(m)个未知数，(m)个状态转移方程看做数学方程，尝试列出增广矩阵。

先看第一类方程：(f[i][1]=frac{1}{3}(f[i][1]+f[i][2]+f[i+1][1])+1)，简单做一下移项分类：

[2f[i][1]-f[i][2]=f[i+1][1]+3 ]

同理，可以化简剩下两个方程：

[-f[i][j-1]+3f[i][j]-f[i][j+1]=f[i+1][j]+4 ]

[-f[i][m-1]+2f[i][m]=f[i+1][m]+3 ]

那就可以写出系数矩阵了：

[left[ egin{array}{ccccccc|c} 2 & -1 & 0 & cdots & 0 & 0 & 0 & f[i+1][1]+3 \ -1 & 3 & -1 & cdots & 0 & 0 & 0 & f[i+1][2]+4 \ vdots & vdots & vdots & ddots & vdots & vdots & vdots & vdots \ 0 & 0 & 0 & cdots& -1 & 3 & -1 & f[i+1][m-1]+4 \ 0 & 0 & 0 & cdots& 0 & -1 & 2 & f[i+1][m]+3 \ end{array} ight] ]

我们发现这个矩阵很特殊，总共只有三条斜列有系数，所以我们可以线性直接消元。具体的说，我们可以从上往下消一遍，将第一斜列的系数消去，同时正确地处理第二第三斜列。再从下往上消一遍，将第三斜列的系数消去，这样就可以直接计算答案了。

利用如上的(dp)以及高斯消元算法，时间复杂度为(O(nm))。

值得注意的是，这三个方程在(m=1)是会出现边界问题，简单推导可知(m-1)时答案就是((n-x)*2)。

(Code:)

#include <bits/stdc++.h>
using namespace std;
const int N = 1020;
int n,m,x,y;
double f[N][N],a[N][N],b[N];
inline void input(void)
{
    scanf("%d%d
%d%d",&n,&m,&x,&y);
}
inline void init(void)
{
    a[1][1] = 2.0 , a[1][2] = -1.0;
    for (int i=2;i<m;i++)
        a[i][i-1] = -1.0 , a[i][i] = 3.0 , a[i][i+1] = -1.0;
    a[m][m-1] = -1.0 , a[m][m] = 2.0;
}
inline void gauss(void)
{
    double rate = a[2][1] / a[1][1];
    a[2][1] = 0.0;
    a[2][2] -= a[1][2] * rate , b[2] -= b[1] * rate;
    for (int i=2;i<m;i++)
    {
        rate = a[i+1][i] / a[i][i];
        a[i+1][i] = 0.0;
        a[i+1][i+1] -= a[i][i+1] * rate , b[i+1] -= b[i] * rate;
    }
    for (int i=m-1;i>=1;i--)
    {
        rate = a[i][i+1] / a[i+1][i+1];
        a[i][i+1] = 0.0 , b[i] -= b[i+1] * rate; 
    }
}
inline void dp(void)
{
    for (int i=1;i<=m;i++)
        f[n][i] = 0;
    for (int i=n-1;i>=1;i--)
    {
        init();
        b[1] = f[i+1][1] + 3.0;
        for (int j=2;j<m;j++) 
            b[j] = f[i+1][j] + 4.0;
        b[m] = f[i+1][m] + 3.0;
        gauss();
        for (int j=1;j<=m;j++) 
            f[i][j] = b[j] / a[j][j];
    }
}
int main(void)
{
    freopen("in.in","r",stdin);
    freopen("out.out","w",stdout);
    input();
    if ( m != 1 ) dp();
    else f[x][y] = ( n - x ) * 2.0;
    printf("%.10lf
",f[x][y]);
    return 0;
}

---

<后记>