luoguP3600 随机数生成器概率与期望Dp

luoguP3600 随机数生成器

题目描述

sol研发了一个神奇的随机数系统，可以自动按照环境噪音生成真·随机数。
现在sol打算生成n个[1,x]的整数a_1…a_n，然后进行一些询问。
q次询问，每次询问i有两个参数li和ri，sol会计算min⁡li<j<riaj\min_{li<j<ri}ajminli<j<riaj
（a数组中下标在li、ri之间的数的最小值）。

最后测试结果会是这些询问得到的结果的最大值。

sol进行了很多次实验，现在他想问问你测试结果的期望大小是多少，对666623333取模。

分析

神仙期望Dp
运用整数概率公式E(x)=∑iinf⁡P(i≤x)E(x)=\sum_i^{\inf}P(i\leq x)E(x)=∑iinfP(i≤x)
可以把题目转化为求maxi{min⁡j=liriaj}≥smax_i\{\min_{j=li}^{ri}a_j\} \ge smaxi{minj=liriaj}≥s的概率
技巧：任意性问题一般比较好入手。
于是转化成
maxi{min⁡j=liriaj}≤s−1 ⟺ ∀i.∃j∈[li,ri]aj≤s−1max_i\{\min_{j=li}^{ri}a_j\} \le s-1\iff \forall i. \exists j\in[l_i,r_i]a_j\le s-1maxi{minj=liriaj}≤s−1⟺∀i.∃j∈[li,ri]aj≤s−1
如此看上去就变成了一个比较可做的问题。
我们发现互相包含的区间可以直接忽略大的那个。
于是询问区间左右端点递增
用f[i][j]f[i][j]f[i][j]表示前iii个数满足了j个区间的概念。转移随便转转可以有O(n3)O(n^3)O(n3)
考虑这个问题的对偶问题，就是每个点可以满足多少个区间。
这个时候每个点满足的区间的编号一定也是一个连续的区间，而且也满足左右端点单调。
点和区间互化之后。问题转化为若干个区间，每个区间选择的概率是p(p=s−1x)p(p=\frac{s-1}{x})p(p=xs−1)，求覆盖所有点的概率。
方程搞一搞，用f[i]f[i]f[i]表示iii是最后一个选的区间，覆盖前面所有点的概率
f[i]=p(∑l[j]≥r[i]−1f[j]⋅(1−p)i−j−1+[l[i]=1](1−p)i−1)f[i]=p(\sum_{l[j]\ge r[i]-1} f[j]\cdot (1-p)^{i-j-1}+[l[i]=1](1-p)^{i-1})f[i]=p(∑l[j]≥r[i]−1f[j]⋅(1−p)i−j−1+[l[i]=1](1−p)i−1)
答案是∑r[i]=nf[i]⋅(1−p)n−i\sum_{r[i]=n} f[i]\cdot (1-p)^{n-i}∑r[i]=nf[i]⋅(1−p)n−i
优化一下
f[i]=p((1−p)i−1∑l[j]≥r[i]−1f[j]⋅(1−p)−j+[l[i]=1](1−p)i−1)f[i]=p((1-p)^{i-1}\sum_{l[j]\ge r[i]-1} f[j]\cdot (1-p)^{-j}+[l[i]=1](1-p)^{i-1})f[i]=p((1−p)i−1∑l[j]≥r[i]−1f[j]⋅(1−p)−j+[l[i]=1](1−p)i−1)
双指针扫一扫就行了。
前半段套路，后半段区间-点转化模型要记得

代码

#include<cstdio>
#include<algorithm>
int ri() {char c = getchar(); int x = 0, f = 1; for(;c < '0' || c > '9'; c = getchar()) if(c == '-') f  = -1;for(;c >= '0' && c <= '9'; c = getchar()) x = (x << 1) + (x << 3) - '0' + c; return x * f;
}
const int P = 666623333, N = 5000;
struct Sec {int l, r;}x[N], y[N];
bool cmp(Sec a, Sec b) {return a.l == b.l ? a.r < b.r : a.l < b.l;}
int f[N], pw[N], ia, n, tp, m;
int Pow(int x, int k) {int r = 1;for(;k; x = 1LL * x * x % P, k >>= 1) if(k & 1) r = 1LL * r * x % P;return r;
}
int Cal(int x) {if(!x) return 0;int p = 1LL * x * ia % P, g = (1 + P - p) % P;pw[tp] = 1; int iv = Pow(g, P - 2), A = 0;for(int i = 1;i <= tp; ++i) pw[tp + i] = 1LL * pw[tp + i - 1] * g % P, pw[tp - i] = 1LL * pw[tp - i + 1] * iv % P;int L = 1, R = 0, s = 0;for(int i = 1;i <= tp; ++i) {for(;L <= R && y[L].r < y[i].l - 1; ++L) s = (s - 1LL * f[L] * pw[tp - L] % P + P) % P;f[i] = 1LL * s * pw[tp + i - 1] % P;if(y[i].l == 1) f[i] = (f[i] + pw[tp + i - 1]) % P;f[i] = 1LL * f[i] * p % P; s = (s + 1LL * f[i] * pw[tp - i]) % P; ++R;if(y[i].r == m) A = (A + 1LL * f[i] * pw[tp + tp - i]) % P;}return A;
}
int main() {n = ri(); int a = ri(); ia = Pow(a, P - 2); int q = ri(); for(int i = 1;i <= q; ++i) x[i].l = ri(), x[i].r = ri();std::sort(x + 1, x + q + 1, cmp);m = 1;for(int i = 2;i <= q; ++i) {for(;x[m].r >= x[i].r;) --m;if(x[m].l < x[i].l) x[++m] = x[i];}int L = 1, R = 0;for(int i = x[1].l;i <= x[m].r; ++i) {for(;L <= R && x[L].r < i; ++L) ;for(;x[R + 1].l <= i && R < m; ++R) ;if(L <= R) y[++tp].l = L; y[tp].r = R;}int A = 0;for(int i = 1;i <= a; ++i) A = (A + 1 + P - Cal(i - 1)) % P;printf("%d\n", A);return 0;
}