BZOJ4555[HEOI2016/TJOI2016]求和
题目链接
前置知识
第二类斯特林数
含义
\(S_n^m\)表示将\(n\)个互不相同的元素划分为\(m\)个非空集合的方案数
递推式
\[S_n^m = S_{n - 1}^{m - 1} + m \cdot S_{n - 1}^m
\]
通项公式
\[S_n^m = \frac{1}{m!} \cdot \sum_{k = 0}^m (-1)^k \cdot {m \choose k} \cdot (m - k)^n
\]
NTT
类似于\(FFT\)，在固定模数的情况下，可以将单位根用模数\(P\)的原根的\(\frac{P - 1}{n}\)次幂代替，避免复数运算掉精度
还是放个链接吧：NTT详解
解析
首先\(j\)的上界可以扩展到\(n\)，因为\(j > i\)的时候\(S(i,j) = 0\)
后两项都只和\(j\)有关，所以换一下枚举顺序，把这两项提出来：
\[f(n) = \sum_{j = 0}^n 2^j \cdot (j!) \cdot \sum_{i = 0}^n S(i, j)
\]
你看那个\(S(i,j)\)就很不友好，我们直接把它干掉，用通项替换：
\[f(n) = \sum_{j = 0}^n 2^j \cdot (j!) \cdot \sum_{i = 0}^n \frac{1}{j!} \cdot \sum_{k = 0}^j (-1)^k \cdot \frac{j!}{k! \cdot (j - k)!} \cdot (j - k)^i
\]
然后一波化简可得：
\[f(n) = \sum_{j = 0}^n 2^j \cdot (j!) \cdot \sum_{i = 0}^n \sum_{k = 0}^j (-1)^k \cdot \frac{1}{k! \cdot (j - k)!} \cdot (j - k)^i
\]
我们把后面那一坨提出来，令:
\[g(j) = \sum_{i = 0}^n \sum_{k = 0}^j (-1)^k \cdot \frac{1}{k! \cdot (j - k)!} \cdot (j - k)^i \\
= \sum_{k = 0}^j \frac{(-1)^k}{k!} \cdot \frac{\sum_{i = 0}^n (j - k)^i}{(j - k)!}
\]
再令：
\[a(i) = \frac{(-1)^i}{i!}，b(i) = \frac{\sum_{k = 0}^n i^k}{i!} = \frac{i^{k + 1} - 1}{(i - 1) \cdot i!}
\]
容易发现\(g\)就是\(a\)和\(b\)的卷积
特别地，\(b(0) = 1, b(1) = n + 1\)
那么用\(NTT\)或者\(FFT\)什么的求出\(g\)后，答案就是\(f(n) = \sum_{i = 0}^n 2^i \cdot (i!) \cdot g(i)\)
代码
我写的\(NTT\)，然而忘了开\(4n\)的数组\(RE\)了\(3\)发。。。。。
#include <cstdio>
#include <cstring>
#include <iostream>
#define MAXN 100005
typedef long long LL;
const LL mod = 998244353ll;
LL N, g[MAXN << 2], a[MAXN << 2], b[MAXN << 2], fact[MAXN];
LL G = 3, ans;
LL qpower(LL, LL);
void NTT(LL *, int, int);
inline LL inverse(LL x) { return qpower(x, mod - 2); }
int main() {
std::ios::sync_with_stdio(false);
std::cin >> N;
fact[0] = 1;
for (int i = 1; i <= N; ++i) fact[i] = fact[i - 1] * i % mod;
for (int i = 0; i <= N; ++i) a[i] = ((i & 1) ? -1 : 1) * inverse(fact[i]);
b[0] = 1, b[1] = N + 1;
for (int i = 2; i <= N; ++i) b[i] = (qpower(i, N + 1) - 1) * inverse(fact[i]) % mod * inverse(i - 1) % mod;
int sz = 1;
while (sz <= (N << 1)) sz <<= 1;
NTT(a, sz, 1);
NTT(b, sz, 1);
for (int i = 0; i < sz; ++i) g[i] = a[i] * b[i] % mod;
NTT(g, sz, -1);
LL pow2 = 1;
for (int i = 0; i <= N; ++i, pow2 = (pow2 << 1) % mod)
ans = (ans + pow2 * fact[i] % mod * g[i] % mod) % mod;
std::cout << ans << std::endl;
return 0;
}
LL qpower(LL x, LL y) {
LL res = 1;
while (y) {
if (y & 1) (res *= x) %= mod;
(x *= x) %= mod;
y >>= 1;
}
return res;
}
void NTT(LL *arr, int sz, int type) {
int rev[MAXN << 2]; rev[0] = 0;
for (int i = 1; i < sz; ++i) {
rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (sz >> 1) : 0);
if (rev[i] > i) std::swap(arr[i], arr[rev[i]]);
}
for (int len = 2, half = 1; len <= sz; len <<= 1, half <<= 1) {
LL wn = qpower(G, (mod - 1) / len);
if (type == -1) wn = inverse(wn);
for (int i = 0; i < sz; i += len) {
LL w = 1;
for (int j = 0; j < half; ++j, w = w * wn % mod) {
LL tmp1 = arr[i + j], tmp2 = arr[i + half + j] * w % mod;
arr[i + j] = tmp1 + tmp2; if (arr[i + j] >= mod) arr[i + j] -= mod;
arr[i + half + j] = tmp1 - tmp2; if (arr[i + half + j] < 0) arr[i + half + j] += mod;
}
}
}
if (type == -1) {
LL inv = inverse(sz);
for (int i = 0; i < sz; ++i) arr[i] = arr[i] * inv % mod;
}
}
//Rhein_E
