Atcoder Regular Contest 113F. Social Distance
给定一个长度为 \(n+1\) 的非负整数序列 \(X_0,X_1,\ldots,X_{n}\)，满足 \(0=X_0<X_1<\ldots<X_{n}\)。
在每个区间 \([X_0,X_1],[X_1,X_2],\ldots,[X_{n-1},X_{n}]\) 中均匀随机一个点，求相邻两个点之间的最短距离的期望。
答案对 \(998244353\) 取模。
\(2\le n\le 20,0\le X_i\le 10^6\)。
首先，此题中的连续随机变量的期望可以用下式计算
考虑一般的离散型连续变量 \(X\)，它有 \(n+1\) 个取值 \(0,\Delta x,2\Delta x,\ldots,n\Delta x\)，\(\Delta x\) 满足 \(\Delta x=\frac{V}{n}\)，期望为
\[\begin{aligned}\text{E}(X)&=\sum_{i=0}^{n}i\Delta x\cdot \text{P}(X=i\Delta x)\\ &=\sum_{i=0}^{n}\text{P}(X=i\Delta x)\sum_{j=0}^{i-1}\Delta x\\&=\sum_{j=0}^{n}\Delta x\sum_{i=j+1}^{n}\text{P}(X=j\Delta x)\\ &=\sum_{j=0}^{n}\Delta x\cdot\text{P}(j\Delta x<X)\end{aligned} \]令 \(n\rightarrow \infty\)，可得
\[\text{E}(X)=\int_{0}^{V}\text{P}(t<X)\text{d}t \]
考虑固定 \(t\) 后计算 \(\text{P}(t<X)\)。设第 \(i\) 个人点的位置为 \(p_i\)，那么有
令 \(p'_{i}=p_{i}-it\)，则有
根据 \(X\) 序列可以求出 \(p'\) 的取值范围，进而问题转化成
有 \(n\) 个连续随机变量，第 \(i\) 个连续随机变量 \(x_i\) 的取值范围为 \([l_i,r_i]\)，求 \(x_1\le x_2\le \ldots\le x_{n}\) 的概率。
先根据所有取值范围的左右边界放到坐标轴上，此时形成 \(2n-1\) 个区间，设第 \(i\) 个区间为 \([L_i,R_i]\)。
设 \(dp_{i,j}\) 为第 \(i\) 个连续随机变量在第 \(j\) 个区间且 \(x_1\le x_2\le\ldots\le x_{i}\) 的概率，则有转移
即枚举最后一个在第 \(j-1\) 个区间的随机变量 \(x_k\)，于是 \(x_{k+1},x_{k+2},\ldots,x_i\) 都在第 \(j\) 个区间，满足 \(x_{k+1},x_{k+2},\ldots,x_i\) 单调不减的概率显然是 \(\frac{1}{(i-k)!}\)。
维护 \(dp_{i,j}\) 关于 \(j\) 的前缀和以及上式的乘积式就可以做到 \(O(n^3)\) 求出答案。
但是 \(t\) 是变量，而不是固定的。考虑到 \(p'_{i}=-i\cdot t+p_i\) 是关于 \(t\) 的一次函数，\(\text{dp}\) 的转移方式只与 \(L_i,R_i\) 的相对顺序有关，而在若干个段之间的 \(t\) 形成的 \(L_i,R_i\) 的相对顺序是相同的，可以发现不同的相对顺序是 \(O(n^2)\) 级别的。维护关于 \(t\) 的多项式进行转移，最后再对结果多项式进行积分即可。
总时间复杂度 \(O(n^6)\)。
#include<bits/stdc++.h>
#define fi first
#define se second
using namespace std;
typedef long long ll;
const double eps=1e-6;
const int N=45,M=5005,mod=998244353;
struct Poly{
int n;ll coe[N];
Poly(){n=0;memset(coe,0,sizeof coe);}
Poly operator+(const Poly&p)const{
Poly ret;ret.n=max(n,p.n);
for(int k=0;k<=ret.n;++k)ret.coe[k]=(coe[k]+p.coe[k])%mod;
while(ret.n&&!ret.coe[ret.n])--ret.n;
return ret;
}
Poly operator*(const Poly&p)const{
Poly ret;ret.n=n+p.n;
for(int j=0;j<=n;++j)for(int k=0;k<=p.n;++k)
(ret.coe[j+k]+=coe[j]*p.coe[k]%mod)%=mod;
while(ret.n&&!ret.coe[ret.n])--ret.n;
return ret;
}
Poly operator*(const int c)const{
Poly ret;ret.n=n;
for(int k=0;k<=n;++k)ret.coe[k]=coe[k]*c%mod;
return ret;
}
}Len[M],dp[N][N];
struct Frac{
ll x,y;
bool operator<(const Frac&a){return x*a.y<a.x*y;}
}seq[M];
int n,X[N],cnt,idL[N],idR[N];ll ans,fac[M],inv[M],Inv[M];
pair<double,int>p[N];
inline ll qpow(ll a,ll b){ll ans=1;for(;b;b>>=1,a=a*a%mod)if(b&1)ans=ans*a%mod;return ans;}
inline void precalc(int n){
fac[0]=inv[0]=Inv[0]=fac[1]=inv[1]=Inv[1]=1;
for(int i=2;i<=n;++i)
fac[i]=fac[i-1]*i%mod,
Inv[i]=(mod-mod/i)*Inv[mod%i]%mod,
inv[i]=inv[i-1]*Inv[i]%mod;
}
inline ll Integral(Poly f,ll x){
ll ans=0,pw=x;
for(int i=0;i<=f.n;++i,pw=pw*x%mod)(ans+=pw*Inv[i+1]%mod*f.coe[i]%mod)%=mod;
return ans;
}
inline ll calc(Frac sL,Frac sR){
int tot=0;ll FR=sL.x*qpow(sL.y,mod-2)%mod,RF=sR.x*qpow(sR.y,mod-2)%mod;
double tL=1.0*sL.x/sL.y+eps;
for(int i=0;i<n;++i)
p[++tot]=make_pair(X[i]-i*tL,i+1),p[++tot]=make_pair(X[i+1]-i*tL,-i-1);
sort(p+1,p+1+tot);
for(int i=1;i<tot;++i)
Len[i].coe[0]=(X[p[i+1].se>0?p[i+1].se-1:-p[i+1].se]-X[p[i].se>0?p[i].se-1:-p[i].se]+mod)%mod,
Len[i].coe[1]=(abs(p[i].se)-abs(p[i+1].se)+mod)%mod,Len[i].n=1;
for(int i=1;i<=tot;++i){
if(p[i].se>0)idL[p[i].se-1]=i;
else idR[-p[i].se-1]=i-1;
}
for(int i=0;i<n;++i)for(int j=0;j<tot;++j)dp[i][j]=Poly();
for(int i=idL[0];i<=idR[0];++i)dp[0][i]=Len[i]*qpow(X[1],mod-2);
for(int i=1;i<tot;++i)dp[0][i]=dp[0][i]+dp[0][i-1];
for(int i=1;i<n;++i){
for(int j=idL[i];j<=idR[i];++j){
Poly sufp;sufp.coe[0]=1;
for(int k=i;~k&&idL[k]<=j&&j<=idR[k];--k){
Poly tmp;
if(!k)tmp.coe[0]=1;else tmp=dp[k-1][j-1];
sufp=sufp*Len[j]*qpow(X[k+1]-X[k],mod-2);
dp[i][j]=dp[i][j]+sufp*tmp*inv[i-k+1];
}
}
for(int j=1;j<tot;++j)dp[i][j]=dp[i][j]+dp[i][j-1];
}
ll Ans=(Integral(dp[n-1][tot-1],RF)-Integral(dp[n-1][tot-1],FR)+mod)%mod;
return Ans;
}
int main(){
precalc(5000);
scanf("%d",&n);
for(int i=0;i<=n;++i)scanf("%d",X+i);
for(int i=0;i<n;++i)
for(int j=i+1;j<n;++j)
seq[++cnt]=(Frac){X[j]-X[i],j-i},seq[++cnt]=(Frac){X[j+1]-X[i+1],j-i},
seq[++cnt]=(Frac){X[j]-X[i+1],j-i},seq[++cnt]=(Frac){X[j+1]-X[i],j-i};
sort(seq+1,seq+1+cnt);
for(int i=1;i<=cnt;++i)(ans+=calc(seq[i],seq[i+1]))%=mod;
return printf("%lld\n",ans),0;
}
