题解 P3803 【【模板】多项式乘法（NTT）】2

## [更好的阅读体验点这里](http://www.cnblogs.com/zwfymqz/p/8980809.html) 很久以前写过一篇FFT的题解 今天有幸学到了NTT—— 一种比FFT更好用的算法，来跟大家分享一下 ## NTT 在FFT中，我们需要用到复数，复数虽然很神奇，但是它也有自己的局限性——需要用double类型计算，精度太低 那有没有什么东西能够代替复数且解决精度问题呢？ 这个东西，叫**原根** ------------ 原根的定义 > 设$m$是正整数，$a$是整数，若$a$模$m$的阶等于$\phi(m)$，则称$a$为模$m$的一个原根 定义中用到了群论的一些知识，不过不会也没关系，不影响接下来的学习 我们定义$P$为素数，$g$为$P$的原根 考虑原根为什么能代替单位根进行运算，(这部分可以跳过) 原因很简单，因为它具有和单位根相同的性质 在FFT中，我们用到了单位根的四条性质，而原根也满足这四条性质 1. 对于所有$\omega_n ^ t (0 \leq t \leq n - 1)$均不相同 不加证明的扔出一个定理 若$P$为素数，假设一个数$g$是$P$的原根，那么$g^i \mod P (1<g<P,0<i<P)$的结果两两不同 不要问我为什么，因为我也不知道。。 2. $\omega_{2n} ^ {2k} = \omega_n ^ k$ 通过代换可以得到 3. $\omega_n ^ { k + \frac{n}{2} } = -\omega_n ^ k$ 根据费马小定理和性质1可以得到 4. $1 + \omega_n ^ k + (\omega_n ^ k) ^ 2 + \dots + (\omega_n ^ k) ^ {n - 1} = 0 $ 由性质3和FFT中傅里叶逆变换的定理可以得到 ------------ 这样我们最终可以得到一个结论 $$\omega_n \equiv g^\frac{p-1}{n} \mod p$$ 然后把FFT中的$\omega_n$都替换掉就好了 $p$建议取$998244353$，它的原根为$3$。 如何求任意一个质数的原根呢？ 可以证明满足$g^r \equiv 1(\mod p)$的最小的$r$一定是$p-1$的约数 对于质数$p$，质因子分解$p-1$，若$g^{\frac{p-1}{p_i}} \neq 1 \pmod p$恒成立，$g$为$p$的原根 ------------ NTT求卷积代码: 确实比FFT快了不少 ```cpp #include<cstdio> #define getchar() (p1 == p2 && (p2 = (p1 = buf) + fread(buf, 1, 1<<21, stdin), p1 == p2) ? EOF : *p1++) #define swap(x,y) x ^= y, y ^= x, x ^= y #define LL long long const int MAXN = 3 * 1e6 + 10, P = 998244353, G = 3, Gi = 332748118; char buf[1<<21], *p1 = buf, *p2 = buf; inline int read() { char c = getchar(); int x = 0, f = 1; while(c < '0' || c > '9') {if(c == '-') f = -1; c = getchar();} while(c >= '0' && c <= '9') x = x * 10 + c - '0', c = getchar(); return x * f; } int N, M, limit = 1, L, r[MAXN]; LL a[MAXN], b[MAXN]; inline LL fastpow(LL a, LL k) { LL base = 1; while(k) { if(k & 1) base = (base * a ) % P; a = (a * a) % P; k >>= 1; } return base % P; } inline void NTT(LL *A, int type) { for(int i = 0; i < limit; i++) if(i < r[i]) swap(A[i], A[r[i]]); for(int mid = 1; mid < limit; mid <<= 1) { LL Wn = fastpow( type == 1 ? G : Gi , (P - 1) / (mid << 1)); for(int j = 0; j < limit; j += (mid << 1)) { LL w = 1; for(int k = 0; k < mid; k++, w = (w * Wn) % P) { int x = A[j + k], y = w * A[j + k + mid] % P; A[j + k] = (x + y) % P, A[j + k + mid] = (x - y + P) % P; } } } } int main() { N = read(); M = read(); for(int i = 0; i <= N; i++) a[i] = (read() + P) % P; for(int i = 0; i <= M; i++) b[i] = (read() + P) % P; while(limit <= N + M) limit <<= 1, L++; for(int i = 0; i < limit; i++) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (L - 1)); NTT(a, 1);NTT(b, 1); for(int i = 0; i < limit; i++) a[i] = (a[i] * b[i]) % P; NTT(a, -1); LL inv = fastpow(limit, P - 2); for(int i = 0; i <= N + M; i++) printf("%d ", (a[i] * inv) % P); return 0; } ``` ## 参考资料 [Menci的博客](https://oi.men.ci/fft-to-ntt/) [Candy?的博客](http://www.cnblogs.com/candy99/p/6641972.html)