XVI.[BJOI2018]治疗之雨
一眼能看出这是道高斯消元题。
我们设\(f_i\)表示当前英雄血量为\(i\)时期望多少次死掉。
则我们有
\[f_i=\dfrac{1}{m+1}\times\Big(\sum\limits_{j=0}^iq_jf_{i+1-j}\Big)+\dfrac{m}{m+1}\times\Big(\sum\limits_{j=0}^{i-1}q_jf_{i-j}\Big)+1
\]
其中\(q_j\)是打掉\(j\)点血的概率，有
\[q_j=\Big(\dfrac{1}{m+1}\Big)^j\Big(\dfrac{m}{m+1}\Big)^{k-j}\dbinom{k}{j}
\]
特殊地，当满血（即\(i=n\)）时，因为必定奶不到英雄，故直接有
\[f_i=\sum\limits_{j=0}^{i-1}q_jf_{i-j}+1
\]
明显这个DP有后效性，必须得高斯消元；但是\(O(n^3)\)你告诉我能跑\(1500\)？
抱歉，还真能。
我们列出转移矩阵，发现它大概长这样：
\[\begin{bmatrix}&?,?,0,0,0,\dots,0\\&?,?,?,0,0,\dots,0\\&?,?,?,?,0,\dots,0\\&\dots\\&?,?,?,?,?,\dots,?\end{bmatrix}
\]
发现它十分接近倒着的上三角矩阵；故我们只需要消掉对角线上面的一行，就能把它消成一个真正的倒置上三角矩阵，然后就能\(O(n^2)\)代回解出了。
明显只消掉对角线上面一行的复杂度也是\(O(n^2)\)的；故总复杂度即为\(O(n^2)\)。
代码：
#include<bits/stdc++.h>
using namespace std;
const int mod=1e9+7;
int T,n,m,S,p,g[2010][2010],q[2010],inv[2010],f[2010];
int ksm(int x,int y){
int z=1;
for(;y;y>>=1,x=1ll*x*x%mod)if(y&1)z=1ll*x*z%mod;
return z;
}
void Gauss(){
for(int i=n;i>1;i--){
if(!g[i][i])swap(g[i],g[i-1]);
int tmp=1ll*g[i-1][i]*ksm(g[i][i],mod-2)%mod;
for(int j=1;j<=i;j++)(g[i-1][j]+=mod-1ll*g[i][j]*tmp%mod)%=mod;
(g[i-1][n+1]+=mod-1ll*g[i][n+1]*tmp%mod)%=mod;
}
// for(int i=1;i<=n;i++){for(int j=1;j<=n+1;j++)printf("%d ",g[i][j]);puts("");}
f[1]=1ll*g[1][n+1]*ksm(g[1][1],mod-2)%mod;
for(int i=2;i<=n;i++){
f[i]=g[i][n+1];
for(int j=1;j<i;j++)(f[i]+=mod-1ll*g[i][j]*f[j]%mod)%=mod;
f[i]=1ll*f[i]*ksm(g[i][i],mod-2)%mod;
}
}
int main(){
scanf("%d",&T);
for(int i=1;i<=2000;i++)inv[i]=ksm(i,mod-2);
while(T--){
scanf("%d%d%d%d",&n,&S,&m,&p);
if(p==0||p==1&&m==0&&(S!=n||n!=1)){puts("-1");continue;}
int P=ksm(m+1,mod-2),Q=1ll*m*P%mod;
for(int i=0;i<=n;i++)q[i]=0;
for(int i=0,j=1;i<=min(n,p);j=1ll*j*(p-i)%mod,i++,j=1ll*j*inv[i]%mod)q[i]=1ll*ksm(P,i)*ksm(Q,p-i)%mod*j%mod;
for(int i=1;i<=n;i++)for(int j=1;j<=n+1;j++)g[i][j]=0;
for(int i=1;i<=n;i++){
for(int j=0;j<i;j++)g[i][i-j]=1ll*q[j]*(i<n?mod+1-P:1)%mod;
if(i<n)for(int j=0;j<i+1;j++)(g[i][i+1-j]+=1ll*q[j]*P%mod)%=mod;
(g[i][i]+=mod-1)%=mod;
g[i][n+1]=mod-1;
}
// for(int i=1;i<=n;i++){for(int j=1;j<=n+1;j++)printf("%d ",g[i][j]);puts("");}
Gauss();
printf("%d\n",f[S]);
}
return 0;
}