voidInv(int N, int *a, int *b){ if (N == 1) { b[0] = Pow(a[0], mod - 2, mod); return; } // 清空 Inv((N + 1) >> 1, a, b);
int lim, l, *c = T; for (lim = 1, l = -1; lim <= N << 1; lim <<= 1, ++l); for (int i = 0; i < lim; ++i) to[i] = (to[i >> 1] >> 1) | ((i & 1) << l);
memcpy(c, a, sizeof(int) * N); memset(c + N, 0, sizeof(int) * (lim - N)); NTT(lim, b, 1), NTT(lim, c, 1); for (int i = 0; i < lim; ++i) b[i] = (2ll - 1ll * b[i] * c[i] % mod + mod) % mod * b[i] % mod; NTT(lim, b, -1); memset(b + N, 0, sizeof(int) * (lim - N)); }
code - 迭代实现
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
voidInv(int N, int *a, int *b){ int *c = T1; memset(b, 0, sizeof(int) * (N << 1)); b[0] = Pow(a[0], mod - 2, mod); for (int m = 2, lim = 4, l = 1; m <= N << 1; m <<= 1, lim <<= 1, ++l) { for (int i = 0; i < lim; ++i) to[i] = (to[i >> 1] >> 1) | ((i & 1) << l);
memcpy(c, a, sizeof(int) * m); memset(c + m, 0, sizeof(int) * (lim - m)); NTT(lim, b, 1), NTT(lim, c, 1); for (int i = 0; i < lim; ++i) b[i] = (2ll - 1ll * b[i] * c[i] % mod + mod) % mod * b[i] % mod; NTT(lim, b, -1); memset(b + m, 0, sizeof(int) * (lim - m)); } }
voidDiv(int N, int M, int *f, int *g, int *Q, int *R){ int lim, m = N - M + 1, *t = T1, *s = T2;
for (lim = 1; lim < m << 1; lim <<= 1); std::reverse_copy(g, g + M, t); // t <- Rev(g) Inv(m, t, s); // s <- Inv(Rev(g)) memset(s + m, 0, sizeof(int) * (lim - m)); // s <- Rev(g) mod x^(deg f - deg g + 1) std::reverse_copy(f, f + N, t); // t <- Rev(f) memset(t + m, 0, sizeof(int) * (lim - m)); // t <- Rev(f) mod x^(deg f - deg g + 1)
NTT(lim, t, 1), NTT(lim, s, 1); for (int i = 0; i < lim; ++i) t[i] = 1ll * t[i] * s[i] % mod; // t <- Rev(Q) NTT(lim, t, 0);
for (lim = 1; lim < N; lim <<= 1); memset(t + m, 0, sizeof(int) * (lim - m)); std::reverse(t, t + m); // t <- Q memcpy(Q, t, sizeof(int) * m); memcpy(s, g, sizeof(int) * M); // s <- g memset(s + M, 0, sizeof(int) * (lim - M));
NTT(lim, t, 1), NTT(lim, s, 1); for (int i = 0; i < lim; ++i) t[i] = 1ll * t[i] * s[i] % mod; // t <- Q * s NTT(lim, t, 0); for (int i = 0; i < M; ++i) R[i] = (f[i] - t[i] + mod) % mod; }
多项式求导 & 积分
设
则
1 2 3 4 5 6 7 8 9 10 11
voidDeri(int N, int *f, int *g){ // [0, N) -> [0, N - 1) for (int i = 1; i < N; ++i) g[i - 1] = 1ll * i * f[i] % mod; g[N - 1] = 0; }
voidInt(int N, int *f, int *g){ // [0, N) -> [0, N + 1) for (int i = N; i; --i) g[i] = 1ll * f[i - 1] * inv[i] % mod; g[0] = 0; }
voidExp(int N, int *f, int *g){ int *t = T3; g[0] = 1; for (int m = 2, lim = 4; m < N << 1; m <<= 1, lim <<= 1) { Ln(m, g, t); // t <- Ln(f_0) for (int i = 0; i < m; ++i) t[i] = (f[i] - t[i] + mod) % mod; t[0] = (t[0] + 1) % mod; // t <- 1 - Ln(f_0) + h
NTT(lim, g, 1), NTT(lim, t, 1); for (int i = 0; i < lim; ++i) g[i] = 1ll * g[i] * t[i] % mod; NTT(lim, g, 0);
memset(g + m, 0, sizeof(int) * (lim - m)); } }
多项式幂函数
problem
给定多项式 $f(x)$ 和整数 $k$,求 $f^k(x)$。
solution
显然可以直接套快速幂,复杂度 $O(n \log n \log k)$。有时候 $k$ 会非常大,需要一个更优秀的算法。