[NOI2017]泳池——概率DP+线性递推

时间:2022-10-14 16:10:27

[NOI2017]泳池

实在没有思路啊~~~

luogu题解

1.差分,转化成至多k的概率减去至多k-1的概率。这样就不用记录“有没有出现k”这个信息了

2.n是1e9,感觉要递推然后利用数列的加速技巧

f[n]表示宽度为n的值,然后枚举最后一个连续高度至少为1的块,dp数组辅助

神仙dp:dp[i][j]表示宽度为i,j的高度出现限制,任意矩形不大于k的概率

设计确实巧妙:宽度利于转移给f,高度利于自己的转移

dp数组转移:枚举第一个到达j的限制的位置,这样,前面部分限制至少是j+1,后面至少是j,中间概率是P^(j-1)*(1-P)

    并且,若i*j<=k,那么新增的就是中间的i*j的矩形,两边都满足的情况下,现在也一定满足

  后缀和优化,限制j的范围,调和级数可以证明状态实际上只有klogk个,转移O(k)

边界:f[0]=1,dp[0][*]=1,注意,*是0~k+2

3.我们现在得到了f的递推式

先递推找到f的前几项

k范围限制不能矩乘

利用[学习笔记]Cayley-Hilmiton

可以倍增加多项式取模

k<=1000可以暴力取模

发现,模数的最高项是1,所以不用逆元!

O(k^2(logn+logk))

代码:

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
#define int long long
using namespace std;
typedef long long ll;
il void rd(int &x){
char ch;x=;bool fl=false;
while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
for(x=numb;isdigit(ch=getchar());x=x*+numb);
(fl==true)&&(x=-x);
}
namespace Miracle{
const int N=;
const int mod=;
ll a[*N],b[*N],f[*N],yu[*N],ret[*N];
ll dp[*N][*N];
ll mi[*N];
ll n,k,x,y;
ll P;
ll ans;
ll qm(ll x,ll y){
ll ret=;
while(y){
if(y&) ret=ret*x%mod;
x=x*x%mod;
y>>=;
}
return ret;
}
void mul(ll *f,ll *g,ll *mo,int n){
memset(b,,sizeof b);
for(reg i=;i<=n;++i){
for(reg j=;j<=n;++j){
b[i+j]=(b[i+j]+f[i]*g[j]%mod)%mod;
}
}
for(reg i=*n;i>=n;--i){
if(b[i])
for(reg j=;j<=n;++j) b[i+j-n]=(b[i+j-n]+mod-b[i]*mo[j]%mod)%mod;
}
for(reg i=;i<=n;++i) f[i]=b[i],b[i]=;
for(reg i=n+;i<=*n;++i) f[i]=,b[i]=;
}
ll calc(int k){
memset(dp,,sizeof dp);
for(reg i=;i<=k+;++i) dp[][i]=;
for(reg i=;i<=k;++i){
for(reg j=k/i+;j>=;--j){
for(reg l=;l<=i;++l){
dp[i][j]=(dp[i][j]+dp[l-][j+]*dp[i-l][j]%mod*(+mod-P)%mod*mi[j-]%mod)%mod;
}
dp[i][j]=(dp[i][j+]+dp[i][j])%mod;
}
}
for(reg i=;i<=k;++i) a[k-i]=(mod-dp[i][]*(+mod-P)%mod)%mod;
a[k+]=;
memset(f,,sizeof f);
f[]=;
for(reg l=;l<=*k;++l){
for(reg i=;i<=l;++i){
if(l-i->=) f[l]=(f[l]+f[l-i-]*dp[i][]%mod*(+mod-P)%mod)%mod;
else f[l]=(f[l]+dp[i][])%mod;
}
}
if(n<=*k) return f[n];
int y=n;
memset(ret,,sizeof ret);
memset(yu,,sizeof yu);
ret[]=;
yu[]=;
while(y){
if(y&) mul(ret,yu,a,k+);
mul(yu,yu,a,k+);
y>>=;
}
ll bac=;
for(reg i=;i<=k;++i){
bac=(bac+ret[i]*f[i]%mod)%mod;
}
return bac;
}
int main(){
scanf("%lld%lld%lld%lld",&n,&k,&x,&y);
P=x*qm(y,mod-)%mod;
int lp=k;
mi[]=;
for(reg i=;i<=k;++i) mi[i]=mi[i-]*P%mod;
printf("%lld\n",(calc(lp)-calc(lp-)+mod)%mod);
return ;
} }
signed main(){
Miracle::main();
return ;
} /*
Author: *Miracle*
Date: 2019/2/25 16:31:18
*/

总结:
1.DP还是有一点套路可循的,f是枚举最后连续的段,dp是枚举第一个限制到j的位置

2.线性递推的优化方法。