洛谷P8967 追寻 | Pursuit of Dream 题解

题意：

在 \(n\) 维空间中有一个梦想。这梦想坐落在 \((d_1, d_2, \ldots, d_n)\) 的地方。而你从 \((0, 0, \ldots, 0)\) 开始，开启寻梦的旅程。

你的步伐轻缓，每一步只能走一个单位长度。你并不知道你的梦想位于哪里，所以你只能随机选择 \(n\) 个正方向中的一个，然后向这个方向走一步。也就是说，在 \([1, n]\) 中均匀随机选择一个正整数 \(h\)，然后，使你在第 \(h\) 维的坐标变成原来的坐标加一。

然而，天有不测风云。在你走每一步的过程中，你会有 \(p = \sum_{i = 1}^k p_i\) 的概率散入天际，并开始一段新的旅程。你会在 \(k\) 个地点中的一个重新开始这段旅程，其中第 \(i\) 个地点的坐标是 \((a_{i,1}, a_{i,2}, \ldots, a_{i,n})\)，从这里重新开始的概率为 \(p_i\)。

那么，期望下，你离到达这个梦想还需要多少步呢？

输入格式：

第一行，两个正整数 \(n,k\)。

第二行，\(n\) 个非负整数 \(d_1, d_2, \ldots, d_n\)。

接下来 \(k\) 行，第 \(i\) 行 \(n + 1\) 个整数 \(a_{i, 1}, a_{i, 2}, \ldots, a_{i, n}, x_i\)，每行最后一个整数 \(x_i\) 表示 \(p_i=x_i\times 10^{-8}\)。

输入的 \(x_i\) 保证了 \(p_i > 0\) 且 \(p < 1\)。

保证每个 \(x_i\) 在所有可能的组合中等概率随机生成。

输出格式：

一行，一个整数，表示答案对 \(998244353\) 取模的结果。

如果你不知道如何进行实数取模：可以说明答案一定是有限的，且是有理数，设它的最简分数形式为 \(\frac{p}{q}\)。如果存在一个整数 \(x\) 满足 \(x \cdot q \equiv p \pmod{998244353}\) 且 \(0 \le x < 998244353\)，那么你只需输出 \(x\) 的值即可。

由于保证了 \(x_i\) 是随机生成的，可以说明以接近 \(1\) 的概率答案在模意义下存在。事实上，一个当 \(x_i\) 尚不确定时以合理地高的概率给出正确答案的算法足以通过本题，考察复杂的模意义下的有理数的处理不是我们的本意。

数据范围：

对于 \(100 \%\) 的数据：

\(1 \le n \le 100\)，\(1 \le k \le 10000\)。

\(d_i \ge 0\)，\(\sum_i d_i \le 10^7\)。

\(0 \le a_{i, j} \le {10}^7\)。

\(x_i \ge 1\)，\(\sum_i x_i < {10}^8\)。此即保证了 \(p_i > 0\) 和 \(p < 1\)。

保证存在一个 \(i \in [1, k]\) 使得对于每个 \(j \in [1, n]\) 均有 \(a_{i,j} \le d_j\)。

保证每个 \((a_{i, 1}, a_{i, 2}, \ldots, a_{i, n})\) 作为空间中的点互不相同。

保证每个 \(x_i\) 在所有可能的组合中等概率随机生成。

由于保证了 \(x_i\) 是随机生成的，可以说明以接近 \(1\) 的概率答案在模意义下存在。事实上，一个当 \(x_i\) 尚不确定时以合理地高的概率给出正确答案的算法足以通过本题，考察复杂的模意义下的有理数的处理不是我们的本意。

样例中的 \(x_i\) 不是随机生成的，仅为理解题意所用。

考虑设 \(f_i\) 表示从第 \(i\) 个点走到终点的期望步数。

显然 \(f_i\) 由两部分构成：一是直接走到终点的期望，另一是散入天际后再走到终点的期望。

先考虑如何计算第一部分。记 \(q_i\) 表示从第 \(i\) 个点直接走到终点的概率，并记 \(P=\sum_{i=1}^{k} p_i\) ，则

若存在 \(j\) 使得 \(a_{i,j} > d_j\) ，则第 \(i\) 个点永远无法直接走到终点，即 \(q_i = 0\) 。
否则第 \(i\) 个点走到终点需要走 \(s_i = \sum_{j = 1}^{n}(d_j - a_{i,j})\) 。
- 每步不散入天际的概率是 \(1-P\) ，总概率为 \((1-P)^{s_i}\) 。
- 走 \(s_i\) 步有 \(n^{s_i}\) 种走法，其中只有 \(\frac{s_{i} !}{\prod_{j=1}^n\left(d_j-a_{i, j}\right) !}\) 种走法能到达终点
因此有 \[ q_i=(1-P)^{s_i} \times \frac{s_{i} !}{n^{s_i}\times\prod_{j=1}^n\left(d_j-a_{i, j}\right) !} \]

那么这部分的答案就是 \(q_i\cdot s_i\) 。

接着考虑较为复杂的第二部分。如果我们直接用 \(f_j\) 来表示 \(f_i\) ，则 \(f\) 就成了一个 \(k\) 元线性方程组

如果使用高斯消元，复杂度就到了 \(\mathcal{O}(k^3)\) ，因此我们考虑其他方法。

注意到无论从哪里散入天际，接着走到终点的期望都是一定的，不妨设为 \(r\) ，显然 \[ r = \sum_{i=1}^{k}\frac{p_i\cdot f_i}{P} \] 则我们只需考虑从第 \(i\) 个点走，期望散入天际的步数 \(g_i\) ，并和 \(r\) 加起来就是第二部分了。

不考虑走到终点停下，期望散入天际的步数有 \[ g_{i1}=\sum_{j=1}^{\infty} j \times(1-P)^{j-1} \times P=\frac{1}{P} \] 而刚才未考虑的这一部分的期望显然为「走到终点的期望」加上「走到终点后的期望」，有 \[ g_{i2}=q_i\times\left(s_i + \frac{1}{P}\right) \] 则 \[ g_i=g_{i1}-g_{i2} = \frac{1}{P}-q_i \times\left(s_i+\frac{1}{P}\right)=\left(1-q_i\right) \times \frac{1}{P}-q_i \cdot s_i \] 由此，我们得知 \(f_i\) 有 \[ f_i =q_i \cdot s_i+\left(1-q_i\right) \cdot r+g_i=\left(1-q_i\right) \times\left(r+\frac{1}{P}\right) \] 根据 \(r = \sum_{i=1}^{k}\frac{p_i\cdot f_i}{P}\) ，将 \(f_i\) 代入，得 \[ r=\sum_{i=1}^k \frac{p_i \cdot\left(1-q_i\right) \cdot\left(r+\frac{1}{P}\right)}{P} \] 整理一下 \[ r=\frac{\sum_{i=1}^k p_i \cdot\left(1-q_i\right)}{P^2-P \times \sum_{i=1}^k p_i \cdot\left(1-q_i\right)} \] 最后把 \(r\) 代回，得到 \[ f_0=(1-q_0)\times\left(r + \frac{1}{P}\right) \] 时间复杂度 \(\mathcal{O}(kn\log V + V)\)

代码：

#include <bits/stdc++.h>
using namespace std;
#define int long long
#define INF 0x3f3f3f3f3f3f3f3f
const int mod = 998244353;
void up(int &x,int y) { x < y ? x = y : 0; }
void down(int &x,int y) { x > y ? x = y : 0; }
#define N ((int)(115))
#define K ((int)(1e4 + 15))
#define M ((int)(1e7 + 15))

const int i1e8 = 205817851;
const int V = 1e7;

int qpow(int a,int b)
{
    int ans = 1;
    for(; b; b >>= 1) {
        if(b & 1) { ans = ans * a % mod; }
        a = a * a % mod;
    }
    return ans;
}
int n,k,q0,r,sump,d[N],a[K][N],p[K],flg[K],fac[M];
signed main()
{
    ios::sync_with_stdio(0);
    cin.tie(0);cout.tie(0);
    // freopen("check.in","r",stdin);
    // freopen("check.out","w",stdout);
    fac[0] = 1;
    for(int i = 1; i <= V; i++) fac[i] = fac[i - 1] * i % mod;
    cin >> n >> k;
    for(int i = 1; i <= n; i++) cin >> d[i];
    for(int i = 1; i <= k; i++)
    {
        flg[i] = 1;
        for(int j = 1; j <= n; j++)
        {
            cin >> a[i][j];
            flg[i] &= (d[j] >= a[i][j]);
        }
        cin >> p[i];
        p[i] = p[i] * i1e8 % mod; sump = (sump + p[i]) % mod;
    }
    for(int i = 1; i <= k; i++)
    {
        if(flg[i])
        {
            int q = 1, s = 0;
            for(int j = 1; j <= n; j++)
            {
                s += d[j] - a[i][j];
                q = q * qpow(fac[d[j] - a[i][j]], mod - 2) % mod;
            }
            q = q * qpow(((1 - sump) % mod + mod) % mod, s) % mod * fac[s] % mod * qpow(qpow(n, s), mod - 2) % mod;
            r = (r + p[i] * (((1 - q) % mod + mod) % mod) % mod) % mod;
        }else { r = (r + p[i]) % mod; }
    }
    r = (r == 0) ? 0 : qpow(((sump * sump % mod * qpow(r, mod - 2) % mod - sump) % mod + mod) % mod, mod - 2);
    q0 = 1; int s = 0;
    for(int i = 1; i <= n; i++) {
        s += d[i], q0 = q0 * qpow(fac[d[i]], mod - 2) % mod;
    }
    q0 = q0 * qpow(((1 - sump) % mod + mod) % mod, s) % mod * fac[s] % mod * qpow(qpow(n, s), mod - 2) % mod;
    cout << ((1 - q0) % mod + mod) % mod * (r + qpow(sump, mod - 2)) % mod << '\n';
    return 0;
}