> > I uploaded the files. If there are any questions, post them here
> > or to comp.lang.fortran, or decode my email address...
> Nifty code! My site started as a way to show my face to a a stranger I
> had agreed to meet at the airport. The rest grew by accretion. CuteFTP
> and CuteHTML made it easy.
Well, I got back from my walk in the park and concluded that maybe
I should put some descriptions on my website of what the files
are and also, to be topical to the original thread (now that the
subject line has changed, of course :) ) I zipped up my code for all
prime power orders <= 64, so it includes the factors 25 and 49,
although as we have seen, this old code is now known not to be
optimal.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Reply by Jerry Avins●January 24, 20042004-01-24
James Van Buskirk wrote:
> ... cool website. I finally got around to figuring
> out how to create a personal web page through my ISP.
>
> http://home.comcast.net/~kmbtib/index.html
>
> I uploaded the files. If there are any questions, post them here
> or to comp.lang.fortran, or decode my email address...
Nifty code! My site started as a way to show my face to a a stranger I
had agreed to meet at the airport. The rest grew by accretion. CuteFTP
and CuteHTML made it easy.
Jerry
--
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
Reply by James Van Buskirk●January 24, 20042004-01-24
"Jerry Avins" <jya@ieee.org> wrote in message
news:4011d483$0$12745$61fed72c@news.rcn.com...
> I'll be happy to host it. Even with my slow dial-up connection, it won't
> be much trouble. Check out http://users.rcn.com/jyavins to see if it
> would be appropriate. It could replace Jeff Briden's simulation that is
> old hat by now, have a line of its own, or be anonymous like
> http://users.rcn.com/jyavins/walala-lena.jpg.
Thanks for the offer; cool website. I finally got around to figuring
out how to create a personal web page through my ISP.
http://home.comcast.net/~kmbtib/index.html
I uploaded the files. If there are any questions, post them here
or to comp.lang.fortran, or decode my email address...
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Reply by Jerry Avins●January 23, 20042004-01-23
James Van Buskirk wrote:
> "Steven G. Johnson" <stevenj@alum.mit.edu> wrote in message
> news:bus6sq$8ib$1@news.fas.harvard.edu...
>
>
>> I hope
>>he publishes his algorithm in the near future, and encourage him to post
>>his code (gzipped it is < 10k).
>
>
> Perhaps if you have a little corner on your website reserved for
> miscellaneous stuff, you could place it there and post a link. Seems
> to me to be the preferred method for communicating largish chunks
> of data in newsgroups, if it's not too much trouble from your
> viewpoint.
I'll be happy to host it. Even with my slow dial-up connection, it won't
be much trouble. Check out http://users.rcn.com/jyavins to see if it
would be appropriate. It could replace Jeff Briden's simulation that is
old hat by now, have a line of its own, or be anonymous like
http://users.rcn.com/jyavins/walala-lena.jpg.
Jerry
P.S. If you look at the site, please tell me what you think of my new
favicon. It's a reproduction of the stamp I put on the jewelry I make.
--
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
Reply by James Van Buskirk●January 23, 20042004-01-23
"Steven G. Johnson" <stevenj@alum.mit.edu> wrote in message
news:bus6sq$8ib$1@news.fas.harvard.edu...
> I hope
> he publishes his algorithm in the near future, and encourage him to post
> his code (gzipped it is < 10k).
Perhaps if you have a little corner on your website reserved for
miscellaneous stuff, you could place it there and post a link. Seems
to me to be the preferred method for communicating largish chunks
of data in newsgroups, if it's not too much trouble from your
viewpoint.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Reply by Steven G. Johnson●January 23, 20042004-01-23
PS. James graciously sent me a copy of his code privately, and I can
verify that it is doing what he claims. This looks like a very
impressive result, and an important advance for the complexity theory of
FFT algorithms: breaking the 36-year-old record for the minimal flop
count in power-of-two sizes, from 1160 to 1152 flops for n=64. I hope
he publishes his algorithm in the near future, and encourage him to post
his code (gzipped it is < 10k).
PPS. I should caution other readers that a small change in the flop
count alone is unlikely to significantly impact FFT speed, which is
dominated these days by cache and pipeline utilization. But the
question of minimal arithmetic complexity is still important on
theoretical grounds.
Reply by James Van Buskirk●January 22, 20042004-01-22
"Steven G. Johnson" <stevenj@alum.mit.edu> wrote in message
news:4010424d$0$564$b45e6eb0@senator-bedfellow.mit.edu...
> Please do let me know when you publish your results, as I try to keep
> abreast of the latest FFT developments.
Well, I would just post the code the n = 64, but it's a little long --
35937 bytes. Is this objectionably long for this group? The trick
that made the difference for n = 64 would also improve the number of
multiplies for some of my smaller modules. I'm not sure at this point
whether I want to rewrite them all (many take a day each) or try to
create a code generator that can do what I do in the short time
available to me just now to ponder this problem again.
> Actually, it has supposedly been rigorously proved [1] that split-radix
> achieves the minimal number of complex additions.
My code is structured for minimum-add so I am not challenging the
above assertion.
> [1] P. Duhamel and M. Vetterli, "Fast Fourier transforms: A tutorial
> review and a state of the art," Signal Processing 19, 259-299 (1990).
> [2] R. Yavne, "An economical method for calculating the discrete Fourier
> transform," Proc. AFIPS 33, 115-125 (1968).
> [3] H. V. Sorensen, D. L. Jones, M. T. Heideman, and C. S. Burrus,
> "Real-valued fast Fourier transform algorithms," IEEE Trans. Acoust.
> Speech Sig. Processing ASSP-35, 849�863 (1987).
Thanks for the comments and the references. Lets see... I've read [1],
and it is highly recommended. [2] just happened to be sitting on my
desk as I read this, but I'm not sure if I've seen [3].
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Reply by Steven G. Johnson●January 22, 20042004-01-22
James Van Buskirk wrote:
>>James Van Buskirk wrote:
>>>But it seems to me that split-radix requires 912 additions and 248
>>>multiplications for n = 64, and 1160 isn't the minimal number of
>>>operations.
>
> The best I know of is 912 additions and 240 multiplications; just
> figured out how to achieve this result yesterday.
Please do let me know when you publish your results, as I try to keep
abreast of the latest FFT developments.
The arithmetic complexity of FFT algorithms was extensively studied in
the 1980's, when it was still considered to be a good proxy for speed
(sadly no longer the case). As of 1991 [1], the lowest-known operation
count for an n=64 complex DFT was 1160 flops, a record that had been set
in 1968 by Yavne [2] and was later matched using a simpler algorithm by
the split-radix technique. For N=2^m, this record is:
flops = 4 N log2(N) - 6N + 8
(This is asymptotically 20% fewer flops than for vanilla radix 2. Note
that you can trade off some multiplies for additions by using the 3 adds
+ 3 multiplies complex-multiply trick.)
Actually, it has supposedly been rigorously proved [1] that split-radix
achieves the minimal number of complex additions. There is also a
rigorous lower bound [1] on the number of non-trivial complex (or real)
multiplications, which split-radix can achieve up to N=16 (and it is
close for larger N). For larger N, there are known algorithms that
achieve the minimal number of multiplies, but they do not have fewer
total flops than split-radix, and have more flops for N > 64. So, the
achievable total flop count has yet not been rigorously proven.
> Algorithm is
> the same as used for all my computational kernals (heavy use of
> real-half-complex DFTs,) just became aware of another possible trick.
(There is the known technique of taking the real-input FFTs of the real
and imaginary parts, and then combining them with (2N-4) additions.
However, if you use the lowest-known operation counts for real-input
FFTs [3], this gives exactly the same known flops for the complex DFT.
e.g. 1160flops for n=64 from two 518flops n=64 real-input FFTs.)
Cordially,
Steven G. Johnson
[1] P. Duhamel and M. Vetterli, "Fast Fourier transforms: A tutorial
review and a state of the art," Signal Processing 19, 259-299 (1990).
[2] R. Yavne, "An economical method for calculating the discrete Fourier
transform," Proc. AFIPS 33, 115-125 (1968).
[3] H. V. Sorensen, D. L. Jones, M. T. Heideman, and C. S. Burrus,
"Real-valued fast Fourier transform algorithms," IEEE Trans. Acoust.
Speech Sig. Processing ASSP-35, 849�863 (1987).
Reply by James Van Buskirk●January 22, 20042004-01-22
"Steven G. Johnson" <stevenj@alum.mit.edu> wrote in message
news:bunj9a$hts$1@news.fas.harvard.edu...
> James Van Buskirk wrote:
> > But it seems to me that split-radix requires 912 additions and 248
> > multiplications for n = 64, and 1160 isn't the minimal number of
> > operations.
> What do you think is the lowest known number of operations for an n=64
> complex DFT, using what algorithm?
The best I know of is 912 additions and 240 multiplications; just
figured out how to achieve this result yesterday. Algorithm is
the same as used for all my computational kernals (heavy use of
real-half-complex DFTs,) just became aware of another possible trick.
The lowest known by anyone is a different question. The more you
look at DFTs, the more you become uncertain that what you or anyone
else does is optimal, even if you get to write the criteria of
optimality yourself.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
Reply by Steven G. Johnson●January 22, 20042004-01-22
James Van Buskirk wrote:
> But it seems to me that split-radix requires 912 additions and 248
> multiplications for n = 64, and 1160 isn't the minimal number of
> operations.
What do you think is the lowest known number of operations for an n=64
complex DFT, using what algorithm?