bzoj4872 [Shoi2017]分手是祝愿 （期望概率DP）

bzoj4872 [Shoi2017]分手是祝愿

原题地址：http://www.lydsy.com/JudgeOnline/problem.php?id=4872

题意：
给出n个灯和它们初始的开关状态，每次操作若选择灯i，则所有i的约数都改变状态（包括i和1）。我们要通过操作把全部灯关掉。每次随机选择一个灯，如果当前最优策略操作数≤k直接用最优策略。问期望操作数*n! %100003的值。

数据范围
1 ≤ n ≤ 100000, 0 ≤ k ≤ n
%50的数据满足k==n。

题解：
因为只想到需要表示状态的DP没想到还有这种操作…
这个k==n的部分数据是一个提示：
在这种情况下，最优方案就是从大到小，如果是亮的就对他操作。
为什么？ 最大的亮的灯如果不操作它本身，只能被其倍数操作，而其倍数就是多按的还要弄灭，以此类推，只会增加多余操作，最后还是要按这个灯本身。
故，对于一个开关状态，其最优方案是固定的唯一的。

于是产生了我们的状态定义：
f[i]表示最优方案是按i个灯 到 最优方案是按i-1个灯的期望操作数。
因为这i个灯是固定的，所以有i/n的概率按对，就直接转过来了，
(n-i)/n的概率按错，赚到了f[i+1]
于是：
$f[i]=\frac{i}{n}+\frac{n-i}{n}(1+f[i+1]+f[i])$
移项：
$f[i]=\frac{(n+(n-i)*f[i+1])}{i}$
发现f[n]=1，于是可以直接推下来，
若给出的状态本来的最优方案是$cnt$，
答案就是，$k+\sum_{i=k+1}^{cnt} f[i]$

---

其实一般不好想到DP维护这个差值，
一般的想法$dp[i]$表示从 最优方案是按i个灯 到 没有灯的期望操作，
$dp[i]=\frac{i}{n}*dp[i-1]+\frac{n-i}{n}dp[i+1]+1$
这东西…
初始化就是$f[k]=k，f[i]=f[i-1]+1$
但是发现移项可以得到更简便的东西，
于是才想到这个差分数组f。

---

猜原来可能是模任意数，所以才要乘n!，刚好对应式子里的分母，这样就不用担心逆元的问题。
结果100003是个质数。

---

代码：

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
using namespace std;
const int mod=100003;
const int N=100005;
int n,k,inv[N],a[N],ans=0;
int f[N],cnt=0;
int main()
{
    scanf("%d%d",&n,&k);
    for(int i=1;i<=n;i++) scanf("%d",&a[i]);
    inv[0]=1; inv[1]=1; 
    for(int i=2;i<=n;i++) inv[i]=1LL*(mod-(mod/i))*inv[mod%i]%mod;
    for(int i=n;i>=1;i--) if(a[i])
    {
        cnt++;
        for(int j=1;j*j<=i;j++) 
        if(i%j==0) { a[j]^=1; if(i/j!=j) a[i/j]^=1; }
    }
    if(cnt<=k) ans=cnt;
    else
    {
        f[n+1]=0; f[n]=1;
        for(int i=n-1;i>k;i--) f[i]=1LL*inv[i%mod]*((n+1LL*(n-i)*f[i+1]%mod)%mod)%mod;
        for(int i=cnt;i>k;i--) ans=(ans+f[i])%mod;
        ans=(ans+k)%mod;
    }
    for(int i=1;i<=n;i++) ans=(1LL*ans*i)%mod;
    printf("%d\n",ans);
    return 0;
}