「NOI2017」泳池

时间:2023-03-08 23:46:33
「NOI2017」泳池

DP式子比后面的东西难推多了

LOJ2304

Luogu P3824

UOJ #316


题意

给定一个长度为$ n$高为$ \infty$的矩形

每个点有$ 1-P$的概率不可被选择

求最大的和底边重合的不包含不可选点的矩形的面积为$ K$的概率

$ n \leq 10^9 k \leq 10^3$


题解

K可以出到50000的

首先考虑DP

面积恰好为$ K$的概率可以差分为不高于$ K$的概率减去不高于$ K-1$的概率

设$ f[i][j]$表示长度为$ i$的矩形,从底边起$ j$行都可选,最大面积不大于$ K$的概率

边界为$ f[0][j]=1,f[i][j]=0 当且仅当i*j>k$

考虑转移,要么第$ j+1$行也都可选,要么第$ j+1$行有不可选的位置

对于第二种情况我们枚举从左到右第一个不可选的位置

有转移方程式

$$f[i][j]=f[i][j+1]*P^i+\sum_{k=1}^iP^{k-1}(1-P)f[k-1][j+1]·f[i-k][j]$$

我们要求的是$ f[n][0]$

由于$ i*j \leq K$因此复杂度大致是$ n·k \log k$的

可以得$ 70$分

容易发现当$ n$远大于$ k$的时候,每连续$ k$列必然有一列最低端有不可选点

令$ F[i]$表示当$ i>k$时,长度为$ i$的矩形的答案

枚举从右往左第一个不可选点,有转移方程式

$$ F[i]=f[i-k]*(1-P)*f[k-1][1]*P^{k-1}$$

这是一个线性递推的标准形式,可以用特征多项式优化到$ O(k^2 \log n)$甚至$ O(k \log k \log n)$

如果采用后一种的话复杂度的瓶颈在于前面的$ O(k^2)DP$,这部分可以用分治$ NTT$优化

然而我并没有写


代码

#include<ctime>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<queue>
#include<vector>
#define p 998244353
#define rt register int
#define ll long long
using namespace std;
inline ll read(){
ll x=;char zf=;char ch=getchar();
while(ch!='-'&&!isdigit(ch))ch=getchar();
if(ch=='-')zf=-,ch=getchar();
while(isdigit(ch))x=x*+ch-'',ch=getchar();return x*zf;
}
void write(ll y){if(y<)putchar('-'),y=-y;if(y>)write(y/);putchar(y%+);}
void writeln(const ll y){write(y);putchar('\n');}
namespace poly{
vector<int>R;
int ksm(int x,int y=p-){
int ans=;
for(rt i=y;i;i>>=,x=1ll*x*x%p)if(i&)ans=1ll*ans*x%p;
return ans;
}
void NTT(int n,vector<int>&A,int fla){
A.resize(n);
for(rt i=;i<n;i++)if(i>R[i])swap(A[i],A[R[i]]);
for(rt i=;i<n;i<<=){
int w=ksm(,(p-)//i);
for(rt j=;j<n;j+=i<<){
int K=;
for(rt k=;k<i;k++,K=1ll*K*w%p){
int x=A[j+k],y=1ll*K*A[i+j+k]%p;
A[j+k]=(x+y)%p,A[i+j+k]=(x-y)%p;
}
}
}
if(fla==-){
reverse(A.begin()+,A.end());
int invn=ksm(n);
for(rt i=;i<n;i++)A[i]=1ll*A[i]*invn%p;
}
}
vector<int>Mul(vector<int>x,vector<int>y){
int lim=,sz=x.size()+y.size()-;
while(lim<=sz)lim<<=;R.resize(lim);
for(rt i=;i<lim;i++)R[i]=(R[i>>]>>)|(i&)*(lim>>);
NTT(lim,x,);NTT(lim,y,);
for(rt i=;i<lim;i++)x[i]=1ll*x[i]*y[i]%p;
NTT(lim,x,-);x.resize(sz);
return x;
}
vector<int>sqr(vector<int>x){
int lim=,sz=x.size()*-;
while(lim<=sz)lim<<=;R.resize(lim);
for(rt i=;i<lim;i++)R[i]=(R[i>>]>>)|(i&)*(lim>>);
NTT(lim,x,);for(rt i=;i<lim;i++)x[i]=1ll*x[i]*x[i]%p;
NTT(lim,x,-);x.resize(sz);
return x;
}
vector<int>Inv(vector<int>A,int n=-){
if(n==-)n=A.size();
if(n==)return vector<int>(,ksm(A[]));
vector<int>b=Inv(A,(n+)/);
int lim=;while(lim<=n+n)lim<<=;R.resize(lim);
for(rt i=;i<lim;i++)R[i]=(R[i>>]>>)|(i&)*(lim>>);
A.resize(n);NTT(lim,A,);NTT(lim,b,);
for(rt i=;i<lim;i++)A[i]=1ll*b[i]*(2ll-1ll*A[i]*b[i]%p)%p;
NTT(lim,A,-);A.resize(n);
return A;
}
vector<int>Div(vector<int>A,vector<int>B){
int n=A.size(),m=B.size();
reverse(A.begin(),A.end());
reverse(B.begin(),B.end());
A.resize(n-m+),B.resize(n-m+);
int lim=;while(lim<=*(n-m+))lim<<=;R.resize(lim);
for(rt i=;i<lim;i++)R[i]=(R[i>>]>>)|(i&)*(lim>>);
vector<int>ans=Mul(A,Inv(B));ans.resize(n-m+);
reverse(ans.begin(),ans.end());
return ans;
}
vector<int>add(vector<int>A,vector<int>B){
int len=max(A.size(),B.size());A.resize(len+);
for(rt i=;i<=len;i++)(A[i]+=B[i])%=p;
return A;
}
vector<int>sub(vector<int>A,vector<int>B){
int len=max(A.size(),B.size());A.resize(len+);
for(rt i=;i<=len;i++)(A[i]-=B[i])%=p;
return A;
}
vector<int>Mod(vector<int>x,vector<int>y){
if(x.size()<=y.size())return x;
vector<int>ans=Div(x,y);
ans=sub(x,Mul(y,ans));
while(!ans[ans.size()-])ans.pop_back();
if(ans.size()>y.size())ans.resize(y.size());
return ans;
}
}
using namespace poly;
int a[];
vector<int>fmo;
vector<int>ksm(vector<int>x,int y){
if(y==)return x;
vector<int>ans=Mod(sqr(ksm(x,y>>)),fmo);
if(y&){
ans.push_back();
for(rt i=ans.size()-;i>=;i--)ans[i+]=ans[i],ans[i]=;
}
return ans;
}
using namespace poly;
int k,m,n,x,y,z,cnt,ans,K,P;
int f[][],mi[];
//f[i][j]长度为i合法高度最低至少为j的合法概率
void calc(int n,int K,int fla){
memset(f,,sizeof(f));
for(rt i=;i<=K+;i++)f[][i]=;
for(rt i=;i<=K+;i++)
for(rt j=K/i;j>=;j--){
f[i][j]=1ll*f[i][j+]*mi[i]%p;
for(rt k=;k<=i;k++)(f[i][j]+=1ll*f[k-][j+]*mi[k-]%p*(+p-P)%p*f[i-k][j]%p)%=p;
}
fmo.resize(K+);
for(rt j=;j<=K+;j++)fmo[K+-j]=-1ll*mi[j-]*f[j-][]%p*(p+-P)%p;
fmo[K+]=;int ret=;
if(n<=K+)ret=f[n][];else {
vector<int>x;x.push_back();x.push_back();
x=ksm(x,n);
for(rt i=;i<=K+;i++)(ret+=1ll*f[i][]*x[i]%p)%=p;
}
(ans+=ret*fla)%=p;
}
int main(){
// file("pool");
n=read();K=read();x=read();y=read();P=1ll*x*ksm(y)%p;
mi[]=;for(rt i=;i<=K+;i++)mi[i]=1ll*mi[i-]*P%p;
calc(n,K,);calc(n,K-,-);cout<<(ans+p)%p;
return ;
}