Fast Fourier Transform

Introduction

The DFT algorithm works very well, but it has one major problem. Calculating the DFT of a sample set of size n takes time proportional to n². While this is not a problem for small sample sets, larger ones start to take rather a long time.

Enter the FFT, or fast Fourier transform. This is nothing more than a divide-and-conquer algorithm for doing Fourier transforms - but it only takes time proportional to n log n, which is a lot better than n²! As such, it provides an interesting and relevant diversion from all the maths we've covered so far.

In order to maximise the speed of the FFT, most sample sets are powers of 2 long. I will, therefore, not cover the adjustment required to cope with arbitrarily sized sets.

What we want is a function FFT(x,m) that returns X[m] for a given x[0..2N-1] (as before, except now we have a sample set size 2N). (For clarity, I have omitted the 1/2N factor present in the FFT.)

The divide-and-conquer algorithm now divides the array x[n] into odd and even parts:

If I call the odd parts x_ODD and the even x_EVEN , we get

Now, the two summations are both FFTs in themselves, so we get a recursive definition of FFT(x, m):

Now, this is clearly of order log N, so the net result is that calculating the full spectrum of FFT(x) is of order Nlog N, as requred. Using a similar derivation, the inverse FFT is given by:

The FFT made a lot of problems prevously thought to be impossible suddenly solvable, and useful. For example, a 4,096-point DFT takes over 300 times as long as a 4,096-point FFT, and the more points you have, the better the FFT is.

What's it For?

Multi-word Integer Multiplication

Let's take a simple example. There are several ways of performing multiplication on big (really big) numbers. One method is the classical long multiplication algorithm (as you should have learned in primary school). Another method is divide-and-conquer, where we can express a multiplication of two 2n-word integers in terms of the muliplication of three n-word integers.

Neither of these methods is really very quick, however. A better method of multiplication is to express each number as a polynomial:

e.g. 57829 = 5x⁴ + 7x³ + 8x² + 2x + 9

where x = 10. We can then perform polynomial multiplication using the FFT algorithm. This will therefore give us a time proportional to N log N, where N is the number of digits in the numbers.

The FFT Polynomial Multiplication Algorithm (coming soon)

During polynomial multiplication, each resulting coefficient is given by:

 

Now, if we compare that to our convolution integral, from before:

 

we can see that polynomial multiplication is actually a discrete form of convolution. It follows that we can solve this set of discrete convolutions by doing a FFT on each source set (a N log N time algorithm), multiplying each individual coefficient (a N time algorithm) and doing an inverse FFT (N log N). This system is therefore O(N log N).

So, to do this, we can perform the following:

• Generate two vectors, of length at least the sum of the two polynomial lengths, and a power of 2 long.
• In the lower halves of these vectors, place the polynomial coefficients.
• Fourier transform these vectors.
• Muliply the corresponding coefficients of these vectors to create a third vector.
• Inverse Fourier transform these vectors.

This algorithm eases considerably the task of large-number multiplication. This sort of thing crops up a lot, in fact, including the task of computing high-precision values of constants like π and e.