Prime-factor FFT algorithm

﻿
Prime-factor FFT algorithm

The Prime-factor algorithm (PFA), also called the Good-Thomas algorithm (1958/1963), is a fast Fourier transform (FFT) algorithm that re-expresses the discrete Fourier transform (DFT) of a size "N" = "N"1"N"2 as a two-dimensional "N"1&times;"N"2 DFT, but "only" for the case where "N"1 and "N"2 are relatively prime. These smaller transforms of size "N"1 and "N"2 can then be evaluated by applying PFA recursively or by using some other FFT algorithm.

PFA should not be confused with the mixed-radix generalization of the popular Cooley-Tukey algorithm, which also subdivides a DFT of size "N" = "N"1"N"2 into smaller transforms of size "N"1 and "N"2. The latter algorithm can use "any" factors (not necessarily relatively prime), but it has the disadvantage that it also requires extra multiplications by roots of unity called twiddle factors, in addition to the smaller transforms. On the other hand, PFA has the disadvantages that it only works for relatively prime factors (e.g. it is useless for power-of-two sizes) and that it requires a more complicated re-indexing of the data based on the Chinese Remainder Theorem (CRT). Note, however, that PFA can be combined with mixed-radix Cooley-Tukey, with the former factorizing "N" into relatively prime components and the latter handling repeated factors.

PFA is also closely related to the nested Winograd FFT algorithm, where the latter performs the decomposed "N"1 by "N"2 transform via more sophisticated two-dimensional convolution techniques. Some older papers therefore also call Winograd's algorithm a PFA FFT.

(Although the PFA is distinct from the Cooley-Tukey algorithm, it is interesting to note that Good's 1958 work on the PFA was cited as inspiration by Cooley and Tukey in their famous 1965 paper, and there was initially some confusion about whether the two algorithms were different. In fact, it was the only prior FFT work cited by them, as they were not then aware of the earlier research by Gauss and others.)

Algorithm

Recall that the DFT is defined by the formula:

:$X_k = sum_\left\{n=0\right\}^\left\{N-1\right\} x_n e^\left\{-frac\left\{2pi i\right\}\left\{N\right\} nk \right\}qquadk = 0,dots,N-1.$

The PFA involves a re-indexing of the input and output arrays, which when substituted into the DFT formula transforms it into two nested DFTs (a two-dimensional DFT).

Re-indexing

Suppose that "N" = "N"1"N"2, where "N"1 and "N"2 are relatively prime. In this case, we can define a bijective re-indexing of the input "n" and output "k" by:

:$n = n_1 N_2 + n_2 N_1 mod N,$:$k = k_1 N_2^\left\{-1\right\} N_2 + k_2 N_1^\left\{-1\right\} N_1 mod N,$

where "N"1-1 denotes the modular multiplicative inverse of "N"1 modulo "N"2 and vice-versa for "N"2-1; the indices "k""a" and "n""a" run from 0,...,"N""a"-1 (for "a" = 1, 2). These inverses only exist for relatively prime "N"1 and "N"2, and that condition is also required for the first mapping to be bijective.

This re-indexing of "n" is called the "Ruritanian" mapping (also "Good's" mapping), while this re-indexing of "k" is called the "CRT" mapping. The latter refers to the fact that "k" is the solution to the Chinese remainder problem "k" = "k"1 mod "N"1 and "k" = "k"2 mod "N"2.

(One could instead use the Ruritanian mapping for the output "k" and the CRT mapping for the input "n", or various intermediate choices.)

A great deal of research has been devoted to schemes for evaluating this re-indexing efficiently, ideally in-place, while minimizing the number of costly modulo (remainder) operations (Chan, 1991, and references).

DFT re-expression

The above re-indexing is then substituted into the formula for the DFT, and in particular into the product "nk" in the exponent. Because "e"2&pi;"i" = 1, this exponent is evaluated modulo "N": any "N"1"N"2 = "N" cross term in the "nk" product can be set to zero. (Similarly, "X""k" and "x""n" are implicitly periodic in "N", so their subscripts are evaluated modulo "N".) The remaining terms give:

:$X_\left\{k_1 N_2^\left\{-1\right\} N_2 + k_2 N_1^\left\{-1\right\} N_1\right\} = sum_\left\{n_1=0\right\}^\left\{N_1-1\right\} left\left( sum_\left\{n_2=0\right\}^\left\{N_2-1\right\} x_\left\{n_1 N_2 + n_2 N_1\right\} e^\left\{-frac\left\{2pi i\right\}\left\{N_2\right\} n_2 k_2 \right\} ight\right) e^\left\{-frac\left\{2pi i\right\}\left\{N_1\right\} n_1 k_1 \right\}.$

The inner and outer sums are simply DFTs of size "N"2 and "N"1, respectively.

(Here, we have used the fact that "N"1-1"N"1 is unity when evaluated modulo "N"2 in the inner sum's exponent, and vice-versa for the outer sum's exponent.)

References

* I. J. Good, "The interaction algorithm and practical Fourier analysis," "J. R. Statist. Soc. B" 20 (2), 361-372 (1958). Addendum, "ibid." 22 (2), 373-375 (1960).
* L. H. Thomas, "Using a computer to solve problems in physics," in "Applications of Digital Computers" (Ginn: Boston, 1963).
* P. Duhamel and M. Vetterli, "Fast Fourier transforms: a tutorial review and a state of the art," "Signal Processing" 19, 259-299 (1990).
* S. C. Chan and K. L. Ho, "On indexing the prime-factor fast Fourier transform algorithm," "IEEE Trans. Circuits and Systems" 38 (8), 951-953 (1991).

Wikimedia Foundation. 2010.

Look at other dictionaries:

• Cooley–Tukey FFT algorithm — The Cooley–Tukey algorithm, named after J.W. Cooley and John Tukey, is the most common fast Fourier transform (FFT) algorithm. It re expresses the discrete Fourier transform (DFT) of an arbitrary composite size N = N1N2 in terms of smaller DFTs… …   Wikipedia

• Cooley-Tukey FFT algorithm — The Cooley Tukey algorithm, named after J.W. Cooley and John Tukey, is the most common fast Fourier transform (FFT) algorithm. It re expresses the discrete Fourier transform (DFT) of an arbitrary composite size N = N 1 N 2 in terms of smaller… …   Wikipedia

• Rader's FFT algorithm — Rader s algorithm (1968) is a fast Fourier transform (FFT) algorithm that computes the discrete Fourier transform (DFT) of prime sizes by re expressing the DFT as a cyclic convolution. (The other algorithm for FFTs of prime sizes, Bluestein s… …   Wikipedia

• Bruun's FFT algorithm — Bruun s algorithm is a fast Fourier transform (FFT) algorithm based on an unusual recursive polynomial factorization approach, proposed for powers of two by G. Bruun in 1978 and generalized to arbitrary even composite sizes by H. Murakami in 1996 …   Wikipedia

• Twiddle factor — A twiddle factor, in fast Fourier transform (FFT) algorithms, is any of the trigonometric constant coefficients that are multiplied by the data in the course of the algorithm. This term was apparently coined by Gentleman amp; Sande in 1966, and… …   Wikipedia

• Multiplication algorithm — A multiplication algorithm is an algorithm (or method) to multiply two numbers. Depending on the size of the numbers, different algorithms are in use. Efficient multiplication algorithms have existed since the advent of the decimal system.… …   Wikipedia

• Fast Fourier transform — A fast Fourier transform (FFT) is an efficient algorithm to compute the discrete Fourier transform (DFT) and its inverse. There are many distinct FFT algorithms involving a wide range of mathematics, from simple complex number arithmetic to group …   Wikipedia

• List of algorithms — The following is a list of the algorithms described in Wikipedia. See also the list of data structures, list of algorithm general topics and list of terms relating to algorithms and data structures.If you intend to describe a new algorithm,… …   Wikipedia

• List of numerical analysis topics — This is a list of numerical analysis topics, by Wikipedia page. Contents 1 General 2 Error 3 Elementary and special functions 4 Numerical linear algebra …   Wikipedia

• List of mathematics articles (P) — NOTOC P P = NP problem P adic analysis P adic number P adic order P compact group P group P² irreducible P Laplacian P matrix P rep P value P vector P y method Pacific Journal of Mathematics Package merge algorithm Packed storage matrix Packing… …   Wikipedia