设第 \(i\) 步后第 \(j\) 维位移量的值域为 \([l_{i, j}, r_{i, j}]\)。
每天死亡的点应该有 \((r_{i, j} - l_{i, j})\) 个，因为 \([1, -l_{i, j}]\) 和 \([n - r_{i, j} + 1, n]\) 中的节点死了，故：
\[Ans = \sum_{i = 0}^T\prod_{j = 1}^k w_j - (r_{i, j} - l_{i, j})
\]
其中 \(\sup(T) = \max\limits_{i = 1}^k\{w_i\}\)。
时间复杂度 \(\mathcal O(nkT)\)。
称 \(n\) 步为一轮，则从第二轮开始，每一轮的答案开始呈现出规律性的变化。
具体地，记 \(\Delta p_{i, j}\) 表示每个周期中走了 \(i\) 步后第 \(j\) 维的死亡点数的变化值，则第 \(x(x \ge 2)\) 轮的答案 \(f(x)\) 为：
\[f(x) = \sum_{i = 1}^n\prod_{j = 1}^kw_j - (r_{i, j} - l_{i, j}) - [(x - 2)\Delta p_{n, j} + \Delta p_{i, j}]
\]
做了前缀和后的式子是一个关于 \(x\) 的 \(k + 1\) 次多项式，拉格朗日插值维护即可，时间复杂度 \(\mathcal O(nk^2)\)。
总时间复杂度为 \(\mathcal O(nk^2)\)，注意步数并不整，最后有余量要单独计算。
代码：
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
constexpr int N = 5e5 + 10, K = 20, MOD = 1e9 + 7;
int n, k, w[N], c[N], d[N], now[K], l[N][K], r[N][K], tl[N][K], tr[N][K], p[N][K];
ll y[K];
ll inv(ll base, int e = MOD - 2) {
ll res = 1;
while (e) {
if (e & 1) res = res * base % MOD;
base = base * base % MOD;
e >>= 1;
}
return res;
}
ll lag(int n, int k) {
if (k <= n) return y[k];
ll res = 0;
for (int i = 2; i <= n; i++) {
ll a = y[i], b = 1;
for (int j = 2; j <= n; j++) {
if (j != i) a = a * (k - j) % MOD, b = b * (i - j) % MOD;
}
res = (res + a * inv(b)) % MOD;
}
return res;
}
int main() {
ios_base::sync_with_stdio(0); cin.tie(nullptr), cout.tie(nullptr);
cin >> n >> k; y[0] = 1;
for (int i = 1; i <= k; i++) cin >> w[i], y[0] = y[0] * w[i] % MOD;
for (int i = 1; i <= n; i++) {
memcpy(l[i], l[i - 1], sizeof(l[i])), memcpy(r[i], r[i - 1], sizeof(r[i]));
cin >> c[i] >> d[i]; now[c[i]] += d[i];
l[i][c[i]] = min(l[i][c[i]], now[c[i]]), r[i][c[i]] = max(r[i][c[i]], now[c[i]]);
}
// -1
int cnt = 0;
for (int i = 1; i <= k; i++) cnt += (r[n][i] - l[n][i] < w[i] && !now[i]);
if (cnt == k) {
cout << -1; return 0;
}
// Round 1
for (int i = 1; i <= n; i++) {
ll prod = 1;
for (int j = 1; j <= k; j++) prod = prod * max(w[j] - (r[i][j] - l[i][j]), 0) % MOD;
y[1] += prod;
}
y[1] %= MOD;
// Calc p
memcpy(tl[1], l[n], sizeof(tl[1])), memcpy(tr[1], r[n], sizeof(tr[1]));
for (int i = 1; i <= n; i++) {
now[c[i]] += d[i];
tl[i][c[i]] = min(tl[i][c[i]], now[c[i]]), tr[i][c[i]] = max(tr[i][c[i]], now[c[i]]);
for (int j = 1; j <= k; j++) p[i][j] = (tr[i][j] - tl[i][j]) - (r[n][j] - l[n][j]);
memcpy(tl[i + 1], tl[i], sizeof(tl[i + 1])), memcpy(tr[i + 1], tr[i], sizeof(tr[i + 1]));
}
// Calc maxr
int maxr = 1e9;
for (int i = 1; i <= k; i++) {
if (p[n][i]) maxr = min(maxr, (w[i] - (r[n][i] - l[n][i])) / p[n][i] + 1);
}
// Calc Round x in [2, k + 3]
for (int x = 2; x <= k + 3; x++) {
for (int i = 1; i <= n; i++) {
ll prod = 1;
for (int j = 1; j <= k; j++) prod = prod * (w[j] - (r[n][j] - l[n][j]) - ((x - 2) * p[n][j] + p[i][j])) % MOD;
y[x] += prod;
}
y[x] %= MOD;
}
for (int i = 1; i <= k + 3; i++) y[i] = (y[i] + y[i - 1]) % MOD;
ll ans = lag(k + 3, maxr ? maxr : 1);
// Clac the rest
for (int i = 1; i <= n; i++) {
ll prod = 1;
for (int j = 1; j <= k; j++) prod = prod * max((w[j] - (r[n][j] - l[n][j]) - ((maxr - 1) * p[n][j] + p[i][j])), 0) % MOD;
ans += prod;
}
cout << ans % MOD;
return 0;
}