洛谷P4457/loj#2513 [BJOI2018]治疗之雨(高斯消元+概率期望)

题面

传送门(loj)

传送门(洛谷)

题解

模拟赛的时候只想出了高斯消元然后死活不知道怎么继续……结果正解居然就是高斯消元卡常?

首先有个比较难受的地方是它一个回合可能不止扣一滴血……我们得算出\(P_i\)表示一回合扣\(i\)滴血的概率,为

\[P_i={{k\choose i}m^{k-i}\over (m+1)^k}
\]

所以这个柿子啥意思?

我们可以把\(k\)次扣血看成一个长度为\(k\)的序列,每个序列有\(m+1\)种选择方法,于是总的选法就是\((m+1)^k\)。我们要钦定它扣\(i\)滴血,就是令其中\(i\)个数强制选择\(1\),方案数为\({k\choose i}\)。剩下的数都不能选\(1\),所以方案数为\(m^{k-i}\)

有了这个\(P_i\)我们就可以考虑找关系了

然而这里还有一个很讨厌的地方就是它时不时会被奶一口……

我们记\(f_i\)表示当血量为\(i\)时被干掉的期望回合数,\(g_i\)表示一回合内打出伤害大于等于\(i\)的概率,然后考虑这东西该怎么转移

\[f_i={m\over m+1}\left(g_i+\sum_{k=0}^{i-1}P_k(f_{i-k}+1)\right)+{1\over m+1}\left(g_{i+1}+\sum_{k=0}^{i}P_k(f_{i+1-k}+1)\right)
\]

边界条件为

\[f_n=g_n+\sum_{k=0}^{n-1}P_k(f_{n-k}+1)
\]

所以这柿子是啥?

\({m\over m+1}\)表示没有被奶到的概率,那么我们枚举它被\(A\)了几下。如果它一回合被干掉了,那么期望局数为\(1\),否则\(A\)了\(k\)下之后血量会到\(i-k\),这一部分的期望步数就是\(f_{i-k}+1\)。后面那个就是如果被奶了之后的情况。顺便因为满血的时候是不可能被奶的,所以\(f_n\)要特殊考虑

然而这柿子一点都不好看,特别是\(g_i\)很不爽,那就继续推倒

\[f_i={m\over m+1}\left(g_i+\sum_{k=0}^{i-1}P_k+\sum_{k=0}^{i-1}P_kf_{i-k}\right)+{1\over m+1}\left(g_{i+1}+\sum_{k=0}^{i}P_k+\sum_{k=0}^{i}P_kf_{i+1-k}\right)
\]

\[f_i={m\over m+1}\left(1+\sum_{k=0}^{i-1}P_kf_{i-k}\right)+{1\over m+1}\left(1+\sum_{k=0}^{i}P_kf_{i+1-k}\right)
\]

同理有\(f_n=1+\sum_{k=0}^{n-1}P_kf_{n-k}\)

那么我们就可以愉快地递推……等会儿这咋推啊……转移好像都成环了啊……

那么我们就把它看成一个方程组来高斯消元求解吧

哈?\(n=1500\)你让我高斯消元?

我们来好好观察一下这个方程组,\(f_i\)所对应的方程只与\(f_1,...f_{i+1}\)的值有关,也就是说每一行的对角线右边只有在它右边一位的那个系数不为\(0\)

因为我们高斯消元的时候是拿自己这行去减下面的,所以每一行中只有\(2\)个系数要去和下面的相减,复杂度就能化到\(O(n^2)\)了

虽然我还是不知道这个复杂度是怎么过去的

//minamoto
#include<bits/stdc++.h>
#define R register
#define fp(i,a,b) for(R int i=a,I=b+1;i<I;++i)
#define fd(i,a,b) for(R int i=a,I=b-1;i>I;--i)
#define go(u) for(int i=head[u],v=e[i].v;i;i=e[i].nx,v=e[i].v)
template<class T>inline bool cmin(T&a,const T&b){return a>b?a=b,1:0;}
using namespace std;
char buf[1<<21],*p1=buf,*p2=buf;
inline char getc(){return p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++;}
int read(){
R int res,f=1;R char ch;
while((ch=getc())>'9'||ch<'0')(ch=='-')&&(f=-1);
for(res=ch-'0';(ch=getc())>='0'&&ch<='9';res=res*10+ch-'0');
return res*f;
}
char sr[1<<21],z[20];int K=-1,Z=0;
inline void Ot(){fwrite(sr,1,K+1,stdout),K=-1;}
void print(R int x){
if(K>1<<20)Ot();if(x<0)sr[++K]='-',x=-x;
while(z[++Z]=x%10+48,x/=10);
while(sr[++K]=z[Z],--Z);sr[++K]='\n';
}
const int N=1505,P=1e9+7;
inline int add(R int x,R int y){return x+y>=P?x+y-P:x+y;}
inline int dec(R int x,R int y){return x-y<0?x-y+P:x-y;}
inline int mul(R int x,R int y){return 1ll*x*y-1ll*x*y/P*P;}
int ksm(R int x,R int y){
R int res=1;
for(;y;y>>=1,x=mul(x,x))if(y&1)res=mul(res,x);
return res;
}
int a[N][N],ans[N],s[N],inv[N],g[N];
int n,m,p,k,T,qwq,aa,bb,tmp,res;
void init(int n){
inv[0]=inv[1]=1;
fp(i,2,n)inv[i]=mul(P-P/i,inv[P%i]);
}
void Gauss(int n){
fp(i,1,n-1){
int t=ksm(a[i][i],P-2);
a[i][i]=1,a[i][n]=mul(a[i][n],t);
if(i!=n-1)a[i][i+1]=mul(a[i][i+1],t);
fp(j,i+1,n-1){
t=a[j][i],a[j][i]=0;
a[j][i+1]=dec(a[j][i+1],mul(t,a[i][i+1])),
a[j][n]=dec(a[j][n],mul(t,a[i][n]));
}
}
ans[n-1]=a[n-1][n];
fd(i,n-2,1)ans[i]=dec(a[i][n],mul(a[i][i+1],ans[i+1]));
}
int main(){
// freopen("testdata.in","r",stdin);
T=read();
init(N-5);
while(T--){
n=read(),p=read(),m=read(),k=read();
if(!k||(!m&&k==1)){puts("-1");continue;}
if(!m){
while(p>0){if(p<n)++p;p-=k,++res;}
printf("%d\n",res),res=0;continue;
}
qwq=ksm(m+1,k),qwq=ksm(qwq,P-2)%P,tmp=1;
fp(i,0,n)s[i]=0;
fp(i,0,min(n,k)){
s[i]=1ll*tmp*ksm(m,k-i)%P*qwq%P;
tmp=1ll*tmp*inv[i+1]%P*(k-i)%P;
}
memset(a,0,sizeof(a));
bb=ksm(m+1,P-2),aa=mul(m,bb);
fp(i,1,n-1){
++a[i][i],++a[i][n+1];
fp(j,0,i-1)a[i][i-j]=dec(a[i][i-j],mul(s[j],aa));
fp(j,0,i)a[i][i+1-j]=dec(a[i][i+1-j],mul(s[j],bb));
}
++a[n][n],++a[n][n+1];
fp(j,0,n-1)a[n][n-j]=dec(a[n][n-j],s[j]);
Gauss(n+1);
printf("%d\n",ans[p]);
}
return 0;
}
上一篇:shell脚本例子集锦(习题总结)


下一篇:SSL证书制作