[BJOI2018]治疗之雨
参考博客
https://www.luogu.org/blog/ShadowassIIXVIIIIV/solution-p4457
题目大意
(T) 组数据
场上有 (m) 个生命值无上下限的随从与 (1) 个生命值上限为 (n) 下限为 (0) 当前为 (p) 的英雄，现在每次随机选择一个生命值不为下限的人物，将其生命值 (+1) ，然后重复 (k) 次，每次将随机选择一个生命值不为下限的人物将其生命值 (-1) ，问期望多少次后生命值归 (0) ，如果不行则输出 (-1) ，答案对 (1000000007) 取模，保证分母不为模数的倍数，不存在 (n=p=k=1,m=0) 的情况
数据范围
(1 le T le 100, 1 le p le n le 1500, 0 le m,k le 1000000000)
时空限制
4000ms, 512MB
分析
我研究了一个 (n^2) 的迭代方法，但是它让我没有兴趣写出来......
我们设 (P(y)) 为在 (k) 次减 (1) 之后，将英雄血量从 (x) 变成 (x - y > 0) 的概率
[P(y) = (dfrac 1{1 + m})^y(dfrac m{1 + m})^{k-y} inom ky
]
设 (f(x)) 为 (p = x) 时的答案
[egin{cases}
f(0) = 0 \
f(n) = sum P(y) f(n - y) \
f(x) = dfrac 1{1 + m} P(0) f(x + 1) + sum left[ dfrac 1{1 + m}P(y + 1) + dfrac m{1 + m }P(y)
ight] f(x - y)
end{cases}
]
那么我们高斯消元即可。。。然而 (1500) 怎么 (n^3)
那么 (O(n^2)) 该怎么办呢（虽然我觉得 (n^2) 能过也很神奇）
我们迭代就可以了
观察一下高斯消元矩阵
[left[
egin{matrix}
1 & 1 & 0 & 0 & 0 & cdots & 0 & 0 \
1 & 1 & 1 & 0 & 0 & cdots & 0 & 0\
vdots & vdots & vdots & vdots & vdots & ddots & vdots & vdots \
1 & 1 & 1 & 1 & 1 & cdots & 1 & 1
end{matrix}
ight]
]
它很像一个三角矩阵不是吗......
所以我们消一消
[left[
egin{matrix}
1 & 1 & 0 & 0 & 0 & cdots & 0 & 0 \
0 & 1 & 1 & 0 & 0 & cdots & 0 & 0\
vdots & vdots & vdots & vdots & vdots & ddots & vdots & vdots \
0 & 0 & 0 & 0 & 0 & cdots & 0 & 1
end{matrix}
ight]
]
于是就是 (n^2) 了
Code
#include <cstdio>
#include <cstring>
#include <iostream>
using namespace std;
inline char nc() {
static char buf[100000], *l = buf, *r = buf;
return l==r&&(r=(l=buf)+fread(buf,1,100000,stdin),l==r)?EOF:*l++;
}
template<class T> void read(T &x) {
x = 0; int f = 1, ch = nc();
while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=nc();}
while(ch>='0'&&ch<='9'){x=x*10-'0'+ch;ch=nc();}
x *= f;
}
typedef long long ll;
const int mod = 1000000007;
const int maxn = 1500 + 5;
int T, n, p, m, k, p1, p2;
int P[maxn], e[maxn][maxn];
inline void sub(int &x, int y) {
x -= y; if(x < 0) x += mod;
}
inline void add(int &x, int y) {
x += y; if(x >= mod) x -= mod;
}
ll qpow(ll x, ll y) {
ll re = 1;
while(y) {
if(y & 1) re = re * x % mod;
x = x * x % mod, y >>= 1;
}
return re;
}
void init() {
memset(P, 0, sizeof(P));
p1 = 1 * qpow(1 + m, mod - 2);
p2 = m * qpow(1 + m, mod - 2) % mod;
int C = 1;
for(int i = 0; i <= min(n, k); ++i) {
P[i] = qpow(p1, i) * qpow(p2, k - i) % mod * C % mod;
C = (ll)C * (k - i) % mod * qpow(i + 1, mod - 2) % mod;
}
}
void Gauss() {
for(int i = 1; i < n; ++i) {
int r = qpow(e[i][i], mod - 2);
for(int j = i + 1; j <= n; ++j) {
int t = (ll)e[j][i] * r % mod; e[j][i] = 0;
sub(e[j][i + 1], (ll)e[i][i + 1] * t % mod);
sub(e[j][n + 1], (ll)e[i][n + 1] * t % mod);
}
}
e[n][n + 1] = (ll)e[n][n + 1] * qpow(e[n][n], mod - 2) % mod, e[n][n] = 1;
for(int i = n - 1; i >= 1; --i) {
sub(e[i][n + 1], (ll)e[i][i + 1] * e[i + 1][n + 1] % mod), e[i][i + 1] = 0;
e[i][n + 1] = (ll)e[i][n + 1] * qpow(e[i][i], mod - 2) % mod, e[i][i] = 1;
}
}
int main() {
// freopen("testdata.in", "r", stdin);
read(T);
for(int kase = 1; kase <= T; ++kase) {
read(n), read(p), read(m), read(k);
if(k == 0 || (m == 0 && k == 1)) {
puts("-1");
continue;
}
init();
memset(e, 0, sizeof(e));
for(int x = 1; x < n; ++x) {
for(int y = 0; y < x; ++y) {
e[x][x - y] = mod - ((ll)p1 * P[y + 1] % mod + (ll)p2 * P[y] % mod) % mod;
}
e[x][x + 1] = mod - (ll)p1 * P[0] % mod;
add(e[x][x], 1), e[x][n + 1] = 1;
}
for(int y = 1; y < n; ++y) {
e[n][n - y] = mod - P[y];
}
e[n][n] = ((1 - P[0]) + mod) % mod, e[n][n + 1] = 1;
Gauss(); printf("%d
", e[p][n + 1]);
}
return 0;
}
总结
特殊矩阵的高斯消元可以优化