【HDOJ】4579 Random Walk

1. 题目描述
一个人沿着一条长度为n个链行走，给出了每秒钟由i到j的概率（$i,j in [1,n]$）。求从1开始走到n个时间的期望。

2. 基本思路
显然是个DP。公式推导也相当容易。不妨设$dp[i], i in [1,n]$表示由i到n的期望时间。
egin{align}
    dp[i] &= Sigma_{j=1}^{n} p(i, j) (dp[j] + 1),    &j<n\
    dp[i] &= 0 &i=n
end{align}
显然这是一个n元方程组，可以高斯消元解。但是因为n很大，因此不能直接套用高斯消元。但是通过观察系数矩阵可以发现规律。
以$n=7, m=2$为例，$ imes$表示$A_{ij}$不为0。

发现每个行向量最多包含$2m+1$个非零向量，即$[i-m, i+m]$。因此，在高斯消元的过程中，实际每次需要减掉的系数最多也就$2m+1$个。
因为$m in [1,5]$，可以直接模拟非零系数的消元。因为$i-m$有可能小于0，因此，将$[i-m,i+m]$映射到$[0,2m]$的区间内，$A_{ii}$恰好映射到$P_{im}$。
最后可以生成几个小规模的n，与高斯消元对拍一下。

3. 代码

/* 4579 */
#include <iostream>
#include <sstream>
#include <string>
#include <map>
#include <queue>
#include <set>
#include <stack>
#include <vector>
#include <deque>
#include <algorithm>
#include <cstdio>
#include <cmath>
#include <ctime>
#include <cstring>
#include <climits>
#include <cctype>
#include <cassert>
#include <functional>
#include <iterator>
#include <iomanip>
using namespace std;
//#pragma comment(linker,"/STACK:102400000,1024000")

#define sti                set<int>
#define stpii            set<pair<int, int> >
#define mpii            map<int,int>
#define vi                vector<int>
#define pii                pair<int,int>
#define vpii            vector<pair<int,int> >
#define rep(i, a, n)     for (int i=a;i<n;++i)
#define per(i, a, n)     for (int i=n-1;i>=a;--i)
#define clr                clear
#define pb                 push_back
#define mp                 make_pair
#define fir                first
#define sec                second
#define all(x)             (x).begin(),(x).end()
#define SZ(x)             ((int)(x).size())
#define lson            l, mid, rt<<1
#define rson            mid+1, r, rt<<1|1

const double eps = 1e-8;
const int maxn = 50005;
const int maxm = 12;
double g[maxn][maxm], p[maxn][maxm];
double v[maxn], x[maxn];
int C[maxn][6];
int n, m;

void solve() {
    int i, j, k;
    
    memset(p, 0, sizeof(p));
    for (i=1; i<n; ++i) {
        int tot = 1;
        double tmp = 0.0;
        for (j=1; j<=m; ++j)
            tot += C[i][j];
        
        for (j=1; j<=m; ++j) {
            if (i-j >= 1) {
                p[i][m-j] = 0.3 * C[i][j] / tot; 
                tmp += p[i][m-j];
            }
            if (i+j <= n) {
                p[i][m+j] = 0.7 * C[i][j] / tot; 
                tmp += p[i][m+j];
            }
        }
        p[i][m] = -tmp;
        v[i] = -1;
    }
    p[n][m] = 1;
    v[n] = 0;
    
    memcpy(g[1], p[1], sizeof(p[1]));
    for (i=2,k=1; i<=n; ++i,++k) {
        int l = max(k-m, 1);
        int r = min(k+m, n);
        for (j=i; j<=n&&j-k<=m; ++j) {
            if (fabs(p[k][m]) < eps)
                continue;
            double t = p[j][k-j+m] / p[k][m];
            for (int kk=k+1; kk<=n&&kk-k<=m; ++kk)
                p[j][kk-j+m] -= t * p[k][kk-k+m];
            v[j] -= t * v[k];
        }
        
        l = max(i-m, 1);
        r = min(i+m, n);
        for (j=l; j<=r; ++j)
            g[i][j-i+m] = p[i][j-i+m];
    }
    
    x[n] = 0;
    for (i=n-1,k=n; i>0; --i,--k){
        for (j=i; j>0&&k-j<=m; --j)
            v[j] -= x[k] * g[j][k-j+m];
        x[i] = v[i] / g[i][m];
    }
    
    printf("%.2lf
", x[1]);
}

int main() {
    ios::sync_with_stdio(false);
    #ifndef ONLINE_JUDGE
        freopen("data.in", "r", stdin);
        freopen("data.out", "w", stdout);
    #endif
    
    while (scanf("%d%d",&n,&m)!=EOF && (n||m)) {
        rep(i, 1, n+1)
            rep(j, 1, m+1)
                scanf("%d", &C[i][j]);
        solve();
    }
    
    #ifndef ONLINE_JUDGE
        printf("time = %d.
", (int)clock());
    #endif
    
    return 0;
}

4. 数据生成器

import sys
import string
from random import randint

    
def GenData(fileName):
    with open(fileName, "w") as fout:
        t = 10
        for tt in xrange(t):
            n = randint(1, 200)
            m = randint(1, 5)
            fout.write("%d %d
" % (n, m))
            L = [0] * m
            for i in xrange(n):
                for j in xrange(m):
                    L[j] = randint(1, 9)
                fout.write(" ".join(map(str, L)) + "
")    
        fout.write("0 0
")
                
        
def MovData(srcFileName, desFileName):
    with open(srcFileName, "r") as fin:
        lines = fin.readlines()
    with open(desFileName, "w") as fout:
        fout.write("".join(lines))

        
def CompData():
    print "comp"
    srcFileName = "F:Qt_prjhdojdata.out"
    desFileName = "F:workspacecpp_hdojdata.out"
    srcLines = []
    desLines = []
    with open(srcFileName, "r") as fin:
        srcLines = fin.readlines()
    with open(desFileName, "r") as fin:
        desLines = fin.readlines()
    n = min(len(srcLines), len(desLines))-1
    for i in xrange(n):
        ans2 = int(desLines[i])
        ans1 = int(srcLines[i])
        if ans1 > ans2:
            print "%d: wrong" % i

            
if __name__ == "__main__":
    srcFileName = "F:Qt_prjhdojdata.in"
    desFileName = "F:workspacecpp_hdojdata.in"
    GenData(srcFileName)
    MovData(srcFileName, desFileName)