Random Ranking

· · 题解

:::info[Statement]{open}

有 n 个实数 a_1,a_2,\dots,a_n，对于每个 1 \le i \le n，已知 a_i \in [l_i,r_i]，其中 l_i,r_i 为整数。对于每个 1 \le i \le n，1 \le j \le n，求 a_i 在所有数中为第 j 小的概率。

原题 n \le 80，但可加强到 n \le 120。

:::

先将区间离散化，拆成若干个小区间 [c_i,c_{i+1})。对于每个小区间，若钦定有 x 个数落在该区间内，则这 x 个数得到任意一个排名的概率都是 \frac{1}{x}。进一步地，若钦定有 y 个数落在该区间左侧，则这 x 个数的排名在 [y+1,y+x] 范围内。

据此，可定义 f_{i,j,u,v} 为：考虑区间 [c_i,c_{i+1}) 及 a_1 \sim a_j，该区间左侧有 u 个数，内部有 v 个数的概率，则：

f_{i,j,u,v}=f_{i,j-1,u,v} \cdot P(a_j \ge c_{i+1})+f_{i,j-1,u-1,v} \cdot P(a_j<c_i)+f_{i,j-1,u,v-1} \cdot P(c_i \le a_j<c_{i+1})

需对每个 a_i 跑一遍上面的 dp，并跳过 j=i。统计答案时，枚举 u,v，将贡献加到 [u+1,u+v+1] 内的每个排名上。复杂度 O(n^5)，足以通过本题。

可用退背包进一步优化。发现每一次只是撤销了 j=i 这一步的转移，且这个转移是具有交换律的，故可先对全部的 n 个数跑一遍 dp，然后倒推一步即可。

具体地，记 P_R=P(a_j \ge c_{i+1})，P_L=P(a_j<c_i)，P_M=P(c_i \le a_j<c_{i+1})，则有

f'_{u,v}=f_{u,v} \cdot P_R+f_{u-1,v} \cdot P_L+f_{u,v-1} \cdot P_M

现在要倒推，那么移一下项，变为

f_{u,v}=\frac{f'_{u,v}-f_{u-1,v} \cdot P_L-f_{u,v-1} \cdot P_M}{P_R}

据此转移即可。统计答案是区间加的形式，可用差分优化，复杂度 O(n^4)。

实现上有个小问题，就是上式中 P_R 可能 =0。转移肯定是要确保 [l_i,r_i] 覆盖了 [c_j,c_{j+1}) 这个区间的，也就是说 P_M 一定 \ne 0，那我们可以把 f_{u,v-1} 扔到左边去，这样分母就是 P_M，就没有问题了。式子是相似的，这里不再列出。

模拟赛把这题改成了取模，我懒得改回浮点数了，将就着看吧。

:::success[O(n^5) Code]

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int MAXN = 150;
const int MOD = 1e9 + 9;
int l[MAXN], r[MAXN], c[MAXN], f[MAXN][MAXN], g[MAXN][MAXN], ans[MAXN][MAXN], inv[MAXN], invi[MAXN], n, tot;
int qpow(int x, int y){
    int res = 1;
    while (y){
        if (y & 1){
            res = res * x % MOD;
        }
        x = x * x % MOD;
        y >>= 1;
    }
    return res;
}
int p(int i, int x){
    if (l[i] >= x){
        return 1;
    }
    if (r[i] <= x){
        return 0;
    }
    return (r[i] - x) * inv[i] % MOD;
}
signed main(){
//  freopen("q.in", "r", stdin);
//  freopen("q.out", "w", stdout);
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    cin >> n;
    for (int i = 1; i <= n; i++){
        cin >> l[i] >> r[i];
        c[++tot] = l[i];
        c[++tot] = r[i];
        inv[i] = qpow(r[i] - l[i], MOD - 2);
    }
    for (int i = 1; i <= n; i++){
        invi[i] = qpow(i, MOD - 2);
    }
    sort(c + 1, c + tot + 1);
    tot = unique(c + 1, c + tot + 1) - c - 1;
    for (int i = 1; i <= n; i++){
        int li = lower_bound(c + 1, c + tot + 1, l[i]) - c;
        int ri = lower_bound(c + 1, c + tot + 1, r[i]) - c;
        for (int j = li; j < ri; j++){
            int lj = c[j], rj = c[j + 1];
            for (int u = 0; u <= n; u++){
                for (int v = 0; v <= n; v++){
                    f[u][v] = 0;
                }
            }
            f[0][0] = 1;
            for (int k = 1; k <= n; k++){
                if (i == k){
                    continue;
                }
                for (int u = 0; u <= n; u++){
                    for (int v = 0; v <= n; v++){
                        g[u][v] = f[u][v];
                        f[u][v] = 0;
                    }
                }
                for (int u = 0; u <= k; u++){
                    for (int v = 0; v <= k; v++){
                        f[u][v] = (f[u][v] + g[u][v] * p(k, rj) % MOD) % MOD;
                        if (u >= 1){
                            f[u][v] = (f[u][v] + g[u - 1][v] * (1 - p(k, lj) + MOD) % MOD) % MOD;
                        }
                        if (v >= 1){
                            f[u][v] = (f[u][v] + g[u][v - 1] * (p(k, lj) - p(k, rj) + MOD) % MOD) % MOD;
                        }
                    }
                }
            }
            for (int u = 0; u < n; u++){
                for (int v = 0; v < n; v++){
                    for (int k = u + 1; k <= min(u + v + 1, n); k++){
                        ans[i][k] = (ans[i][k] + (rj - lj) * inv[i] % MOD * f[u][v] % MOD * invi[v + 1] % MOD) % MOD;
                    }
                }
            }
        }
    }
    for (int i = 1; i <= n; i++){
        for (int j = 1; j <= n; j++){
            cout << ans[i][j] << " ";
        }
        cout << "\n";
    }
    return 0;
}

:::

:::success[O(n^4) Code]{open}

#include <bits/stdc++.h>
#define int long long
using namespace std;
const int MAXN = 400;
const int MOD = 1e9 + 9;
int l[MAXN], r[MAXN], c[MAXN], f[MAXN][MAXN], g[MAXN][MAXN], ans[MAXN][MAXN], inv[MAXN], invi[MAXN], n, tot;
int qpow(int x, int y){
    int res = 1;
    while (y){
        if (y & 1){
            res = res * x % MOD;
        }
        x = x * x % MOD;
        y >>= 1;
    }
    return res;
}
int p(int i, int x){
    if (l[i] >= x){
        return 1;
    }
    if (r[i] <= x){
        return 0;
    }
    return (r[i] - x) * inv[i] % MOD;
}
signed main(){
//  freopen("q.in", "r", stdin);
//  freopen("q.out", "w", stdout);
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    cin >> n;
    for (int i = 1; i <= n; i++){
        cin >> l[i] >> r[i];
        c[++tot] = l[i];
        c[++tot] = r[i];
        inv[i] = qpow(r[i] - l[i], MOD - 2);
    }
    for (int i = 1; i <= n; i++){
        invi[i] = qpow(i, MOD - 2);
    }
    sort(c + 1, c + tot + 1);
    tot = unique(c + 1, c + tot + 1) - c - 1;
    for (int i = 1; i < tot; i++){
        int lj = c[i], rj = c[i + 1];
        vector <int> act;
        int cntl = 0;
        for (int j = 1; j <= n; j++){
            if (r[j] <= lj){
                cntl++;
            }
            if (l[j] <= lj && r[j] >= rj){
                act.push_back(j);
            }
        }
        int len = (int)act.size();
        if (!len){
            continue;
        }
        for (int u = 0; u <= len; u++){
            for (int v = 0; v <= len; v++){
                f[u][v] = 0;
            }
        }
        f[0][0] = 1;
        for (int k = 0; k < len; k++){
            int j = act[k];
            for (int u = 0; u <= len; u++){
                for (int v = 0; v <= len; v++){
                    g[u][v] = f[u][v];
                    f[u][v] = 0;
                }
            }
            for (int u = 0; u <= k + 1; u++){
                for (int v = 0; v <= k + 1; v++){
                    f[u][v] = (f[u][v] + g[u][v] * p(j, rj) % MOD) % MOD;
                    if (u >= 1){
                        f[u][v] = (f[u][v] + g[u - 1][v] * ((1 - p(j, lj) + MOD) % MOD) % MOD) % MOD;
                    }
                    if (v >= 1){
                        f[u][v] = (f[u][v] + g[u][v - 1] * ((p(j, lj) - p(j, rj) + MOD) % MOD) % MOD) % MOD;
                    }
                }
            }
        }
        for (int j : act){
            for (int u = 0; u <= len; u++){
                for (int v = 0; v <= len; v++){
                    g[u][v] = 0;
                }
            }
            int pm = qpow((p(j, lj) - p(j, rj) + MOD) % MOD, MOD - 2);
            for (int u = 0; u < len; u++){
                for (int v = len - 1; v >= 0; v--){
                    g[u][v] = f[u][v + 1];
                    g[u][v] = (g[u][v] - g[u][v + 1] * p(j, rj) % MOD + MOD) % MOD;
                    if (u >= 1){
                        g[u][v] = (g[u][v] - g[u - 1][v + 1] * ((1 - p(j, lj) + MOD) % MOD) % MOD + MOD) % MOD;
                    }
                    g[u][v] = g[u][v] * pm % MOD;
                }
            }
            for (int u = 0; u < len; u++){
                for (int v = 0; v < len; v++){
                    int li = cntl + u + 1, ri = min(cntl + u + v + 1, n);
                    int val = (rj - lj) * inv[j] % MOD * g[u][v] % MOD * invi[v + 1] % MOD;
                    if (li <= ri){
                        ans[j][ri + 1] = (ans[j][ri + 1] - val + MOD) % MOD;
                        ans[j][li] = (ans[j][li] + val) % MOD;
                    }
                }
            }
        }
    }
    for (int i = 1; i <= n; i++){
        for (int j = 1; j <= n; j++){
            ans[i][j] = (ans[i][j] + ans[i][j - 1]) % MOD;
            cout << ans[i][j] << " ";
        }
        cout << "\n";
    }
    return 0;
}

:::