Luogu 4705 玩遊戲
看見這個題依稀想起了$5$月月賽時候的事情，到現在仍然它感覺非常神仙。
遊戲$k$次價值的期望答案
$$ans_k = \frac{1}{nm}\sum_{i = 1}^{n}\sum_{j = 1}^{n}(a_i + b_j)^k$$
二項式定理展開
$$ans_k=\frac{1}{nm}\sum_{i = 1}^{n}\sum_{j = 1}^{m}\sum_{t = 0}^{k}\binom{k}{t}a_i^tb_j^{k - t}$$
$$ = \frac{1}{nm}\sum_{t = 0}^{k}\binom{k}{t}\sum_{i = 1}^{n}a_i^t\sum_{j = 1}^{m}b_j^{k - t}$$
$$= \frac{k!}{nm}\frac{\sum_{i = 1}^{n}a_i^t}{t!}\frac{\sum_{j = 1}^{m}b_j^{k - t}}{(k - t)!}$$
容易發現這後面是一個卷積的形式。
我們設$A(x) = \frac{\sum_{i = 1}^{n}a_i^x}{x!}$，$B(x) = \frac{\sum_{i = 1}^{n}b_i^x}{x!}$，答案就變成了
$$\frac{k!}{nm}(A * B)(k)$$
現在考慮怎麽算這個$\sum_{i = 1}^{n}a_i^k$。
不會了，點開題解
我們構造一個多項式$F(x) = \prod_{i = 1}^{n}(1 - a_ix)$，這個$F(x)$可以通過分治算出來，時間復雜度$T(n) = 2T(\frac{n}{2}) + O(nlogn)$趨向$O(nlogn)$？？？。
接下來要開始變魔術了，
設$G(x) = lnF(x)$，
$$G(x) = \sum_{i = 1}^{n}ln(1 - a_ix)$$
對裏面的$ln$在$x_0 = 0$處泰勒展開。
我們知道
$$ln(1 - x) = 0 - \frac{x}{1} - \frac{x^2}{2} - \frac{x^3}{3} - \cdots = -\sum_{i = 1}^{\infty}\frac{x^i}{i}$$
所以
$$G(x) = \sum_{i = 1}^{n}\sum_{j = 1}^{k}-\frac{a_i^j}{j}x^j$$
因為這個題對$x^{k + 1}$次取模了，所以後面的項可以不用再寫了。
變形一下
$$G(x) = \sum_{j = 1}^{k}(-\frac{1}{j}\sum_{i = 1}^{n}a_i^j)x^j$$
發現$1$到$k$的$sum_{i = 1}^{n}a_i^k$直接算出來了，而$i = 0$時候這東西顯然為$n$。
鼓掌~~~
現在回過頭來考慮一下這個魔術是怎麽變出來的。
考慮構造每一個$a_i$的生成函數
$$1 + a_ix + a_i^2x^2 + a_i^3x^3 + \cdots$$
把每一個$a_i$的生成函數都構造出來然後加起來的第$i$項系數就是$k = i$時候的答案了。
$$F(x) = \sum_{i = 1}^{n}\frac{1}{1 - a_ix}$$
發現這個$F(x)$並不是很好算，而
$$ln‘(1 - a_ix) = \frac{1}{1 - a_ix}$$
$$(ln(1 - a_ix))‘ = \frac{-a_i}{1 - a_ix}$$
所以先計算
$$G(x) = \sum_{i = 1}^{n}\frac{-a_i}{1 - a_ix} = (ln\prod_{i = 1}^{n}(1 - a_ix))‘$$
把這兩個式子寫回到生成函數的形式，發現
$$F(x) = -x * G(x) + n$$
於是大功告成。
用$vector$寫多項式真爽。
#include <cstdio> #include <cstring> #include <algorithm> #include <vector> using namespace std; typedef long long ll; typedef vector <ll> poly; const int N = 1e5 + 5; const int Maxn = 1e5; int n, m, K; ll a[N], b[N], fac[N], ifac[N], inv[N]; template <typename T> inline void read(T &X) { X = 0; char ch = 0; T op = 1; for (; ch > ‘9‘ || ch < ‘0‘; ch = getchar()) if (ch == ‘-‘) op = -1; for (; ch >= ‘0‘ && ch <= ‘9‘; ch = getchar()) X = (X << 3) + (X << 1) + ch - 48; X *= op; } inline void deb(poly c) { for (int i = 0; i < (int)c.size(); i++) printf("%lld%c", c[i], " \n"[i == (int)c.size() - 1]); } namespace Poly { const int L = 1 << 20; const ll P = 998244353LL; int lim, pos[L]; inline ll fpow(ll x, ll y) { ll res = 1; for (; y > 0; y >>= 1) { if (y & 1) res = res * x % P; x = x * x % P; } return res; } template <typename T> inline void inc(T &x, T y) { x += y; if (x >= P) x -= P; } template <typename T> inline void sub(T &x, T y) { x -= y; if (x < 0) x += P; } inline void prework(int len) { int l = 0; for (lim = 1; lim < len; lim <<= 1, ++l); for (int i = 0; i < lim; i++) pos[i] = (pos[i >> 1] >> 1) | ((i & 1) << (l - 1)); } inline void ntt(poly &c, int opt) { c.resize(lim, 0); for (int i = 0; i < lim; i++) if (i < pos[i]) swap(c[i], c[pos[i]]); for (int i = 1; i < lim; i <<= 1) { ll wn = fpow(3, (P - 1) / (i << 1)); if (opt == -1) wn = fpow(wn, P - 2); for (int len = i << 1, j = 0; j < lim; j += len) { ll w = 1; for (int k = 0; k < i; k++, w = w * wn % P) { ll x = c[j + k], y = w * c[j + k + i] % P; c[j + k] = (x + y) % P, c[j + k + i] = (x - y + P) % P; } } } if (opt == -1) { ll inv = fpow(lim, P - 2); for (int i = 0; i < lim; i++) c[i] = c[i] * inv % P; } } inline poly operator * (const poly &x, const poly &y) { poly res, u = x, v = y; prework(u.size() + v.size() - 1); ntt(u, 1), ntt(v, 1); for (int i = 0; i < lim; i++) res.push_back(v[i] * u[i] % P); ntt(res, -1); res.resize(u.size() + v.size() - 1); return res; } poly getInv(poly x, int len) { x.resize(len); if (len == 1) { poly res; res.push_back(x[0]); return res; } poly y = getInv(x, (len + 1) >> 1); prework(len << 1); poly u = x, v = y, res; ntt(u, 1), ntt(v, 1); for (int i = 0; i < lim; i++) res.push_back(v[i] * (2LL - u[i] * v[i] % P + P) % P); ntt(res, -1); res.resize(len); return res; } inline void direv(poly &c) { for (int i = 0; i < (int)c.size() - 1; i++) c[i] = c[i + 1] * (i + 1) % P; c[c.size() - 1] = 0; } inline void integ(poly &c) { for (int i = (int)c.size() - 1; i > 0; i--) c[i] = c[i - 1] * inv[i] % P; c[0] = 0; } inline poly getLn(poly c) { poly a = getInv(c, (int)c.size()); poly b = c; direv(b); poly res = b * a; res.resize(c.size()); integ(res); return res; } inline poly solve(int l, int r, ll *c) { if (l == r) { poly res; res.push_back(1), res.push_back((P - c[l]) % P); return res; } int mid = ((l + r) >> 1); poly u = solve(l, mid, c), v = solve(mid + 1, r, c); poly res = u * v; res.resize(u.size() + v.size() - 1); return res; } inline poly mul(const poly &x, const poly &y) { return x * y; } } using Poly::P; using Poly::fpow; using Poly::inc; using Poly::sub; using Poly::solve; using Poly::getLn; inline void prework() { fac[0] = inv[0] = inv[1] = 1; for (int i = 1; i <= Maxn; i++) { fac[i] = fac[i - 1] * i % P; if (i > 1) inv[i] = (P - P / i) * inv[P % i] % P; } ifac[Maxn] = fpow(fac[Maxn], P - 2); for (int i = Maxn - 1; i >= 0; i--) ifac[i] = ifac[i + 1] * (i + 1) % P; } int main() { /* #ifndef ONLINE_JUDGE freopen("Sample.txt", "r", stdin); #endif */ freopen("input.txt", "r", stdin); freopen("my.out", "w", stdout); prework(); read(n), read(m); for (int i = 1; i <= n; i++) read(a[i]); for (int i = 1; i <= m; i++) read(b[i]); read(K); ++K; poly Fa = solve(1, n, a), Fb = solve(1, m, b); // deb(Fa), deb(Fb); Fa.resize(K, 0), Fb.resize(K, 0); poly Ga = getLn(Fa), Gb = getLn(Fb); Ga.resize(K, 0), Gb.resize(K, 0); // deb(Ga), deb(Gb); poly A, B, C; A.push_back(n), B.push_back(m); for (int i = 1; i < K; i++) { A.push_back((P - Ga[i]) % P * i % P * ifac[i] % P); B.push_back((P - Gb[i]) % P * i % P * ifac[i] % P); } // deb(A), deb(B); C = Poly::mul(A, B); C.resize(K); // deb(C); // printf("\n"); ll invn = fpow(n, P - 2), invm = fpow(m, P - 2); for (int i = 1; i < K; i++) printf("%lld\n", invn * invm % P * fac[i] % P * C[i] % P); return 0; }View Code
Luogu 4705 玩遊戲