HDOJ-4089 Activation(概率DP)

 题意：

n个人在服务器上排队激活仙剑5游戏，番茄排在第m位。服务器给每个人激活时，可能会发生一下4种情况：

事件1：激活失败，队伍保持不变，概率为p1；

事件2：连接失败：队首的人回到队尾，概率为p2；

事件3：激活成功：队首的人出队，其他人前移，概率为p3；

事件4：服务器故障：事件发生后所有未激活的人不能激活，概率为p4；

规定事件A：番茄排到小于等于第K个位置时候，出现服务器故障的概率。

做法：

自己入门概率DP求概率的第一题，参考l其他博客的思路。菜鸟一只，如果下面说的不对，还请多多指教。

设f[i][j]为队伍一共有i个人番茄排在第j位时能发生事件A点概率，我们要求的结果就为f[n][m]，注意不是恰好是在第j位时遇到服务器故障的概率，而是在第j个位置排队时，能到达目标状态的概率，比如f[10][j],可能是番茄在第j位置是发生事件A，也可能在1、2、3、4、5、6......j-1位置，发生的事件A（可能我表述的不太清楚，请读者结合下面给出的方程再理解一下）

转移方程：

j=1时：f[i][1]=p1*f[i][1]+p2*f[i][i]+p4;

可以这么理解：番茄排在第一位时，如果发生事件1，那么他要继续排在队首等待；如果发生事件2，那么他要回到队尾（到队尾也是可能发生事件A的，队首激活成功或者链接失败，队伍会移动）；如果发生事件3，那么就不再可能发生事件A了，所以没有这一项；如果发生事件4，直接达成事件A；

2<=j<=k时：f[i][j]=p1*f[i][j]+p2*f[i][j-1]+p3*f[i-1][j-1]+p4;

同样：番茄排在这个区间时，队首如果发生事件1，则要继续等待；如果发生事件2，队伍总人数不变，番茄前进1名；发生事件3，队伍少1人，番茄也前进1名；发生事件4，直接达成事件A；

k<j<=i时：f[i][j]=p1*f[i][j]+p2*f[i][j-1]+p3*f[i-1][j-1]；

和第二种情况类似，只是番茄如果排在k位之后，队首不能发生服务器故障，否则A发生的概率就为0了，所以转移的时候没有p4.

这样的方程是没法转移的，首先要进行移项。

令k2=p2/(1-p1)    k3=p3/(1-p1)     k4=p4/(1-p1)

j=1:           f[i][1]=k2*f[i][i] + | k4

2<=j<=k:  f[i][j]=k2*f[i][j-1]+ | k3*f[i-1][j-1]+k4

k<j<=i:      f[i][j]=k2*f[i][j-1]+ | k3*f[i-1][j-1]

当求到f的第i行时，第i-1行是已知数，所以分割线|右边的可以提前预处理出来当做常数c[i]

得到i元1次方程组(简写，略去f的第一维坐标)：

f[1]=k2*f[i]+c[1]

f[2]=k2*f[1]+c[2]

f[3]=k2*f[2]+c[3]

......

f[n]=k2*f[n-1]+c[n]

把第一个式子代入到第二个，得到f[2]关于f[i]的表达式，再把这个表达式代入到第三个方程，得到f[3] 关于f[i]的表达式，这样一步步迭代，最后得到：

f[n]=(k2^n)f[n]+(k2^(n-1))c[1]+(k2^(n-2))c[2]+......+(k2^0)c[n]

这样就能求出f[n],于是就可以根据j=1时候的方程得出f[1],下面进行dp就好了,不要用解方程组的方法再去求前面各项（除法比较慢的，会超时）。

做之前要判断下特殊情况，比如服务器上来就会崩溃或者服务器永远不能激活完所有的人，避免分母出现0或者加快一下速度，减少不必要运算；

一开始RE了——在函数里开了2000×2000的数组，后来开到全局变量里，交上去之后——超内存了（32MB），而且是刚刚超，看到f[i]只用到了f[i-1]，干脆就用滚动数组了，也挺好改的，所有的i项都与1就好了，数组在0行与1行之间来回转移。

#include <cmath>
#include <cstdio>
#include <cstdlib>
#include <cstring>
using namespace std;
const int N = 2000 + 5;


void work(int n, int m, int k, double * p)
{
    int    i, j;
    double temp;
    double k2, k3, k4;//1/(1-pi)
    double f[2][N];   //DP数组,开成[N][N]超32MB内存，用滚动数组
    double c[N];      //常数数组
    double q[N];      //预处理系数的幂

    //特判1，永远不可能把这队人激活完，而且服务器也不会瘫痪
    if (fabs(1 - p[1] - p[2])<1e-9)
    {
        printf("0.00000\n");
        return;
    }

    //特判2，服务器一开始就瘫痪
    if (fabs(1 - p[4])<1e-9)
    {
        if (m <= k)
            printf("1.00000\n");
        else
            printf("0.00000\n");
        return;
    }

    //计算用到的系数
    p[1] = 1 - p[1];
    k2 = p[2] / p[1];
    k3 = p[3] / p[1];
    k4 = p[4] / p[1];

    //计算解方程时用到的k2的i次幂
    q[0] = 1;
    q[1] = k2;
    for (i = 2; i <= n; i++)
        q[i] = q[i - 1] * k2;


    memset(f, 0, sizeof(0));
    f[1][1] = p[4] / (p[3] + p[4]);
    c[1] = k4;

    for (i = 2; i <= n; i++)
    {
        //计算f[i]中涉及到的常数项
        for (j = 2; j <= k; j++)
            c[j] = k3*f[(i - 1)&1][j - 1] + k4;
        for (j = k + 1; j <= i; j++)
            c[j] = k3*f[(i - 1)&1][j - 1];

        //解方程求出f[i][i]
        temp = 0;
        for (j = 1; j <= i; j++)
            temp += c[j] * q[i - j];
        f[i&1][i] = temp / (1 - q[i]);

        //终于能愉快地DP了,常数数组再次派上用场
        f[i&1][1] = k2*f[i&1][i] + k4;
        for (j = 2; j < i; j++)
            f[i&1][j] = k2*f[i&1][j - 1] + c[j];

    }

    printf("%0.5lf\n", f[n&1][m]);

}

int main()
{
    int i, n, m, k;
    double p[10];
    while (~scanf("%d%d%d", &n, &m, &k))
    {
        for (i = 1; i <= 4; i++) scanf("%lf", &p[i]);
        work(n, m, k, p);
    }
    return 0;
}