--- title: Fast Fourier Transform module: Mathematical Algorithms moduleNumber: 10 lessonNumber: 6 order: 1006 summary: | Multiplying two degree-$n$ polynomials by the schoolbook method costs $\Theta(n^2)$. Evaluating them at the **$n$-th roots of unity** turns multiplication into pointwise products, and the **Cooley–Tukey FFT** computes all those evaluations in $\Theta(n\log n)$ by splitting even and odd coefficients. The inverse FFT interpolates back, giving $\Theta(n\log n)$ polynomial and big-integer multiplication. topics: [Number Theory] sources: - book: CLRS ref: "Ch. 30 — Polynomials and the FFT" - book: Erickson ref: "Ch. — (the Fast Fourier Transform)" - book: Skiena ref: "§ — Divide and Conquer / Numerical Algorithms" practice: - title: 'Multiply Strings' slug: multiply-strings difficulty: Medium - title: 'Add Two Numbers' slug: add-two-numbers difficulty: Medium - title: 'Maximum Score Of Spliced Array' slug: maximum-score-of-spliced-array difficulty: Hard --- [Karatsuba](/algorithms/divide-and-conquer/fast-multiplication) cut integer multiplication below quadratic by trading one multiplication for a few additions, reaching $\Theta(n^{\log_2 3})\approx\Theta(n^{1.585})$. This lesson goes further, to $\Theta(n\log n)$, with the **Fast Fourier Transform**. The route is indirect. The obstacle to fast multiplication is the **coefficient representation** of a polynomial, in which the product is a convolution costing $\Theta(n^2)$. If instead we represent a polynomial by its **values** at enough points, the product is computed pointwise in linear time. The FFT converts between the two representations in $\Theta(n\log n)$, and the points it chooses — the complex **roots of unity** — are what make that conversion fast. ## The problem: polynomial multiplication is convolution A polynomial of degree $n-1$ is given by its coefficient vector $a = (a_0,\dots,a_{n-1})$, meaning $A(x) = \sum_{j=0}^{n-1} a_j x^j$. The product $C(x) = A(x)\,B(x)$ has coefficients $$ c_k = \sum_{j=0}^{k} a_j\, b_{k-j}, \qquad k = 0,1,\dots,2n-2. $$ This sum is the **convolution** $a \ast b$. Computed directly, each of the $\sim 2n$ output coefficients is a sum over up to $n$ products, so the cost is $\Theta(n^2)$ — the schoolbook method, and exactly the cost of grade-school long multiplication when $a$ and $b$ are the digit vectors of two integers. Convolution is everywhere (signal filtering, string matching, probability of sums), so a faster algorithm pays off far beyond polynomials. ## A polynomial is its values A degree-$(n-1)$ polynomial is **uniquely determined** by its values at any $n$ distinct points — the **interpolation** fact that $n$ points fix a degree-$(n-1)$ curve (two points fix a line, three fix a parabola). So besides the coefficient vector there is a second, equally faithful description: the **point-value representation**, a list of pairs $\{(x_0, A(x_0)), \dots, (x_{n-1}, A(x_{n-1}))\}$. > **Definition (Point-value representation).** A point-value representation of a degree > $\le n-1$ polynomial $A$ is a set of $n$ pairs $\{(x_i, y_i)\}_{i=0}^{n-1}$ with the > $x_i$ distinct and $y_i = A(x_i)$. By the uniqueness of interpolation, it determines > $A$ exactly. In this representation, multiplication is cheap. If $A$ and $B$ are both sampled at the **same** points $x_i$, then the product $C = AB$ is sampled there by a single multiply per point: $$ C(x_i) = A(x_i)\, B(x_i). $$ That is $\Theta(n)$ work, not $\Theta(n^2)$. The catch is that the product has degree up to $2n-2$, so $2n-1$ points are needed to determine it; we must therefore sample $A$ and $B$ at (at least) $2n$ points up front, padding their coefficient vectors with zeros to length $2n$. The strategy is now a three-step detour: $$ % caption: Multiply via point-value form: evaluate both polynomials (DFT), multiply % pointwise, interpolate back (inverse DFT) \begin{tikzpicture}[ every node/.style={font=\small}, box/.style={draw, minimum width=24mm, minimum height=9mm, inner sep=2pt, align=center}, >=stealth, x=1cm, y=1cm] \definecolor{acc}{HTML}{2348F2} \node[box] (coef) at (0,0) {input vectors\\ $a$ and $b$}; \node[box] (val) at (5.8,0) {values\\ $A(x_i)$ and $B(x_i)$}; \node[box] (cval) at (5.8,-2.4) {product values\\ $A(x_i)\,B(x_i)$}; \node[box] (ccoef) at (0,-2.4) {output vector\\ $c$}; \draw[->, acc, thick] (coef) -- node[above, font=\footnotesize, acc] {evaluate (FFT)} (val); \draw[->, acc, thick] (val) -- node[right, font=\footnotesize, acc] {multiply} (cval); \draw[->, acc, thick] (cval) -- node[below, font=\footnotesize, acc] {interp\/olate (inverse FFT)} (ccoef); \end{tikzpicture} $$ Evaluate, multiply pointwise, interpolate. The middle step is linear, so the total cost depends on the two end steps: evaluating a polynomial at $2n$ points, and recovering coefficients from $2n$ values. Naively each is $\Theta(n^2)$ (evaluating at one point by Horner's rule is $\Theta(n)$, times $2n$ points), so a poor choice of points gains nothing. For **one special set of points**, both steps drop to $\Theta(n\log n)$. ## The roots of unity The special points are the **$n$-th roots of unity**: the $n$ complex numbers $\omega$ with $\omega^n = 1$. They are equally spaced around the unit circle in the complex plane, $$ \omega_n = e^{2\pi i / n}, \qquad \omega_n^k = e^{2\pi i k / n}\ \ (k = 0,1,\dots,n-1), $$ where $\omega_n$ is the **principal** $n$-th root and its powers $\omega_n^0,\omega_n^1, \dots,\omega_n^{n-1}$ enumerate all $n$ of them. Geometrically, multiplying by $\omega_n$ rotates a point by $2\pi/n$ radians; $n$ such rotations come full circle back to the start, which is the statement $\omega_n^n=1$. $$ % caption: The eight $8$-th roots of unity on the unit circle, labelled by their exponent % $k$ in $\omega_8^k=e^{2\pi i k/8}$, equally spaced by $\pi/4$; squaring % doubles the angle, folding the eight points onto the four $4$-th roots \begin{tikzpicture}[ every node/.style={font=\small}, >=stealth, x=1cm, y=1cm] \definecolor{acc}{HTML}{2348F2} \draw[black] (-2.2,0) -- (2.2,0); \draw[black] (0,-2.2) -- (0,2.2); \node[font=\footnotesize, anchor=north west] at (2.05,-0.05) {Re}; \node[font=\footnotesize, anchor=south west] at (0.12,2.05) {Im}; \draw[thick] (0,0) circle (1.8cm); \foreach \k in {0,...,7}{ \fill[acc] (\k*45:1.8) circle (2.2pt); } \node[acc, anchor=west, font=\footnotesize] at (1.95,0.05) {$k{=}0$}; \node[acc, anchor=south west, font=\footnotesize] at (1.4,1.4) {$k{=}1$}; \node[acc, anchor=east, font=\footnotesize] at (-0.18,1.95) {$k{=}2$}; \node[acc, anchor=south east, font=\footnotesize] at (-1.4,1.4) {$k{=}3$}; \node[acc, anchor=east, font=\footnotesize] at (-1.95,0.05) {$k{=}4$}; \draw[acc, thick, ->] (0.6,0.55) arc (42:88:0.82); \node[acc, font=\footnotesize, anchor=south west] at (0.6,0.5) {square}; \end{tikzpicture} $$ Two algebraic properties of these points make the FFT work. > **Lemma (Halving / squaring).** If $n$ is even, the squares of the $n$-th roots of unity > are the $\tfrac{n}{2}$-th roots of unity, each appearing **twice**: > $$\parens{\omega_n^{k}}^2 = \omega_n^{2k} = \omega_{n/2}^{k}, \qquad \omega_n^{k + n/2} = -\,\omega_n^{k}.$$ > **Proof.** Since $\omega_n = e^{2\pi i/n}$, we have $\omega_n^{2} = e^{4\pi i/n} = > e^{2\pi i/(n/2)} = \omega_{n/2}$, so squaring a root halves the order. And > $\omega_n^{n/2} = e^{\pi i} = -1$, giving $\omega_n^{k+n/2} = \omega_n^k\cdot\omega_n^{n/2} > = -\omega_n^k$: the second half of the roots is the **negation** of the first half. $\qed$ The $n$ evaluation points, when squared, become only $n/2$ **distinct** points, and they pair up as $\pm$ of each other. These two facts turn one size-$n$ evaluation into two size-$(n/2)$ evaluations. ## The DFT and the Cooley–Tukey split Evaluating $A$ at all $n$ roots of unity is, by definition, the **Discrete Fourier Transform** of its coefficient vector. > **Definition (DFT).** The DFT of $a = (a_0,\dots,a_{n-1})$ is the vector > $\hat a = (\hat a_0,\dots,\hat a_{n-1})$ with > $$\hat a_k = A(\omega_n^{k}) = \sum_{j=0}^{n-1} a_j\,\omega_n^{kj}.$$ Computed as written, that is $n$ evaluations of an $n$-term sum: $\Theta(n^2)$. The FFT computes the same vector in $\Theta(n\log n)$ by a divide-and-conquer split on the **parity** of the coefficient index. Separate $A$ into its even- and odd-indexed coefficients and factor an $x$ out of the odd part: $$ A(x) = \underbrace{\parens{a_0 + a_2 x^2 + a_4 x^4 + \cdots}}_{A_{\text{even}}(x^2)} + x\,\underbrace{\parens{a_1 + a_3 x^2 + a_5 x^4 + \cdots}}_{A_{\text{odd}}(x^2)}, $$ where $A_{\text{even}}$ and $A_{\text{odd}}$ are each polynomials of half the degree, in the variable $x^2$. Repeating the split recursively builds a tree: at each level the input is sorted by parity into two halves, dark for the **even** indices and light for the **odd**, until single coefficients sit at the leaves. $$ % caption: The even/odd recursive split for $n=8$. Each level partitions the indices by % parity — dark tint = even-indexed coefficients, light tint = odd-indexed — so % after $\log_2 n$ levels every leaf is one coefficient $a_j$ \begin{tikzpicture}[ every node/.style={font=\footnotesize}, ev/.style={draw=acc, fill=acc!28, minimum width=30mm, minimum height=7mm, inner sep=2pt}, od/.style={draw=acc, fill=acc!8, minimum width=30mm, minimum height=7mm, inner sep=2pt}, lf/.style={draw=acc, fill=acc!14, minimum width=7mm, minimum height=7mm, inner sep=1pt}, >=stealth, x=1cm, y=1cm] \definecolor{acc}{HTML}{2348F2} \node[draw=acc, fill=acc!14, minimum width=62mm, minimum height=7mm] (root) at (4,3) {a0 a1 a2 a3 a4 a5 a6 a7}; \node[ev] (L) at (1.6,1.6) {a0 a2 a4 a6}; \node[od] (R) at (6.4,1.6) {a1 a3 a5 a7}; \node[ev, minimum width=14mm] (LL) at (0.4,0.3) {a0 a4}; \node[od, minimum width=14mm] (LR) at (2.8,0.3) {a2 a6}; \node[ev, minimum width=14mm] (RL) at (5.2,0.3) {a1 a5}; \node[od, minimum width=14mm] (RR) at (7.6,0.3) {a3 a7}; \node[lf] at (-0.1,-1) {a0}; \node[lf] at (0.9,-1) {a4}; \node[lf] at (2.3,-1) {a2}; \node[lf] at (3.3,-1) {a6}; \node[lf] at (4.7,-1) {a1}; \node[lf] at (5.7,-1) {a5}; \node[lf] at (7.1,-1) {a3}; \node[lf] at (8.1,-1) {a7}; \draw[->, acc] (root) -- (L); \draw[->, acc] (root) -- (R); \draw[->, acc] (L) -- (LL); \draw[->, acc] (L) -- (LR); \draw[->, acc] (R) -- (RL); \draw[->, acc] (R) -- (RR); \draw[->, acc] (LL) -- (-0.1,-0.6); \draw[->, acc] (LL) -- (0.9,-0.6); \draw[->, acc] (LR) -- (2.3,-0.6); \draw[->, acc] (LR) -- (3.3,-0.6); \draw[->, acc] (RL) -- (4.7,-0.6); \draw[->, acc] (RL) -- (5.7,-0.6); \draw[->, acc] (RR) -- (7.1,-0.6); \draw[->, acc] (RR) -- (8.1,-0.6); \end{tikzpicture} $$ Now evaluate at a root $\omega_n^k$. By the halving lemma, $(\omega_n^k)^2 = \omega_{n/2}^{k}$, so $$ A(\omega_n^{k}) = A_{\text{even}}(\omega_{n/2}^{k}) + \omega_n^{k}\,A_{\text{odd}}(\omega_{n/2}^{k}). $$ The right-hand side asks for $A_{\text{even}}$ and $A_{\text{odd}}$ **evaluated at the $\tfrac n2$-th roots of unity** — two DFTs of size $n/2$. And here the $\pm$ pairing pays off: the point $\omega_n^{k+n/2} = -\omega_n^k$ squares to the **same** $\omega_{n/2}^k$, so the two half-size results serve **both** $k$ and $k+n/2$ at once, with only a sign flip on the odd part: $$ \begin{aligned} A(\omega_n^{k}) &= A_{\text{even}}(\omega_{n/2}^{k}) + \omega_n^{k}\,A_{\text{odd}}(\omega_{n/2}^{k}),\\ A(\omega_n^{k+n/2}) &= A_{\text{even}}(\omega_{n/2}^{k}) - \omega_n^{k}\,A_{\text{odd}}(\omega_{n/2}^{k}), \end{aligned} \qquad 0\le k < \tfrac n2. $$ These paired equations are the **butterfly**: two outputs assembled from two inputs by one multiply ($\omega_n^k$ times the odd value, the **twiddle factor**) and one add/subtract. A full level of $n/2$ butterflies combines two half-DFTs into one full DFT in $\Theta(n)$ work. $$ % caption: One butterfly. The two half-DFT outputs $E=A_{\text{even}}(\omega_{n/2}^k)$ and % $O=A_{\text{odd}}(\omega_{n/2}^k)$ form the twiddle $t=\omega_n^k\cdot O$; the top % output is $A(\omega_n^k)=E+t$, the bottom $A(\omega_n^{k+n/2})=E-t$ \begin{tikzpicture}[ every node/.style={font=\small}, dot/.style={circle, fill=acc, minimum size=4pt, inner sep=0}, >=stealth, x=1cm, y=1cm] \definecolor{acc}{HTML}{2348F2} % inputs \node[anchor=east] (e) at (0,0.9) {even half $E$}; \node[anchor=east] (o) at (0,-0.9) {odd half $O$}; \node[dot] (ed) at (0.25,0.9) {}; \node[dot] (od) at (0.25,-0.9) {}; % outputs \node[dot, fill=acc] (top) at (4.6,0.9) {}; \node[dot, fill=acc, draw=acc] (bot) at (4.6,-0.9) {}; \node[anchor=west, acc] at (4.85,0.9) {top output: sum $E$ plus $t$}; \node[anchor=west, acc] at (4.85,-0.9) {next output: $E$ min\/us $t$}; % butterfly edges: solid = the plus path, dashed = the minus path \draw[->, acc, thick] (ed) -- (top); \draw[->, acc, thick, dashed] (ed) -- (bot); \draw[->, acc, thick] (od) -- (top); \draw[->, acc, thick, dashed] (od) -- (bot); \node[acc, font=\footnotesize, anchor=south] at (1.45,1.0) {scale $O$ by the root: t\/widdle $t$}; \node[acc, font=\footnotesize, anchor=north, align=center] at (2.3,-1.7) {the two outputs are $E$ plus and\\ min\/us the same t\/widdle $t$}; \end{tikzpicture} $$ Recursing this split bottoms out at size-$1$ DFTs (a single coefficient evaluates to itself), and at every level above it does $\Theta(n)$ combining work, giving the recurrence $$ T(n) = 2\,T\!\parens{\tfrac n2} + \Theta(n). $$ This is the **mergesort recurrence**. By the [Master Theorem](/algorithms/foundations/recurrences) ($a=2$, $b=2$, $f(n)=\Theta(n)$, the balanced case $\log_b a = 1$), it solves to $$ T(n) = \Theta(n\log n). $$ ```algorithm caption: $\textsc{FFT}(a)$ — recursive Cooley–Tukey DFT, $n$ a power of two number: 1 $n \gets \text{length}(a)$ if $n = 1$ then return $a$ // size-1 DFT is the identity $a_{\text{even}} \gets (a_0, a_2, \dots, a_{n-2})$ // even-indexed $a_{\text{odd}} \gets (a_1, a_3, \dots, a_{n-1})$ // odd-indexed $E \gets \textsc{FFT}(a_{\text{even}})$ // half-size DFT $O \gets \textsc{FFT}(a_{\text{odd}})$ // half-size DFT $\omega \gets 1,\quad \omega_n \gets e^{2\pi i / n}$ for $k \gets 0$ to $n/2 - 1$ do $t \gets \omega \cdot O_k$ // twiddle the odd part $\hat a_k \gets E_k + t$ // butterfly: top output $\hat a_{k + n/2} \gets E_k - t$ // butterfly: bottom output $\omega \gets \omega \cdot \omega_n$ // advance to next root return $\hat a$ ``` ::impl{algo="fft"} The algorithm requires $n$ to be a **power of two** so the halving is exact; pad the coefficient vector with zeros up to the next power of two (and to at least $2n$, for the product's degree) before transforming. A non-recursive, in-place variant first permutes the inputs into **bit-reversed** order and then runs the butterflies bottom-up, which is what production libraries implement, but the recursion above is the cleanest statement of why it is $\Theta(n\log n)$. ## Inverting the transform Evaluation is half the round trip; we still must turn the $2n$ **product values** back into coefficients — the **inverse DFT**. The DFT is a linear map, multiplication of $a$ by the matrix $V$ with $V_{kj} = \omega_n^{kj}$ (a **Vandermonde** matrix at the roots of unity). Its inverse has a simple form. > **Lemma (Inverse DFT).** The inverse of the DFT matrix is its own conjugate, scaled: > $$\parens{V^{-1}}_{kj} = \tfrac{1}{n}\,\omega_n^{-kj}, \qquad\text{so}\qquad a_j = \frac{1}{n}\sum_{k=0}^{n-1} \hat a_k\,\omega_n^{-kj}.$$ > **Proof.** It suffices to show $V\overline{V} = nI$, where $\overline V_{jk} = > \omega_n^{-jk}$. The $(k,k')$ entry of $V\overline V$ is > $\sum_{j} \omega_n^{kj}\omega_n^{-jk'} = \sum_{j}\omega_n^{(k-k')j}$. If $k=k'$ this sums > $n$ ones, giving $n$. If $k\ne k'$ it is a **geometric series** with ratio > $\omega_n^{k-k'}\ne 1$, summing to $\frac{\omega_n^{(k-k')n}-1}{\omega_n^{k-k'}-1}=0$ > because $\omega_n^{(k-k')n} = (\omega_n^n)^{k-k'} = 1$. So $V\overline V = nI$. $\qed$ The consequence is that the inverse transform is **the same algorithm**: replace $\omega_n$ by its conjugate $\omega_n^{-1} = e^{-2\pi i/n}$, run $\textsc{FFT}$ unchanged, and divide every output by $n$. We never invert a matrix or solve a linear system — interpolation at the roots of unity is just another FFT. $$ % caption: Forward and inverse FFT are the same machine. The forward pass evaluates at the % roots $\omega_n^k$; the inverse runs the identical butterflies on the conjugate % roots $\omega_n^{-k}$ and scales the result by $1/n$ \begin{tikzpicture}[ every node/.style={font=\small}, box/.style={draw=acc, fill=acc!8, minimum width=26mm, minimum height=11mm, inner sep=2pt, align=center}, >=stealth, x=1cm, y=1cm] \definecolor{acc}{HTML}{2348F2} \node[box] (coef) at (0,0) {input list\\ of numbers}; \node[box] (val) at (6.4,0) {list of values\\ at the roots}; \draw[->, acc, thick] (coef.20) to[bend left=18] node[above, font=\footnotesize, acc, align=center] {forw\/ard FFT: evaluate at the roots} (val.160); \draw[->, acc, thick, dashed] (val.200) to[bend left=18] node[below, font=\footnotesize, acc, align=center] {inverse FFT: conjugate roots, then scale} (coef.340); \end{tikzpicture} $$ ## Putting it together: fast multiplication The three-step plan now has an $\Theta(n\log n)$ cost on every leg. Pad $a$ and $b$ to length $N$, the next power of two $\ge 2n$; transform both; multiply the two value vectors pointwise; inverse-transform. Three FFTs plus one linear pointwise multiply: ```algorithm caption: $\textsc{Poly-Multiply}(a, b)$ — convolution in $\Theta(n\log n)$ number: 2 $N \gets$ next power of two $\ge 2n$ pad $a, b$ with zeros to length $N$ $\hat a \gets \textsc{FFT}(a)$ // evaluate at the $N$-th roots $\hat b \gets \textsc{FFT}(b)$ for $k \gets 0$ to $N - 1$ do $\hat c_k \gets \hat a_k \cdot \hat b_k$ // pointwise product $c \gets \tfrac{1}{N}\,\textsc{FFT-conj}(\hat c)$ // inverse FFT (conjugate roots, /N) return $c$ // the convolution $a \ast b$ ``` ::impl{algo="poly_multiply"} > **Theorem (Fast convolution).** The convolution of two length-$n$ vectors, equivalently > the product of two degree-$(n-1)$ polynomials, can be computed in $\Theta(n\log n)$ > arithmetic operations.[^clrs-fft] **Worked example (a size-$4$ transform).** Take $A(x) = 1 + 2x$, coefficient vector $a = (1, 2, 0, 0)$ padded to length $4$. The $4$-th roots of unity are $\omega_4^0 = 1,\ \omega_4^1 = i,\ \omega_4^2 = -1,\ \omega_4^3 = -i$. Splitting by parity, $a_{\text{even}} = (1,0)$ gives $E = (1, 1)$ (a constant evaluates to itself at both $2$nd roots), and $a_{\text{odd}} = (2,0)$ gives $O = (2, 2)$. The butterflies with twiddles $\omega_4^0 = 1$ and $\omega_4^1 = i$ assemble $$ \hat a_0 = E_0 + 1\cdot O_0 = 3,\quad \hat a_1 = E_1 + i\cdot O_1 = 1 + 2i,\quad \hat a_2 = E_0 - 1\cdot O_0 = -1,\quad \hat a_3 = E_1 - i\cdot O_1 = 1 - 2i. $$ Check directly: $A(1) = 3$, $A(i) = 1 + 2i$, $A(-1) = -1$, $A(-i) = 1 - 2i$ — the transform is exactly $A$ evaluated at the four roots. To square $A$ (compute $(1+2x)^2 = 1 + 4x + 4x^2$), multiply the transform pointwise — $\hat c = (9,\ (1+2i)^2,\ 1,\ (1-2i)^2) = (9,\ -3+4i,\ 1,\ -3-4i)$ — then run the inverse FFT (conjugate roots, divide by $4$), recovering $c = (1, 4, 4, 0)$, the coefficients of $1 + 4x + 4x^2$. This is the asymptotic payoff that beats both schoolbook $\Theta(n^2)$ and Karatsuba's $\Theta(n^{1.585})$. The application that started us off, **big-integer multiplication**, is the same algorithm: an $n$-digit integer in base $b$ is the polynomial $\sum d_j b^j$ evaluated at $x=b$, so multiplying two big integers is multiplying their digit polynomials, then **carrying** the over-large coefficients — the problem behind _Multiply Strings_, and at large enough sizes it is what bignum libraries do; the same FFT pipeline (with carries) is the route from Karatsuba's $\Theta(n^{\log_2 3})$ down toward the near-linear $\Theta(n\log n)$ that the [divide-and-conquer lesson](/algorithms/divide-and-conquer/mergesort) flagged as the next rung. _Add Two Numbers_ shares the carry mechanics on a per-digit scale. > **Remark (Floating point vs. exact: the NTT).** The FFT works over the **complex** > numbers, so its arithmetic is floating-point and accumulates rounding error; for exact > integer results one rounds the outputs and must keep $N$ and the coefficient sizes small > enough that the error stays below $1/2$. When the answer is required **modulo a prime** > (the competitive-programming norm), one instead runs the **Number-Theoretic Transform**: > replace the complex root of unity $e^{2\pi i/n}$ by a primitive $n$-th root of unity in > $\mathbb Z_p$ for a prime $p$ with $n \mid p-1$ (such as $998244353$). Every step of the > FFT carries over verbatim — the halving lemma, the butterflies, the inverse-by-conjugate > argument all hold in $\mathbb Z_p$ — but the arithmetic is now **exact**, using the > [modular inverse](/algorithms/mathematical-algorithms/number-theory-basics) of $n$ for > the final scaling. ::impl{algo="ntt,big_integer_multiply"} ## The FFT beyond powers of two Cooley and Tukey's 1965 paper is one of the most-cited in all of computing, and the ideas around it stretch from signal processing to the fastest multiplication algorithm known. **Not just powers of two.** The recursion above needs $n$ to be a power of two. Real FFT libraries handle any length by **mixed-radix** splitting (factor $n = n_1 n_2$ and do $n_1$ transforms of size $n_2$ then $n_2$ of size $n_1$), and for a **prime** length they use **Bluestein's algorithm**, which re-expresses the size-$n$ DFT as a convolution and pads _that_ to a convenient power of two.[^bluestein] So the $\Theta(n\log n)$ bound holds for every $n$, not only the padded powers of two — a detail that matters when you cannot afford to round the transform size up. **The FFT was not new.** Cooley and Tukey rediscovered and popularized the algorithm, but Gauss had written down the same even/odd recursion around 1805 to interpolate asteroid orbits — before Fourier's own work was published.[^heideman] The 1965 paper's real impact was timing: it arrived just as digital computers made $\Theta(n\log n)$ spectral analysis suddenly practical, launching modern digital signal processing. **Toward optimal multiplication.** FFT-based multiplication leads to the asymptotically fastest integer multiplication known. **Schönhage–Strassen** (1971) runs the transform over a ring of the form $\mathbb Z/(2^m+1)\mathbb Z$ — no floating-point, no NTT-prime constraints — to multiply $n$-bit integers in $O(n\log n\log\log n)$.[^schonhage] The $\log\log n$ factor stood for nearly fifty years until Harvey and van der Hoeven (2021) reached the conjectured optimum $O(n\log n)$.[^harvey] These are galactic algorithms — the crossover point is enormous — but Schönhage–Strassen is genuinely used by GMP for numbers beyond a few tens of thousands of digits. **Beyond multiplication.** The same transform accelerates string matching (match a pattern against a text by convolving indicator vectors, catching wildcards), the computation of the autocorrelation and power spectrum of a signal, large-scale polynomial interpolation, and the multiplication step inside the Fiduccia/Kitamasa recurrence solver from the [matrix-exponentiation lesson](/algorithms/mathematical-algorithms/matrix-exponentiation). All of these reduce to convolution, which the FFT makes cheap. ## Takeaways - Multiplying polynomials (or big integers) is **convolution**, $\Theta(n^2)$ in the coefficient representation; but in the **point-value representation** the product is a $\Theta(n)$ **pointwise** multiply at shared sample points. - The bottleneck is converting between the two forms. Sampling at the **$n$-th roots of unity** makes both conversions cheap because squaring **halves** the root set and the second half is the **negation** of the first. - The **DFT** evaluates $A$ at all roots of unity; the **Cooley–Tukey FFT** splits on even/odd coefficient parity into two half-size DFTs combined by **butterflies**, giving $T(n)=2T(n/2)+\Theta(n)=\Theta(n\log n)$ by the **Master Theorem**. - The **inverse FFT** is the same algorithm with conjugate roots and a $1/n$ scaling, because the DFT matrix's inverse is $\tfrac1n\,\overline V$ (a geometric-series argument). - Three FFTs plus a pointwise multiply give $\Theta(n\log n)$ **polynomial and big-integer multiplication**; the **NTT** does the same exactly modulo a prime. [^clrs-fft]: **CLRS**, Ch. 30 — Polynomials and the FFT: the coefficient and point-value representations, the recursive Cooley–Tukey FFT with the $2T(n/2)+\Theta(n)$ recurrence, the inverse DFT via $V^{-1}=\tfrac1n\overline V$, and $\Theta(n\log n)$ polynomial multiplication. [^bluestein]: L. Bluestein, "A linear filtering approach to the computation of discrete Fourier transform," _IEEE Trans. Audio and Electroacoustics_ **18**(4), 1970; the chirp-z transform for arbitrary and prime lengths. [^heideman]: M. T. Heideman, D. H. Johnson, C. S. Burrus, "Gauss and the history of the fast Fourier transform," _IEEE ASSP Magazine_ **1**(4), 1984 — tracing the algorithm to Gauss (c. 1805), before Cooley–Tukey (1965). [^schonhage]: A. Schönhage and V. Strassen, "Schnelle Multiplikation großer Zahlen," _Computing_ **7**, 1971 — integer multiplication in $O(n\log n\log\log n)$. [^harvey]: D. Harvey and J. van der Hoeven, "Integer multiplication in time $O(n\log n)$," _Annals of Mathematics_ **193**(2), 2021.