题目描述

scx 在考试中取得了好成绩,因此赢得了去 Alushta 夏令营的机会,那是一个在海边的小岛,但是不幸地,正值 Alushta 的飓风季,因此帐篷有一定概率倒塌!

营地建筑可以看成一个高 $n + 2$,宽 $m$ 的砖块构成的矩形,每个白天,从海边会吹来一阵风,除了最上、最下那两行,其余每行最左边的砖块都有独立的 $p = \dfrac ab$ 的概率被吹走,类似地,每个晚上,会有另一阵风吹向海边,除了最上、最下那两行,其余每行最右边的砖块也各有独立的 $p$ 的概率吹走,记住,最上面的一行和最下面的一行的砖块特别坚固,因此只有 $n \cdot m$ 的砖块可能被吹走。

飓风季一共会持续 $k$ 天 $k$ 夜,如果在某一天,营地建筑被分成两个及以上连通分量,它就会被倒塌,scx 不得不寻找另一个地方度假。

求 scx 不需要寻找其它地方,即她能这个营地度假的概率。

输入格式

第一行包含两个正整数 $n, m$ ($n, m \leq 1500$),表示建筑物可能被倒塌部分的高度和宽度。

第二行包含两个正整数 $a, b$ ($a < b \leq 10^9$),表示每个砖块被吹走的概率为 $p = \dfrac ab$。

第三行包含一个正整数 $k$ ($k \leq 10^5$),表示飓风季持续的天数 (夜数)。

输出格式

输入一行一个整数,表示 scx 能在这个营地度假的概率在模 $10^9 + 7$ 意义下的取值。

题解

依旧考虑 DP。

因为所有的行是独立的,所以记 $pl_i$ (pl[i]) 表示某一行中,$k$ 个白天后最左边的砖块为第 $i$ 块的概率,类似地,记 $pr_i$ (pr[i]) 表示某一行中,$k$ 个晚上后最右边的砖块为第 $i$ 块的概率。

可以看出,如果将 $pl_i$ 中的 $i$ 看作一个随机变量,即 $k$ 个白天后吹掉了多少个砖块,则它应该服从二项分布,即,记 $q = 1 - p$ 是一个砖块不被吹走的概率,则有

$$ pl_{i+1} = \binom ki q^i p^{k-i} $$

由对称性,应有 $pr_{c-i} = pl_{i+1}$。再记 $f_{i, l, r}$ 表示 $k$ 天后,从下往上数第 $i$ 行仅剩下砖块 $l, l + 1, \cdots, r$ (且下面 $i$ 行不倒塌) 的概率。

可以看出,如果营地建筑最后没有倒塌 (只有一个连通分量),那么任意相邻的两行一定是有交集的。

因此,我们只要枚举有交集的部分转移即可,时间复杂度 $O(n^5)$,然而这太慢啦。

撕烤一下可以得到,有交集的部分不怎么好枚举,但是没交集的部分,可以看出,都是在两边,还是比较好枚举的。

因此,记 $F_{i, r}$ (F[i][r]) 表示 $k$ 天后,下数第 $i$ 行最右边的砖块为第 $r$ 块的概率,即

$$ F_{i, r} = \sum_{l=1}^r f_{i, l, r} $$

再记 $\sigma r_{i, j}$ (sr[i][j]) 表示前面条件相同,最右边的砖块为第 $j$ 块或更靠左的概率,即

$$ \sigma r_{i, j} = \sum_{r=1}^j F_{i, r} = \sum_{r=1}^j \sum_{l=1}^r f_{i, l, r} = \sum_{1 \leq l \leq r \leq j} f_{i, l, r} $$

同理,记 $\sigma l_{i, j}$ (sl[i][j]) 表示前面条件相同,最左边的砖块为第 $j$ 块或更靠右的概率,即

$$ \sigma l_{i, j} = \sum_{j \leq l \leq r \leq m} f_{i, l, r} $$

可以由对称性由 $\sigma r_{i, c-j+1}$ 的值推得。于是转移方程就是

$$ f_{i, l, r} = pl_l \cdot pr_r \cdot \sum_{[j, k] \cap [l, r] \neq \varnothing} f_{i-1, j, k} = pl_l \cdot pr_r \cdot \left( \sum_{1 \leq j \leq k \leq m} f_{i-1, j, k} - \sum_{1 \leq j \leq k < l} f_{i-1, j, k} - \sum_{r < j \leq k \leq m}f_{i-1, j, k} \right) = pl_l \cdot pr_r \cdot (\sigma r_{i-1, m} - \sigma r_{i-1, l-1} - \sigma l_{i-1, r+1}) $$

对其求前缀和,再和式变换,得

$$ F_{i, r} = \sum_{l=1}^r f_{i, l, r} = \sum_{l=1}^r \left( pl_l \cdot pr_r \cdot \left( \sigma r_{i-1, m} - \sigma r_{i-1, l-1} - \sigma l_{i-1, r+1} \right) \right) = pr_r \cdot \left( \left( \sum_{l=1}^r pl_l \right) \cdot (\sigma r_{i-1, m} - \sigma l_{i-1, r+1}) - \sum_{l=1}^r (pl_l \cdot \sigma r_{i-1, l-1}) \right) $$

于是记 $\sigma pl_i = \sum\limits_{k=1}^i pl_k$ (spl[i]),$\sigma pl \sigma r_i = \sum\limits_{k=1}^i pl_k \cdot \sigma r_{i-1, k-1}$ (splsr[i]),于是

$$ F_{i, r} = pr_r \cdot \left( \sigma pl_r \cdot (\sigma r_{i-1, m} - \sigma l_{i-1, r+1}) - \sigma pl \sigma r_r \right) $$

这样就能做到总时间复杂度 $O(nm)$,最终答案就是 $\sigma r_{n, m}$。

代码

#include <bits/stdc++.h>
#define N 2034
#define K 100034
using namespace std;

typedef long long ll;
const ll mod = 1000000007;

int r, c, k, i, j;
ll p, q; // p: not collapse, q: collapse
ll P[K], Q[K], C[N], inv[N]; // p^i, q^i, C(k, i), i^-1
ll pl[N], pr[N], spl[N]; // probability of a level to i
ll F[N][N], sl[N][N], sr[N][N], splsr[N]; // Sum[f[i][j]], Sum[pl[i] * sr[i - 1]]

ll PowerMod(ll a, ll n, ll m = mod){
    if(!n || a == 1) return 1ll;
    ll s = PowerMod(a, n >> 1, m);
    (s *= s) %= m;
    return n & 1 ? s * a % m : s;
}

ll Mul(int count, ...){
    va_list ptr;
    va_start(ptr, count);
    ll res = 1;
    for(int i = 0; i < count; ++i) (res *= va_arg(ptr, ll)) %= mod;
    va_end(ptr);
    return res;
}

int main(){
    scanf("%d%d%d%d%d", &r, &c, &i, &j, &k);
    q = (ll)i * PowerMod(j, mod - 2) % mod;
    (p = 1 - q) < 0 ? (p += mod) : p;
    inv[1] = C[0] = P[0] = Q[0] = 1; C[1] = k;
    for(i = 1; i <= k; ++i) {P[i] = P[i - 1] * p % mod; Q[i] = Q[i - 1] * q % mod;} // power of p(q)
    for(i = 2; i <= c; ++i){
        inv[i] = (mod - mod / i) * inv[mod % i] % mod;
        C[i] = Mul(3, C[i - 1], (ll)(k - i + 1), inv[i]); // binomial coefficients C(k, i)
    }
    for(i = 0; i <= k && i < c; ++i) pr[c - i] = pl[i + 1] = Mul(3, C[i], Q[i], P[k - i]); // binomial distribution
    for(i = 1; i <= c; ++i) (spl[i] = spl[i - 1] + pl[i]) >= mod ? spl[i] -= mod : 0; // Sum[pl[i]]
    sr[0][c] = 1;
    for(i = 1; i <= r; ++i){
        for(j = 1; j <= c; ++j)
            splsr[j] = (splsr[j - 1] + pl[j] * sr[i - 1][j - 1]) % mod; // Sum[pl[j] * sr[j - 1]]
        for(j = 1; j <= c; ++j){
            F[i][j] = pr[j] * ((sr[i - 1][c] - sl[i - 1][j + 1]) * spl[j] % mod - splsr[j]) % mod; // DP
            F[i][j] < 0 ? F[i][j] += mod : 0;
        }
        for(j = 1; j <= c; ++j){
            (sr[i][j] = sr[i][j - 1] + F[i][j]) >= mod ? sr[i][j] -= mod : 0; // Sum[F[i][j]]
            sl[i][c - j + 1] = sr[i][j];
        }
    }
    printf("%d\n", (int)sr[r][c]);
    return 0;
}

坑1:各个数组的大小不要开错,特别是计算 pl[i] 时,值大于 0 需要同时满足 $i \leq k \wedge i \leq c$。

坑2:在转移到的时候,记得及时取模,不要到时候爆 long long 了。