Luogu4457「BJOI2018」治疗之雨
Description
有 \(m+1\)个数：第一个为 \(p\)，最小值为 \(0\)，最大值为 \(n\)；剩下 \(m\) 个都是无穷，没有最小值或最大值。你可以进行任意多轮操作，每轮操作如下：
在不为最大值的数中等概率随机选择一个（如果没有则不操作），把它加一；
进行 \(k\) 次这个步骤：在不为最小值的数中等概率随机选择一个（如果没有则不操作），把它减一。
现在问期望进行多少轮操作以后第一个数会变为最小值 \(0\)。
Solution
按照分手是祝愿的做法：定义 \(f_i\) 为剩下 \(i\) 滴血的期望轮数
转移考虑每次的增减方式，因为只需要考虑第一个数的值，所以
这东西可以用高斯消元来做，把式子一推就是个没有分数的式子了
然后发现上面的式子是有问题的：没有考虑可以是任意 \(k\) 个地方掉血，所以正经的式子可以写成：
这题目只想到了这里，然后发现这个式子高斯消元完和初始量就没关系了，初始化的细节还是不太会想
发现这是个概率期望题目，所以状态定义是有问题的：应该考虑逆着做，就是定义 \(f_i\) 为有 \(i\) 滴血被杀死的期望轮数,\(g_i\) 表示一轮里面打出伤害大于等于 \(i\) 的概率,式子大概差不多：
注意 \(f_n\) 就不考虑后面一半了
然后发现这个式子有个有趣的地方：\(f_{i-j}+1\)
那么后面的概率部分可以拿出来，就能跟前面的 \(g_i\) 合并了
然后长成这个样子：
其中 \(P_j\) 表示一轮内打掉 \(k\) 滴血的概率
然后直接答案就是 \(f_p\) 就做完了
时间复杂度 \(O(n^3)\) 期望得分：\(60\to 80\ pts\)
下面是卡过的方法：
观察矩阵的形式：每行都是一些连续有值的和一些连续没有值的
然后我们发现这东西比较好的在于我们可以每次直接把 \(a_{i,i},a_{i,i+1}\) 消成 \(1\)
顺着仅更新下面行的 \(a_{j,i+1},a_{j,n}\) 就好了
\(upd:orz\ skyh\)
因为这个矩阵有值的部分是一个直角梯形，那么可以考虑让每个 \(x_i=a\times x_1 +b\)
然后全怼到第 \(n\) 行里面消元，再反解其他的就完成了
这样得到的矩阵是一个类似于这样的
最后一行没有多出来的一个有数的值，所以可以直接计算出答案，最后返回去带入就做完了
Code
#include<bits/stdc++.h>
using namespace std;
#define reg register
#define For(i,a,b) for(reg int i=a;i<=b;++i)
#define Down(i,a,b) for(reg int i=a;i>=b;--i)
namespace yspm{
inline int read()
{
int res=0,f=1; char k;
while(!isdigit(k=getchar())) if(k=='-') f=-1;
while(isdigit(k)) res=res*10+k-'0',k=getchar();
return res*f;
}
const int mod=1e9+7;
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
inline int del(int x,int y){return x-y<0?x-y+mod:x-y;}
inline int mul(int x,int y){return 1ll*x*y-1ll*x*y/mod*mod;}
inline int min(int x,int y){return x<y?x:y;}
inline int max(int x,int y){return x>y?x:y;}
inline int ksm(int x,int y)
{
int res=1;
for(;y;y>>=1,x=mul(x,x)) if(y&1) res=mul(res,x);
return res;
}
const int N=1510;
int inv[N],ans[N],n,p,m,k,P[N],a[N][N];
inline int calc(int mx,int now,int d)
{
int ans=0;
while(now>0)
{
if(now<mx) now++;
now-=d; ++ans;
}return ans;
}
int cnt;
inline void work()
{
n=read(); p=read(); m=read(); k=read();
if(!k) return puts("-1"),void();
if(k==1&&!m) return puts("-1"),void();
if(!m) return printf("%d\n",calc(n,p,k)),void();
memset(a,0,sizeof(a));
int now=1,t1=ksm(m+1,mod-2);
For(i,0,n) P[i]=0,ans[i]=0;
For(i,0,min(n,k)) P[i]=mul(now,mul(ksm(t1,k),ksm(m,k-i))),now=mul(now,mul(inv[i+1],k-i));
For(i,1,n-1)
{
a[i][i]=1; a[i][n+1]=1;
For(j,0,i-1) a[i][i-j]=del(a[i][i-j],mul(P[j],mul(m,t1)));
For(j,0,i) a[i][i-j+1]=del(a[i][i-j+1],mul(P[j],t1));
}
a[n][n]=1; a[n][n+1]=1;
For(i,0,n-1) a[n][n-i]=del(a[n][n-i],P[i]);
For(i,1,n)
{
int tmp=ksm(a[i][i],mod-2); a[i][i]=1; a[i][n+1]=mul(a[i][n+1],tmp);
if(i!=n) a[i][i+1]=mul(a[i][i+1],tmp);
For(j,i+1,n)
{
tmp=a[j][i],a[j][i]=0;
a[j][i+1]=del(a[j][i+1],mul(tmp,a[i][i+1]));
a[j][n+1]=del(a[j][n+1],mul(tmp,a[i][n+1]));
}
}ans[n]=a[n][n+1];
Down(i,n-1,1) ans[i]=del(a[i][n+1],mul(a[i][i+1],ans[i+1]));
return printf("%d\n",ans[p]),void();
}
signed main()
{
inv[0]=inv[1]=1;
For(i,2,N-1) inv[i]=del(mod,mul(mod/i,inv[mod%i]));
int T=read(); while(T--) work();
return 0;
}
}
signed main(){return yspm::main();}
