Fast Fourier Transform

Overview

The Fast Fourier Transform (FFT) algorithm is a divide and conquer algorithm that computes the Discrete Fourier Transform (DFT) of a sequence, or its inverse (IDFT). The DFT is a linear transformation that maps a vector of complex numbers to another vector of complex numbers. The DFT is widely used in signal processing and related fields to analyze frequencies contained in a signal.

The important observation behind the FFT is that the DFT contains a great deal of repeated structure. Rather than computing each output independently, we can divide the original sequence into its even- and odd-indexed elements, recursively compute the DFTs of those smaller sequences, and combine their results.

This turns what would otherwise be a $\theta(n^2)$ computation into a $\theta(n\log n)$ algorithm.

The DFT

Let our input sequence be $x_0,x_1,\ldots,x_{n-1}$. The DFT produces another sequence $X_0,X_1,\ldots,X_{n-1}$ according to

$$X_k = \sum_{j=0}^{n-1}x_j\omega_n^{jk}, \qquad k=0,1,\ldots,n-1$$

where

$$\omega_n = e^{-2\pi i/n}.$$

A straightforward implementation evaluates this sum separately for every $k$. Since there are $n$ values of $k$ and each sum contains $n$ terms, this requires $\theta(n^2)$ operations.

The FFT gets around this by noticing that we do not actually need to treat all $n$ terms as unrelated.

An Observation

Suppose, for simplicity, that $n$ is even. We can separate the input into its even- and odd-indexed elements:

$$X_k = \sum_{j=0}^{n/2-1}x_{2j}\omega_n^{2jk} + \sum_{j=0}^{n/2-1}x_{2j+1}\omega_n^{(2j+1)k}.$$

Since $\omega_n^{2}=\omega_{n/2}$, this becomes

$$X_k = \sum_{j=0}^{n/2-1}x_{2j}\omega_{n/2}^{jk} + \omega_n^k \sum_{j=0}^{n/2-1}x_{2j+1}\omega_{n/2}^{jk}.$$

The first sum is simply the DFT of the even-indexed elements, while the second sum is the DFT of the odd-indexed elements. Let these two transforms be $E_k$ and $O_k$:

$$E_k = \sum_{j=0}^{n/2-1}x_{2j}\omega_{n/2}^{jk}$$

$$O_k = \sum_{j=0}^{n/2-1}x_{2j+1}\omega_{n/2}^{jk}.$$

Then the original DFT can be written as

$$X_k = E_k + \omega_n^kO_k.$$

This is the key observation. An $n$-point DFT can be reduced to two $(n/2)$-point DFTs, followed by a relatively small amount of work to combine their results.

The Butterfly

There is one more useful property. Because $\omega_n^{k+n/2}=-\omega_n^k$, we can obtain the second half of the output without performing another pair of DFTs:

$$X_{k+n/2}=E_k-\omega_n^kO_k.$$

Thus, for each $k$ in the first half, we compute two outputs using the same two values $E_k$ and $O_k$:

$$X_k=E_k+\omega_n^kO_k$$

$$X_{k+n/2}=E_k-\omega_n^kO_k.$$

This pair of operations is commonly called a butterfly. The same pattern is repeated at every level of the recursion.

The Algorithm

Roughly,

def FFT(x):
    n = len(x)

    # Base case
    if n == 1:
        return x

    # Divide
    even = FFT(x[0::2])
    odd  = FFT(x[1::2])

    # Combine
    X = [0] * n
    for k in range(n // 2):
        omega = exp(-2*pi*i*k/n)

        X[k] = even[k] + omega * odd[k]
        X[k + n//2] = even[k] - omega * odd[k]

    return X
  

At each recursive level, the problem is divided into two equally sized problems. Once the smaller transforms have been computed, the butterfly operations combine them into the transform of the original sequence.

Why Is It Fast?

The recurrence for the algorithm is

$$T(n)=2T(n/2)+\theta(n).$$

The two recursive calls account for the $2T(n/2)$ term, while the $\theta(n)$ term comes from combining the two smaller transforms. By the Master Theorem,

$$T(n)=\theta(n\log n).$$

This is the central improvement over the straightforward DFT algorithm, whose runtime is $\theta(n^2)$. The FFT therefore allows us to compute the same transformation while using dramatically less computation as the size of the input grows.

A Broader Connection

The usefulness of the FFT comes from more than simply making the DFT faster. By transforming a signal from the time or spatial domain into the frequency domain, we can expose information about the frequencies that make up the signal.

The same mathematical structure also appears in problems far removed from traditional signal processing. In particular, the FFT can be used to accelerate polynomial multiplication by transforming polynomials into a representation where multiplication can be performed efficiently, and then transforming the result back.

The FFT is therefore a particularly nice example of a broader algorithmic idea: a transformation can reveal structure that is difficult to exploit in the original representation, and a carefully chosen divide-and-conquer strategy can make that transformation computationally practical.