卷積導向快速傅立葉變換（FFT/NTT）教程

1 Forewords

卷積，但不止卷積 - FFT 漫談

• 先有 FT，再有 DFT，才有 FFT
• 時頻轉換是最初的用途
• 發現單位根優秀性質，James Cooley, John Tukey 發明現代 FFT 加速 DFT，但此前相似的發現早已有之
• 後來將 DFT 與卷積定理聯絡，FFT 才被用於計算多項式乘法
• 複數運算精度誤差推動了 NTT 的發展
• 應用：任何需要頻率和卷積的地方．頻譜、濾波器、音樂、雷達、影像處理……
• OI/XCPC 中主要關心卷積

推薦食用方法

• 初步要求

  • 知道 DFT、FFT 可用於快速計算多項式卷積
  • 掌握 FFT 加速 DFT 計算的原理和實現
  • 會應用結論改動 FFT 加速 NTT 計算
  • 見識一些卷積解決的基本問題，初步瞭解生成函式在組合計數中的應用
  • 題目可選擇性完成，請多花時間消化原理和思想

• 學有餘力 / 集訓後繼續消化

  • 系統學習生成函式
  • 實現多項式全家桶
  • 對原理感興趣的同學可對數學部分做進一步研究．本講內容是線性代數、抽象代數、數論等多領域的綜合應用．歡迎討論．
  • 學習集合冪級數相關知識點（FMT，FWT，……），體會其思想與 FFT 的同與異
  • 學習 FFT 在訊號、頻譜等非演演算法競賽向實際問題中的應用

• Learn for fun :)

記號說明

• 記 $[n]=\{0,1,\dots ,n-1\}$，此時可用 $k\in [n]$ 代替下標取值範圍 $k=0,1,\dots ,n-1$ 的記法．

  • 集合論中已定義 $n=\{0,1,\dots ,n-1\}$，這裡的中括號是為了強調其集合含義．

• 使用 Iverson 括號約定：設 $P$ 是一個命題，記 $$[P]:=\{\begin{array}{ll}1& P\text{ is true}\\ 0& \text{otherwise}\end{array}$$

• 多項式的規模定義為多項式的次數加一．特別的，零多項式的規模為 $0$．

  • 以後會混用 $n-1$ 次多項式和規模為 $n$ 的多項式的說法．

2 FFT/NTT in a nutshell

2.1 FFT

多項式卷積

給定兩個至多 $n-1$ 次的多項式 $$A(x)=\sum _{k=0}^{n-1}{a}_k{x}^k,B(x)=\sum _{k=0}^{n-1}{b}_k{x}^k$$ 如何快速計算兩者的卷積，即它們相乘得到的多項式的係數？ $$\begin{array}{rl}(A\ast B)(x)=A(x)B(x)& =\sum _{i=0}^{n-1}{a}_i{x}^i\sum _{j=0}^{n-1}{b}_j{x}^j\\ & =\sum _{k=0}^{2n-2}{x}^k\sum _{i+j=k}{a}_i{b}_j\end{array}$$ 上式給出了 $O(n^2)$ 的樸素做法．

係數 - 點值 - 係數

• 代入任意 $x$ 可得到多項式在 $x$ 處的點值

• 點值意義下的多項式乘法是 $O(n)$ 的 $$(A\ast B)(x)=A(x)B(x)$$

• $n$ 點確定一個至多 $n-1$ 次的多項式¹

若計算至多 $n-1$ 次的多項式的某 $n$ 個點值存在快速演演算法，而透過多項式的某 $n$ 個點值確定原多項式係數亦存在快速演演算法，就有透過 $O(n)$ 的點值乘法加速多項式乘法計算的可能．

係數 - 點值 - 係數 - 快速轉換？

• 樸素計算任意指定 $n$ 個位置點值需要 $O(n^2)$．

• Lagrange[1] 插值給出了 $O(n^2)$ 將任意位置 $n$ 個點值還原為多項式係數的演演算法．

• 能否選取 $n$ 個特殊的點值使係數 - 點值、點值 - 係數的變換支援快速計算？

Discrete Fourier Transform

離散傅立葉變換（Discrete Fourier Transform, DFT）接受一個至多 $n-1$ 次的多項式的 $n$ 個係數，將複數域上的 $n$ 個 $n$ 次單位根代入係數表示式以得到給定多項式的 $n$ 個點值． $$a_0,a_1,\dots ,a_{n-1}⟶A(1),A({\omega }_n),\dots ,A({\omega }_n^{n-1})$$

得益於單位根的特殊運算性質，二者均有被稱為快速傅立葉變換（Fast Fourier Transform, FFT）的快速演演算法．

複數域單位根

複數域上的 $n$ 個 $n$ 次單位根² $${\omega }_n^k:=e^{\frac{2\pi k}{n}i}=\cos \frac{2\pi k}{n}+i\sin \frac{2\pi k}{n},k\in [n]$$ 是複平面單位圓上的 $n$ 等分點，易驗證它們是複數域中唯一滿足方程 $z^n=1$ 的解．

所有單位根模長均為 $1$．第 $k$ 個單位根的輻角為 $\frac{2\pi k}{n}$．複數乘法”模長相乘，輻角相加”的性質告訴我們 ${\omega }_n^i{\omega }_n^j={\omega }_n^{i+j}$．

Remark (Euler’s formula for nerds). Euler 公式 $e^{it}=\cos t+i\sin t$ 的一種證明可用指數函式的另一定義 $\exp z=\underset{n\to \mathrm{\infty }}{lim}(1+\frac{z}{n}{)}^n$ 分析複數處極座標的極限 [2]．嚴格來講，如果三角函式和指數函式都由級數定義 [3]，Euler 公式幾乎是顯然的．總之，這裡僅將 Euler 公式作為一種緊湊的記號使用，細節不做要求．

複數域單位根 - 三個重要性質³

Theorem 1 (消去引理) $${\omega }_{dn}^{dk}={\omega }_n^k,\mathrm{\forall }d\in {\mathbb{N}}_{+}$$

Theorem 2 (折半引理) $$\{{\omega }_{2n}^{2k}:k\in [2n]\}=\{{\omega }_n^k:k\in [n]\}$$

消去 / 折半引理將在 FFT 的推導中使用．

複數域單位根 - 三個重要性質

Theorem 3 (求和引理) $$\frac{1}{n}\sum _{k=0}^{n-1}{\omega }_n^{ik}=[n\mid i]$$

求和引理的證明使用了等比數列求和公式．將在 IDFT 的推導中用到．

Fast Fourier Transform

考慮將至多 $2n-1$ 次的待變換多項式 $A(x)=\sum _{k=0}^{2n-1}{a}_k{x}^k$ 奇偶分項兩個至多 $n-1$ 次的多項式 $A(x)=A_0(x^2)+x{A}_1(x^2)$，其中 $$A_0(x)=\sum _{k=0}^{n-1}{a}_{2k}{x}^k,A_1(x)=\sum _{k=0}^{n-1}{a}_{2k+1}{x}^k$$ 代入 $x={\omega }_{2n}^k,k\in [2n]$，用單位根消去或折半引理（${\omega }_{2n}^{2k}={\omega }_n^k$）得 $$A({\omega }_{2n}^k)=A_0({\omega }_n^k)+{\omega }_{2n}^k{A}_1({\omega }_n^k)$$ 再用 ${\omega }_{2n}^n=-1$ 得 $$\begin{array}{rl}A({\omega }_{2n}^k)& =A_0({\omega }_n^k)+{\omega }_{2n}^k{A}_1({\omega }_n^k)\\ A({\omega }_{2n}^{n+k})& =A_0({\omega }_n^k)-{\omega }_{2n}^k{A}_1({\omega }_n^k)\end{array},k\in [n]$$ 原多項式 $A(x)$ 規模為 $2n$ 的 DFT 轉化為規模為 $n$ 的 $A_0(x)$ 和 $A_1(x)$ 的 DFT．遞迴計算就可得到 $O(n\log n)$ 的演演算法．

DFT 的矩陣表示

記 $$\begin{array}{c}\mathit{a}=\left(\begin{array}{c}{a}_0\\ a_1\\ ⋮\\ a_{n-1}\end{array}\right),F={\left({\omega }_n^{ij}\right)}_{(i,j)\in n\times n}=\left(\begin{array}{cccc}1& 1& \dots & 1\\ 1& {\omega }_n& \dots & {\omega }_n^{n-1}\\ ⋮& ⋮& \ddots & ⋮\\ 1& {\omega }_n^{n-1}& \dots & {\omega }_n^{(n-1)(n-1)}\end{array}\right)\end{array}$$ 則 DFT 的變換結果（給定係數 $a_0,a_1,\dots ,a_{n-1}$ 的多項式在 $n$ 個 $n$ 次單位根處的點值）可表示為 $$\hat{\mathit{a}}=\left(\begin{array}{c}A(1)\\ A({\omega }_n)\\ ⋮\\ A({\omega }_n^{n-1})\end{array}\right)=F\mathit{a}$$

IDFT

由單位根的消去引理可證，DFT 矩陣 $F$ 的逆矩陣為⁴ $$F^{-1}=\frac{1}{n}{\left({\omega }_n^{-ij}\right)}_{(i,j)\in n\times n}=\frac{1}{n}\left(\begin{array}{cccc}1& 1& \dots & 1\\ 1& {\omega }_n^{-1}& \dots & {\omega }_n^{-(n-1)}\\ ⋮& ⋮& \ddots & ⋮\\ 1& {\omega }_n^{-(n-1)}& \dots & {\omega }_n^{-(n-1)(n-1)}\end{array}\right)$$ 於是 IDFT 的過程可表示為 $$\mathit{a}=F^{-1}\hat{\mathit{a}}$$

故快速計算 IDFT 的方法與 FFT 幾乎一致，只需將計算 DFT 時使用的本原單位根 ${\omega }_n$ 替換為 ${\omega }_n^{-1}$ 並對最終結果除以 $n$ 即可．

梳理

• 怎麼計算卷積？

  • 把至多 $n-1$ 次的多項式 $A(x)$ 和至多 $m-1$ 次的多項式 $B(x)$ 寫成至多 $n+m-2$ 次的多項式（高位補 $0$）．為計算 FFT 方便，還要繼續補 $0$ 至一個大於其次數的 $2$ 的冪．
  • 對 $A(x)$ 和 $B(x)$ 多點求值．
  • 把兩個多項式的點值逐點相乘．
  • 多點插值還原 $(A\ast B)(x)$ 的係數．

• 怎麼快速求值？

  • 選點選單位根就是 DFT．
  • 消去引理和折半引理使我們可以用 FFT 演演算法遞迴地計算 DFT．
  • 推導已經給出了遞迴的寫法，之後還會介紹常數更優的迭代實現．

• 怎麼快速插值？

  • 求和引理給出了 DFT 矩陣的逆矩陣．
  • 計算方法很相似，最後逐項除掉一個規模．

FFT 遞迴實現 - DFT 部分

方便起見，我們只處理 $n$ 為 $2$ 的冪的情形．以下 C 風格的程式碼實現了遞迴的 DFT 和 IDFT．

#include<bits/stdc++.h>
#include<complex>
using namespace std;
typedef long long ll; typedef complex<double> CP;
const ll MXN=4E6+5; const double PI=3.14159265358979323846l;
CP tmp[MXN];
void _DFT(CP A[],ll n,ll typ){
    n/=2; if(n==0) return;
    for(ll k=0;k<n;k++) tmp[k]=A[2*k],tmp[n+k]=A[2*k+1];
    for(ll k=0;k<2*n;k++) A[k]=tmp[k];
    _DFT(A,n,typ); _DFT(A+n,n,typ);
    CP w(cos(2*PI/(2*n)),typ*sin(2*PI/(2*n))), wk=1;
    for(ll k=0;k<n;k++){
        tmp[  k]=A[k]+wk*A[n+k];
        tmp[n+k]=A[k]-wk*A[n+k];
        wk*=w;
    } for(ll k=0;k<2*n;k++) A[k]=tmp[k];
} void DFT(CP A[],ll n,ll typ){
    _DFT(A,n,typ); if(typ==-1) for(ll i=0;i<n;i++) A[i]*=1.0l/n;
}

FFT 遞迴實現 - 卷積部分

// alternatively, use std::__lg() in GCC
ll log2ceil(ll n){ll cnt=0; for(ll t=1;t<n;t<<=1) cnt++; return cnt;} 
CP A[MXN],B[MXN],C[MXN]; ll outC[MXN];
ll* conv(ll inA[],ll aN,ll inB[],ll bN){
    ll n=1LL<<log2ceil(aN+bN-1);
    for(ll i=0;i<aN;i++) A[i]=inA[i];
    for(ll i=0;i<bN;i++) B[i]=inB[i];
    DFT(A,n,1); DFT(B,n,1);
    for(ll i=0;i<n;i++) C[i]=A[i]*B[i];
    DFT(C,n,-1); for(ll i=0;i<n;i++) outC[i]=round(C[i].real());
    return outC;
}

Drawbacks?

• 遞迴實現慢
• 臨時陣列醜
• 封裝性為零

FFT 迭代

迭代地實現 FFT 不僅在常數上更加優秀，亦更便於使用 C++ 的容器進行封裝．這並不困難，只需自底向上模擬 FFT 遞迴過程即可．

唯一的問題——最底層的順序？

來觀察一輪 $2^3$-FFT 自頂向下的置換過程 $$\begin{array}{ccccccccc}{2}^3& 0/000& 1/001& 2/010& 3/011& 4/100& 5/101& 6/110& 7/111\\ 2^2& 0/000& 2/010& 4/100& 6/110& 1/001& 3/011& 5/101& 7/111\\ 2^1& 0/000& 4/100& 2/010& 6/110& 1/001& 5/101& 3/011& 7/111\\ 2^0& 0/000& 4/100& 2/010& 6/110& 1/001& 5/101& 3/011& 7/111\end{array}$$

你發現了什麼？

蝶形運算

在計算點值前，$2^n$-FFT 事實上完成了一次 $n$-位逆序置換．分解來看，規模為 $2^k$ 的層的置換完成了對 $2^{n-k}$ 對應二進位制位的分類．

我們有 $O(n)$ 的遞推方法獲得這一置換．

void spawnrev(ll n){
    rev[0]=0;
    for(ll i=1;i<(1<<n);i++)
        rev[i]=(rev[i>>1]>>1)+((i&1)<<(n-1));
}

FFT 迭代實現

void DFT(CP A[],ll n,ll typ){ // rev[] should be spawned in advance
    for(ll i=0;i<n;i++) if(i<rev[i]) swap(A[i],A[rev[i]]); // a one-to-one permutation
    for(ll hf=1;hf<n;hf*=2){
        CP w(cos(2*PI/(2*n)),typ*sin(2*PI/(2*n))), wk=1;
        for(ll i=0;i<n;i+=hf*2){
            CP wk=1;
            for(ll k=0;k<hf;k++){
                CP x=A[i+k],y=wk*A[i+hf+k];
                A[i+k]=x+y; A[i+hf+k]=x-y;
                wk=wk*w;
            }
        }
    }
    if(typ==-1) for(ll i=0;i<n;i++) A[i]*=1.0l/n;
}

請自行實現更易用的容器封裝版本．

2.2 NTT

NTT 速成⁵

FFT 的缺點？浮點數帶來的大常數與精度問題．

$998244353=2^{23}\times 7\times 17+1$，同時是一個質數．

我們指出，在係數和點值模 $p=998244353$ 的意義下，當規模 $n\le p$ 時⁶，至多 $n-1$ 次的多項式仍可由其 $n$ 個點值唯一確定，故仍可使用係數-點值-係數的方法求得多項式卷積．

我們指出，對於滿足 $n\mid p-1$ 的 $n$，依 ${\omega }_n:=3^{\frac{p-1}{n}}$ 定義的 ${\omega }_n$ 在模 $p$ 意義下與複數域中定義的 ${\omega }_n$ 發揮相同的作用，仍可進行規模至多為 $2^{23}$ 的 FFT 作為模 $p=998244353$ 意義下 NTT 的快速演演算法．

只需修改單位根定義，把複數運算改為整數取模，就得到了能算 NTT 的 FFT 的實現．

const PR=3,MOD=998244353;
ll w=qpow(PR,(MOD-1)/(hf*2)); if(typ==-1) w=inv(w);

FFT/NTT in a nutshell - 小結：概念區分

• 關於 DFT

  • Discrete Fourier Transform, DFT, 離散傅立葉變換
  • Fast Fourier Transform, FFT, 快速傅立葉變換
  • FFT 是計算 DFT 的快速演演算法

• 關於 NTT

  • Number Theoretic Transform, NTT, 數論變換
  • FFT 在複數域上的多項式環 $\mathbb{C}[x]$ 中進行，而 NTT 在模 $p$ 剩餘類域 ${\mathbb{Z}}_p[x]$ 上進行
  • 快速計算 DFT / IDFT, NTT 都用 FFT，故一般不使用”FNTT”的說法

3 Applications

3.1 基本應用

基本應用

大整數乘法

十進位制數可拆解為多項式表示，計算卷積後處理進位即可．由於數字最大隻是 $9$，合理資料範圍下捲起來不會爆模數，直接用 NTT 實現即可．

揹包計數

兩個揹包的合併就是多項式卷積．

濾波器

反轉多項式的係數陣列再做卷積，可以快速得到兩個多項式滑動視窗式的內積．

位運算

有些位運算可以寫成卷積的形式．模 $2$ 意義下”異或”是加法，“與”是乘法，“或”可以透過取反轉化為”與”．

字串

透過巧妙設定字串距離函式，FFT 可解決更廣泛的字串匹配問題．

值域次數化

當值域較小時，將待計算的值放在多項式次數上統計貢獻次數，可以繞開某些極難求解的數值問題．例如 Vandermonde 行列式的快速計算．

基本應用 - 分治 FFT

對多個長度相同的多項式的卷積，分治地捲起來可降低時間複雜度．長度不一時，挑小的先卷也可減小常數（用堆維護）．

另有一種 CDQ 風格的分治 FFT．CDQ 長於處理帶偏序的二元點對貢獻，在處理形如 $c_k=\sum _{i>j}[i+j=k]a_i{b}_j$ 的帶偏序卷積時可以應用．

當卷積的前後項存在依賴關係時，也可使用此法保證處理順序恰當．然而此類依賴問題往往也可透過解生成函式方程的方法求得封閉形式．

時間複雜度均為 $$T(n)=2T\left(\frac{n}{2}\right)+O(n\log n)=O(n{\log }^{2}n)$$

3.2 生成函式初步

生成函式初步 - 導言

生成函式是一種對數列的操作技巧．透過將數列表示為多項式或形式冪級數，數列間複雜的和式操作可用簡單函式的乘法、複合等運算進行表示，從而大大降低了數列變換技巧的使用門檻．

生成函式在組合數學中應用廣泛，且生成函式的部分操作在組合意義下也有較為直觀的理解．本節將帶大家初窺其中的奧妙．限於篇幅和主講人能力，我們僅以題帶點地講解，期冀為大家建立構造生成函式的直覺．請感興趣的同學下來做進一步研究．

熟悉 Taylor 展開的同學或能較快上手此部分內容．

Ordinary Generating Function

序列 $a_n$ 的普通生成函式（Ordinary Generating Function, OGF）定義為其對應的多項式（形式冪級數）$A(x)=\sum _{n=0}^{+\mathrm{\infty }}{a}_n{x}^n$．這樣書寫的目的是為了便於進行數列間的卷積操作．特別的，Taylor 展開（或廣義二項式定理）風格的 $\frac{1}{1-x}=1+x+x^2+\dots$ 也是常用技巧之一⁷．

Exercise 1 寫出下列數列的 OGF．下標從 $0$ 開始．

• $1,0,0,0\dots$
• $1,1,1,1,\dots$
• $1,0,1,0,\dots$
• $1,-1,1,-1,\dots$
• $1,2,3,4,\dots$
• $1,2,4,8,\dots$

OGF 組合意義

OGF 相乘，是揹包，是卷積，是兩塊無標號組合物件的有序拼接．

Exercise 2 寫出下列計數問題的 OGF，均以 $n$ 作為數列的下標．

• 在 $m$ 個物品中選出 $n$ 個的方案數．
• 容量為 $n$ 的揹包裝下體積分別為 $a_i$ 的 $m$ 個物品的方案數．
• 同上，但每個物品有無限個．
• $n$ 個無標號球放入 $r$ 個有標號盒的方案數，要求盒非空．
• 同上，但盒可空．
• 將整數 $n$ 分拆為若干正整數之和的方案數．

Exponential Generating Function

序列 $a_n$ 的指數生成函式（Exponential Generating Function, EGF）定義為 $$A(x)=\sum _{n=0}^{+\mathrm{\infty }}{a}_n\frac{x^n}{n!}$$

Exercise 3 寫出下列序列的 EGF，下標從 $0$ 開始．

• $1,0,0,0,\dots$
• $1,1,1,1,\dots$
• $1,2,4,8,\dots$
• $0,1,0,-1,0,1,0,-1,\dots$
• $1,0,-1,0,1,0,-1,0,\dots$
• $0,1,-1,2!,-3!,4!,\dots$

EGF 組合意義

觀察兩個 EGF 的乘積 $$\begin{array}{rl}A(x)B(x)& =\sum _{i=0}^{+\mathrm{\infty }}{a}_i\frac{x^i}{i!}\sum _{j=0}^{+\mathrm{\infty }}{b}_j\frac{x^j}{j!}\\ & =\sum _{i=0}^{+\mathrm{\infty }}\sum _{j=0}^{+\mathrm{\infty }}{a}_i{b}_j\frac{x^{i+j}}{i!j!}\\ & =\sum _{n=0}^{+\mathrm{\infty }}\sum _{i=0}^n{a}_i{b}_{n-i}\frac{x^n}{i!(n-i)!}\\ & =\sum _{n=0}^{+\mathrm{\infty }}\sum _{i=0}^n{a}_i{b}_{n-i}(\genfrac{}{}{0ex}{}{n}{i})\frac{x^n}{n!}\end{array}$$ EGF 的乘積，相當於兩塊有標號組合物件的有序拼接．

EGF 組合意義

$$A(x)B(x)=\sum _{i=0}^{+\mathrm{\infty }}{a}_i\frac{x^i}{i!}\sum _{j=0}^{+\mathrm{\infty }}{b}_j\frac{x^j}{j!}=\sum _{n=0}^{+\mathrm{\infty }}\sum _{i=0}^n{a}_i{b}_{n-i}(\genfrac{}{}{0ex}{}{n}{i})\frac{x^n}{n!}$$

Exercise 4 寫出下列計數問題的 EGF，均以 $n$ 作為數列的下標．

• 長度為 $n$ 的排列的構型數．
• 長度為 $n$ 的圓排列的構型數．
• 將 $n$ 個有標號球放入 $r$ 個有標號盒的方案數，要求盒非空．
• 將 $n$ 個有標號球放入 $r$ 個無標號盒的方案數，要求盒非空．
• 將 $n$ 元集合劃分為 $r$ 個等價類的方案數．
• 劃分 $n$ 元集合的方案數．

4 Mathematics behind

4.1 NTT 原理

NTT 原理 - 導言

本節介紹 NTT 的原理．

FFT 加速卷積演演算法的核心，一是多項式的求值插值原理，二是單位根帶來的分治快速演演算法．我們將在本節中證明，模 $p$ 剩餘類域 ${\mathbb{Z}}_p$ 中，只要 $n\le p$，仍可透過 $n$ 個點值唯一確定至多 $n-1$ 次的 ${\mathbb{Z}}_p$ 上的多項式，這保證了係數-點值-係數方法的正確性；此外，基於原根構造出的 ${\mathbb{Z}}_p$ 上的本原單位根，使得使用與 FFT 相同的分治演演算法加速 NTT 的計算成為可能．

NTT 原理涉及原根等數論內容．本講的目標是建立 DFT 變換和 FFT 演演算法的通用數學框架，而非具體研究其某一特例．故我們只講解 NTT 所需的基礎數論知識，無關的細節則略過處理．對數論感興趣的同學可前往 OI Wiki 學習．

模 $p$ 剩餘類域 ${\mathbb{Z}}_p$

模 $p$ 剩餘類域 ${\mathbb{Z}}_p$ 是刻畫取模運算下整數加法、乘法運算規律的代數結構，其良定義性由 $$\{\begin{array}{l}{a}_1\equiv a_2(\mathrm{mod}p)\\ b_1\equiv b_2(\mathrm{mod}p)\end{array}⟹\{\begin{array}{l}{a}_1+b_1\equiv a_2+b_2(\mathrm{mod}p)\\ a_1{b}_1\equiv a_2{b}_2(\mathrm{mod}p)\end{array}$$ 保證．

模 $p$ 剩餘類域的最大特點是其中每個元素都存在（唯一）逆元．這一點是數論中 Bézout 定理的直接應用．

模 $p$ 剩餘類域的另一優良性質是其滿足消去律，即 $$ab\equiv 0(\mathrm{mod}p)⟹a\equiv 0\vee b\equiv 0(\mathrm{mod}p)$$ 而模一般的數 $m$ 形成的代數結構則不滿足此性質．

${\mathbb{Z}}_p[x]$ 上的多點插值

之後記係數均在 ${\mathbb{Z}}_p$ 上的多項式組成的集合為 ${\mathbb{Z}}_p[x]$．

Lemma 1 (${\mathbb{Z}}_p[x]$ 上的多項式餘式定理) 設 $A(x)\in {\mathbb{Z}}_p[x]$，用 $Z_p[x]$ 中的一次多項式 $(x-a)$ 對 $A(x)$ 進行帶餘除法，得到的餘式為至多零次的多項式 $A(a)$．

${\mathbb{Z}}_p[x]$ 上的多點插值

Theorem 4 (Lagrange 定理) 設 $p$ 是一個質數，設 $A(x)\in Z_p[x]$．同餘方程 $A(x)\equiv 0(\mathrm{mod}p)$ 只有至多 $\mathrm{deg}A(x)$ 個模 $p$ 意義下不同的整數解，除非這多項式的係數全為零．

Corollary 1 設 $A(x),B(x)\in {\mathbb{Z}}_p[x]$ 是至多 $n-1$ 次的兩個多項式．若同餘方程 $A(x)\equiv B(x)(\mathrm{mod}p)$ 有 $n$ 個不同的模 $p$ 意義下的整數解，則兩多項式一定相等．

推論告訴我們，欲確定 ${\mathbb{Z}}_p$ 上的一個至多 $n-1$ 次的多項式，只需找到其 $n$ 個模 $p$ 意義下的點值即可．因此，使用係數 - 點值 - 係數方法加速 ${\mathbb{Z}}_p[x]$ 上的卷積是可行的．

本原單位根

回顧 DFT 中複數域 $\mathbb{C}$ 上單位根的定義 $${\omega }_n^k:=e^{\frac{2\pi k}{n}i}=\cos \frac{2\pi k}{n}+i\sin \frac{2\pi k}{n},k\in [n]$$ 其中 ${\omega }_n$ 是最重要的單位根，因為它生成了所有其它的單位根．

事實上所有的 ${\omega }_n^k:gcd(n,k)=1$ 也都是重要的，它們的 $0$ 至 $n-1$ 次方也能生成所有的單位根．這也是數論中 Bézout 定理的直接應用．

我們把這一類重要的單位根稱為本原單位根．抽象的來說，$n$ 次本原單位根 ${\omega }_n$ 可定義如下 $${\omega }_n^n=1;{\omega }_n^k\ne 1,k=1,2,\dots ,n-1$$

${\mathbb{Z}}_m$ 上的階

如何定義 ${\mathbb{Z}}_m$ 上的“本原單位根”呢？

稱在模 $m$ 意義下使得 $a^n=1$ 的最小正整數 $n$ 為 $a$ 的階，記為 ${\mathrm{ord}}_m(a)$．在不至混淆的情況下，也簡記為 $\mathrm{ord}(a)$．

• 階最高有多高？

Theorem 5 (Fermat 小定理) 若 $p$ 是質數，$a$ 是非零數，則 $a^{p-1}\equiv 1(\mathrm{mod}p)$．

Theorem 6 (Euler 定理) 若 $gcd(a,m)=1$，則 $a^{\phi (m)}\equiv 1(\mathrm{mod}m)$，其中 $\phi (m)$ 是數論中的 Euler 函式⁸．

Euler 定理的證明 - 簡化剩餘系

${\mathbb{Z}}_m^R:=\{a\in {\mathbb{Z}}_m:gcd(a,m)=1\}$ 被稱為 ${\mathbb{Z}}_m$ 的簡化剩餘系．簡化剩餘系對模 $m$ 乘法封閉．簡化剩餘系中的每個元素都具有模 $m$ 意義下的唯一逆元，且其也在 ${\mathbb{Z}}_m^R$ 中．此二性質易由 $gcd$ 的線性組合意義證明．

簡化剩餘系的大小即 $\phi (m)$．這是 Euler 函式的其中一種定義．

設 $a\in {\mathbb{Z}}_m^R$，其在 ${\mathbb{Z}}_m$ 中存在逆元，故我們有 $$\mathrm{\forall }x,y\in {\mathbb{Z}}_m,x\equiv y(\mathrm{mod}m)⟺ax\equiv ay(\mathrm{mod}m)$$ 當限定 $x,y\in {\mathbb{Z}}_m^R$ 時，結合 ${\mathbb{Z}}_m^R$ 上述兩個性質，就有 $${\mathbb{Z}}_m^R=a{\mathbb{Z}}_m^R=\{axmodm:x\in {\mathbb{Z}}_m^R\}$$ 故 $$a^{\phi (m)}\prod _{x\in {\mathbb{Z}}_m^R}x\equiv \prod _{x\in {\mathbb{Z}}_m^R}ax\equiv \prod _{x\in a{\mathbb{Z}}_m^R}x\equiv \prod _{x\in {\mathbb{Z}}_m^R}x⟹a^{\phi (m)}\equiv 1(\mathrm{mod}m)$$

$a^n$ 的取值 - $a$, $m$ 互質時

補充討論 $a^{n}modm,n\in {\mathbb{Z}}_{+}$ 可取得哪些 ${\mathbb{Z}}_m$ 中的元素．

當 $gcd(a,m)=1$ 時，用 $gcd$ 的線性組合含義，有 $$\begin{array}{rl}gcd(a,m)=1& ⟹\mathrm{\exists }{x}_0,a{x}_0=1(\mathrm{mod}m)\\ & ⟹\mathrm{\exists }{x}_0,a^n{x}_0^n=1(\mathrm{mod}m)\\ & ⟹gcd(a^n,m)=1\end{array}$$ 故此時 $a^{n}modm$ 只可能取得 ${\mathbb{Z}}_m$ 的簡化剩餘系 ${\mathbb{Z}}_m^R$ 中的元素．

$a^n$ 的取值 - 一般情況⁹

當 $d_1=gcd(a,m)\ne 1$ 時，設 $a=d_1{a}_1,m=d_1{m}_1$，當 $n\ge 1$ 時， $$\begin{array}{rl}{a}^{n}modm& =d_1{a}_1{a}^{n-1}mod{d}_1{m}_1=d_1(a_1{a}^{n-1}mod{m}_1)\end{array}$$ 此時 $gcd(a_1,m_1)=1$，根據線性同餘方程理論，括號內表示式的取值集合僅與 $a^{n-1}mod{m}_1$ 有關．此時若 $d_2=gcd(a,m_1)$ 仍不等於 $1$，則令 $a_1=d_2{a}_2,m_1=d_2{m}_2$．類似地繼續化歸，最終有 $d_{k+1}=gcd(a,m_k)=1$， $$a^{n}modm=d_1(a_1{d}_2(\dots a_{k-1}{d}_k(a_k{a}^{n-k}mod{m}_k)\dots )mod{m}_1)$$ 故當 $n\ge k$ 時，$a^{n}modm$ 的取值集合為 $$\{d_1{d}_2\dots d_k(a^{t}mod{m}_k):t\in {\mathbb{Z}}_{+}\}\subset d_1{d}_2\dots d_k{\mathbb{Z}}_{m_k}^R$$ 事實上，$\prod _{i=1}^k{d}_i$ 恰為 $gcd(a,m)$ 的各素因子在 $m$ 中全部成分的乘積．

原根¹⁰

需要注意的是，Euler 定理只給出了 ${\mathbb{Z}}_m$ 上元素階的一個上界．換言之，其說明 $$\mathrm{\forall }a\in {\mathbb{Z}}_m^R⟹{\mathrm{ord}}_m(a)\mid \phi (m)$$ 而對於那些使得 $gcd(a,m)\ne 1$ 的元素 $a$，用 $gcd$ 的線性組合性質容易證明其永遠不可能透過自乘變為 $1$．我們設定這類元素的階為 $0$．

稱 ${\mathbb{Z}}_m$ 上那些確能達到階數上界 $\phi (m)$ 的元素是 ${\mathbb{Z}}_m$ 上的原根，用符號 $g$ 表示．

Theorem 7 (原根存在定理) ${\mathbb{Z}}_m$ 上存在原根，當且僅當 $m=2,4,p^{\alpha },2p^{\alpha }$ ，其中 $p$ 是一奇質數．

求階和原根

• 求階，用定理 ${\mathrm{ord}}_m(a)\mid \phi (m)$，求因子 + 快速冪即可 $O(\sqrt{m}\log m)$．
• 找原根，從小到大用上述求階方法暴力即可．最小原根一般不會太大．
• 更快的方法請參考洛谷求原根模板題題解．

下面再介紹原根的兩個定理．

Theorem 8 (原根判定定理) 若 $gcd(g,m)=1$，則 $g$ 是 ${\mathbb{Z}}_m$ 上原根的充要條件是，對 $\phi (m)$ 的每個素因數 $p$，都有 $g^{\frac{\phi (m)}{p}}\not\equiv 1(\mathrm{mod}m)$.

必要性顯然．充分性，反證出所有 $\phi (m)$ 的非平凡因子均不是 $g$ 的階即可．

Theorem 9 (原根個數定理) 若 ${\mathbb{Z}}_m$ 上存在原根，則它原根的個數為 $\phi (\phi (m))$.

${\mathbb{Z}}_m$ 上的本原單位根構造

設 $a$ 是 ${\mathbb{Z}}_m$ 上一元素，當 $n\mid {\mathrm{ord}}_m(a)$，令 ${\omega }_n=a^{\frac{{\mathrm{ord}}_m(a)}{n}}$，容易發現 ${\omega }_n$ 就是 ${\mathbb{Z}}_m$ 上的 $n$ 次本原單位根．

使用數論中 Bézout 定理，${\omega }_n^k$ 是 ${\mathbb{Z}}_m$ 上的 $\frac{n}{gcd(n,k)}$ 次本原單位根．這是 ${\mathbb{Z}}_m$ 上本原單位根的消去引理．

遺憾的是，一般的 ${\mathbb{Z}}_m$ 中，本原單位根不一定滿足求和引理¹¹ $\frac{1}{n}\sum _{k=0}^{n-1}{\omega }_n^{ik}=[n\mid i]$．這是因為，求和引理的證明中使用了等比數列求和和消去律，但一般的 ${\mathbb{Z}}_m$ 中消去律並不成立．

${\mathbb{Z}}_m$ 上的 NTT - $m$ 是質數

求和引理是保障 NTT 逆變換對應矩陣確為 $n^{-1}{\left(\begin{array}{c}{\omega }_n^{-ij}\end{array}\right)}_{(i,j)\in n\times n}$ 的關鍵性質．如果我們還想保留這一點，就必須要求 $m$ 是質數．

為保證 $n$ 在 ${\mathbb{Z}}_m$ 中存在逆元，必須要求 $gcd(n,m)=1$．已經要求 $m$ 是質數，故無需做出額外要求．

綜上，逆變換矩陣確為 $n^{-1}{\left(\begin{array}{c}{\omega }_n^{-ij}\end{array}\right)}_{(i,j)\in n\times n}$ 的 $n$ 點 NTT 要求 $m$ 是質數，且 $n\mid {\mathrm{ord}}_m(a)$．已經要求 $m$ 是質數，故 ${\mathbb{Z}}_m$ 中必存在原根，不妨直接取 ${\mathbb{Z}}_m$ 的一個原根 $g$ 構造 $n$ 次本原單位根 ${\omega }_n=g^{\frac{\phi (m)}{n}}=g^{\frac{m-1}{n}}$，這樣的構造支援至多 $m-1$ 次單位根的存在，顯然是最優的選擇．

${\mathbb{Z}}_m$ 上的 FFT 與卷積加速 - 998244353

為滿足 FFT 對 $2^m$ 次本原單位根的需求，只需選擇 $p=k{2}^m+1$ 型的奇質數 $p$，就可以在 ${\mathbb{Z}}_p$ 上支援規模至多為 $2^m$ 的 NTT/FFT 及其逆變換執行．

$998244353=2^{23}\times 7\times 17+1$，同時是一個質數，在 int 型中的單次加減操作不會溢位，是 OI/XCPC 計數題中不可多得的優秀模數¹²．

結合前述關於 ${\mathbb{Z}}_p[x]$ 上多點插值的討論，透過 NTT/FFT 加速卷積運算的正確性得到完整證明．

${\mathbb{Z}}_m$ 上的 NTT - 對一般的 $m$

當 $m$ 不一定是質數時，NTT 逆變換不能表示為 $n^{-1}{\left(\begin{array}{c}{\omega }_n^{-ij}\end{array}\right)}_{(i,j)\in n\times n}$ 的形式．但這並不代表 NTT 不可逆．因此，儘管不再實用，一般 ${\mathbb{Z}}_m$ 上的 NTT 變換仍有討論價值，即研究 ${\left(\begin{array}{c}{\omega }_n^{ij}\end{array}\right)}_{(i,j)\in n\times n}$ 的可逆性．

我們找到一篇有關該問題的參考文獻[7, Secs. 3 and appendix B]，但尚不確定其證明的正確性．友情提示讀者：${\mathbb{Z}}_m$ 不是域，甚至不是整環，故線性空間中矩陣的性質不能直接應用於 ${\mathbb{Z}}_m$ 上的矩陣，討論可逆性時還需小心謹慎．

從分析到代數

後續數學內容導讀 NTT 原理雖已非常”數學”，但也只是 DFT 在有限域上的一個例項．本節往後，我們要嘗試為多項式係數位於複數域 $\mathbb{C}$ 上的 DFT、位於模 $p$ 剩餘類域 ${\mathbb{Z}}_p$ 上的 NTT 及它們的快速演演算法 FFT 建立一個統一的數學框架．這需要我們剖析求值插值的基本原理，提煉出 FFT 演演算法成立的根本要求．

後續數學內容不再要求掌握．望同學們在紛繁的定理定義中抓住要旨，窺見抽象數學背後蘊藏的規律．熟悉高等代數和抽象代數的同學或會對某些內容感到熟悉．抽象地討論 FFT 的資料並不多見，後續內容多為主講人的新進探索，或有謬誤，望不吝指正．

4.2 求值與插值

求值與插值 - 導言

本節將重新審視已經熟知的多項式，把抽象的、代數的多項式和具體的、分析的多項式函式區分開來．我們指出，多項式和多項式函式不同但關聯緊密，形式冪級數與冪級數亦有此類聯絡．這些抽象的討論將幫助我們剖析多項式求值插值的基本原理．

除常見代數書目（如[8]），也推薦參考 OI Wiki 的多項式基礎[9]和 Wikipedia 的形式冪級數[10]．

多項式

設（無窮）序列 $\{a_n\}$ 是一個只有有限個非零項的序列，其元素均在環 $R$ 上．環 $R$ 上的多項式環[^14] $R[x]$ 是所有滿足上述條件的序列構成的集合以及在它們之間定義的兩種運算 $+,\times$ 的合稱．此語境下，我們也將構成 $R[x]$ 的序列稱為 $R$ 上的多項式．序列中的元素被稱為多項式的係數．序列間定義的加法和乘法是普通序列升級為多項式的關鍵．

多項式環上的加法、乘法的定義已經為大家所熟知．係數所處的環保證了多項式加法和乘法的良定義，而在這兩種運算下，$R$ 上多項式的集合也構成一個環的結構．

Remark (群、環、域). 群、環、域都是常見的代數結構，其中的元素在給定運算下封閉，並滿足特定的運算性質．簡單來說，環¹³上定義了加法和可能不可逆的、不一定交換的乘法，域上定義了加減乘除所有四則運算．交換環中的乘法滿足交換律．除環中的所有元素有乘法逆元．域是交換除環．

多項式

習慣上也會將多項式 $A$ 寫為 $A(x)=\sum _{k=0}^{n-1}{a}_k{x}^k$ 的形式．其中未定元 $x$ 只是一個符號，僅為方便加、乘法的理解而使用，並無任何含義．

多項式 $A$ 的次數 $\mathrm{deg}A$ 定義為其最高非零項所處的位置下標．特別的，定義零多項式的次數為 $-\mathrm{\infty }$．對任意多項式 $A,B\in R[x]$，$\mathrm{deg}(A+B)\le max\{\mathrm{deg}A,\mathrm{deg}B\}$．當 $R$ 是整環¹⁴時，$A,B$ 的首項乘積非零，故 $\mathrm{deg}(A\ast B)=\mathrm{deg}A+\mathrm{deg}B$．

下面額外為多項式定義一種新的運算．多項式 $A$ 和多項式 $B$ 的複合 $A\circ B$ 定義為 $$A\circ B:=\sum _{k=0}^{+\mathrm{\infty }}{a}_k{B}^k$$ 其中 $B^k$ 代表 $k$ 個多項式 $B$ 的乘積．習慣上也會將多項式 $A$ 和多項式 $B$ 的複合寫為 $A(B(x))=\sum _{k=0}^{+\mathrm{\infty }}{a}_k{B}^k(x)$ 的形式．由於多項式的次數有限，複合運算中僅包含了有限次多項式加法和乘法．

形式冪級數

形式冪級數定義與多項式的唯一區別是其不要求 $\{a_n\}$ 只有有限項非零．類似的，$R$ 上的形式冪級數環記為 $R[[x]]$．

由於涉及無限次運算，形式冪級數的複合運算需考慮環 $R$ 上的收斂問題．為迴避此問題，一般規定複合右側函式的常數項須為零．

Remark. DFT, NTT 與多項式環 DFT/FFT 加速的多項式卷積在複數域 $\mathbb{C}$ 上的多項式環 $\mathbb{C}[x]$ 中進行，NTT/FFT 加速的多項式卷積在模 $p$ 剩餘類域 ${\mathbb{Z}}_p$ 上的多項式環 ${\mathbb{Z}}_p[x]$ 中進行．

帶餘除法

整環是無零因子的交換麼環．所謂無零因子，即環中任意元素 $a,b$ 滿足 $ab=0⟹a=0\vee b=0$．所有域都是整環．當 $R$ 是整環時，$R[x]$ 也是整環¹⁵．整環上可以定義整除相關理論．

定義有帶餘除法的環被稱為 Euclid 整環．域上的多項式環都是 Euclid 整環[8, 第 7 章第 2 節定理 3, p. 11]．值得注意的是，$R$ 的 Euclid 性不可傳遞至 $R[x]$¹⁶．

如在帶餘除法中保證除式是首一多項式¹⁷，則帶餘除法的良定義和進行過程也均可在整環上實現．

多項式函式 / 冪級數

剛剛強調，多項式 / 形式冪級數只是定義了加法和乘法的序列．現在介紹多項式函式和冪級數．它們不是序列，而是對映¹⁸．

多項式 $A$ 對應的環 $R$ 上的多項式函式 $A(x)$ 定義為對映¹⁹ $$A:R\to R,x\mapsto \sum _{k=0}^{+\mathrm{\infty }}{a}_k{x}^k$$ 其中 $\{a_n\}$ 只有有限項非零．環 $R$ 上的冪級數則無此限制²⁰．

與多項式 / 形式冪級數不同，這裡的 $x\in R$ 代表函式的自變數，是會按環 $R$ 上運演演算法則參與運算的有意義的變數．

多項式函式的加法和乘法定義為函式的加法和乘法，即 $$\begin{array}{rl}A+B:x& \mapsto A(x)+B(x)\\ A\otimes B:x& \mapsto A(x)B(x)\end{array}$$

求值

多項式和多項式函式似乎在許多情況下有著平行的關係．下面介紹一個較直觀的結論．

若環 $R$ 是交換環，則環 $R$ 上每個多項式都唯一確定一個 $R$ 上的多項式函式，且該對映保持加法和乘法運算²¹．形式化地可以記為 $$\begin{array}{rl}(A+B)(x)& =A(x)+B(x)\\ (A\ast B)(x)& =A(x)B(x)\end{array}$$ 可用於證明上述結論的定理證明見[8, 第 7 章第 1 節定理 4, p. 7]²²．由於多項式的複合由多項式加法和乘法組成，而多項式函式的複合（定義為函式的複合）亦可在環上平行地展開，故上述對映也保持了複合運算．

多項式和多項式函式的這一關係為多項式在任意點的求值操作提供了理論基礎．在多項式函式的 $t$ 處求值也被稱為用 $t$ 代入多項式 $A$，記為 $A(t)$．

插值

需要注意的是，前述結論的逆命題不一定成立，即環 $R$ 上的一個多項式函式不一定與 $R$ 上多項式一一對應．一個經典反例是有限域上多項式環 ${\mathbb{Z}}_p[x]$ 上的 $(x+1{)}^p\equiv x^p+1(\mathrm{mod}p)$．我們指出其成立的一個充分條件是 $R$ 是無限整環²³．

該逆命題的本質是透過多項式函式的所有函式值反過來確定多項式（係數）的過程．而如果在確定時只使用一部分函式值，就是所謂的多點插值過程．很多時候待求多項式的次數是已知的，這在相當程度上縮小了待定多項式的範圍．我們指出，只要 $R$ 是整環，就可以透過多項式任意 $n$ 個不同位置的點值確定一個至多 $n-1$ 次的多項式．

證明的關鍵是用帶餘除法²⁴討論多項式根與其一次因式的關係，即多項式餘式定理或小 Bézout 定理[8, 第 7 章第 6 節定理 6, p. 35]．

需要強調，上述結論只能說明，只有那些確實可透過多項式 $n$ 點求值得到的點值組，才能與多項式建立保持加法和乘法運算的雙射．對於任取的 $R$ 中的一個 $n$ 元組，則無法確定是否存在與之對應的多項式．

求值與插值的線性表示

到這裡，多項式的求值與插值的線性表示已經呼之欲出了．若將交換環 $R$ 上的至多 $n-1$ 次的多項式看做模²⁵ $R^n$ 中的一個向量，則多項式多點求值的過程就是 $R^n$ 上的一個線性變換（記為 $V$），其矩陣表示正是 Vandermonde 矩陣．由先前對求值的討論，我們還知道其保持模 $R^n$ 上的加法運算，並將結果不超過 $n-1$ 次的卷積運算轉化為 $R^n$ 上的逐項乘法．

若將對 $R$ 的要求加強到域，注意到 Vandermonde 行列式在引數互不相同時非零，故此時 $V$ 可逆．事實上，Vandermonde 矩陣求逆的過程就是 Lagrange 插值，從而給出了將任意 $R^n$ 中向量還原為其對應多項式的構造方法²⁶．這與數論中中國剩餘定理的構造有異曲同工之妙[11, section: A perspective from linear algebra]．

求值與插值 - 小結

本節的核心是多項式和多項式函式的區別與聯絡，兩個方向的”確定”分別給出了多項式多點求值和多點插值的理論基礎．

由於多項式和多項式函式這種若即若離的關係，往往在記號上也有意無意地混淆了它們，某些情況下加大了區分的難度．本篇使用的記號體系將盡量用單個字母 $A$ 表示多項式，用 $A(x)$ 表示多項式 $A$ 對應的多項式函式，用 $A(t)$ 表示用 $t$ 代入多項式 $A$，同時在需要區分的符號前加上適當的名詞．

4.3 環上的 DFT

環上的 DFT - 導言

前面討論了在多項式任意點處求值插值的基本原理，但 DFT/FFT 的執行只需在單位根處求值和插值．本節將進一步放寬對多項式環的限制，介紹定義在有主要單位根的環上的一般的 DFT 及其快速演演算法 FFT．

本節內容主要參考了 [12, Secs. 2, pp. 983–984] 和 [13]．

環上的單位根 - 兩個定義

定義環 $R$ 上的 $n$ 次主要單位根（principal n-th root of unity）${\omega }_n$ 是滿足如下條件的環 $R$ 中的一個元素： $${\omega }_n^n=1;\sum _{k=0}^{n-1}{\omega }_n^{ik}=0,i=1,2,\dots ,n-1$$ 由定義立得 $\sum _{k=0}^{n-1}{\omega }_n^{ik}=[n\mid i]n$，可見主要單位根的定義即求和引理本身．

定義環 $R$ 上的 $n$ 次本原單位根（primitive n-th root of unity）${\omega }_n$ 是滿足如下條件的環 $R$ 中的一個元素： $${\omega }_n^n=1;{\omega }_n^k\ne 1,k=1,2,\dots ,n-1$$

本原單位根在許多情況下與主要單位根等價，但亦非完全相同．

環上的單位根 - 區別與聯絡

Proposition 1 若環 $R$ 是整環，則環上的本原單位根也是一個主要單位根．

Proof. 對任意正整數 $i\in [n]$，令 $S=\sum _{k=0}^n{\omega }_n^{ik}$，則 $(1-{\omega }_n^i)S=S-{\omega }_n^{i}S=1-{\omega }_n^n=0$，由 $R$ 是整環及 ${\omega }_n^i=0$ 得 $S=0$．

環上的單位根 - 區別與聯絡

Proposition 2 若環 $R$ 的特徵²⁷ $\mathrm{char}R$ 滿足 $\mathrm{char}R\nmid n$，則環上的主要單位根也是一個本原單位根．

Proof. 反證．若存在一正整數 $i\in [n]$ 使得 $n$ 次主要單位根 ${\omega }_n$ 的某一冪次 ${\omega }_n^i=1$，則 $\sum _{k=0}^{n-1}{\omega }_n^{ik}=\sum _{k=0}^{n-1}{1}^k=n1\ne 0$，與主要單位根定義矛盾．

環上的單位根 - 其它性質

分別根據定義和數論中的 Bézout 定理，容易證明主要單位根和本原單位根的消去引理：設 ${\omega }_n$ 是環上的 $n$ 次主要（本原）單位根，則 ${\omega }_n^k$ 是環上的 $\frac{n}{gcd(n,k)}$ 次主要（本原）單位根．

注意到 ${\left({\omega }_{2n}^n\right)}^2=1⟹({\omega }_{2n}^n+1)({\omega }_{2n}^n-1)=0$，故整環上的 $2n$ 次主要單位根或本原單位根 ${\omega }_{2n}$ 均滿足 ${\omega }_{2n}^n=-1$．

環上的 $n$ 次主要單位根或本原單位根 ${\omega }_n$ 都存在逆元 ${\omega }_n^{-1}={\omega }_n^{n-1}$．

環上的 DFT

設 ${\omega }_n$ 是環 $R$ 上的一個 $n$ 次主要單位根，其對應的 $R^n$ 上的 $n$ 點 DFT 定義為線性對映（或矩陣） $$F={\left({\omega }_n^{ij}\right)}_{(i,j)\in n\times n}=\left(\begin{array}{cccc}1& 1& \dots & 1\\ 1& {\omega }_n& \dots & {\omega }_n^{n-1}\\ ⋮& ⋮& \ddots & ⋮\\ 1& {\omega }_n^{n-1}& \dots & {\omega }_n^{(n-1)(n-1)}\end{array}\right)$$ 若 $n1$ 在環 $R$ 中存在乘法逆元，則 $R^n$ 上的 DFT 可逆，其逆對映為²⁸ $$F^{-1}=(n1{)}^{-1}{\left({\omega }_n^{-ij}\right)}_{(i,j)\in n\times n}=(n1{)}^{-1}\left(\begin{array}{cccc}1& 1& \dots & 1\\ 1& {\omega }_n^{-1}& \dots & {\omega }_n^{-(n-1)}\\ ⋮& ⋮& \ddots & ⋮\\ 1& {\omega }_n^{-(n-1)}& \dots & {\omega }_n^{-(n-1)(n-1)}\end{array}\right)$$

環上的 FFT

除無法在任意環上使用 ${\omega }_{2n}^n=-1$，現有的抽象已足夠讓我們寫出與之前類似的推導 FFT 的過程．

設 $R$ 是任意環，${\omega }_{2n}$ 是環 $R$ 上的一個 $2n$ 次主要單位根．由主要單位根的消去引理，${\omega }_{2n}^2$ 是環 $R$ 上的 $\frac{2n}{gcd(2n,2)}=n$ 次主要單位根．方便起見，記 ${\omega }_n={\omega }_{2n}^2$．

考慮對模 $R^{2n}$ 中向量 $\mathit{a}={(a_0,a_1,\dots ,a_{2n-1})}^T$ 做 ${\omega }_{2n}$ 對應的 $2n$ 點 DFT 變換，得到向量 $\hat{\mathit{a}}={({\hat{a}}_0,{\hat{a}}_1,\dots ,{\hat{a}}_{2n-1})}^T$，其中 ${\hat{a}}_i=\sum _{k=0}^{2n-1}{a}_k{\omega }_{2n}^{ik}$．將其奇偶分項為 ${\hat{a}}_i=A_0(i)+A_1(i){\omega }_{2n}^i$，其中 $$\begin{array}{rl}{A}_0(i)& :=\sum _{k=0}^{n-1}{a}_{2k}{\omega }_n^{ik}=\sum _{k=0}^{n-1}{a}_{2k}{\omega }_{2n}^{2ik}\\ A_1(i)& :=\sum _{k=0}^{n-1}{a}_{2k+1}{\omega }_n^{ik}=\sum _{k=0}^{n-1}{a}_{2k+1}{\omega }_{2n}^{2ik}\end{array}$$

環上的 FFT

注意到 $A_0(i)=A_0(n+i)$，$A_1(i)=A_1(n+i)$，故 $$\begin{array}{rl}{\hat{a}}_i& =A_0(i)+A_1(i){\omega }_{2n}^i\\ {\hat{a}}_{n+i}& =A_0(i)+A_1(i){\omega }_{2n}^{n+i}\end{array},i\in [n]$$ 而 $A_0(i),i\in [n]$ 的計算即計算向量 ${\mathit{a}}_{\mathbf{0}}=(a_0,a_2,\dots ,a_{2n-2})$ 對應於 ${\omega }_n$ 的 $n$ 點 DFT 的過程；$A_1(i),i\in [n]$ 的計算即計算向量 ${\mathit{a}}_{\mathbf{1}}=(a_1,a_3,\dots ,a_{2n-1})$ 對應於 ${\omega }_n$ 的 $n$ 點 DFT 的過程．

因此，令 $n=2^m$，則只要環 $R$ 上存在 $2^m$ 次主要單位根，如上形式的 FFT 遞迴演演算法就可在 $O(n\log n)$ 的時間複雜度內快速計算 $R^n$ 上的 $n$ 點 DFT．而若 $2^{k}1,0\le k\le m$ 在 $R$ 內均存在逆元，$R^n$ 上的 $n$ 點 IDFT 也可類似地快速計算．

4.4 迴圈卷積與卷積定理

迴圈卷積與卷積定理 - 導言

上一節中，我們建立了在有主要單位根的環上的 DFT 及其快速演演算法 FFT 的相關理論，但由於放寬了對多項式環 $R[x]$ 的限制，多項式的求值插值理論無法在此直接得到應用，環上 DFT 加速多項式卷積的理論尚需重新構建．

在求值與插值部分已經介紹求值變換 $V$ 在 $R^n$ 上的線性表示，其保持 $R^n$ 上的加法運算，並將規模不超過 $n$ 的兩向量的卷積轉化為 $R^n$ 上的逐項乘法．本節中，我們來證明 DFT 變換 $F$ 也滿足同樣的性質，且將這一結果推廣到卷積規模超過 $n$ 時的情況．

迴圈卷積

設 $R$ 是交換環，我們記 $A\otimes B$ 代表模 $R^n$ 上兩個向量的逐項乘法．

設 $F$ 是 $R^n$ 上主要單位根 ${\omega }_n$ 對應的 $n$ 點可逆 DFT 變換，$A,B$ 是 $R^n$ 中兩個任意取定的向量．顯然 $F(A+B)=FA+FB$，故 DFT 變換保持 $R^n$ 上的加法運算．現在的主要問題是，求值與插值部分計算卷積的方法 $F^{-1}(FA\otimes FB)$ 在放寬環限制的 $R^n$ 上進行時，將會得到什麼結果？

令 $$\begin{array}{c}\hat{\mathit{a}}=F\mathit{a},\hat{\mathit{b}}=F\mathit{b}\\ \hat{\mathit{c}}=\hat{\mathit{a}}\otimes \hat{\mathit{b}},\mathit{c}=F^{-1}\hat{\mathit{c}}\end{array}$$ 則 $${\hat{a}}_i=\sum _{j=0}^{n-1}{a}_j{\omega }_n^{ij},{\hat{b}}_i=\sum _{k=0}^{n-1}{b}_k{\omega }_n^{ik}$$

迴圈卷積

$${\hat{c}}_i={\hat{a}}_i{\hat{b}}_i=\sum _{j=0}^{n-1}{a}_j{\omega }_n^{ij}\sum _{k=0}^{n-1}{b}_k{\omega }_n^{ik}=\sum _{j=0}^{n-1}\sum _{k=0}^{n-1}{a}_j{b}_k{\omega }_n^{i(j+k)}$$ $$c_i=\sum _{t=0}^{n-1}{\omega }_n^{-it}{\hat{c}}_t=\sum _{t=0}^{n-1}{\omega }_n^{-it}\sum _{j=0}^{n-1}\sum _{k=0}^{n-1}{a}_j{b}_k{\omega }_n^{t(j+k)}=\sum _{j=0}^{n-1}\sum _{k=0}^{n-1}{a}_j{b}_k\sum _{t=0}^{n-1}{\omega }_n^{t(j+k-i)}$$ 對最裡側的求和使用主要單位根的定義（求和引理），就有 $$c_i=\sum _{j=0}^{n-1}\sum _{k=0}^{n-1}{a}_j{b}_k[n\mid j+k-i]=\sum _{j+k\equiv imodn}{a}_j{b}_k$$ 這便是所謂的迴圈卷積．以後將 $R^n$ 中向量 $\mathit{a}$ 和 $\mathit{b}$ 的迴圈卷積記為 $\mathit{a}\ast \mathit{b}$．形象地來看，迴圈卷積中，次數超過 $n-1$ 的卷積項被模意義地地疊加到了以次數為 $0$ 為始的項上．

卷積定理

前述討論已經證明²⁹ $$F(A\ast B)=FA\otimes FB$$ 這便是一般交換環上 DFT 變換的（迴圈）卷積定理．

我們指出，當要求所作變換可逆時，卷積定理反過來也要求所作變換是一類似 DFT 對映的變換³⁰．

Exercise 5 嘗試將任意序列 DFT 兩次，觀察結果．證明你的結論．

Exercise 6 求 DFT 矩陣的行列式．儘可能縮小可行解範圍．

5 Advanced Operations

多項式全家桶 - 序言

• 有哪些？

  • 求逆、開根、對數、指數、快速冪、複合等

• 咋推的？

  • 求解思路幾乎都是倍增，時間複雜度幾乎都是大常數 $O(n\log n)$
  • Newton 迭代法是推導全家桶迭代公式的通法
  • 嚴格化需要進一步的形式冪級數理論，主講人不會

• 有啥用？

  • 當你一波操作化出生成函式發現不會求係數

• 怎麼講？

  • 受限於篇幅和主講人能力，我們講不完
  • 只講求逆和對數，其餘請左轉 OI Wiki

多項式求逆

給定一多項式 $A(x)$，求解滿足 $A(x)B(x)=1(\mathrm{mod}{x}^n)$ 的多項式 $B(x)$．

多項式逆元存在的充分必要條件是其常數項非零（這是因為邊界條件 $b_0=\frac{1}{a_0}$），若存在則在模意義下一定唯一．這結論可直接由下述求解方法得到．不失一般性，只研究 $A(x)$ 的次數至多為奇數 $2n-1$ 的情況．設 $$\begin{array}{rl}A(x)& =A_0(x)+x^n{A}_1(x)\\ B(x)& =B_0(x)+x^n{B}_1(x)\end{array}$$ 考慮遞迴地在已知 $$A(x)B_0(x)=A_0(x)B_0(x)=1(\mathrm{mod}{x}^n)$$ 的基礎上求解 $A(x)$ 的逆元．先指出該方法的時間複雜度為 $$T(n)=T(\frac{n}{2})+O(n\log n)=O(n\log n)$$

多項式求逆

以下簡記 $A(x)$ 為 $A$，其它多項式同理．注意到 $$\{\begin{array}{l}A{B}_0=1(\mathrm{mod}{x}^n)\\ AB=1(\mathrm{mod}{x}^n)\end{array}⟹A(B-B_0)=0(\mathrm{mod}{x}^n)$$ 由於 $Amod{x}^n$ 非零，故 $$B-B_0=0(\mathrm{mod}{x}^n)$$ 這也說明，多項式逆元在模不同 $x^n$ 下的字首保持一致．

兩邊平方得 $$B^2-2B{B}_0+B_0^2=0(\mathrm{mod}{x}^{2n})$$ 兩側同乘 $A$ 並移項得 $$B=2B_0-A{B}_0^2(\mathrm{mod}{x}^{2n})$$

多項式求逆 - 實現

Poly inv(const Poly &A){
    ll n=A.len(); Poly B(1); B[0]=inv(A[0]);
    for(ll hf=1;hf<n;hf<<=1){
        B=B*2-B*B*A.subpoly(0,hf*2); B.resize(hf*2);
    } B.resize(n);
    return B;
}

• 常數巨大的寫法，僅作演示，請勿學習．
• 正確的寫法是在 DFT 後的點值上操作多項式，請小心實現封裝．

多項式 $\ln$

給定一多項式 $A(x)$，求解滿足 $B(x)=\ln A(x)(\mathrm{mod}{x}^n)$ 的多項式 $B(x)$．

次數為 $+\mathrm{\infty }$ 的形式冪級數的 $\ln$ 存在的充分必要條件為其常數項非零（這是因為邊界條件 $b_0=\ln a_0$），同樣一旦存在則唯一．注意到僅整數 $a_0=1$ 時，$\ln a_0$ 可取得整數，故合理的 $a_0$ 只能是 $1$．

推導是容易的．方程兩側同時求導得 $$B^{\prime }(x)=\frac{A^{\prime }(x)}{A(x)}(\mathrm{mod}{x}^{n-1})$$ 兩側再積分得 $$B(x)=\int \frac{A^{\prime }(x)}{A(x)}\mathrm{d}x+C(\mathrm{mod}{x}^n)$$ 其中 $C=\ln a_0$．多項式求逆、求導、積分即可．時間複雜度 $O(n\log n)$．

多項式 $\ln$ - 實現

Poly drv(Poly A){ // derivative
    for(ll i=0;i<A.len();i++) A[i]=(i+1)*A[i+1]_;
    A.pop_back(); return A;
}
Poly itg(Poly A,ll c){ // integral
    A.push_back(0); for(ll i=A.len();i>=1;i--) A[i]=A[i-1]*inv(i)_;
    A[0]=c; return A;
}
Poly ln(const Poly &A){
    return itg((drv(A)*inv(A)).subpoly(0,A.len()-1),0/*log(A[0])*/);
}

Acknowledgements

• 感謝 keke_046 學長教授 FFT、集合冪級數與生成函式．微言大義，博大精深，至今仍在消化．

• 感謝隊友 ItzDesert 提供位運算典題一道並提供內容編排建議．

題單

主講人練題少，僅供參考．

模板

• 洛谷 P3803 【模板】多項式乘法（FFT）

• 洛谷 P6091 【模板】原根

• 洛谷 P4238 【模板】多項式乘法逆

• 洛谷 P5205 【模板】多項式開根

• 洛谷 P4725 【模板】多項式對數函式（多項式 ln）

• 洛谷 P4726 【模板】多項式指數函式（多項式 exp）

• 洛谷 P5245 【模板】多項式快速冪

大整數乘法

• 洛谷 P1919 【模板】A*B Problem 升級版（FFT 快速傅立葉變換）

基礎計數

• SPOJ-TSUM Triple Sums

• BZOJ3513-MUTC2013 Idiots

  上面兩道題都是 OGF 消序，較 EGF 消序困難．一般的方法是使用 Polya 計數原理．

位運算

• ABC291G OR Sum

  也是濾波器的應用．

分治 FFT

• 百度之星 2023 初賽第二場 T8

  容斥後需要計算若干一次多項式乘積，分治 NTT 即可．

• 百度之星 2023 初賽第一場 T6

  求 $\sum _{i>j}\frac{a_i}{a_i+a_j}$．值域次數化後 CDQ 處理偏序．

• 洛谷 P4721 【模板】分治 FFT

  CDQ 偏序化處理前後項依賴．也可解生成函式方程再多項式求逆．

字串

• ABC196F Substring 2

• 洛谷 P4173 殘缺的字串

  帶單字元萬用字元的字串匹配．體會設定字串距離函式的方法．

全家桶

• 洛谷 P4721 【模板】分治 FFT

  體驗生成函式簡化複雜數列遞推的威力．

• 洛谷 P4389 付公主的揹包

  解決此問題的方法也可用於分拆數計算．

• 洛谷 P4841 【集訓隊作業2013】城市規劃

  EGF 消序劃分典題．體驗生成函式風格的計數方法．

• LOJ6538 烷基計數 - 加強版 - 加強版

  對生成函式使用的 Polya 計數．

  Polya ex. 苯環碳接 $-{\mathrm{C}}_4{\mathrm{H}}_9$ 同分異構體計數

其它

• QOJ5748-UCUP2023-Stage7-K Determinant, or...?

  看似人畜無害的 ${\left(a_{i\mathrm{or}j}\right)}_{(i,j)\in n\times n}$ 行列式求值，解法的背後卻潛藏 FMT 的思想．此類行列式與各類反演、卷積有密切聯絡，我們給出一篇研究了其與偏序集上反演的關係的參考文獻[15]，歡迎討論．

References

[1]

OI-Wiki, “拉格朗日插值.” https://oi-wiki.org/math/numerical/lagrange/.

[2]

張築生, “數學分析新講（重排本）（第二冊）,” 2nd ed.北京: 北京大學出版社, 2021, pp. 256–262.

[3]

楊樹森, “三角函式的嚴格定義.” https://zhuanlan.zhihu.com/p/58814328/, 2023.

[4]

T. H. Cormen, C. E. Leiserson, R. L. Rivest, and C. Stein, “演演算法導論（原書第三版）,” 北京: 機械工業出版社, 2013.

[5]

OI-Wiki, “離散對數.” https://oi-wiki.org/math/number-theory/discrete-logarithm/.

[6]

OI-Wiki, “原根.” https://oi-wiki.org/math/number-theory/primitive-root/.

[7]

R. Agarwal and C. Burrus, “Fast convolution using fermat number transforms with applications to digital filtering,” IEEE Transactions on Acoustics, Speech, and Signal Processing, vol. 22, no. 2, pp. 87–97, 1974, doi: 10.1109/TASSP.1974.1162555.

[8]

丘維聲, “高等代數 下冊,” 3rd ed.北京: 高等教育出版社, 2015.

[9]

OI-Wiki, “多項式與生成函式簡介.” https://oi-wiki.org/math/poly/intro/.

[10]

Wikipedia, “Formal power series.” https://en.wikipedia.org/wiki/Formal_power_series.

[11]

Wikipedia, “Lagrange polynomial.” https://en.wikipedia.org/wiki/Lagrange_polynomial.

[12]

M. Fürer, “Faster integer multiplication,” SIAM Journal on Computing, vol. 39, no. 3, pp. 979–1005, 2009, doi: 10.1137/070711761.

[13]

Wikipedia, “Discrete fourier transform over a ring.” https://en.wikipedia.org/wiki/Discrete_Fourier_transform_over_a_ring.

[14]

I. Baraquin and N. Ratier, “Uniqueness of the discrete fourier transform,” Signal Processing, vol. 209, p. 109041, 2023.

[15]

H. S. Wilf, “Hadamard determinants möbius functions, and the chromatic number of a graph,” 1968.

Footnotes

1. 一種常見的證法是使用 Vandermonde 行列式證明矩陣可逆．後面會介紹多項式環風格的證明．↩︎

2. 有的文獻定義 ${\omega }_n:=e^{-\frac{2\pi }{n}i}$，或是因為訊號處理領域常用 IDFT 將訊號時域取樣資料變為頻域資訊．事實上，DFT/IDFT 的說法也常有反轉，但這只是形式問題．↩︎

3. 這三個引理是《演演算法導論》[4] 引入的．↩︎

4. 這也表明適當歸一化後的 DFT 矩陣是一個酉矩陣．↩︎

5. NTT 原理需較多筆墨，稍後介紹．↩︎

6. 否則只有 $p$ 個不同元素的 ${\mathbb{Z}}_p$ 中根本取不到 $n$ 個不同位置的點值．後面會深入討論．↩︎

7. 多項式求逆等多項式進階操作，我們後續講解．↩︎

8. 對質數 $p$，$\phi (p)=p-1$．故 Euler 定理是 Fermat 小定理的一個推廣．↩︎

9. ExBSGS 求解離散對數 [5] 的推導與此相似．↩︎

10. 部分證明稍複雜，我們略過處理．感興趣的同學請參考 OI Wiki [6]．↩︎

11. 這是後文所述定理“整環上的本原單位根也是主要單位根”在一般環上的一個反例．↩︎

12. 網傳此模數由 UOJ 站長 vfleaking 提出並推廣．在所有需要取模的題目中使用該模數，可使選手無法透過模數判斷題目的做法．↩︎

13. 本篇中環的定義包含乘法單位元，即麼環．↩︎

14. 無零因子的交換麼環，稍後介紹．↩︎

15. 可用前述多項式乘積次數公式證明．↩︎

16. 試試在 $\mathbb{Z}[x]$ 上用 $2x+1$ 去除 $x+1$．↩︎

17. 首項為 $1$ 的多項式．↩︎

18. 函式和對映幾乎是等價名詞．有時函式特指值域包含於複數域 $\mathbb{C}$ 的對映．↩︎

19. 為良定義 $x^0$，環 $R$ 必須有單位元．↩︎

20. 這裡再次涉及環 $R$ 上的收斂問題．由於實踐中只關心形式冪級數的前有限項，後續討論係數-點值-係數法轉換卷積時不需要用到冪級數理論，可以迴避．↩︎

21. 這種保持結構不變的對映被稱為同態（homomorphism）．↩︎

22. 其證明瞭域上一元多項式環的通用性質．仿照該證明應可證明環上的版本，從而證明這一同態關係．↩︎

23. 該定理是下方高亮定理的一個自然的推論．↩︎

24. 由於一次因式均為首一多項式，可以在整環上對其使用帶餘除法．↩︎

25. 模是定義在環上的“線性空間”．↩︎

26. Lagrange 插值的構造用到了除法，且行列式非零推出矩陣可逆僅在域上的線性空間中適用，因此必須要求 $R$ 是域．↩︎

27. 使得 $m1=\sum _{k=0}^{m-1}1=0$ 的最小正整數 $m$．不存在則記為 $0$．$\mathrm{char}\mathbb{C}=0$，$\mathrm{char}{\mathbb{Z}}_p=p$．可以證明域的特徵一定是 $0$ 或一質數[8, 第 7 章第 11 節定理 3，p. 70]．↩︎

28. 證明使用主要單位根的定義（求和引理）即可．↩︎

29. 雖然前面用到了 $F$ 可逆的要求，但該定理在 $F$ 不可逆時也成立．只需類似地驗證兩邊相等即可．↩︎

30. 具體來說，該變換隻能是 DFT 矩陣的某個行置換．在 ${\mathbb{C}}^n$ 上的證明可參見[14]，主講人目前正在研究整環上的版本，歡迎討論．↩︎