【期望】【P5081】Tweetuzki 爱取球

Description

Tweetuzki 有一个袋子，袋子中有 \(N\) 个无差别的球。Tweetuzki 每次随机取出一个球后放回。求取遍所有球的期望次数。

取遍是指，袋子中所有球都被取出来过至少一次。

Input

一行一个整数 \(N\)

Output

一行一个整数，表示期望次数对 \(200403013\) 取模的结果

Hint

\(1~\leq~N~\leq~10^7\)

Solution

期望嘛，就是你xjb乱推好几个看起来很靠谱事实上最多有一个是对的甚至可能一个也不对的式子，然后用蒙特卡洛方法拍一拍，挑一个最靠谱的交上去

感谢@_rqy提供指导

先给出结论：

\[E(n)~=~n~\times~(\sum_{i=1}^{n}~\frac{1}{i}) \]

证明如下：

引理：对于一次试验，有 \(p~(0~\leq~p~\leq~1)\) 的概率成功，有 \(1-p\) 的概率不成功，如果一次试验后不成功则重复进行试验，否则停止。则期望在 \(\frac{1}{p}\) 步后停止试验

证明：

设期望 \(x\) 步后停止试验，根据已知以及期望的线性相加性，可得方程：

\[x~=~1~+~p~\times~0~+~(1~-~p)~\times~x \]

上述方程的直观解释是当前进行完一步后有 \(p\) 的概率再进行 \(0\) 次（停止），\(1~-~p\) 的概率继续进行，由于下次试验和这次试验完全相同，所以有 \(1~-~p\) 的概率期望再进行 \(x\) 次。

解方程，得

\[x~=~\frac{1}{p} \]

引理证毕。

于是考虑摸球，首先有 \(n\) 个球没有被摸到的时候，摸一次摸到新球的概率为 \(1\)，根据引理，期望摸 \(1\) 次。

然后剩下 \(n~-~1\) 个球没有被摸到，摸一次摸到新球的概率为 \(\frac{n-1}{n}\)，根据引理，期望摸 \(\frac{n}{n-1}\) 次。

然后剩下 \(n~-~2\) 个球没有摸到，摸一次摸到新球的概率为 \(\frac{n-2}{n}\) ，根据引理，期望摸 \(\frac{n}{n-2}\) 次。

……（一位读者及时砸坏了复读机，停止了复读）

按照球的个数进行数学归纳，根据期望的可加性（非线性相加性），可以得到摸遍 \(n\) 个球的期望

\[E(n)~=~n~+~\frac{n}{n-1}~+~\frac{n}{n-1}~+~\dots~+~\frac{n}{1} \]

提取公因数 \(n\)，可得

\[E(n)~=~n~\times~(\sum_{i=1}^{n}~\frac{1}{i}) \]

证毕。

然鹅最后一个式子只是为了让结论变得优美一点，事实上用倒数第二个式子就可以做题了233

于是先筛一遍逆元，然后求答案即可。

Code

#include <cstdio>
#ifdef ONLINE_JUDGE
#define freopen(a, b, c)
#endif
#define rg register
#define ci const int
#define cl const long long

typedef long long int ll;

namespace IPT {
	const int L = 1000000;
	char buf[L], *front=buf, *end=buf;
	char GetChar() {
		if (front == end) {
			end = buf + fread(front = buf, 1, L, stdin);
			if (front == end) return -1;
		}
		return *(front++);
	}
}

template <typename T>
inline void qr(T &x) {
	rg char ch = IPT::GetChar(), lst = ' ';
	while ((ch > '9') || (ch < '0')) lst = ch, ch=IPT::GetChar();
	while ((ch >= '0') && (ch <= '9')) x = (x << 1) + (x << 3) + (ch ^ 48), ch = IPT::GetChar();
	if (lst == '-') x = -x;
}

template <typename T>
inline void ReadDb(T &x) {
	rg char ch = IPT::GetChar(), lst = ' ';
	while ((ch > '9') || (ch < '0')) lst = ch, ch = IPT::GetChar();
	while ((ch >= '0') && (ch <= '9')) x = x * 10 + (ch ^ 48), ch = IPT::GetChar();
	if (ch == '.') {
		ch = IPT::GetChar();
		double base = 1;
		while ((ch >= '0') && (ch <= '9')) x += (ch ^ 48) * ((base *= 0.1)), ch = IPT::GetChar();
	}
	if (lst == '-') x = -x;
}

namespace OPT {
	char buf[120];
}

template <typename T>
inline void qw(T x, const char aft, const bool pt) {
	if (x < 0) {x = -x, putchar('-');}
	rg int top=0;
	do {OPT::buf[++top] = x % 10 + '0';} while (x /= 10);
	while (top) putchar(OPT::buf[top--]);
	if (pt) putchar(aft);
}

const int maxn = 10000010;
const int MOD = 20040313;

int n, ans;
int inv[maxn];

void GetInv();

int main() {
	freopen("1.in", "r", stdin);
	qr(n);
	GetInv();
	for (rg int i = 0; i < n; ++i) {
		ans = (ans + 1ll * n * inv[n - i]) % MOD;
	}
	qw(ans, '\n', true);
}

void GetInv() {
	inv[1] = 1;
	for(int i = 2; i <= n; ++i) inv[i] = 1ll * (MOD - MOD / i) * inv[MOD % i] % MOD;
}

Summary

对于一次试验，有 \(p~(0~\leq~p~\leq~1)\) 的概率成功，有 \(1-p\) 的概率不成功，如果一次试验后不成功则重复进行试验，否则停止。则期望在 \(\frac{1}{p}\) 步后停止试验