[ARC113F]Social Distance

题目

点这里看题目。

分析

首先需要知道，在此题中连续随机变量的期望可以如下计算：

[E(X)=int_{0}^{+infty} P(t< X)mathop{}!mathrm d t ]

关于这个东西的理解：

首先考虑一般的离散随机变量 (X) ，它有 (n) 个取值 (0,delta,2delta,3delta,dots,(n-1)delta) ，其中 (delta) 满足 (delta=frac{M}{n}) ，且 (M) 为常数。那么 它的期望就是：

[egin{aligned} E(x)=&sum_{k=0}^{n-1}kdelta imes P(X=kdelta)\ =&sum_{k=0}^{n-1}P(X=kdelta)sum_{j=0}^{k-1}delta\ =&sum_{j=0}^{n-1}deltasum_{k=j+1}^{n-1}P(X=kdelta)\ =&sum_{j=0}^{n-1}delta imes P(jdelta<X) end{aligned} ]

考虑 (n ightarrow +infty) 的时候，可以发现：

[E(X)=int_{0}^{M}mathop{}!mathrm dt imes P(t< X) ]

由于 (X) 的取值被限制在 ([0,M)) 里面，所以积分上界也可以被修改。

在此题中，我们直接设 (X) 为任意两个人之间距离的最小值，那么考虑 (t) 确定的时候该如何计算 (P(t<X)) 。很显然，如果设第 (k) 个人的位置为 (p_k) （其中 (k) 为 (0sim N-1) ），那么必有：

[forall 0le k<N-1,p_{k+1}-p_{k}>t ]

适当地做一下变形，让变量独立：设 (p'_{k}=p_k-kt) ，那么就有：

[forall 0le k<N-1,p'_{k+1}>p'_k ]

此时每个 (p') 也有其对应的取值区间，那么问题就变成了：

给定 (n) 个区间，有 (n) 个连续随机变量，第 (k) 个的取值范围为 ([l_k,r_k]) 。求这 (n) 个随机变量组成的序列为单调递增的概率。

这是经典的问题。我们可以根据所有的 ([l_k,r_k]) 的端点，将所有可用值划分为若干段区间，其中第 (j) 段为 ([L_j,R_j)) ，于是可以设计 DP ：

(f_{i,j}) 表示前 (i) 个随机变量组成的序列递增，且第 (i) 个随机变量落在第 (j) 个区间里的概率。

不难写出转移：

[f_{i,j}=sum_{k=1}^{i}left(left(sum_{l=1}^{j-1}f_{k-1,l} ight) imes left(prod_{p=k}^i[[L_j,R_j)sube [l_p,r_p]]frac{R_j-L_j}{r_p-l_p} ight) ight) ]

使用前缀和优化，转移可以做到 $O(n^3w) $，其中 (w) 是对元素做单次运算的复杂度。

但是 (t) 实际上是个变量，我们不可能对每一个 (t) 都计算一遍。不过可以注意到，转移的方式和 ([l_k,r_k]) 这些端点的相对位置有关，而转移的结果才和 (t) 有关。根据转移，我们可以发现在转移方式确定的时候，转移结果是 (t) 的多项式，那么这样就很好处理了，对于同样的转移方式，我们可以得出结果多项式；同时不难发现同样转移方式对应的 (t) 必然是一段区间，因此可以直接算出该段的积分结果。

其一：不太清楚这里可不可以代入 (n+1) 个点，从而直接维护多项式的点值，这样加法、乘法都是 (O(n)) 的；

其实这里乘的就是一次函数，本身乘法就是 (O(deg)) 的

其二：可以发现这里的 (P(t<X)) 其实就是一个分段函数，分段函数自然就要分段积分。

这里转移方式对应的就是区间端点相对位置，可以发现区间端点相对位置只会变化 (O(n^2)) 次，因此 (t) 的区间也只有 (O(n^2)) 段。每一段暴力计算，时间复杂度即是 (O(n^6)) 。

小结：

1. 关于连续随机变量，用得虽然不多，但是一定要了解，不然这道题都下不去手；
2. 其中关于区间的转化其实挺巧的，主要还是一个分离变量的思想；
3. 这里将 DP 的含义变为了对多项式进行 DP ，虽然还是在应对 (t) 是连续的限制，不过这也提醒了我们：
   1. DP 的值理论上可以是任何东西，只要能满足一定的运算即可；
   2. DP 的转移方式和转移结果并不是绑在一块的东西——通过向 DP 值中引入变量我们可以让它们两个脱钩。

代码

#include <cstdio>
#include <cstring>
#include <utility>
#include <algorithm>
using namespace std;

#define rep( i, a, b ) for( int i = (a) ; i <= (b) ; i ++ )
#define per( i, a, b ) for( int i = (a) ; i >= (b) ; i -- )

typedef long long LL;
typedef pair<double, int> Point;

const double eps = 1e-6;
const int mod = 998244353, inv2 = mod + 1 >> 1;
const int MAXN = 35;

template<typename _T>
void read( _T &x )
{
	x = 0; char s = getchar(); int f = 1;
	while( s < '0' || '9' < s ) { f = 1; if( s == '-' ) f = -1; s = getchar(); }
	while( '0' <= s && s <= '9' ) { x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar(); }
	x *= f;
}

template<typename _T>
void write( _T x )
{
	if( x < 0 ) putchar( '-' ), x = - x;
	if( 9 < x ) write( x / 10 );
	putchar( x % 10 + '0' );
}

template<typename _T>
_T ABS( const _T a )
{
	return a < 0 ? -a : a;
}

inline LL Gcd( LL, LL );
inline int Mul( int, int );
inline int Sub( int, int );
inline int Add( int, int );
inline int Inv( const int );

struct Fraction
{
	LL a, b;
	Fraction() { a = 0, b = 1; }
	Fraction( int V ) { a = V, b = 1; }
	int operator () () const { return Mul( a, Inv( b ) ); }
	Fraction( int A, int B ) { LL d = Gcd( A, B ); a = A / d, b = B / d; }
	bool operator < ( const Fraction &g ) const { return a * g.b < b * g.a; }
	bool operator != ( const Fraction &g ) const { return a ^ g.a || b ^ g.b; }
	bool operator == ( const Fraction &g ) const { return a == g.a && b == g.b; }
};

int smallInv[MAXN];

struct Poly
{
	int coe[MAXN], n;
	Poly() { n = 0, memset( coe, 0, sizeof coe ); }
	Poly( int B ) { memset( coe, 0, sizeof coe ), coe[n = 0] = B; }
	Poly( int K, int B ) { memset( coe, 0, sizeof coe ), coe[n = 1] = K, coe[0] = B; }
	Poly& operator /= ( const int &g ) { return *this = *this / g; }
	Poly& operator *= ( const int &g ) { return *this = *this * g; }
	Poly& operator *= ( const Poly &g ) { return *this = *this * g; }
	Poly& operator += ( const Poly &g ) { return *this = *this + g; }
	Poly& operator -= ( const Poly &g ) { return *this = *this - g; }
	
	int operator () ( const int x ) const
	{
		int ret = 0;
		for( int i = 0, pw = 1 ; i <= n ; i ++, pw = Mul( pw, x ) )
			ret = Add( ret, Mul( coe[i], pw ) );
		return ret;
	}
	
	Poly operator + ( const Poly &g ) const
	{
		Poly ret; 
		for( int &k = ret.n ; k <= n || k <= g.n ; k ++ ) ret.coe[k] = Add( coe[k], g.coe[k] );
		for( int &k = ret.n ; k && ! ret.coe[k] ; k -- );
		return ret;
	}
	
	Poly operator - ( const Poly &g ) const 
	{
		Poly ret; 
		for( int &k = ret.n ; k <= n || k <= g.n ; k ++ ) ret.coe[k] = Sub( coe[k], g.coe[k] );
		for( int &k = ret.n ; k && ! ret.coe[k] ; k -- );
		return ret;
	}
	
	Poly operator * ( const int d ) const
	{
		Poly ret;
		for( int &k = ret.n ; k <= n ; k ++ ) ret.coe[k] = Mul( coe[k], d );
		for( int &k = ret.n ; k && ! ret.coe[k] ; k -- );
		return ret;
	}
	
	Poly operator / ( const int d ) const
	{
		Poly ret; int v = Inv( d );
		for( int &k = ret.n ; k <= n ; k ++ ) ret.coe[k] = Mul( coe[k], v );
		for( int &k = ret.n ; k && ! ret.coe[k] ; k -- );
		return ret;
	}

	Poly Integral() const
	{
		Poly ret; ret.n = n + 1;
		rep( i, 0, n ) ret.coe[i + 1] = Mul( coe[i], smallInv[i + 1] );
		return ret;
	}

	Poly operator * ( const Poly &g ) const 
	{
		Poly ret; ret.n = n + g.n;
		rep( i, 0, n ) rep( j, 0, g.n )
			ret.coe[i + j] = Add( ret.coe[i + j], Mul( coe[i], g.coe[j] ) );
		for( int &k = ret.n ; k && ! ret.coe[k] ; k -- );
		return ret;
	}
};

Poly dp[MAXN][MAXN << 1];

Point pnt[MAXN << 1];

Fraction seq[MAXN * MAXN * 4];
int fac[MAXN], ifac[MAXN];
int X[MAXN], fir[MAXN], lst[MAXN];
int N, tot = 0;

inline int Qkpow( int, int );
inline int Mul( int x, int v ) { return 1ll * x * v % mod; }
inline int Inv( const int a ) { return Qkpow( a, mod - 2 ); }
inline LL Gcd( LL a, LL b ) { return b ? Gcd( b, a % b ) : a; }
inline int Sub( int x, int v ) { return ( x -= v ) < 0 ? x + mod : x; }
inline int Add( int x, int v ) { return ( x += v ) >= mod ? x - mod : x; }

inline int Qkpow( int base, int indx )
{
	int ret = 1;
	while( indx )
	{
		if( indx & 1 ) ret = Mul( ret, base );
		base = Mul( base, base ), indx >>= 1;
	}
	return ret;
}

int Integral( const Poly &f, const int L, const int R ) 
{ return Sub( f.Integral()( R ), f.Integral()( L ) ); }

Poly Length( const int id ) 
{
	Poly l, r;
	int idL = pnt[id].second, idR = pnt[id + 1].second;
	l = Poly( - ABS( idL ), idL < 0 ? X[- idL] : X[idL - 1] );
	r = Poly( - ABS( idR ), idR < 0 ? X[- idR] : X[idR - 1] );
	return r - l; 
}

int Calc( const Fraction lef, const Fraction rig )
{
	double dif = 1.0 * lef.a / lef.b + eps;
	int len = 0;
	rep( i, 0, N - 1 ) 
		pnt[++ len] = Point( X[i] - i * dif, i + 1 ), 
		pnt[++ len] = Point( X[i + 1] - i * dif, - i - 1 );
	sort( pnt + 1, pnt + 1 + len );
	rep( i, 1, len )
	{
		if( pnt[i].second > 0 ) fir[pnt[i].second - 1] = i;
		else lst[- pnt[i].second - 1] = i - 1;
	}
	rep( i, 0, N - 1 ) rep( j, 0, len - 1 ) dp[i][j] = Poly();
	rep( i, fir[0], lst[0] ) dp[0][i] = Length( i ) / ( X[1] - X[0] );
	rep( i, 1, len - 1 ) dp[0][i] += dp[0][i - 1];
	Poly tmp, poss;
	rep( i, 1, N - 1 ) 
	{
		rep( j, fir[i], lst[i] )
		{
			poss = 1;
			for( int k = i ; ~ k && fir[k] <= j && j <= lst[k] ; k -- )
			{
				if( ! k ) tmp = 1;
				else tmp = dp[k - 1][j - 1];
				poss *= Length( j ) / ( X[k + 1] - X[k] );
				dp[i][j] += poss * tmp * ifac[i - k + 1];
			}
		}
		rep( j, 1, len - 1 ) dp[i][j] += dp[i][j - 1];
	}
	int L = lef(), R = rig();
	return Integral( dp[N - 1][len - 1], L, R );
}

void Init( const int n )
{
	fac[0] = 1; rep( i, 1, n ) fac[i] = Mul( fac[i - 1], i );
	ifac[n] = Inv( fac[n] ); per( i, n - 1, 0 ) ifac[i] = Mul( ifac[i + 1], i + 1 );
	smallInv[1] = 1; rep( i, 2, n ) smallInv[i] = Mul( smallInv[mod % i], mod - mod / i );
}

int main()
{
	read( N ), Init( N + 1 );
	rep( i, 0, N ) read( X[i] );
	rep( i, 0, N - 1 )
		rep( j, i + 1, N - 1 )
		{
			seq[++ tot] = Fraction( X[j] - X[i], j - i );
			seq[++ tot] = Fraction( X[j + 1] - X[i + 1], j - i );
			seq[++ tot] = Fraction( X[j + 1] - X[i], j - i );
			seq[++ tot] = Fraction( X[j] - X[i + 1], j - i );
		}
	sort( seq + 1, seq + 1 + tot );
	tot = unique( seq + 1, seq + 1 + tot ) - seq - 1;
	
	seq[++ tot] = 1e7; int ans = 0;
	rep( i, 1, tot ) 
		ans = Add( ans, Calc( seq[i], seq[i + 1] ) );
	write( ans ), putchar( '
' );
	return 0;
}