hdu4579 Random Walk

题目：

Yuanfang is walking on a chain. The chain has n nodes numbered from 1 to n. Every second, he can move from node i to node j with probability: 

c(i,j) is an element in a given parameter matrix which is n×m. (1 <= c(i, j) <= 9) 
Yuanfang wants to know the expectation time for him to walk from node 1 to node n.

Input

There are no more than 10 test cases. 
In each case, there are two integers n (2 <= n <= 50000), m (1 <= m <= 5), in the first line, meaning that there are n nodes and the parameter matrix is n×m . There are m integers in each of the next n lines which describe the parameter matrix . 
The input ends with 0 0. 

Output

For each case, output the expectation time for Yuanfang to walk from node 1 to node n in one line. The answer should be rounded to 2 digits after decimal point.

Sample Input

3 1
1
1
1
5 2
1 2
2 1
3 2
2 3
1 3
0 0

Sample Output

6.94
8.75

思路：设dp[i]表示从i到n的步数期望值

那么dp[n]=0

列出n个变元的n个方程，这样就能用高斯消元解出来了。

代码：

#pragma comment(linker, "/STACK:1024000000,1024000000")
#include<iostream>
#include<algorithm>
#include<ctime>
#include<cstdio>
#include<cmath>
#include<cstring>
#include<string>
#include<vector>
#include<map>
#include<set>
#include<queue>
#include<stack>
#include<list>
#include<numeric>
using namespace std;
#define LL long long
#define ULL unsigned long long
#define INF 0x3f3f3f3f
#define mm(a,b) memset(a,b,sizeof(a))
#define PP puts("*********************");
template<class T> T f_abs(T a){ return a > 0 ? a : -a; }
template<class T> T gcd(T a, T b){ return b ? gcd(b, a%b) : a; }
template<class T> T lcm(T a,T b){return a/gcd(a,b)*b;}
// 0x3f3f3f3f3f3f3f3f
// 0x3f3f3f3f

const int maxn=5e4+50;
const double eps=1e-6;
double c[maxn][6],p[maxn][11];
double a[maxn][11],b[maxn];
double dp[maxn];
int main(){

    int n,m;
    while(~scanf("%d%d",&n,&m)){
        if(!n&&!m) break;
        for(int i=1;i<=n;i++)
            for(int j=1;j<=m;j++)
                scanf("%lf",&c[i][j]);
        for(int i=1;i<n;i++){
            double sum=1,s=0;
            for(int j=1;j<=m;j++)
                sum+=c[i][j];
            for(int j=1;j<=m&&i-j>=1;j++){
                p[i][m-j]=0.3*c[i][j]/sum;
                s+=p[i][m-j];
            }
            for(int j=1;j<=m&&i+j<=n;j++){
                p[i][m+j]=0.7*c[i][j]/sum;
                s+=p[i][m+j];
            }
            p[i][m]=-s;
            b[i]=-1;
        }
        for(int i=1;i<n;i++){
            int L=max(1,i-m);
            int R=min(n,i+m);
            for(int j=L;j<i;j++){
                if(f_abs(p[i][m+j-i])<=eps) continue;
                double t=p[i][m+j-i]/a[j][1];
                for(int k=1;k<=m+1&&j+k-1<=n;k++){
                    p[i][m+j-i+k-1]-=t*a[j][k];
                }
                b[i]-=t*b[j];
            }
            for(int j=1;j<=R-i+1;j++)
                a[i][j]=p[i][m+j-1];
        }
        dp[n]=0;
        for(int i=n-1;i>=1;i--){
            for(int j=2;j<=m+1&&i+j-1<=n;j++)
                b[i]-=a[i][j]*dp[i+j-1];
            dp[i]=b[i]/a[i][1];
        }
        printf("%.2f\n",dp[1]);
    }
    return 0;
}