luoguP4457 [BJOI2018]治疗之雨 概率期望 + 高斯消元

应该是最后一道紫色的概率了....然而颜色啥也代表不了....

首先看懂题意：

你现在有$p$点体力，你的体力上限为$n$

在一轮中，

1.如果你的体力没有满，你有$frac{1}{m + 1}$的几率回复一点体力

2.紧接着有$k$轮攻击，每轮攻击都有$frac{1}{m + 1}$的几率使你掉一点体力

如果一轮后，你的体力$ leq 0$，那么游戏结束

询问游戏结束的期望轮数

看懂题应该就懂了什么吧....

设状态$f[i]$表示生命值为$i$游戏结束的期望轮数

那么

$$f[i] = egin{cases}
& 0;;(i = 0)\
& 1 + sumlimits_{j = 1}^i f[j] * P[i - j];;(i = n)\
& 1 + sumlimits_{j = 1}^i f[j] * (frac{P[i - j + 1]}{m + 1} + frac{m * P[i - j]}{m + 1}) + f[i + 1]*frac{P[0]}{m + 1} ;;(else)
end{cases}$$

其中，$P[i]$表示在一轮攻击中受到$i$点攻击的概率

（$f[0]$结束游戏，因此不能向别的状态转移）

考虑怎么求$P[i]$

由于确定了$i$种要攻击，其余的不攻击，那么一种攻击方式的概率为$(frac{1}{m + 1})^i * (frac{m}{m + 1})^{k - i}$

对于攻击$i$次而言，总共有$C(k, i)$种攻击方式

因此$P[i] = C(k, i) *(frac{1}{m + 1})^i * (frac{m}{m + 1})^{k - i} $

直接$O(n^2)$暴力全部求出来...

那么，有了$P[i]$后，高斯消元的复杂度过高

但是，注意到方程中$0$的数量非常的多，因此在消元的时候把$0$项跳过

酱紫，复杂度就到$O($玄学$)$了，仔细优化一下就$O(n^2)$了，具体看代码吧....

注：记得判无解

注2：不知道为什么$P[]$数组莫名的要多预处理一些.....

注（注2）：仿佛是$n = 0$的时候挂了.....出题人有猫病啊.....

复杂度$O(Tn^2)$，没怎么卡常

#include <cstdio>
using namespace std;

extern inline char gc() {
    static char RR[23456], *S = RR + 23333, *T = RR + 23333;
    if(S == T) fread(RR, 1, 23333, stdin), S = RR;
    return *S ++;
}
inline int read() {
    int p = 0, w = 1; char c = gc();
    while(c > '9' || c < '0') { if(c == '-') w = -1; c = gc(); }
    while(c >= '0' && c <= '9') p = p * 10 + c - '0', c = gc();
    return p * w;
}

#define sid 1600
#define eid 200050
#define ri register int
const int mod = 1000000007;

int n, p, m, k;
int inv[eid], P[sid], f[sid][sid];

void Init_Inv() {
    inv[0] = inv[1] = 1;
    for(ri i = 2; i <= 200000; i ++)
    inv[i] = 1ll * (mod - mod / i) * inv[mod % i] % mod;
    for(ri i = 2; i <= 200000; i ++)
    inv[i] = 1ll * inv[i - 1] * inv[i] % mod;
}

int fp(int a, int k) {
    int ret = 1;
    for( ; k; k >>= 1, a = 1ll * a * a % mod)
    if(k & 1) ret = 1ll * ret * a % mod;
    return ret;
}

void Init_P() {
    int invm = fp(fp(m + 1, mod - 2), k);
    for(ri i = 0; i <= n + 5; i ++) {
        if(i > k) { P[i] = 0; continue; }
        int C = 1;
        for(ri j = k - i + 1; j <= k; j ++) C = 1ll * C * j % mod;
        C = 1ll * C * inv[i] % mod;
        P[i] = 1ll * C * fp(m, k - i) % mod * invm % mod;
    }
}

void Init_Guass() {
    int m11 = fp(m + 1, mod - 2);
    int mm1 = 1ll * m * m11 % mod;
    for(ri i = 1; i <= n + 1; i ++)
    for(ri j = 1; j <= n + 1; j ++)
    f[i][j] = 0;
    for(ri i = 1; i < n; i ++) {
        f[i][i] = 1; f[i][n + 1] = 1;
        f[i][i + 1] = (mod - 1ll * P[0] * m11 % mod);
        for(ri j = 1; j <= i; j ++) {
            f[i][j] = (f[i][j] - 1ll * m11 * P[i - j + 1] % mod + mod) % mod;
            f[i][j] = (f[i][j] - 1ll * mm1 * P[i - j] % mod + mod) % mod;
        }
    }
    f[n][n] = 1; f[n][n + 1] = 1;
    for(ri i = 1; i <= n; i ++) f[n][i] = (f[n][i] - P[n - i] + mod) % mod;
}

void Guass() {
    for(ri i = 1; i <= n; i ++) {
        int inv = fp(f[i][i], mod - 2);
        for(ri j = i + 1; j <= n; j ++) {
            int t = 1ll * f[j][i] * inv % mod;
            for(ri k = i; k <= i + 1; k ++) 
            f[j][k] = (f[j][k] - 1ll * f[i][k] * t % mod + mod) % mod;
            f[j][n + 1] = (f[j][n + 1] - 1ll * f[i][n + 1] * t % mod + mod) % mod;
        }
    }
    for(ri i = n; i >= 1; i --) {
        f[i][n + 1] = 1ll * f[i][n + 1] * fp(f[i][i], mod - 2) % mod;
        f[i - 1][n + 1] = (f[i - 1][n + 1] - 1ll * f[i - 1][i] * f[i][n + 1] % mod + mod) % mod; 
    }
}

int main() {
    int Tt = read();
    Init_Inv();
    while(Tt --) {
        n = read(); p = read(); m = read(); k = read();
        if(k == 0) { printf("-1
"); continue; }
        if(m == 0 && k == 1) { printf("-1
"); continue; }
        Init_P(); Init_Guass(); Guass();
　　　　　int ans = f[p][n + 1];
        printf("%d
", ans);
    }
    return 0;
}