OpenCL Fast Fourier Transform

Eric Bainville - May 2010

The Discrete Fourier Transform

I give here a brief introduction on the subject, in order to make this article self-contained. The interested reader will find a lot of material on the Internet. As a starting point, refer to the Wikipedia DFT page. The DFT is defined as follows.

Given a sequence of N complex numbers x_j, j=0..N-1, the DFT of x is a sequence of N complex numbers y=DFT_N(θ,x) defined by:
y_k = Σ θ^j.k x_j, for k=0..N-1.

θ is a primitive N-th root of 1, usually taken as θ = e^-2iπ/N.

In other terms, we associate a polynomial of degree N-1 to a sequence of length N: P(U)=Σ x_j U^j, and the DFT is the vector of all values taken by this polynomial for all successive powers of θ:
y = (P(θ⁰),P(θ¹),...,P(θ^N-1).

This transformation is invertible, and the inverse is another DFT (up to a scalar):
if y = DFT_N(θ,x), then x = (1/N)*DFT_N(θ^-1,y).

The main interest of the DFT is how it operates with cyclic convolutions. Given two vectors u,v of dimension N, the cyclic convolution of u and v is the dimension N vector u×v defined by
(u×v)_k = Σ u_j * v_k-j,
where the index k-j is taken modulo N. Then we have immediately:
P_u×v = P_u * P_v mod (Z^N-1),
and consequently:
DFT(θ,u×v) = DFT(θ,u) * DFT(θ,v).
The product at the right hand side is a coordinate by coordinate product. Computing the convolution using the definition would take O(N²) time. The last equation shows that the time complexity of the convolution is at most equal to the time complexity of the DFT. Fast Fourier Transform (FFT) algorithms compute the DFT in time O(N.log(N)): FFT.