快速數論變換(NTT)小結

NTT

在FFT中，咱們須要用到複數，複數雖然很神奇，可是它也有本身的侷限性——須要用double類型計算，精度過低

那有沒有什麼東西可以代替複數且解決精度問題呢？

這個東西，叫原根

原根

階

若\(a,p\)互素，且\(p>1\)，

對於\(a^n \equiv 1 \pmod{p}\)最小的\(n\)，咱們稱之爲\(a\)模\(p\)的階,記作\(\delta_p(a)\)

例如：

\(\delta_7(2)=3\)，

\(2^1 \equiv 2 \pmod{7}\)

\(2^2 \equiv 4 \pmod{7}\)

\(2^3 \equiv 1 \pmod{7}\)

原根

原根的定義

設\(p\)是正整數，\(a\)是整數，若\(\delta_p(a)\)等於\(\phi(p)\)，則稱\(a\)爲模\(p\)的一個原根

\(\delta_7(3)=6=\phi(7)\)，所以\(3\)是模\(7\)的一個原根

注意原根的個數是不惟一的

若是模數\(p\)有原根，那麼它必定有\(\phi(\phi(p))\)個原根

原根存在的重要條件爲\(m = 2,4,p^a,2p^{a}\)，其中\(p\)爲奇素數\(a \ge 1\)

例如

原根有一個很是重要的定理：

• 若\(P\)爲素數，假設一個數\(g\)是\(P\)的原根，那麼\(g^i \mod P (1<g<P,0<i<P)\)的結果兩兩不一樣

不要問我爲何，由於我也不知道。。

考慮原根爲何能代替單位根進行運算，(這部分能夠跳過)

緣由很簡單，由於它具備和單位根相同的性質

在FFT中，咱們用到了單位根的四條性質，而原根也知足這四條性質

這樣咱們最終能夠獲得一個結論

\[\omega_n \equiv g^\frac{p-1}{n} \mod p\]

而後把FFT中的\(\omega_n\)都替換掉就行了

\(p\)建議取\(998244353\)，它的原根爲\(3\)。

如何求任意一個質數的原根呢？

對於質數\(p\)，質因子分解\(p−1\)，若\(g^{\frac{p-1}{p_i}} \neq 1 \pmod p\)恆成立，\(g\)爲\(p\)的原根

實現

NTT求卷積代碼:

確實比FFT快了很多

#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;
}