【集训队作业2018】喂鸽子
该文不被密码保护。
【集训队作业2018】喂鸽子 题解
算法标签
概率与期望, 生成函数, EGF, 二项式定理, 牛顿二项式定理, 组合计数, 多项式, min-max 容斥
分析题面
有 \(n\) 只鸽子，现在要喂它们玉米粒，每次随机选择一只鸽子喂一粒。一只鸽子吃了至少 \(k\) 粒才能吃饱，问期望多久才能喂饱所有 鸽子。
\(n \le 50, k \le 1000\)
问题分析
法一
因为求的是最晚一只鸽子被喂饱时的期望，不难想到 \(\min - \max\) 容斥
网上题解大多是这种，可以自行查阅
法二
对于这道题，pb 讲了另外一种的推式子方式
先把期望转化为概率
这里有一个结论：
\[\begin{aligned}
E &= \sum_i i \cdot P(X = i)\\
&= P(X = 1) + 2 P(X = 2) + 3 P(X = 3) + \dots\\
&= (P(X = 1) + P(X = 2) + P(X = 3) + \dots) + (P(X = 2) + P(X = 3) + \dots) + \dots\\
&=\sum_i P(X \ge i)
\end{aligned}
\]
设 \(P(X = i)\) 为喂 \(i\) 次喂饱所有鸽子的概率，那么原问题所求即可转化为
\[\begin{aligned}
E &=\sum_t P(X \ge t) \\
\end{aligned}
\]
而 \(P(X \ge t)\) 的意义是喂 \(t - 1\) 次没喂饱的概率，整体平移一下，原问题即求喂 \(t\) 次没喂饱所有鸽子的概率之和
根据一只鸽子喂饱的方案数的 EGF 为
\[\sum_{i = K} ^ {\infty} \frac {x^i}{i!}
\]
所以 \(n\) 只鸽子喂 \(t\) 次喂饱的方案数的 EGF 为
\[(\sum_{i = K} ^ {\infty} \frac {x^i}{i!})^n
\]
可以理解为包含 \(n\) 种元素的多重集的排列数，其中每种元素都有无穷个，但要求每种元素都出现至少 \(K\) 次
由于喂 \(t\) 次的总方案数是 \(n^t\) ，我们就可以得出相应的概率，那我们的答案就是
\[E = \sum_t 1 - \frac{t!}{n^t} [x^t](\sum_{i = K} ^ {\infty} \frac {x^i}{i!})^n
\]
由于这个多项式是无限的，所以直接做是做不了的
但又由于最低项系数的特殊性，我们可以将这个无限项的多项式拆了
令
\[F(x) = \sum_{i = 0} ^ {K - 1} \frac {x^i} {i!}
\]
那么
\[\sum_{i = K} ^ {\infty} \frac {x^i}{i!} = e^x - F(x)
\]
于是
\[\begin{aligned}
E &= \sum_t 1 - \frac{t!}{n^t} [x^t](\sum_{i = K} ^ {\infty} \frac {x^i}{i!})^n\\
&= \sum_t 1 - \frac{t!}{n^t} [x^t](e^x - F(x))^n\\
&= \sum_t 1 - \frac{t!}{n^t} [x^t] \sum_{i = 0} ^ n \binom n i e^{ix} (-1)^{n - i} F^{n - i}(x)\\
&= \sum_t 1 + \sum_{i = 0} ^ n \binom n i (-1)^{n - i + 1} \frac{t!}{n^t} [x^t] [e^{ix} F^{n - i}(x)]\\
\end{aligned}
\]
又因为，当 \(i = n\) 时
\[\begin{aligned}
&\binom n i (-1)^{n - i + 1} e^{ix} F^{n - i}(x) = 1 \cdot (-1) \cdot e^{nx} \cdot 1 = -e^{nx}\\
\\
&\frac {t!}{n^t}[x^t] (-e^{nx}) = -\frac {t!}{n^t} \cdot \frac {n^t}{t!} = -1
\end{aligned}
\]
于是
\[E = \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1}\sum_t \frac{t!}{n^t} [x^t] [e^{ix} F^{n - i}(x)]
\]
注意到 \(F^{n - i}(x)\) 是有限项多项式
设 \(f_{i, j}\) 表示 \([x^j]F^i(x)\)，这个可以通过 \(n\) 次 NTT 求出
于是可以将 \(t\) 次项系数利用加法卷积形式展开
\[\begin{aligned}
E &= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_t \sum_j \frac{t!}{n^t} ([x^{t - j}]e^{ix}) \cdot ([x^j]F^{n - i}(x))\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_t \sum_j \frac{t!}{n^t} \frac {i^{t - j}}{(t - j)!} \cdot f_{n - i,j}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \sum_t \frac{t!}{n^t} \frac {i^{t - j}}{(t - j)!}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot j! \sum_t \frac{t!}{j!(t - j)!}\frac {i^{t - j}}{n^t}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot j! \sum_t \binom{t}{j} \frac {i^{t - j}}{n^{t - j} n^j}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot \frac {j!}{n^j} \sum_t \binom{t}{t - j} (\frac {i}{n})^{t - j}\\
\end{aligned}
\]
根据牛顿二项式定理
\[\frac {1}{(1 - x) ^ n} = \sum_{i = 0} ^ {\infty}\binom {n + i - 1}{i} x^i
\]
发现和答案式子中后面那一堆有点像
\[\begin{aligned}
E &= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot \frac {j!}{n^j} \frac {1}{(1 - \frac i n) ^{j + 1}}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot \frac {j!}{n^j \cdot \frac{(n - i)^{j + 1}}{n^{j + 1}}}\\
&= \sum_{i = 0} ^ {n - 1} \binom n i (-1)^{n - i + 1} \sum_j f_{n - i,j} \cdot \frac {nj!}{(n - i)^{j + 1}}\\
\end{aligned}
\]
由于 \(f_i\) 最高为 \(n(K - 1)\) 次多项式，所以 \(j \le n(K - 1)\)
所以时间复杂度为 \(\mathcal O(n^2K\log(nK) + n^2K)\)
法三
还有 \(\mathcal O(n^2K)\) 的做法，但是我不会/kk
实现细节
NTT 空间开 \(\mathcal O(nK)\) 的，注意阶乘逆元预处理的范围
代码实现
const int N = 52, K = 1e3 + 3, inf = 0x3f3f3f3f, P = 998244353;
int fac[N * K], inv[K];
inline int C(int x, int y) {return 1ll * fac[x] * inv[y] % P * inv[x - y] % P;}
int f[N][N * K];
inline int ksm(int a, int b)
{
int res = 1;
while (b)
{
if (b & 1) res = 1ll * res * a % P;
a = 1ll * a * a % P;
b >>= 1;
}
return res % P;
}
const int g = 3, ig = ksm(3, P - 2);
int rev[66000];
inline void NTT(int *A, const int siz, const bool f)
{
for (re int i = 1; i < siz; ++i) if (i < rev[i]) swap(A[i], A[rev[i]]);
for (re int len = 1; len < siz; len <<= 1)
{
const int rlen = len << 1, wn = ksm(f ? ig : g, (P - 1) / rlen);
for (re int i = 0; i < siz; i += rlen)
{
re int w = 1;
for (re int j = 0; j < len; ++j, w = 1ll * w * wn % P)
{
const int x = A[i + j], y = 1ll * w * A[i + j + len] % P;
A[i + j] = (x + y) % P, A[i + j + len] = (x - y + P) % P;
}
}
}
if (f)
{
const int inv = ksm(siz, P - 2);
for (re int i = 0; i < siz; ++i) A[i] = 1ll * A[i] * inv % P;
}
}
int A[66000], B[66000];
inline void Mul(int *H, int *F, int *G, int n, int m)
{
int siz = 1, L = 0;
while (siz <= n + m) siz <<= 1, ++L;
for (re int i = 1; i < siz; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << L - 1);
for (re int i = 0; i < siz; ++i) A[i] = B[i] = 0;
for (re int i = 0; i <= n; ++i) A[i] = F[i];
for (re int i = 0; i <= m; ++i) B[i] = G[i];
NTT(A, siz, 0), NTT(B, siz, 0);
for (re int i = 0; i < siz; ++i) A[i] = 1ll * A[i] * B[i] % P;
NTT(A, siz, 1);
for (re int i = 0; i <= n + m; ++i) H[i] = A[i];
}
signed main()
{
ios::sync_with_stdio(false);
cin.tie(nullptr);
int n, k; cin >> n >> k;
const int lim = n * (k - 1);
int up = max(lim, max(n, k));
fac[0] = 1;
for (re int i = 1; i <= up; ++i) fac[i] = 1ll * fac[i - 1] * i % P;
up = max(n, k);
inv[up] = ksm(fac[up], P - 2);
for (re int i = up - 1; ~i; --i) inv[i] = 1ll * inv[i + 1] * (i + 1) % P;
for (re int i = 0; i < k; ++i) f[1][i] = inv[i];
for (re int i = 2; i <= n; ++i)
Mul(f[i], f[i - 1], f[1], (k - 1) * (i - 1), k - 1);
for (re int i = 0; i <= lim; ++i) cout << 1ll * f[n][i] * ksm(ksm(n, i), P - 2) % P << ' ';
int ans = 0;
for (re int i = 0; i < n; ++i)
{
int sum = 0;
for (Re int j = 0; j <= lim; ++j)
sum = (sum + 1ll * f[n - i][j] * fac[j] % P * n % P * ksm(ksm(n - i, j + 1), P - 2) % P) % P;
sum = 1ll * sum * C(n, i) % P;
if (n - i + 1 & 1) ans = (ans - sum + P) % P;
else ans = (ans + sum) % P;
}
cout << ans << '\n';
return 0;
}
总结
依旧是利用二项式定理化无穷为 \(\log n\) (快速幂)
浙公网安备 33010602011771号