【bzoj4872】[Shoi2017]分手是祝愿 【概率与期望dp】

题目链接
题解：首先原题有50%的数据n=k。我们只需要从大往小依次关就好。
否则我们先求一遍至少要关多少次，如果<=k就直接输出，否则进行dp。
我们让f[i]$f[i]$表示现在还有i个灯开着，全部关掉的期望次数。
可以得到状态转移方程：
f[i]=in∗f[i−1]+n−in∗f[i+1]+1$f[i]=\frac{i}{n}\ast f[i-1]+\frac{n-i}{n}\ast f[i+1]+1$
然而这个式子不是特别好dp，不过也可以用一个O(n)的神秘算法消一遍元解方程。
我们有一个更简单的做法。
令g[i]=f[i]−f[i−1]$g[i]=f[i]-f[i-1]$
则f[i]=∑ij=1g[j]$f[i]=\sum _{j=1}^{i}g[j]$
把这个式子带入上面的dp方程
∑ij=1g[j]=in∗∑i−1j=1g[j]+n−in∗∑i+1j=1g[j]+1$\sum _{j=1}^{i}g[j]=\frac{i}{n}\ast \sum _{j=1}^{i-1}g[j]+\frac{n-i}{n}\ast \sum _{j=1}^{i+1}g[j]+1$
∑ij=1g[j]=∑i−1j=1g[j]+n−in∗(g[i]+g[i+1])+1$\sum _{j=1}^{i}g[j]=\sum _{j=1}^{i-1}g[j]+\frac{n-i}{n}\ast (g[i]+g[i+1])+1$
g[i]=n−in∗(g[i]+g[i+1])+1$g[i]=\frac{n-i}{n}\ast (g[i]+g[i+1])+1$
ing[i]=n−ing[i+1]+1$\frac{i}{n}g[i]=\frac{n-i}{n}g[i+1]+1$
i×g[i]=(n−i)×g[i+1]+n$i\times g[i]=(n-i)\times g[i+1]+n$
g[i]=(n−i)×g[i+1]+ni$g[i]=\frac{(n-i)\times g[i+1]+n}{i}$
然后g[1..k]=1$g[\mathrm{1..}k]=1$，因为f[i=1..k]=i$f[i=\mathrm{1..}k]=i$
边界：g[n]=1$g[n]=1$，由dp方程易知。
我们从后往前dp一次就行了。
最后不要忘记乘上阶乘。
代码:

  #include<cstdio>
#include<algorithm>
using namespace std;
const int N=100005,mod=100003;
int n,k,tot,ans,jc=1,a[N],f[N];
int fastpow(int a,int x){
    int res=1;
    while(x){
        if(x&1){
            res=1LL*res*a%mod;
        }
        x>>=1;
        a=1LL*a*a%mod;
    }
    return res;
}
int main(){
    scanf("%d%d",&n,&k);
    for(int i=1;i<=n;i++){
        scanf("%d",&a[i]);
        jc=1LL*jc*i%mod;
    }
    for(int i=n;i>=1;i--){
        if(a[i]){
            for(int j=1;j*j<=i;j++){
                if(i%j==0){
                    a[j]^=1;
                    if(j*j!=i){
                        a[i/j]^=1;
                    }
                }
            }
            tot++;
        }
    }
    if(tot<=k){
        printf("%lld
",1LL*tot*jc%mod);
        return 0;
    }
    f[n]=1;
    for(int i=n-1;i>k;i--){
        f[i]=(1LL*(n-i)*f[i+1]%mod+n)%mod*fastpow(i,mod-2)%mod;
    }
    for(int i=1;i<=tot;i++){
        ans+=i<=k?1:f[i];
        ans%=mod;
    }
    printf("%lld
",1LL*ans*jc%mod);
    return 0;
}