NTT
在FFT中，我们需要用到复数，复数虽然很神奇，但是它也有自己的局限性——需要用double类型计算，精度太低
那有没有什么东西能够代替复数且解决精度问题呢？
这个东西，叫原根
原根
原根的定义
设(m)是正整数，(a)是整数，若(a)模(m)的阶等于(phi(m))，则称(a)为模(m)的一个原根
定义中用到了群论的一些知识，不过不会也没关系，不影响接下来的学习
我们定义(P)为素数，(g)为(P)的原根
接下来不加证明的扔出一个很重要定理
- 若(P)为素数，假设一个数(g)是(P)的原根，那么(g^i mod P (1<g<P,0<i<P))的结果两两不同
不要问我为什么，因为我也不知道。。
考虑原根为什么能代替单位根进行运算，(这部分可以跳过)
原因很简单，因为它具有和单位根相同的性质
在FFT中，我们用到了单位根的四条性质，而原根也满足这四条性质
1 . 对于所有(omega_n ^ t (0 leq t leq n - 1))均不相同
这一条可以由上面的定理得到
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快了不少
代码语言：javascript复制#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;
}
