[CTS2019] 珍珠
Problem
有 \(n\) 个范围在 \([1,D]\) 内的整数均匀随机变量，求至少能选出 \(m\) 个瓶子，使得存在一种方案，选择一些变量，并把选出来的每一个变量放到一个瓶子中，满足每个瓶子都恰好装两个值相同的变量的概率。答案乘上 \(D^n\) 对 \(998244353\) 取模。
\(D\le 10^5\)，\(n,m\le 10^9\)。
Sol
首先根据奇偶性，写出来关于每种颜色的多项式（\(x\) 表示选出来的数目，\(y\) 表示能组合的 pair \(\times 2\)）
\[F(x,y)=\frac{e^{xy}+e^{-xy}}2+\frac{e^{xy}-e^{-xy}}{2y}
\]
答案为
\[Ans=n![x^ny^{\ge2m}]F^D(x,y)
\]
\[\Rightarrow n![x^ny^{\ge2m-n}]\left(\frac{e^x+e^{-x}}2+\frac{e^x-e^{-x}}{2y}\right)^D
\]
\[\Rightarrow n![x^ny^{\ge 2m-n+D}]\left(\frac{e^x+e^{-x}}2y+\frac{e^x-e^{-x}}2\right)^D
\]
无必要的 \(y\) 消掉了。
展开右式
\[\sum_{i=0}^D{D\choose i}y^i\left(\frac{e^x+e^{-x}}2\right)^i\left(\frac{e^x-e^{-x}}2\right)^{D-i}
\]
答案变成
\[Ans=n![x^n]\sum_{i=a}^D{D\choose i}\left(\frac{e^x+e^{-x}}2\right)^i\left(\frac{e^x-e^{-x}}2\right)^{D-i}
\]
这里 \(a=2m-n+D\)。
继续考察右式
\[\begin{aligned}
&\sum_{i=a}^D{D\choose i}\left(\frac{e^x+e^{-x}}2\right)^i\left(\frac{e^x-e^{-x}}2\right)^{D-i}\\
=&2^{-D}e^{-Dx}\sum_{i=a}^D{D\choose i}(e^{2x}+1)^i(e^{2x}-1)^{D-i}
\end{aligned}
\]
这里重点来看后半部分。暂时用 \(x\) 替换 \(e^{2x}\)。
\[\begin{aligned}
\sum_{i=a}^D{D\choose i}(x+1)^i(x-1)^{D-i}
\end{aligned}
\]
它一定可以变成
\[\sum_{i=0}^Da_ix^i
\]
这样子就可以解决问题了，现在我们来尝试化简。这里令 \(z=x+1\)，则式子变成
\[\begin{aligned}
&\sum_{i=a}^D{D\choose i}z^i(z-2)^{D-i}\\
=&\sum_{i=a}^D{D\choose i}\sum_{j=0}^{D-i}{D-i\choose j}(-2)^{D-i-j}z^{i+j}\\
=&\sum_{i=a}^D(-2)^{D-i}z^i\sum_{j=0}^{i-a}{D\choose i-j}{D-(i-j)\choose j}\\
=&\sum_{i=a}^D(-2)^{D-i}z^i{D\choose i}\sum_{j=0}^{i-a}{i\choose j}
\end{aligned}
\]
发现 \(\sum_{j=0}^{i-a}{i\choose j}\) 是可以递推的！于是我们可以在 \(\mathcal O(D)\) 的时间求出上面的多项式，对 \(z\) 做平移得到原来的多项式，问题得解。最后得到
\[Ans=n![x^n]\sum_{i=0}^Df_ie^{(2i-D)x}
\]
直接解系数即可。复杂度 \(\mathcal O(D\log n)\)。
附平移多项式推导
\[\begin{aligned}
f(x+c)&=\sum_{i=0}^na_i(x+c)^i\\
&=\sum_{i=0}^na_i\sum_{j=0}^i{i\choose j}x^jc^{i-j}\\
&=\sum_{j=0}^n\frac{x^j}{j!}\sum_{i=j}^ni!a_i\frac{c^{i-j}}{(i-j)!}
\end{aligned}
\]
NTT 即可。
Code
#include <bits/stdc++.h>
using std::max;
const int N = 100005, P = 998244353;
int D, n, m, a[N * 4], b[N * 4], s;
int inc(int a, int b) { return (a += b) >= P ? a - P : a; }
int qpow(int a, int b) {
int t = 1;
for (; b; b >>= 1, a = 1LL * a * a % P)
if (b & 1) t = 1LL * t * a % P;
return t;
}
int W[N * 4], inv[N], fac[N], ifac[N];
int C(int n, int m) {
if (n < 0 || m > n) return 0;
return 1LL * fac[n] * ifac[m] % P * ifac[n - m] % P;
}
void prework(int n) {
for (int i = 1; i < n; i <<= 1)
for (int j = W[i] = 1, Wn = qpow(3, (P - 1) / i >> 1); j < n; j++)
W[i + j] = 1LL * W[i + j - 1] * Wn % P;
}
void fft(int *a, int n, int op) {
static int rev[N * 4] = {0};
for (int i = 0; i < n; i++)
if ((rev[i] = rev[i >> 1] >> 1 | (i & 1 ? n >> 1 : 0)) > i) std::swap(a[i], a[rev[i]]);
for (int q = 1; q < n; q <<= 1)
for (int p = 0; p < n; p += q << 1)
for (int i = 0, t; i < q; i++)
t = 1LL * W[q+i] * a[p+q+i] % P, a[p+q+i] = inc(a[p+i], P-t), a[p+i] = inc(a[p+i], t);
if (op) return;
for (int i = 0, inv = qpow(n, P - 2); i < n; i++)
a[i] = 1LL * a[i] * inv % P;
std::reverse(a + 1, a + n);
}
int getsz(int n) { int x = 1; while (x < n) x <<= 1; return x; }
int main() {
scanf("%d%d%d", &D, &n, &m);
prework((D + 1) * 2);
inv[1] = fac[0] = ifac[0] = 1;
for (int i = 2; i <= D; i++)
inv[i] = 1LL * (P - P / i) * inv[P % i] % P;
for (int i = 1; i <= D; i++)
fac[i] = 1LL * fac[i - 1] * i % P, ifac[i] = 1LL * ifac[i - 1] * inv[i] % P;
s = 2 * m - n + D;
int tmp = 1;
for (int i = max(0, s); i <= D; i++) {
a[i] = 1LL * qpow(P - 2, D - i) * C(D, i) % P * tmp % P * fac[i] % P;
tmp = (2LL * tmp + C(i, i - s + 1)) % P;
}
for (int i = 0; i <= D; i++)
b[i] = ifac[i];
std::reverse(b, b + D + 1);
int l = getsz(2 * (D + 1));
fft(a, l, 1), fft(b, l, 1);
for (int i = 0; i < l; i++) a[i] = 1LL * a[i] * b[i] % P;
fft(a, l, 0);
int ans = 0;
for (int i = 0; i <= D; i++) {
a[i + D] = 1LL * a[i + D] * ifac[i] % P;
ans = (ans + 1LL * a[i + D] * qpow((2 * i - D + P) % P, n)) % P;
}
ans = 1LL * ans * qpow((P + 1) / 2, D) % P;
printf("%d\n", ans);
return 0;
}
