Luogu P3600 随机数生成器（期望+dp）

题意

有一个长度为 (n) 的整数列 (a_1, a_2, cdots, a_n) ,每个元素在 ([1, x]) 中的整数中均匀随机生成。

• 有 (q) 个询问,第 (i) 个询问的结果是下标在 ([l_i , r_i ]) 的元素的最小值。
• 求这 (q) 个询问结果的最大值的期望，(mod 666623333) 。

数据范围

(1 le n, q, x le 2000) 。

题解

参考了 fjzzq2002 的题解。

有个很有用的结论。

对于非负实数变量 (X) 我们有：

[mathbb E(X) = int_{0}^{infty} P(X ge x) mathrm dx ]

对于 (X) 是离散的非负整数变量的情况，我们有：

[mathbb E(X) = sum_{x = 1}^{infty} P(X ge x) ]

如何理解呢？我们可以考虑期望的定义式 (mathbb E(X) = sum_{x} P(X = x) x) 其实是一个数面积的过程，也就是积分。
然后等式两边的其实就是数面积的两种过程，一个是竖着数，另外一个是横着数。

推广一下其实如果是求 (mathbb E(f(X))) 我们考虑求导后积分也就是：

[mathbb E(f(x)) = int_{0}^{infty} f'(x) P(X ge x) mathrm dx ]

我们套用这个后就有

[mathbb E(ans) = sum_{i = 1}^infty P(ans ge i) = sum_{i = 0}^{x - 1} (1 - P(ans le i)) ]

这样就好做多了，考虑枚举 (i) ，然后算 (P(ans le i)) 也就是，每个区间的 (min) 都 (le i) 的概率，假设这个 (i = t) 。

我们发现如果一个区间包含另外一个区间的话，那么它的最小值一定不会比小区间大，所以一定不会贡献到答案，那么我们可以去掉大区间。转化后，所有区间就是互不包含的了。

考虑记 (f_i) 为右端点不超过 (i) 的区间询问结果都 (le t) 的概率。

如果不存在一个 (r = i) ，那么就是 (f_i = f_{i - 1}) 。

否则记右端点为 (i) 的区间左端点为 (j) ，枚举区间内最后一个 (le t) 的数的位置 (l) ，那么我们有

[f_i = sum_{l = j}^{i} f_{l - 1} frac tx (1 - frac tx)^{i - l} ]

然后发现是 (mathcal O(xn^2)) 的，发现它能跑过QAQ

显然我们可以使它更加正确，推导一下式子就有（注意有 (t < x) ）

[f_i = frac tx (1 - frac tx)^i sum_{l = j}^{i} f_{l - 1} (1 - frac tx)^{- l} ]

然后我们就可以用前缀和优化到 (mathcal O(xn)) 。

其实还可以继续优化。。

观察转移方程发现答案其实是一个关于 (displaystyle frac tx) 的 (mathcal O(n)) 次多项式，也就是关于 (i) 的 (mathcal O(n)) 次多项式。

所以它的前缀和仍然是一个 (mathcal O(n)) 次多项式，我们算出前 (mathcal O(n)) 项和前缀和，插值出 (x) 处的前缀和即可。复杂度就优化成 (mathcal O(n^2)) 啦。

代码

#include <bits/stdc++.h>

#define For(i, l, r) for (register int i = (l), i##end = (int)(r); i <= i##end; ++i)
#define Fordown(i, r, l) for (register int i = (r), i##end = (int)(l); i >= i##end; --i)
#define Rep(i, r) for (register int i = (0), i##end = (int)(r); i < i##end; ++i)
#define Set(a, v) memset(a, v, sizeof(a))
#define Cpy(a, b) memcpy(a, b, sizeof(a))
#define debug(x) cout << #x << ": " << (x) << endl
#define l first
#define r second

using namespace std;

typedef pair<int, int> PII;

template<typename T> inline bool chkmin(T &a, T b) { return b < a ? a = b, 1 : 0; }
template<typename T> inline bool chkmax(T &a, T b) { return b > a ? a = b, 1 : 0; }

inline int read() {
    int x(0), sgn(1); char ch(getchar());
    for (; !isdigit(ch); ch = getchar()) if (ch == '-') sgn = -1;
    for (; isdigit(ch); ch = getchar()) x = (x * 10) + (ch ^ 48);
    return x * sgn;
}

void File() {
#ifdef zjp_shadow
    freopen ("P3600.in", "r", stdin);
    freopen ("P3600.out", "w", stdout);
#endif
}

const int N = 2e3 + 1e2, Mod = 666623333;

inline int fpm(int x, int power) {
    int res = 1;
    for (; power; power >>= 1, x = 1ll * x * x % Mod)
        if (power & 1) res = 1ll * res * x % Mod;
    return res;
}

int n, x, q, len, f[N]; PII S[N];

int V[N], powb[N], invpowb[N], sum[N];

int main () {

    File();

    n = read(); x = read(); q = read();

    For (i, 1, q)
        S[i].l = read(), S[i].r = read();
    sort(S + 1, S + q + 1);
    For (i, 1, q) {
        bool flag = true;
        For (j, i + 1, q)
            if (S[i].l <= S[j].l && S[j].r <= S[i].r) flag = false;
        if (flag)
            V[S[i].r] = S[i].l;
    }

    int ans = 0, invx = fpm(x, Mod - 2);
    Rep (v, x) {
        f[0] = powb[0] = invpowb[0] = 1;
        int coef = 1ll * v * invx % Mod;
        powb[1] = (Mod + 1 - coef) % Mod;
        invpowb[1] = fpm(powb[1], Mod - 2);
        For (i, 2, n) {
            powb[i] = 1ll * powb[i - 1] * powb[1] % Mod;
            invpowb[i] = 1ll * invpowb[i - 1] * invpowb[1] % Mod;
        }
        For (i, 1, n) {
            sum[i] = (sum[i - 1] + 1ll * f[i - 1] * invpowb[i]) % Mod;
            if (!V[i]) f[i] = f[i - 1];
            else f[i] = 1ll * (Mod + sum[i] - sum[V[i] - 1]) * coef % Mod * powb[i] % Mod;
        }
        (ans += Mod + 1 - f[n]) %= Mod;
    }
    printf ("%d
", ans);

    return 0;

}