Talk:Fast Fourier transform

From Wikipedia, the free encyclopedia

You have new messages (last change).

Hey, does anyone have an issue number for the citation? The citation reads:

Cooley, James W., and John W. Tukey, 1965, "An algorithm for the machine calculation of complex Fourier series," Math. Comput. 19: 297–301

But there is no issue number.

A scan of the paper (of uncertain legality) can be found here. I believe it was in the April issue, but strangely enough the scanned pages don't list the journal issue or volume. In any case, the volume and page number uniquely identify the article, and should be sufficient to look it up in your library. —Steven G. Johnson 21:59, 2 April 2007 (UTC)

The "Other FFT Algorithms" segment is just an unreadable blob of text; I've separated the algorithms out (I think), which doesn't look as nice but is easier to read. Frankly, however, the whole thing ought to be rewritten for coherency, or expanded and put into its own section. I'm not qualified to do either of those things.

I think breaking this into more than three paragraphs is overkill. As for expanding it, any additional information about the FFT algorithms in question should go in a separate article, like for Prime-factor FFT algorithm and Rader's FFT algorithm. —Steven G. Johnson 23:58, 29 May 2006 (UTC)


More information on 2d or 3d FFT's would be super.

Alex.Szatmary

Done. At least as a start (the fancier multidimensional algorithms such as vector-radix could eventually use their own pages). —Steven G. Johnson 04:41, August 23, 2005 (UTC)

More information about the practical uses of FFT (e.g. how it is used for digital signal processing) would be very useful. If I ever find out what it does, I may write it up myself... Jevon 12:11, 23 November 2005 (UTC)

As is described in the first paragraph, information about applications belongs in discrete Fourier transform. —Steven G. Johnson 16:34, 23 November 2005 (UTC)

[edit] in finite fields

The whole article seems to address floating-point arithmetics as an approximation to the complex field.

However, the same algorithms may be used in any field where there is a nth prime root of the unit, where n is the length of the vector, including finite fields GF(pk). Since the multiplicative group of a finite field is cyclic, the condition for the existence of such a root is n | pk-1. If k=1, such a field is easily implementable in machine by modular arithmetics.

There is an algorithm for multiplying long integers by

  • consider them written in base β: a=\sum a_k.\beta^k and b=\sum b_k.\beta^k, as polynomials
  • transforming them by FFT modulo p
  • writing the convolution product a.b=\sum c_k.\beta^k where c_k=\sum_{j=0}^k a_j.b_{k-j}
  • computing this convolution product by inverse FFT of the product of the two FFTs
  • doing a little carrying.

A little web search brings up: J. M. Pollard, Mathematics of Computation, Vol. 25, No. 114. (Apr., 1971), pp. 365-374. Stable URL: http://links.jstor.org/sici?sici=0025-5718%28197104%2925%3A114%3C365%3ATFFTIA%3E2.0.CO%3B2-U

I'm not an expert in the field, so I don't think comfortable writing about this. I'll look how GNU MP does it. David.Monniaux 22:07, 13 October 2006 (UTC)

Only one section of the article discusses floating-point approximation, so I'm not sure why you're saying that the "whole article" is concerned with this. However, it's true that the article is concerned with algorithms for the discrete Fourier transform over vectors of complex numbers. It's certainly true that many of the algorithms have direct analogues for any finite field (although I can give counter-examples of some that can't, such as the Rader-Brenner algorithm); in particular, those algorithms that only use the fact that ωN is a primitive root of unity (e.g. Cooley-Tukey, PFA, Rader...) are generalizable in this way. But such transforms are generally given different names, like the number theoretic transform. —Steven G. Johnson 22:49, 13 October 2006 (UTC)
I've added a brief mention of such generalizations in the introductions. Any further details belong in articles like number theoretic transform. —Steven G. Johnson 23:02, 13 October 2006 (UTC)
Ah, interesting. I actually didn't know different names were used in English (I learned mathematics in French :-) ). David.Monniaux 08:27, 14 October 2006 (UTC)

[edit] example

Past century, I was frustrated with complexity of FFT progtams suggested in the literature. I believe, beginners can add all the accessories and improvements, if the short algorithm alreaady works. So, I add the short portable example. If you can write even shorter copylefted example, it may be even better. dima 03:23, 24 January 2007 (UTC)

I removed this. First, Wikipedia is not a code repository, and using a particular programming language (as opposed to English pseudo-code) is inappropriate in a mathematical article. Second, it makes the common mistake of confusing an "FFT algorithm" (which is usually considered in the literature to be an abstract decomposition of the DFT) with an "FFT implementation". Third, an example of the Cooley-Tukey algorithm would belong on Cooley-Tukey FFT algorithm, not here. Fourth, that page already gives the radix-2 example, but as math and English, not C++ code. Fifth, your code was far from minimal. Here is a shorter (but still not minimal) example in C (using the C99 complex type) that you would call with rec_fft(-1, N, in, out, 1, 1), for N a power of 2:
void rec_fft(int sign, int N, complex double *x, complex double *y, int is, int os)
{
     if (N == 1) *y = *x;
     else {
          int N2 = N/2;
          rec_fft(sign, N2, x, y, is*2, os);
          rec_fft(sign, N2, x+is, y+os*N2, is*2, os);
          complex double omega = 1.0, omega1 = cexp(sign * I * 6.2831853071795864769 / N);
          for (int i = 0; i < N2; ++i) {
               complex double a = y[os*i], b = omega * y[os*(i + N2)];
               y[os*i] = a + b;
               y[os*(i + N2)] = a - b;
               omega *= omega1;
          }
     }
}
(Making the code work out-of-place, and recursively, allows you to basically copy the two-line mathematical derivation of the radix-2 case, and leads to the simplest implementations because you don't have to worry about bit-reversal etc. It is not nearly as fast as a highly optimized code, but then again no textbook radix-2 FFT is remotely optimal these days.) If you want a code library that is pretty short and simple and of reasonable speed, Google "kissfft".
I would also disagree that "beginners can add all the accessories and improvements". A typical textbook radix-2 FFT implementation is an extremely poor starting point for a full-featured, highly-optimized FFT program.
—Steven G. Johnson 07:00, 25 January 2007 (UTC)

[edit] Example

Dear Steven G. Johnson.

You could declare that any examples of algorithmic codes are prohibited at Wiki, and I would understand;

but you seem to believe, that the kissfft, as well as your example, are simpler than the example you have removed.

I just checked the kissfft; the "minimal" version requires several files insteaad of two;

and the files required are significantly longer than the example you removed;

even the istructiuons of use in the README file are 4792 byte long.

In addition, the use is not straight-forward; it requires skills and/or assistance, because some additional files are supposed to be pre-installed at the user's computer, for example, the bench-user.h header required by the shortest example doit.c

As for the example you suggest, it requires 6 arguments instead of 3, and therefore is not simpler. I would recommend your example to the article recursion.

Sincerely, Dima.

I didn't say that kissfft was simpler than the code you posted; kissfft is more full-featured, faster, and supports non-power-of-two sizes, among other things, and so it is more complicated. But it is still pretty simple and easy to learn from; while still being vastly better than a textbook radix-2 code. (Nor does it require any code to be preinstalled; you're apparently confused because the kissfft author included an optional test program based on the benchfft framework. Kissfft only requires you to compile one .c file, along with two .h files; everything else is optional extras.) I mentioned it to help you since you claimed to be frustrated with the complexity of what you found on the web. (You can find zillions of textbook radix-2 FFT codes on the web and usenet similar to the one you posted, by the way—i.e. nested bit-reversal loops followed or preceded by nested butterfly loops—some of it dating back to the 1960s. Nowadays, however, these are mainly useful as pedagogical exercises and bear little resemblance to the programs actually used for serious numerical calculations.)
Regarding the code I posted above as an example, you have a very odd way of measuring complexity if you measure it by the number of arguments, as opposed to by the number of instructions. And the whole idea of most FFT algorithms is based on recursion; implementations that do not use recursion explicitly and work in-place are certainly possible and even common, but are more complicated (require more instructions and are conceptually less direct) than the simplest non-recursive out-of-place versions. If you want to actually understand FFTs and why they work, as opposed to cutting-and-pasting code, you need to understand recursion. And if you don't need to understand the guts of FFTs, and just want to use them to do something, then it hardly matters whether the code is 5 lines or 500 lines as long as it has the features you require and a reasonable API.
(Not to mention the fact that any user who has trouble with a 5k README file is unlikely to be capable of using a fast Fourier transform in a reasonable way, in my opinion.)
In any case, as I said, Wikipedia is not a code repository, and it doesn't matter whether we agree on the relative merits of different FFT implementations because none of them belongs in the article. Nor are Wikipedia articles on mathematical subjects intended to address the problem of installing and compiling and calling source code on a user's computer. Creating a code repository for numerical code, while a worthwhile mission (and being addressed elsewhere, see e.g. GNU Scientific Library) is quite a different problem from that of creating an encyclopedia. (I didn't say that algorithmic source code is "prohibited" on Wikipedia, just that it's not appropriate here. Code fragments would be appropriate, for example, in an article on programming languages or techniques.) Of course, if you want to create your own numerical code library (somewhere else, perhaps on Wikibooks) for "minimal" simplistic numerical codes, as a pedagogical tool, knock yourself out.
—Steven G. Johnson 04:23, 27 January 2007 (UTC)

Dear Steven. Thank you for the explanaiton. I had to modify your code (the compiler complained about "complex double"), and the modification below seems to work well. Is such modification portable?

using namespace std;
#include <complex>
void rec_fft(int sign, int N, complex<double> *x, complex<double> *y, int is, int os)
{ if (N==1) *y = *x;
  else{ int N2 = N/2;
       rec_fft(sign, N2, x, y, is*2, os);
       rec_fft(sign, N2, x+is, y+os*N2, is*2, os);
       double f=sign*2*M_PI/N;
       complex<double> omega =1., omega1=complex<double>(cos(f), sin(f));
       for(int i=0; i<N2; ++i)
       {       complex<double> a = y[os*i], b = omega * y[os*(i + N2)];
               y[os*i] = a + b;
               y[os*(i + N2)] = a - b;
               omega *= omega1;
       }
     }
}
main(){ int n,N=8,m=1; //play with values of N and m;
       complex<double> x[N],y[N];
       for(n=0;n<N;n++){x[n]=complex<double>(cos(2*M_PI/N*n),sin(-2*M_PI/N*n));}
       rec_fft(1,N,x,y,1,1);
       for(n=0;n<N;n++){printf("%1d %6.3f %6.3f %6.3f %6.3f\n",
                       n,x[n].real(),x[n].imag(),y[n].real(),y[n].imag());}
       }

I agree with your argumentation and I like your code. If it is still too large for Wiki, will you consider to post your example somewhere, and to supply the link at the FFT article? For Wiki, it is better than "kissfft": such example requires neither installation nor instruction. Sincerely, dima 05:23, 27 January 2007 (UTC)

(The reason you had to modify it is that I posted C code, whereas you are using C++. C++ does not support C's complex-number types.)
No, I don't think it's appropriate for Wikipedia. This is an encyclopedia article on a mathematical topic, not an article on C or C++ programming. Readers can follow the links, or use Google, if they want code in their language of choice. As for posting it somewhere else and linking to it, that's a bit too much like sneaking around Wikipedia policy. For Wikipedia to cite or link to a web page, it should generally be notable already. That is, if you create a page with the above code, and it becomes popular (e.g. lots of people link to it, it ranks high in Google, etcetera), then Wikipedia can link to it. Using Wikipedia to try and make a site popular is linkspam and is prohibited by policy.
Also, the code above still requires installation—you still need to compile it. And it still requires instruction: you have to clearly define what its inputs and outputs are, at least (and remember that there are several different conventions for the DFT). Moreover, who is the audience? If you just want to use an FFT code, then I would strongly recommend something more sophisticated. And if the above code is intended for instructional purposes, then it needs a lot more explanation because you have to explain how the code works.
(I hereby release my code above into the public domain, by the way. Feel free to use it wherever you want outside Wikipedia.)
—Steven G. Johnson 07:26, 27 January 2007 (UTC)
Retrieved from "http://en.wikipedia.org../../../f/a/s/Talk%7EFast_Fourier_transform_5c3b.html"