BZOJ4872【SHOI2017】分手是祝愿 <概率DP>

Problem

【SHOI2017】分手是祝愿

$\mathrm{Time\;Limit:\;20\;Sec}$
$\mathrm{Memory\;Limit:\;512\;MB}$

Description

Zeit und Raum trennen dich und mich.
时空将你我分开。$\mathrm{B 君}$在玩一个游戏，这个游戏由 $n$ 个灯和 $n$ 个开关组成，给定这 $n$ 个灯的初始状态，下标为从 $1$ 到 $n$ 的正整数。每个灯有两个状态亮和灭，我们用 $1$ 来表示这个灯是亮的，用 $0$ 表示这个灯是灭的，游戏的目标是使所有灯都灭掉。但是当操作第 $i$ 个开关时，所有编号为 $i$ 的约数（包括 $1$ 和 $i$）的灯的状态都会被改变，即从亮变成灭，或者是从灭变成亮。$\mathrm{B 君}$发现这个游戏很难，于是想到了这样的一个策略，每次等概率随机操作一个开关，直到所有灯都灭掉。这个策略需要的操作次数很多， $\mathrm{B 君}$想到这样的一个优化。如果当前局面，可以通过操作小于等于 $k$ 个开关使所有灯都灭掉，那么他将不再随机，直接选择操作次数最小的操作方法（这个策略显然小于等于 $k$ 步）操作这些开关。$\mathrm{B 君}$想知道按照这个策略（也就是先随机操作，最后小于等于 $k$ 步，使用操作次数最小的操作方法）的操作次数的期望。这个期望可能很大，但是 $\mathrm{B 君}$发现这个期望乘以 $n!$ 一定是整数，所以他只需要知道这个整数对 $100003$ 取模之后的结果。

Input

第一行两个整数 $n$, $k$。
接下来一行 $n$ 个整数，每个整数是 $0$ 或者 $1$，其中第 $i$ 个整数表示第 $i$ 个灯的初始情况。
$1\le n\le 100000$, $0\le k\le n$

Output

一行，为操作次数的期望乘以 $n!$ 对 $100003$ 取模之后的结果。

Sample Input

4 0
0 0 1 1

Sample Output

512

Source

黑吉辽沪冀晋六省联考

标签：概率DP

Solution

显然一种可行的最优解方案是从大往小操作，如果当前灯是亮的，则将其及其约数编号的灯状态翻转。
预处理初始状态通过这种状态操作到全部灭掉的次数，记为$cnt$。
令$f[i]$为还剩$i$步全部灭掉时的期望步数，有

$$f[i]= \begin{cases} i&i\le k\\ \frac{i}{n}f[i-1]+\frac{n-i}{n}f[i+1]+1&i>k\\ \end{cases}$$

发现$f$的转移带环。将$f[i]$进行差分，记为$g[i]$，即$g[i]=f[i]-f[i-1]$。那么有

$$g[i]= \begin{cases} 1&i\le k\\ \frac{n-i}{i}g[i+1]+\frac{n}{i}&k<i<n\\ 1&i=n\\ \end{cases}$$

由此可在$O(n)$时间内算出$g[1]\sim g[n]$，累加可知$f[cnt]$。

Code

#include <bits/stdc++.h>
#define P 100003
#define MAX_N 100000
using namespace std;
typedef long long lnt;
template <class T> inline void read(T &x) {
    x = 0; int c = getchar(), f = 1;
    for (; !isdigit(c); c = getchar()) if (c == 45) f = -1;
    for (; isdigit(c); c = getchar()) (x *= 10) += f*(c-'0');
}
int n, m, cnt; lnt ans;
lnt fac[MAX_N+5], inv[MAX_N+5];
bool a[MAX_N+5]; lnt g[MAX_N+5];
void init() {
    fac[0] = inv[0] = inv[1] = 1;
    for (int i = 1; i <= n; i++) fac[i] = fac[i-1]*i%P;
    for (int i = 2; i <= n; i++) inv[i] = (P-P/i*inv[P%i]%P)%P;
}
int main() {
    read(n), read(m), init();
    for (int i = 1; i <= n; i++) read(a[i]);
    for (int i = n; i; cnt += a[i--])
        for (int j = (i<<1); j <= n; j += i) a[i] ^= a[j];
    g[n] = 1; for (int i = 1; i <= m; i++) g[i] = 1;
    for (int i = n-1; i > m; i--)
        g[i] = (1LL*(n-i)*inv[i]%P*g[i+1]%P+1LL*n*inv[i]%P)%P;
    for (int i = 1; i <= cnt; i++) ans = (ans+g[i])%P;
    return printf("%lld\n", ans*fac[n]%P), 0;
}