Luogu P4705 玩游戏

求自然数等幂和的几乎模板题。

题意：从 $a_1, a_2, \cdots, a_n$ 和 $b_1, b_2\cdots, b_m$ 中任意选数 $x$ 和 $y$，求 $(x + y) ^ k$ 的期望。你需要求出 $k\in [1, t]$ 的所有答案。$n, m, t\leq 10 ^ 5$，3s。

首先大力拆贡献，可以得到：

$$\begin{aligned} ans_k =& \dfrac 1{nm} \sum_{i = 1} ^ n\sum_{j = 1} ^ m (a_i + b_j) ^ k\\ =& \dfrac 1{nm} \sum_{i = 1} ^ n \sum_{j = 1} ^ m \sum_{l = 0} ^ k \binom kl a_i ^ l b_j ^ {k - l}\\ =& \dfrac 1{nm} \sum_{l = 0} ^ k \binom kl \sum_{i = 1} ^ n a_i ^ l\sum_{j = 1} ^ m b_j ^ {k - l}\\ =& \dfrac {k!}{nm} \sum_{l = 0} ^ k \dfrac{\sum_{i = 1} ^ n a_i ^ l}{l!} \dfrac{\sum_{j = 1} ^ m b_j ^ {k - l}}{(k - l)!} \end{aligned}$$

容易发现如果我们已经得到了 $\sum_{i = 1} ^ n a_i ^ l$，卷积一下即可在 $O(n\log n)$ 时间内完成。于是问题变为了快速求出 $\sum_{i = 1} ^ n a_i ^ l$。

考虑生成函数 $F(x) = \prod_{i = 1} ^ n (a_ix + 1)$，容易发现这个可以使用分治在 $O(n\log ^ 2 n)$ 的时间内求出。

考虑对两边求 $\ln$，可以得到：$\ln F(x) = \sum_{i = 1} ^ n \ln (a_ix + 1)$。

考虑对 $\ln(1 + x)$ 泰勒展开，可以得到 $\ln(1 + x) = x - \dfrac{x ^ 2}2 + \dfrac{x ^ 3}3 - \dfrac{x ^ 4}{4}+\cdots$。于是我们带入即可得到：

$$\begin{aligned} \ln F(x) =& -\sum_{i = 1} ^ n \sum_{j = 1}\dfrac{(-a_ix) ^ j}{j}\\ =& -\sum_{j = 1} ^ n \dfrac{(-x) ^ j}{j} \sum_{i = 1} ^ n a_i ^ j \end{aligned}$$

容易发现后半段就是我们要求的东西，于是对 $F(x)$ 取 $\ln$，稍加变换即可得到 $\sum_{i = 1} ^ n a_i ^ j$，于是时间复杂度 $O(n\log ^ 2n)$，稍加卡常即可通过。注意上面的式子无法计算 $\sum_{i = 1} ^ n a_i ^ 0$，需要特殊计算。

poly calcpow(poly a, int m)
{
	auto solve = [&](auto &self, poly &a, int l, int r) -> poly {
		if (l == r) return {1, a[l]};
		int mid = (l + r) >> 1;
		return self(self, a, l, mid) * self(self, a, mid + 1, r);
	};
	auto b = solve(solve, a, 0, a.size() - 1);
	b.resize(m + 1), b = Ln(b);
	for (int i = 0; i <= m; ++ i)
	{
		if (!(i & 1)) adj(b[i] = -b[i]);
		b[i] = (LL) b[i] * i % Mod;
	}
	return b;
}

int main()
{
	init();
	int n, m, t;
	std::cin >> n >> m;
	poly a(n), b(m);
	for (int &x : a) scanf("%d", &x);
	for (int &x : b) scanf("%d", &x);
	std::cin >> t;
	a = calcpow(a, t), a[0] = n;
	for (int i = 0; i <= t; ++ i) a[i] = (LL) a[i] * infact[i] % Mod;
	b = calcpow(b, t), b[0] = m;
	for (int i = 0; i <= t; ++ i) b[i] = (LL) b[i] * infact[i] % Mod;
	a = a * b;
	int Inv = qpow((LL) n * m % Mod);
	for (int i = 0; i <= t; ++ i) a[i] = (LL) a[i] * fact[i] % Mod * Inv % Mod;
	for (int i = 1; i <= t; ++ i) printf("%d\n", a[i]);
	return 0;
}