HDU 4418 Time travel 概率DP 高斯消元

第一次写高斯消元...

题目大意：

就是现在黑衣人中的K探员使用时空穿梭器在几个时间点之间穿梭，由于机器故障现在他从第X个点开始, 每次有p[i] (1 <= i <= M)的几率在经过i个点之后停下, 时间点是N个从左至右， 初始时有一个移动方向D， D = 0表示从左向右, D = 1从右向左, D = -1表示在两端处, 现在每次到两端时会转身向相反方向移动，如此循环, 问从X以方向D出发, 停留到达Y所需要经过的点数的期望

大致思路：

思路见代码注释部分

代码如下：

Result  :  Accepted     Memory  :  1472 KB     Time  :  296 ms

/*
 * Author: Gatevin
 * Created Time:  2014/12/23 15:46:41
 * File Name: Sora_Kasugano.cpp
 */
#include<iostream>
#include<sstream>
#include<fstream>
#include<vector>
#include<list>
#include<deque>
#include<queue>
#include<stack>
#include<map>
#include<set>
#include<bitset>
#include<algorithm>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cctype>
#include<cmath>
#include<ctime>
#include<iomanip>
using namespace std;
const double eps(1e-8);
typedef long long lint;

/*
 * 由于走的路是往返不停, 对于除起点和终点外的点都有两种状态
 * 所以将原序列以终点为中心对称一下, 起点和终点不对称， 比如1234变成123432
 * 这样到达终点之后直接就回到第一个点, 对称出来的点正好是在之前与本身的点相反方向运动的代表
 * 用E[i]表示在第i个点处到达目标点需要的期望经过点数
 * 则E[i] = sigma((E[i + j] + j)*p[i]) 1 <= j <= M, 初始时E[Y] = E[2*N - 2 - Y](对称点) = 0
 * 对于不能到达的点(p[i]可以为0)就不加入方程组, 建立方程是遍历一下看是否能到达即可判断Impossible的情况
 */

double a[210][210], x[210];
/*
 * a是方程左边的矩阵
 * x是等式右边的值
 */
int var, equ;//var变量数, equ方程个数

int Gauss()//Gauss消元求解方程组
{
    //int i, j, k, col, max_r;
    for(int k = 0, col = 0; k < equ && col < var; k++, col++)
    {
        int max_r = k;
        for(int i = k + 1; i < equ; i++)
            if(fabs(a[i][col]) > fabs(a[max_r][col]))
                max_r = i;//max_r表示第k到equ - 1行中第col列绝对值最大的
        if(fabs(a[max_r][col]) < eps) return 0;//col列都是0无解
        if(k != max_r)//col列(k到equ - 1行中) 系数绝对值最大的不是k
        {
            for(int j = col; j < var; j++)//交换第k行和第max_r行
                swap(a[k][j], a[max_r][j]);
            swap(x[k], x[max_r]);
        }
        x[k] /= a[k][col];//将第k行的第col个变量的系数变成1
        for(int j = col + 1; j < var; j++) a[k][j] /= a[k][col];
        a[k][col] = 1;
        for(int i = 0; i < equ; i++)//将第k行以外的所有的方程的第col个系数变成0
            if(i != k)
            {
                x[i] -= x[k]*a[i][k];//x[i] -= x[k]*a[i][k];
                for(int j = col + 1; j < var; j++) a[i][j] -= a[k][j]*a[i][col];
                a[i][col] = 0;
            }
    }//最终第i行i列的系数不是0,其他的都是0, 右边x[i]是对应的第i个变量的解
    return 1;
}

int N, M, X, Y, D, beg, end;
double p[110];
bool vis[210];
int hash[210];
int n;
double cons;

/*
 * hash[i]为点i的对应的在矩阵中的位置, 即第几个变量
 * 因为解方程时不能将不能到的点的加入方程组，否则无解
 */

bool bfs()
{
    queue <int> Q;
    Q.push(beg);
    memset(vis, 0, sizeof(vis));
    memset(hash, -1, sizeof(hash));
    memset(a, 0, sizeof(a));
    hash[beg] = 0;
    int cnt = 0;
    equ = 0;
    while(!Q.empty())
    {
        int now = Q.front();
        Q.pop();
        if(now == end || now == (n - end))
        {
            a[equ][hash[now]] = 1;
            x[equ] = 0;
            for(int i = 1; i <= M; i++)
                if(!vis[(now + i) % n] && fabs(p[i]) > eps)
                {
                    if(hash[(now + i) % n] == -1) hash[(now + i) % n] = ++cnt;
                    vis[(now + i) % n] = 1;
                    Q.push((now + i ) % n);
                }
            equ++;
            continue;
        }
        a[equ][hash[now]] = 1;
        for(int i = 1; i <= M; i++)
        {
            if(!vis[(now + i) % n] && fabs(p[i]) > eps)
            {
                vis[(now + i) % n] = 1;
                Q.push((now + i) % n);
            }
            if(hash[(now + i) % n] == -1 && fabs(p[i]) > eps) hash[(now + i) % n] = ++cnt;
            if(hash[(now + i) % n] != -1)
                a[equ][hash[(now + i) % n]] -= p[i];//考虑到(now + i) % n对于不同的i可能相等, 不能用=
        }
        x[equ] = cons;
        equ++;
    }
    if(!vis[end] && !vis[n - end])//根据对称出来的想法,终点可能有两个
        return false;
    return true;
}

int main()
{
    int t;
    scanf("%d", &t);
    while(t--)
    {
        scanf("%d %d %d %d %d", &N, &M, &Y, &X, &D);
        cons = 0;
        for(int i = 1;  i <= M; i++)
        {
            scanf("%lf", p + i);
            p[i] /= 100;
            cons += i*p[i];
        }
        if(X == Y)//特判一下
        {
            printf("0.00\n");
            continue;
        }
        n = 2*N - 2;
        end = Y;
        if(D == -1 || D == 0) beg = X;
        else beg = n - X;//虚拟的起点编号
        if(!bfs())//没有从期待你到终点的路
        {
            printf("Impossible !\n");
            continue;
        }
        var = equ;
        Gauss();
        printf("%.2f\n", x[hash[beg]]);
    }
    return 0;
}