F - Social Distance¶
题意¶
長さ N+1 の整数列 X_0,X_1,\ldots,X_N が与えられます。ここで、0=X_0<X_1<\ldots<X_N です。
これから、1 から N までの番号のついた N 人の人が、数直線上に現れます。人 i は、区間 [X_{i−1},X_i] の中から一様ランダムに選ばれた実数座標に出現します。
人と人の距離の最小値の期待値を \bmod 998244353 で求めてください。
解析¶
考虑求出人与人距离的最小值的互补累积分布函数 f(z) (即 f(z) 表示人与人的距离的最小值大于等于 z 的概率),这样我们就可以得到答案为 \displaystyle\int_{0}^{+\infty}f(z)\mathrm{d}z。
先尝试对固定的 z 求出 f(z)。设第 i 个人的位置为 x_i,则对于固定的 z,我们需要对 i=1\ldots n-1 满足 x_i+z\le x_{i+1}。设 y_i=x_i-iz,则我们的要求可以改写为 y_1\le y_2\le\ldots\le y_n。于是我们只需求出在 n 个区间均匀随机分布的变量 y_i\in[X_{i-1}-iz,X_i-iz] 满足 y_1\le y_2\le\ldots\le y_n 的概率。
这可以用动态规划 (严格来讲,这不是动态规划,而是递推) 在 O(n^3) 的时间复杂度内求出。将所有区间 [X_{i-1}-iz,X_i-iz] 的左右端点排序,从小到大的第 i 个端点记做 v_i。设 dp_{i,j} 表示 [v_1,v_i] 这段区间内分布了 y_1\le y_2\le\ldots\le y_j 的概率,我们如下转移:
- 若 \forall p\in[j+1,k],[v_{i-1},v_i]\in[X_{p-1}-pz,X_p-pz],使 dp_{i,k}\gets dp_{i,k}+dp_{i-1,j}\dfrac{1}{(k-j+1)!}\prod\limits_{j<p\le k}\dfrac{v_i-v_{i-1}}{X_p-X_{p-1}}
现在考虑变化的 z。由于每个 v_i 都可视作关于 z 的线性函数,在所有区间端点的大小关系不变时,观察递推形式,发现 f(z) 是一个关于 z 的不超过 n 次的多项式函数。我们求出区间内 n+1 个不同的点插值得到 f(z) 在该区间的解析式后,再用微积分基本定理求 f(z) 在这段区间的积分即可。
由于所有区间的端点的大小关系只有 O(n^2) 种,我们可以在 O(n^6) 时间解决本题。
实现¶
#include <bits/stdc++.h>
typedef long long int ll;
const int maxn = 23, mod = 998244353;
ll gcd(ll a, ll b){
return b ? gcd(b, a % b) : a;
}
int qpow(int a, int b){
int res = 1;
while(b){
if(b & 1) res = (ll)res * a % mod;
a = (ll)a * a % mod, b >>= 1;
}
return res;
}
struct frac{
ll x, y;
frac(ll a = 0){x = a, y = 1; }
frac(ll a, ll b){
ll g = gcd(a, b);
a /= g, b /= g;
if(b < 0) x = -a, y = -b;
else x = a, y = b;
}
void out(void){
std::printf("%lld/%lld", x, y);
}
int val(void){
return (x % mod) * qpow(y % mod, mod - 2) % mod;
}
};
frac operator-(const int &a, const frac &b){
return frac(a * b.y - b.x, b.y);
}
frac operator+(const frac &a, const frac &b){
return frac(a.x * b.y + b.x * a.y, a.y * b.y);
}
frac operator*(const frac &a, const int &b){
return frac(a.x * b, a.y);
}
frac operator/(const frac &a, const int &b){
return frac(a.x, a.y * b);
}
bool operator==(const frac &a, const frac &b){
return a.x == b.x && a.y == b.y;
}
bool operator!=(const frac &a, const frac &b){
return a.x != b.x || a.y != b.y;
}
bool operator<(const frac &a, const frac &b){
return a.x * b.y < b.x * a.y;
}
int eval(const std::vector<int> &f, int z){
int res = 0, x = 1;
for(int i = 0; i < (int)f.size(); ++i)
res = (res + (ll)f[i] * x) % mod, x = (ll)x * z % mod;
return res;
}
int inv[maxn];
std::vector<int> integral(const std::vector<int> &f){
std::vector<int> g(f.size() + 1);
for(int i = f.size(); i; --i)
g[i] = (ll)f[i - 1] * inv[i] % mod;
return g;
}
std::vector<int> operator+(const std::vector<int> &f, const std::vector<int> &g){
std::vector<int> h(f);
if(g.size() > h.size()) h.resize(g.size());
for(int i = 0; i < (int)g.size(); ++i)
h[i] = (h[i] + g[i]) % mod;
return h;
}
std::vector<int> operator*(const std::vector<int> &f, const std::vector<int> &g){
std::vector<int> h(f.size() + g.size() - 1);
for(int i = 0; i < (int)f.size(); ++i)
for(int j = 0; j < (int)g.size(); ++j)
h[i + j] = (h[i + j] + (ll)f[i] * g[j]) % mod;
return h;
}
std::vector<int> operator*(const std::vector<int> &f, const int &b){
std::vector<int> h(f);
for(int i = 0; i < (int)h.size(); ++i)
h[i] = (ll)h[i] * b % mod;
return h;
}
int n, X[maxn], l[maxn], r[maxn];
std::pair<int, bool> v[maxn << 1];
void init(const frac &z){
for(int i = 1; i <= n; ++i){
v[i * 2 - 1] = std::make_pair(i, false);
v[i * 2] = std::make_pair(i, true);
}
std::sort(v + 1, v + 1 + n * 2, [&z](const std::pair<int, bool> &a, const std::pair<int, bool> &b){
return X[a.first + a.second - 1] - z * a.first != X[b.first + b.second - 1] - z * b.first ?
X[a.first + a.second - 1] - z * a.first < X[b.first + b.second - 1] - z * b.first :
a.first > b.first;
});
for(int i = 1; i <= n * 2; ++i)
if(v[i].second) r[v[i].first] = i;
else l[v[i].first] = i;
}
frac next(void){
frac res = 1e6;
for(int i = 1; i < n * 2; ++i)
if(v[i].first < v[i + 1].first && v[i].first + v[i].second != v[i + 1].first + v[i + 1].second)
res = std::min(res, frac(X[v[i + 1].first + v[i + 1].second - 1] - X[v[i].first + v[i].second - 1], v[i + 1].first - v[i].first));
return res;
}
bool check(void){
int last = 0;
for(int i = 1; i <= n; ++i)
if(last < r[i]) last = std::max(last, l[i]);
else return false;
return true;
}
int leninv[maxn], dp[maxn << 1][maxn], y[maxn];
int f(int z){
for(int i = 1; i <= n * 2; ++i)
y[i] = (X[v[i].first + v[i].second - 1] - (ll)z * v[i].first) % mod;
dp[1][0] = 1;
for(int i = 2; i <= n * 2; ++i)
for(int j = 0; j <= n; ++j){
dp[i][j] = dp[i - 1][j];
for(int k = j, mul = 1; k >= 1 && l[k] <= i - 1 && r[k] >= i; --k){
mul = (ll)mul * (y[i] - y[i - 1]) % mod * leninv[k] % mod * inv[j - k + 1] % mod;
dp[i][j] = (dp[i][j] + (ll)dp[i - 1][k - 1] * mul) % mod;
}
}
return dp[n * 2][n];
}
int solve(int a, int b){
std::vector<int> F;
for(int i = 0; i <= n; ++i){
std::vector<int> G{1};
int mul = 1;
for(int j = 0; j <= n; ++j)
if(i != j)
G = G * std::vector<int>{-j, 1}, mul = (ll)mul * (i - j) % mod;
mul = (ll)f(i) * qpow(mul, mod - 2) % mod;
F = F + G * mul;
}
F = integral(F);
return (eval(F, b) - eval(F, a)) % mod;
}
int main(){
std::scanf("%d", &n);
for(int i = 0; i <= n; ++i)
std::scanf("%d", X + i);
for(int i = 1; i <= n + 1; ++i)
leninv[i] = qpow(X[i] - X[i - 1], mod - 2), inv[i] = qpow(i, mod - 2);
frac l(0), r; int ans = 0;
while(true){
init(l), r = next();
if(!check()) break;
ans = (ans + solve(l.val(), r.val())) % mod;
l = r;
}
std::printf("%d\n", (ans + mod) % mod);
return 0;
}
来源¶
AtCoder Regular Contest 113 F - Social Distance