题解 AT4690 e
题意简述
- 一张有 \(M\) 个位置的长椅,初始为空。
- 有 \(M\) 个人依次落座,每个人会在长椅中满足 自己与相邻共三个位置均未被占据 的位置中随机选择一个坐下。若不存在满足条件的位置则不落座。
- 随机选择 \(N\) 个连续的位置,求当 \(M\) 趋近于无穷大时,这 \(N\) 个位置的落座位置与给出序列 \(s\) 一致的概率的极限(\(s\) 中
X
表示有人,-
表示无人)。 - 可以证明答案可以表示成 \(p+\frac{q}{e}+\frac{r}{e^2}\),其中 \(p,q,r\) 均为有理数,对 \(10^9+7\) 取模。
- \(1 \leq |S|=N \leq 1000\)。
题解
Part 1
先考虑 \(s=X\)(样例 \(1\)),即求落座比例的期望。
题目中的过程较为繁琐,每个人会随机选择自己的位置。我们不妨令每个人只对应一个位置,这样最后的情况只与落座顺序有关。
具体而言,随机排列 \(P_1,\cdots,P_M\),然后让 \(M\) 个人依次落座 \(P_1,\cdots,P_M\)(如果符合题意的话)。
但是由于 \(M\) 趋近于无穷大,枚举排列不好处理。考虑给每个位置赋一个 \([0,1]\) 中的值 \(X_i\) 来代替排列 \(P\)。让 \(M\) 个人按 \(X_i\) 从小到大的顺序依次落座,仍是等价的。
考虑一个位置 \(i\) 最终是否被占据。
设其向左与向右的极长连续下降长度分别为 \(L\) 与 \(R\)(即 \(X_{i-L-1}<X_{i-L}<\cdots<X_i<\cdots<X_{i+R}>x_{i+R+1}\)),那么 \(i-L\) 与 \(i+R\) 一定会被占据,同理 \(i-L+2\) 与 \(i+R-2\) 也会被占据……
故 \(i\) 被占据当且仅当 \(L\) 与 \(R\) 都是偶数。
设 \(X_i=x\),那么 \(L\) 为偶数的概率为:
故答案为 \(\int_0^1 e^{-2x}=\dfrac{1}{2}-\dfrac{1}{2e^2}\)。
Part 2
考虑原问题,可以考虑 DP。
假设选出的 \(N\) 个位置就是 \(1,2,\cdots,N\),我们来从左向右填 \(X_i\)(\(-\infty\) 到 \(0\) 均已填好了),同时计算与 \(s\) 匹配的概率。
这样会有一个问题:在 DP 过程中,我们只知道 \(L\) 的奇偶性,但 \(R\) 的奇偶性还未确定。
但是 \(R\) 的奇偶性只有两种,不妨倒推。
设 \(dp_{i,id_1,id_2,id_3}\) 为关于 \(x\) 的概率密度函数,四个变量分别表示到第 \(i\) 个位置,\(L\) 的奇偶性,\(R\) 奇偶性为 \(0\) 时是否合法,\(R\) 奇偶性为 \(1\) 时是否合法。
转移时仍然设 \(X_i=x\),考虑 \(X_{i-1}\) 与 \(x\) 的大小关系。
- \(X_{i-1}<x\)。则 \(dp_{i-1,id_1 \oplus 1,1,0/1} \rightarrow dp_{i,id_1,X,X}\),具体而言:
- 如果该位为
-
:\(dp_{i-1,1,1,0/1} \rightarrow dp_{i,0,0,1}\);\(dp_{i-1,0,1,0/1} \rightarrow dp_{i,1,1,1}\)。 - 如果该位为
X
:\(dp_{i-1,1,1,0/1} \rightarrow dp_{i,0,1,0}\)。
- 如果该位为
- \(X_{i-1}>x\)。则 \(dp_{i-1,0/1,id_2,id_3} \rightarrow dp_{i,0,X,X}\),具体而言:
- 如果该位为
-
:\(dp_{i-1,0/1,1,0/1} \rightarrow dp_{i,0,0,1}\)。 - 如果该位为
X
:\(dp_{i-1,0/1,0/1,1} \rightarrow dp_{i,0,1,0}\)。
- 如果该位为
初始状态为 \(dp_{0,0,1,1}=e^{-x}\) 及 \(dp_{0,1,1,1}=1-e^{-x}\)。
实现时需要计算 \(\int_0^x dp_{i-1,X,X,X} dy\) 或是 \(\int_x^1 dp_{i-1,X,X,X} dy\),注意到 \(c \int e^{-x}=-ce^{-x}+C\),故每个 DP 值一定是 \(P(x)+Q(x)e^{-1}+ce^{-x}\),其中 \(P(x)\) 与 \(Q(x)\) 均为关于 \(x\) 的不超过 \(i\) 次多项式。
计算答案时,要将 DP 状态乘上 \(e^{-x}\) 或是 \(1\) 或是 \(1-e^{-x}\) 再求在 \([0,1]\) 上的定积分。
只需要计算两个积分:\(\int_0^1 x^k e^{-x} dx\) 及 \(\int_0^1 e^{-2x}dx=\dfrac{1}{2}-\dfrac{1}{2e^2}\)。
分部积分法:
总时间复杂度 \(O(N^2)\)。
代码
#include<bits/stdc++.h>
using namespace std;
const long long N=1100,mod=1e9+7,inv2=(mod+1)/2;
long long n;
char S[N];
long long invv[N],ans[3];
struct abc{
long long P[N],Q[N],c;
abc(){
memset(P,0,sizeof(P));
memset(Q,0,sizeof(Q));
c=0;
}
}dp[N][2][2][2];
abc operator +(abc X,abc Y){
abc T;
for(long long i=0;i<=n;i++) T.P[i]=(X.P[i]+Y.P[i])%mod,T.Q[i]=(X.Q[i]+Y.Q[i])%mod;
T.c=(X.c+Y.c)%mod;
return T;
}
abc operator -(abc X,abc Y){
abc T;
for(long long i=0;i<=n;i++) T.P[i]=(X.P[i]+mod-Y.P[i])%mod,T.Q[i]=(X.Q[i]+mod-Y.Q[i])%mod;
T.c=(X.c+mod-Y.c)%mod;
return T;
}
abc INTx(abc X){
abc T;
for(long long i=0;i<=n;i++) T.P[i+1]=X.P[i]*invv[i+1]%mod;
for(long long i=0;i<=n;i++) T.Q[i+1]=X.Q[i]*invv[i+1]%mod;
T.c=(mod-X.c)%mod;
return T;
}//求 x 处积分
abc INT0(abc X){
abc T;
T.P[0]=(mod-X.c)%mod;
return T;
}//求 0 处积分
abc INT1(abc X){
abc T;
for(long long i=0;i<=n;i++) T.P[0]=(T.P[0]+X.P[i]*invv[i+1]%mod)%mod;
for(long long i=0;i<=n;i++) T.Q[0]=(T.Q[0]+X.Q[i]*invv[i+1]%mod)%mod;
T.Q[0]=(T.Q[0]+mod-X.c)%mod;
return T;
}//求 1 处积分
abc cal1(abc X){
return INTx(X)-INT0(X);
}
abc cal2(abc X){
return INT1(X)-INTx(X);
}
void solve1(abc X,long long w){
abc T=INT1(X)-INT0(X);
ans[0]=(ans[0]+T.P[0]*w%mod)%mod;
ans[1]=(ans[1]+T.Q[0]*w%mod)%mod;
}
void solvex(abc X,long long w){
ans[0]=(ans[0]+X.c*inv2%mod*w%mod)%mod;
ans[2]=(ans[2]+mod-X.c*inv2%mod*w%mod)%mod;
ans[0]=(ans[0]+X.P[0]*w%mod)%mod;
ans[1]=(ans[1]+mod-X.P[0]*w%mod)%mod;
ans[1]=(ans[1]+X.Q[0]*w%mod)%mod;
ans[2]=(ans[2]+mod-X.Q[0]*w%mod)%mod;
for(long long k=1;k<=n;k++){
long long tot=1;
long long nw=1;
for(long long i=k;i>=1;i--){
nw=nw*i%mod;
tot=(tot+nw)%mod;
}
ans[1]=(ans[1]+mod-X.P[k]*tot%mod*w%mod)%mod;
ans[0]=(ans[0]+X.P[k]*nw%mod*w%mod)%mod;
ans[2]=(ans[2]+mod-X.Q[k]*tot%mod*w%mod)%mod;
ans[1]=(ans[1]+X.Q[k]*nw%mod*w%mod)%mod;
}
}//求答案
int main(){
invv[1]=1;
for(long long i=2;i<N;i++) invv[i]=(mod-mod/i)*invv[mod%i]%mod;
cin>>n;
scanf("\n%s",S+1);
dp[0][0][1][1].c=1;
dp[0][1][1][1].P[0]=1,dp[0][1][1][1].c=mod-1;
for(long long i=1;i<=n;i++){
if(S[i]=='-'){
dp[i][0][0][1]=dp[i][0][0][1]+cal1(dp[i-1][1][1][0]+dp[i-1][1][1][1]);
dp[i][1][1][1]=dp[i][1][1][1]+cal1(dp[i-1][0][1][0]+dp[i-1][0][1][1]);
dp[i][0][0][1]=dp[i][0][0][1]+cal2(dp[i-1][0][1][0]+dp[i-1][1][1][0] +dp[i-1][0][1][1]+dp[i-1][1][1][1]);
}
else{
dp[i][0][1][0]=dp[i][0][1][0]+cal1(dp[i-1][1][1][0]+dp[i-1][1][1][1]);
dp[i][0][1][0]=dp[i][0][1][0]+cal2(dp[i-1][0][0][1]+dp[i-1][1][0][1] +dp[i-1][0][1][1]+dp[i-1][1][1][1]);
}
}
for(long long i=0;i<=1;i++)
for(long long j=0;j<=1;j++)
for(long long k=0;k<=1;k++){
if(j && k) solve1(dp[n][i][j][k],1);
if(j && !k) solvex(dp[n][i][j][k],1);
if(!j && k) solve1(dp[n][i][j][k],1),solvex(dp[n][i][j][k],mod-1);
}
cout<<ans[0]<<" "<<ans[1]<<" "<<ans[2];
}