Fast Fourier transform
These notes derive the radix-2 Fast Fourier Transform (FFT) by viewing the discrete Fourier transform (DFT) as a linear operator on an $N$-dimensional Hilbert space. The key steps are: an even–odd decomposition of the Hilbert space, expressing the $N$-point DFT in terms of two $M$-point DFTs ($N = 2M$), and writing an explicit operator factorization in Dirac notation, followed by a block-matrix form and explicit $N=4$ and $N=8$ examples.
I. Discrete Fourier transform as an operator
Let $H_N$ be an $N$-dimensional Hilbert space with orthonormal basis
\[ \{\,\lvert 0\rangle,\lvert 1\rangle,\dots,\lvert N-1\rangle\,\}. \]
A vector (signal) is written as
\[ \lvert x\rangle = \sum_{n=0}^{N-1} x_n \lvert n\rangle. \]
We define the discrete Fourier transform (DFT) operator $\hat{F}_N$ via
\[ \lvert X\rangle = \hat{F}_N \lvert x\rangle, \qquad X_k = \langle k\vert X\rangle, \]
with matrix elements
\[ \langle k\vert \hat{F}_N \vert n\rangle = \frac{1}{\sqrt{N}}\, \omega_N^{kn}, \qquad \omega_N := e^{-2\pi i / N}. \]
Thus, in the $\{\lvert n\rangle\}$ basis, the DFT matrix $F_N$ has entries
\[ (F_N)_{k n} = \frac{1}{\sqrt{N}}\, \omega_N^{k n}, \]
so that
\[ F_N = \frac{1}{\sqrt{N}} \begin{pmatrix} 1 & 1 & 1 & \cdots & 1 \\ 1 & \omega_N & \omega_N^2 & \cdots & \omega_N^{N-1} \\ 1 & \omega_N^2 & \omega_N^4 & \cdots & \omega_N^{2(N-1)} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & \omega_N^{N-1} & \omega_N^{2(N-1)} & \cdots & \omega_N^{(N-1)(N-1)} \end{pmatrix}. \]
II. Even–odd decomposition of the Hilbert space
Assume $N = 2M$. We identify
\[ H_N \simeq H_2 \otimes H_M, \]
by
\[ \lvert 2m\rangle \leftrightarrow \lvert 0\rangle \otimes \lvert m\rangle, \qquad \lvert 2m+1\rangle \leftrightarrow \lvert 1\rangle \otimes \lvert m\rangle, \qquad m = 0,\dots,M-1. \]
A general input vector is
\[ \lvert x\rangle = \sum_{n=0}^{N-1} x_n \lvert n\rangle = \sum_{m=0}^{M-1} x_{2m}\lvert 2m\rangle + \sum_{m=0}^{M-1} x_{2m+1}\lvert 2m+1\rangle. \]
Using the tensor-product structure
\[ \lvert 2m\rangle = \lvert 0\rangle \otimes \lvert m\rangle, \qquad \lvert 2m+1\rangle = \lvert 1\rangle \otimes \lvert m\rangle, \]
we obtain
\[ \lvert x\rangle = \lvert 0\rangle \otimes \left( \sum_{m=0}^{M-1} x_{2m}\lvert m\rangle \right) + \lvert 1\rangle \otimes \left( \sum_{m=0}^{M-1} x_{2m+1}\lvert m\rangle \right). \]
Define the even and odd subvectors
\[ \lvert x^{(e)}\rangle := \sum_{m=0}^{M-1} x_{2m}\lvert m\rangle, \qquad \lvert x^{(o)}\rangle := \sum_{m=0}^{M-1} x_{2m+1}\lvert m\rangle. \]
Then
\[ \lvert x\rangle = \lvert 0\rangle \otimes \lvert x^{(e)}\rangle + \lvert 1\rangle \otimes \lvert x^{(o)}\rangle. \]
III. Action of $\hat{F}_N$ on the decomposed state
Let
\[ \lvert X\rangle = \hat{F}_N \lvert x\rangle = \sum_{k=0}^{N-1} X_k \lvert k\rangle. \]
By definition,
\[ X_k = \langle k\vert \hat{F}_N \vert x\rangle = \frac{1}{\sqrt{N}} \sum_{n=0}^{N-1} x_n \omega_N^{k n}. \]
A. Even and odd indices
Split the sum over $n$ into even and odd parts:
\[ \begin{aligned} X_k &= \frac{1}{\sqrt{N}} \left( \sum_{m=0}^{M-1} x_{2m}\,\omega_N^{k(2m)} + \sum_{m=0}^{M-1} x_{2m+1}\,\omega_N^{k(2m+1)} \right) \\ &= \frac{1}{\sqrt{N}} \left( \sum_{m=0}^{M-1} x_{2m}\,\omega_N^{2km} + \omega_N^k \sum_{m=0}^{M-1} x_{2m+1}\,\omega_N^{2km} \right). \end{aligned} \]
Define
\[ \tilde{E}_k := \sum_{m=0}^{M-1} x_{2m}\,\omega_N^{2km}, \qquad \tilde{O}_k := \sum_{m=0}^{M-1} x_{2m+1}\,\omega_N^{2km}, \]
so
\[ X_k = \frac{1}{\sqrt{N}} \left( \tilde{E}_k + \omega_N^k \tilde{O}_k \right). \]
Note that
\[ \omega_N^2 = e^{-2\pi i/N \cdot 2} = e^{-2\pi i/(N/2)} = \omega_{N/2}. \]
Hence
\[ \tilde{E}_k = \sum_{m=0}^{M-1} x_{2m}\,\omega_{N/2}^{k m}, \qquad \tilde{O}_k = \sum_{m=0}^{M-1} x_{2m+1}\,\omega_{N/2}^{k m}. \]
B. Relation to the smaller Fourier transform $\hat{F}_M$
Define $\hat{F}_M$ on $H_M$ by
\[ \hat{F}_M\lvert m\rangle = \frac{1}{\sqrt{M}} \sum_{k=0}^{M-1} \omega_M^{km}\,\lvert k\rangle, \qquad \omega_M := e^{-2\pi i/M} = \omega_{N/2}. \]
Let
\[ \lvert E\rangle := \hat{F}_M\lvert x^{(e)}\rangle = \sum_{k=0}^{M-1} E_k \lvert k\rangle, \qquad \lvert O\rangle := \hat{F}_M\lvert x^{(o)}\rangle = \sum_{k=0}^{M-1} O_k \lvert k\rangle. \]
Then
\[ \begin{aligned} E_k &= \langle k\vert \hat{F}_M\lvert x^{(e)}\rangle = \frac{1}{\sqrt{M}} \sum_{m=0}^{M-1} x_{2m}\,\omega_M^{km} = \frac{1}{\sqrt{M}}\,\tilde{E}_k, \\ O_k &= \langle k\vert \hat{F}_M\lvert x^{(o)}\rangle = \frac{1}{\sqrt{M}} \sum_{m=0}^{M-1} x_{2m+1}\,\omega_M^{km} = \frac{1}{\sqrt{M}}\,\tilde{O}_k. \end{aligned} \]
Hence
\[ \tilde{E}_k = \sqrt{M}\,E_k, \qquad \tilde{O}_k = \sqrt{M}\,O_k, \]
and
\[ \begin{aligned} X_k &= \frac{1}{\sqrt{N}} \left( \sqrt{M}\,E_k + \omega_N^k \sqrt{M}\,O_k \right) = \sqrt{\frac{M}{N}} \left( E_k + \omega_N^k O_k \right). \end{aligned} \]
Since $N = 2M$, we have $\sqrt{M/N} = 1/\sqrt{2}$, so
\[ X_k = \frac{1}{\sqrt{2}} \left( E_k + \omega_N^k O_k \right), \qquad k = 0,\dots,M-1. \]
Now compute $X_{k+M}$:
\[ \begin{aligned} X_{k+M} &= \frac{1}{\sqrt{N}} \sum_{n=0}^{N-1} x_n \omega_N^{(k+M)n} = \frac{1}{\sqrt{N}} \sum_{n=0}^{N-1} x_n \omega_N^{kn} \omega_N^{Mn}. \end{aligned} \]
Note that
\[ \omega_N^{M} = e^{-2\pi i M /N} = e^{-2\pi i (N/2)/N} = e^{-\pi i} = -1 \quad\Rightarrow\quad \omega_N^{Mn} = (-1)^n. \]
Therefore
\[ \begin{aligned} X_{k+M} &= \frac{1}{\sqrt{N}} \left( \sum_{m=0}^{M-1} x_{2m}\,\omega_N^{k(2m)}(+1) + \sum_{m=0}^{M-1} x_{2m+1}\,\omega_N^{k(2m+1)}(-1) \right) \\ &= \frac{1}{\sqrt{N}} \left( \sum_{m=0}^{M-1} x_{2m}\,\omega_N^{2km} - \omega_N^k \sum_{m=0}^{M-1} x_{2m+1}\,\omega_N^{2km} \right) \\ &= \frac{1}{\sqrt{N}} \left( \tilde{E}_k - \omega_N^k \tilde{O}_k \right) = \sqrt{\frac{M}{N}} \left( E_k - \omega_N^k O_k \right) \\ &= \frac{1}{\sqrt{2}} \left( E_k - \omega_N^k O_k \right). \end{aligned} \]
Hence the butterfly relations
\[ X_k = \frac{1}{\sqrt{2}} \big( E_k + \omega_N^k O_k \big), \qquad X_{k+M} = \frac{1}{\sqrt{2}} \big( E_k - \omega_N^k O_k \big), \quad k = 0,\dots,M-1. \]
IV. Operator factorization in Dirac notation
We now express the above procedure as a factorization of the operator $\hat{F}_N$ on $H_N \simeq H_2 \otimes H_M$.
A. Two smaller Fourier transforms
Define
\[ \hat{G} := \hat{I}_2 \otimes \hat{F}_M. \]
Acting on the decomposed state,
\[ \begin{aligned} \hat{G}\lvert x\rangle &= (\hat{I}_2 \otimes \hat{F}_M) \left( \lvert 0\rangle \otimes \lvert x^{(e)}\rangle + \lvert 1\rangle \otimes \lvert x^{(o)}\rangle \right) \\ &= \lvert 0\rangle \otimes \hat{F}_M\lvert x^{(e)}\rangle + \lvert 1\rangle \otimes \hat{F}_M\lvert x^{(o)}\rangle \\ &= \lvert 0\rangle \otimes \lvert E\rangle + \lvert 1\rangle \otimes \lvert O\rangle. \end{aligned} \]
B. Twiddle-factor operator
Define
\[ \hat{D} := \sum_{k=0}^{M-1} \omega_N^k \lvert k\rangle\langle k\rvert \]
on $H_M$, and the operator
\[ \hat{T} := \lvert 0\rangle\langle 0\rvert \otimes \hat{I}_M + \lvert 1\rangle\langle 1\rvert \otimes \hat{D}. \]
Then
\[ \hat{T} \left( \lvert 0\rangle \otimes \lvert E\rangle + \lvert 1\rangle \otimes \lvert O\rangle \right) = \lvert 0\rangle \otimes \lvert E\rangle + \lvert 1\rangle \otimes \hat{D}\lvert O\rangle. \]
Let
\[ \lvert O'\rangle := \hat{D}\lvert O\rangle. \]
C. Hadamard on the first qubit
Define the Hadamard operator
\[ \hat{H} = \frac{1}{\sqrt{2}} \big( \lvert 0\rangle\langle 0\rvert + \lvert 0\rangle\langle 1\rvert + \lvert 1\rangle\langle 0\rvert - \lvert 1\rangle\langle 1\rvert \big) = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix}. \]
Consider $\hat{U} := \hat{H} \otimes \hat{I}_M$. Acting on
\[ \lvert \psi\rangle := \lvert 0\rangle \otimes \lvert E\rangle + \lvert 1\rangle \otimes \lvert O'\rangle, \]
\[ \begin{aligned} \hat{U}\lvert \psi\rangle &= (\hat{H} \otimes \hat{I}_M) \left( \lvert 0\rangle \otimes \lvert E\rangle + \lvert 1\rangle \otimes \lvert O'\rangle \right) \\ &= (\hat{H}\lvert 0\rangle) \otimes \lvert E\rangle + (\hat{H}\lvert 1\rangle) \otimes \lvert O'\rangle \\ &= \frac{1}{\sqrt{2}} (\lvert 0\rangle + \lvert 1\rangle) \otimes \lvert E\rangle + \frac{1}{\sqrt{2}} (\lvert 0\rangle - \lvert 1\rangle) \otimes \lvert O'\rangle \\ &= \frac{1}{\sqrt{2}} \left[ \lvert 0\rangle \otimes (\lvert E\rangle + \lvert O'\rangle) + \lvert 1\rangle \otimes (\lvert E\rangle - \lvert O'\rangle) \right]. \end{aligned} \]
Writing
\[ \lvert E\rangle = \sum_{k=0}^{M-1} E_k \lvert k\rangle, \qquad \lvert O'\rangle = \sum_{k=0}^{M-1} \omega_N^k O_k \lvert k\rangle, \]
the coefficient of $\lvert 0\rangle \otimes \lvert k\rangle$ is
\[ \frac{1}{\sqrt{2}} (E_k + \omega_N^k O_k) = X_k, \]
and the coefficient of $\lvert 1\rangle \otimes \lvert k\rangle$ is
\[ \frac{1}{\sqrt{2}} (E_k - \omega_N^k O_k) = X_{k+M}. \]
Identifying $\lvert 0\rangle \otimes \lvert k\rangle$ with $\lvert k\rangle$ and $\lvert 1\rangle \otimes \lvert k\rangle$ with $\lvert k+M\rangle$, we see that $\hat{U}\hat{T}(\hat{I}_2 \otimes \hat{F}_M)$ reproduces the DFT amplitudes.
Up to a permutation $\hat{P}$ that maps the standard basis $\{\lvert n\rangle\}$ to the even–odd tensor basis, we obtain
\[ \hat{F}_N = (\hat{H} \otimes \hat{I}_M)\,\hat{T}\, (\hat{I}_2 \otimes \hat{F}_M)\,\hat{P}. \]
V. Block-matrix form in the even–odd basis
In the even–odd basis, so that $\hat{P} = \hat{I}$, the matrix representations are:
- $\,\hat{I}_2 \otimes \hat{F}_M\,$: \[ I_2 \otimes F_M = \begin{pmatrix} F_M & 0 \\ 0 & F_M \end{pmatrix}. \]
- $\hat{T}$: \[ T = \begin{pmatrix} I_M & 0 \\ 0 & D \end{pmatrix}, \qquad D = \mathrm{diag}(1,\omega_N^1,\dots,\omega_N^{M-1}). \]
- $\hat{H} \otimes \hat{I}_M$: \[ H \otimes I_M = \frac{1}{\sqrt{2}} \begin{pmatrix} I_M & I_M \\ I_M & -I_M \end{pmatrix}. \]
Thus
\[ \begin{aligned} F_N^{(\mathrm{even\text{-}odd})} &= (H \otimes I_M)\,T\,(I_2 \otimes F_M) \\ &= \frac{1}{\sqrt{2}} \begin{pmatrix} I_M & I_M \\ I_M & -I_M \end{pmatrix} \begin{pmatrix} F_M & 0 \\ 0 & D F_M \end{pmatrix} \\ &= \frac{1}{\sqrt{2}} \begin{pmatrix} F_M & D F_M \\ F_M & -D F_M \end{pmatrix}, \qquad M = N/2. \end{aligned} \]
Therefore,
\[ F_N^{(\mathrm{even\text{-}odd})} = \frac{1}{\sqrt{2}} \begin{pmatrix} F_M & D F_M \\ F_M & -D F_M \end{pmatrix}, \qquad M = \frac{N}{2}. \]
To return to the standard ordering $(0,1,\dots,N-1)$, one uses
\[ F_N = F_N^{(\mathrm{even\text{-}odd})} P, \]
where $P$ is the permutation matrix implementing the even–odd reordering.
VI. Explicit examples: $N=4$ and $N=8$
A. The case $N = 4$
For $N = 4$ we have $M = N/2 = 2$ and
\[ \omega_4 = e^{-2\pi i/4} = e^{-i\pi/2} = -i. \]
The $2 \times 2$ Fourier matrix $F_2$ is
\[ F_2 = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix}. \]
The twiddle-factor matrix $D$ (of size $M \times M$) for $N = 4$ is
\[ D_2 = \mathrm{diag}(1,\omega_4^1) = \mathrm{diag}(1,-i) = \begin{pmatrix} 1 & 0 \\ 0 & -i \end{pmatrix}. \]
We can then compute
\[ D_2 F_2 = \begin{pmatrix} 1 & 0 \\ 0 & -i \end{pmatrix} \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ 1 & -1 \end{pmatrix} = \frac{1}{\sqrt{2}} \begin{pmatrix} 1 & 1 \\ -i & i \end{pmatrix}. \]
In the even–odd basis, where the input is ordered as $(x_0,x_2,x_1,x_3)^{T}$, the factorization
\[ F_4^{(\mathrm{even\text{-}odd})} = \frac{1}{\sqrt{2}} \begin{pmatrix} F_2 & D_2 F_2 \\ F_2 & -D_2 F_2 \end{pmatrix} \]
yields, after simplifying, the matrix often written as
\[ F_4^{(\mathrm{even\text{-}odd})} = \frac{1}{2} \begin{pmatrix} 1 & 1 & 1 & 1 \\ 1 & -1 & -i & i \\ 1 & 1 & -1 & -1 \\ 1 & -1 & i & -i \end{pmatrix}. \]
To relate this to the standard 4-point DFT in the natural ordering $(x_0,x_1,x_2,x_3)^{T}$, introduce the permutation matrix
\[ P = \begin{pmatrix} 1 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 1 \end{pmatrix}, \qquad P \begin{pmatrix} x_0 \\ x_1 \\ x_2 \\ x_3 \end{pmatrix} = \begin{pmatrix} x_0 \\ x_2 \\ x_1 \\ x_3 \end{pmatrix}. \]
Then
\[ F_4 = F_4^{(\mathrm{even\text{-}odd})} P = \frac{1}{2} \begin{pmatrix} 1 & 1 & 1 & 1 \\ 1 & -i & -1 & i \\ 1 & -1 & 1 & -1 \\ 1 & i & -1 & -i \end{pmatrix}, \]
which is precisely the 4-point DFT matrix with the normalization $1/\sqrt{N} = 1/2$.
B. The case $N = 8$
For $N = 8$ we have $M = N/2 = 4$ and
\[ \omega_8 = e^{-2\pi i / 8} = e^{-i\pi/4}. \]
The $4 \times 4$ Fourier matrix $F_4$ (in the standard ordering) is
\[ F_4 = \frac{1}{2} \begin{pmatrix} 1 & 1 & 1 & 1 \\ 1 & \omega_4 & \omega_4^2 & \omega_4^3 \\ 1 & \omega_4^2 & \omega_4^4 & \omega_4^6 \\ 1 & \omega_4^3 & \omega_4^6 & \omega_4^9 \end{pmatrix}, \qquad \omega_4 = e^{-2\pi i/4} = -i. \]
Using $\omega_4^2 = -1$ and $\omega_4^3 = i$, this can be written explicitly as
\[ F_4 = \frac{1}{2} \begin{pmatrix} 1 & 1 & 1 & 1 \\ 1 & -i & -1 & i \\ 1 & -1 & 1 & -1 \\ 1 & i & -1 & -i \end{pmatrix}. \]
For $N = 8$, the twiddle-factor matrix $D_4$ is
\[ D_4 = \mathrm{diag} \big( 1,\omega_8^1,\omega_8^2,\omega_8^3 \big) = \mathrm{diag} \big( 1,e^{-i\pi/4},e^{-i\pi/2},e^{-3i\pi/4} \big). \]
In the even–odd basis for $N = 8$, the block factorization reads
\[ F_8^{(\mathrm{even\text{-}odd})} = \frac{1}{\sqrt{2}} \begin{pmatrix} F_4 & D_4 F_4 \\ F_4 & -D_4 F_4 \end{pmatrix}. \]
This $8 \times 8$ matrix acts on the input ordered as
\[ (x_0,x_2,x_4,x_6,x_1,x_3,x_5,x_7)^{T}, \]
and its block structure directly encodes the top-level radix-2 FFT step:
- $F_4$ applied independently to the even and odd subsequences,
- multiplication of the “odd” 4-point transform by $D_4$,
- a final “butterfly” layer given by the block matrix \[ \frac{1}{\sqrt{2}} \begin{pmatrix} I_4 & I_4 \\ I_4 & -I_4 \end{pmatrix}. \]
To obtain the standard 8-point DFT matrix $F_8$ in the natural ordering $(x_0,x_1,\dots,x_7)^{T}$, one composes $F_8^{(\mathrm{even\text{-}odd})}$ with the corresponding even–odd permutation matrix $P_8$:
\[ F_8 = F_8^{(\mathrm{even\text{-}odd})} P_8. \]
The structure of $F_8^{(\mathrm{even\text{-}odd})}$ shows the recursive nature of the FFT, as it is entirely built from smaller Fourier transforms ($F_4$) and diagonal twiddle-factor matrices ($D_4$), combined by a simple $2 \times 2$ block “Hadamard”.