05.3 1965：譜方法、FFT 與指數級精確度#

一段塵封一百六十年的手稿#

1965 年，IBM 的 James Cooley 與貝爾實驗室的 John Tukey 發表了〈An Algorithm for the Machine Calculation of Complex Fourier Series〉，這篇只有三頁的論文，日後被譽為二十世紀最重要的數值演算法。故事卻要從更早說起。

1805 年，Gauss 在處理小行星 Pallas 軌道內插問題時，已經寫下了本質上等同於 FFT 的離散傅立葉變換快速演算法。但這份手稿以拉丁文寫成，未曾發表，直到 1866 年才收入《Gauss 全集》。Gauss 甚至先於 Cooley-Tukey 一百六十年就注意到：把長度為 \(NM\) 的變換拆成 \(N\) 個長度 \(M\) 的子變換，可以把運算量從 \(O((NM)^2)\) 降到 \(O(NM(N+M))\)。只是在那個沒有電子計算機的年代，計算量多十倍與少十倍，不過是計算員多熬幾個夜晚的差別，革命性的想法竟無人聞問。

歷史的弔詭在於：同一個演算法，在 1805 年只是天才的筆記，在 1965 年卻成為改變世界的發明。差別不在數學，而在「計算的成本」發生了革命。

冷戰、地震波與核試驗禁令#

1963 年，美蘇簽署《部分禁止核試驗條約》（Limited Test Ban Treaty），禁止在大氣層、太空與水下核試。此後地下核試的偵測成為談判的核心：美國是否願意簽署更進一步的全面禁試，取決於能否用地震儀可靠地偵測並識別蘇聯境內的地下核爆。地震信號分析的核心工具是功率譜估計，而功率譜估計的核心運算正是離散傅立葉變換。

Tukey 當時是總統科學顧問委員會（PSAC）的成員，所在的委員會立場傾向支持核試禁令，因為他們判斷大氣中放射性落塵的危害被低估了。Tukey 在委員會討論中聽到：若能有效分析地震頻譜，偵測阈值就能降低，禁試條約便可驗證。但當時用電腦算 DFT，\(N\) 個資料點要 \(N^2\) 次複數乘法，取 \(N\) 稍大便不切實際。

傳說中，Tukey 在一次會議上把這個「如何算得更快」的問題帶給了數學社群。1976 年 Tukey 與 Cooley 回憶此事時提到，這個想法的出發點確實與 PSAC 的禁試驗證研究有關。Gordon Welchman 等人當時參與的 VELA 計畫——在衛星與地面布設核爆偵測儀器——提供了制度與經費的動力。Cooley 當時在 IBM 約克敦高地的研究中心任職，他推導並程式化了這個演算法。1965 年論文發表後，Cooley 與 Stockham 又整理了既有技術，確認 Gauss 的先例，並澄清了 Sande 也曾獨立提出分解方式。

演算法的核心想法可以寫成：把 \(N=N_1 N_2\) 點的 DFT 拆開，

先對每個子序列做小變換，再合併。當 \(N = 2^m\) 時，運算量從 \(N^2\) 降到 \(\tfrac{N}{2}\log_2 N\)。取 \(N = 2^{20}\)，這是從「一兆次」到「兩千萬次」的差別——不是改良，是把不可能變成日常。

從震測到整個數位文明#

FFT 問世後，應用爆炸性擴張。1970 年代，德州儀器等公司推出專用 DSP 晶片，其存在的核心理由就是快速做捲積與傅立葉變換；語音的數位化、電話線路的回音消除（AT&T 的迴路）、石油探勘的震測資料處理，全都建立在 FFT 上。經濟上，這是一個由國防預算催生、再被電信產業商業化的典型冷戰技術轉移故事。Cooley 本人在 1990 年代回顧時整理了這段應用史，指出 FFT 之所以能普及，與 1970 年代多位學者——包括他本人與 John Lewis——撰寫的應用綜述文章有關。

譜方法：指數級精確度#

FFT 之所以值得在數值分析的教科書中佔一章，是因為它開啟了譜方法（spectral methods）這一整個領域。

傳統的有限差分法用局部樣本近似導數，收斂階是固定的（一階、二階、四階……）。譜方法則改用全域基底——傅立葉級數或 Chebyshev 多項式——展開解，用 FFT 在物理空間與譜空間之間快速往返。若函數是解析的，譜方法的誤差不是按多項式速度下降，而是按指數速度下降：

這就是「指數級精確度」的含義：多加一個基函數，精確度不是提升一點，而是乘上一個固定倍率。對光滑問題，幾十個自由度往往可以達到傳統方法上千個格點的精確度。

1969 至 1970 年代初，Steven Orszag 將譜方法推廣到流體力學模擬，用偽譜法（pseudospectral method，在物理空間算非線性項、在譜空間算導數）解 Navier–Stokes 與湍流問題；Heinz-Otto Kreiss 與其合作者則建立了譜方法的穩定性理論，說明為什麼這些方法在剛性問題上依然可靠。偽譜法之所以可行，正是因為 FFT 把物理空間與譜空間之間的往返成本降到 \(O(N \log N)\)。沒有 FFT，譜方法只是一個理論優美的想法；有了 FFT，它成為氣候模擬、湍流計算與量子化學的主力。

Chebfun：把「無窮」變成資料結構#

譜方法的思想在二十一世紀演化成一個更激進的形式。2000 年代初，牛津大學的 Trefethen 與其學生開發了 Chebfun：把定義在區間上的函數表示為自動截斷到機器精確度的 Chebyshev 展開——一個函數變成一個長度可變的係數向量。微分變成對係數的矩陣運算，積分、求根、求極值都有對應的高精度演算法。在 Chebfun 中，\([\text{f},\text{g}] = \min(\text{f},\text{g})\) 這種寫法是合法的，因為兩個函數的最小值本身就是一個新的函數。

Chebfun 是譜方法思想的極致體現：當誤差是指數下降的，你可以把「逼近一個函數」做到像「儲存一個數字」一樣可靠。這也呼應了本章的主題：1965 年那篇三頁論文的真正遺產，不只是快了一個數量級的演算法，而是一種新的思維方式——光滑函數的世界裡，存在著比想像中便宜得多的計算。

小結#

FFT 的故事串起了 Gauss 的未發表手稿、冷戰的地震偵測、PSAC 的核試政策辯論、IBM 的研究文化與 DSP 產業的興起。數學史上很少有演算法像 FFT 這樣：思想上它早熟了 160 年，實現上它等待了「計算成本崩跌」這個經濟條件，應用上它又反過來催生了一個產業。譜方法與 Chebfun 則說明，這個演算法最深的價值在於它讓「全域、高階、指數精確」的逼近變得可行。下一章，我們轉向另一個古老的問題：如何解方程。