Bottom     Previous     Contents

An Algorithm for Prime Length Sequences.

We have seen that if the size of a DFT N is a composite (non-prime) number, then it can be recursively decomposed to multiple DFTs, the (minimum) sizes of which are given by the prime factors of N (see the Mixed Radix Section). This still leaves the problem of finding an algorithm for a prime sized FFT. The algorithm described here uses a neat application of number theory to express an FFT of size N (N is prime,>2) as a circular convolution of length (N-1) which can be evaluated using 2 (N-1) point FFT's. (See Annex C if your not familiar with circular convolution). If N is prime (>2) then N-1 can't be prime (at the very least it must be even). It should be emphasised that this algorithm only works if N is prime.

An Example.

To illustrate the algorithm, consider a N=5 point un-scaled transform expressed as a matrix multiplication:

equation

The W factors are periodic:

equation

To simplify the notation we shall write this matrix as (taking account of the periodicity):

equation

Note that all the elements of this matrix (exponents of W) are now modulo N (modulo 5). Note also, that each row and column of the matrix, other than the first, contains some permutation of every possible value in the range 0..4. This is not a special property of the number N=5. If (and only if) N is prime, then each row and column of the matrix, other than the first, contains some permutation of every possible value in the range 0..N-1.

F(0) can simply be evaluated by summing f(n) for n= 0..N-1. The remaining transform can be re-written:

equation

At first glance, this matrix multiplication would appear to require (N-1)2 (=16) complex multiplies. However, the fact that each row of the matrix contains the same values in a different permutation suggests that we might be able to swap rows and columns in such a way that the matrix represents a circular convolution. (In a matrix representing a circular convolution, each row is obtained from the previous row by rotating right one step. A matrix in this form is said to be 'right circulant'. The corresponding result is obtained by convolving the first column of the matrix with the column vector f.

In our example one way (but by no means the only way) to achieve this is by:

  1. Swapping the last 2 matrix rows (and the last 2 elements of the column vector F').
  2. Rotating the last 3 matrix columns left (and the last 3 elements of the column vector f upwards).

If you do this, you end up with:

equation

Note that by doing this we haven't changed the values of F'(k), all we've changed is the order of arithmetic operations. Our matrix is now right circulant, so we can evaluate F'(k) using a fast (FFT based) 4 point convolution with the first column of the matrix:

equation

Clearly finding the appropriate row and column swaps using ad-hoc reasoning like this is going to be very difficult for large primes. If you enjoy puzzles like 'Rubics Cube' this might be an attractive feature of the algorithm. Fortunately (for most of us), there is a very simple method we can use to get the necessary permutations.

A Little of Number Theory.

So it works for N=5, but what we really need is a systematic way of generating the necessary row and column permutations for any prime N (>2). The answer to this problem is provided by a branch of number theory called 'Galois Field Theory' (which is very useful for all sorts of things other than FFTs).

Galois field theory is concerned with the arithmetic properties of numbers modulo N, where N is prime. A Galois field is the set of numbers (0,1,2, .. N-1). Any two values drawn from this set can be added or multiplied, but the result is modulo N (the remainder left after dividing by N). So the result of adding or multiplying any two numbers drawn from this set is also contained in this set (it lies in the range 0 .. N-1).

One useful property of a Galois field is that every number (n) drawn from the set (1,2, .. N-1) has a unique multiplicative inverse (n-1) which is also a member of the set:

equation

This provides a kind of modulo division, which enables us to raise numbers to negative powers.

Fermat's Theorem.

For any number (n) drawn from the set (1,2, .. N-1):

equation

This theorem provides a way to find multiplicative inverses:

equation

Generally, we can see that raising a number to an (integer) power is periodic, with periodicity of N-1:

equation

Generators.

Another useful property of a Galois field is that it always contains at least one number (a) called a generator (or primitive root), which is has the property that the set of numbers (ai | i= 0 .. N-2) is a permutation of the set (1,2, .. N-1). We also know (from Fermat's theorem) that the sequence generated is periododic. So, for any integer (m), the set of numbers (ai | i= m .. m+N-2) is a permutation of the set (1,2, .. N-1).

To test whether or not a number is a generator, simply count the number of times you have to multiply (modulo N) to get a result of 1. Fermat's theorem tells us that this process is guaranteed to terminate after (N-1) multiplies, if not before. If it does terminate early, the number is not a generator.

It is also worth noting that if (a) is a generator, then so is its multiplicative inverse (a-1). It will become apparent that this property is useful for generating matched forward and inverse transform permutation pairs.

A Systematic Method.

Recall the un-scaled DFT definition:

equation

First we break up this up, as we did in the example:

equation

equation

Again, to simplify, we shall use the following notation:

equation

So we are left with the problem of evaluating:

equation

To do this, we need to find a suitable generator (a) (using the method described earlier). With this generator, we can define new sequences h(n') and H(k'), which are simply permutations of our original sequences f(n) and F'(k):

equation

equation

The indices run from 0 .. N-2 instead of 1 .. N-1 for consistency with our usual notation (periodicity means these will still generate valid permutations). Also, in the definition of h(n') we've used the periodicity property to avoid negative powers of a. The expression we're trying to evaluate becomes:

equation

This is the usual form of a circular convolution:

equation

The Example Again.

We can apply all this to our example (N=5), which has a generator a=2:

equation

So..

equation

The corresponding matrix form is:

equation

This is exactly the same as was derived by ad-hoc reasoning earlier. We may note that 3 (=2-1) is also a generator for N=5. If we use 3 instead of 2 we get:

equation

This is would be a good permution choice for an inverse transform to complement the forward transform.

©1999 - Engineering Productivity Tools Ltd.


Next     Top