|Re: [eigen] A complex FFT for Eigen|
[ Thread Index |
| More lists.tuxfamily.org/eigen Archives
- To: eigen@xxxxxxxxxxxxxxxxxxx
- Subject: Re: [eigen] A complex FFT for Eigen
- From: "Benoit Jacob" <jacob.benoit.1@xxxxxxxxx>
- Date: Fri, 28 Nov 2008 20:03:03 +0100
- Dkim-signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=domainkey-signature:received:received:message-id:date:from:to :subject:in-reply-to:mime-version:content-type :content-transfer-encoding:content-disposition:references; bh=AhlNQoWgGed9eWwZMxBbC+3ZD1OncnyQovOzUfh2lJk=; b=Z5kkScqnH8oS83b3IWTeIvF5IzqjqQDAWRTUssvgejTwvb9xZ0XegZ4GNR0J6KSFRX pxPORK6EDXMQAfAm9F6DAO0PH8NyIFLesHbBJKSUxBoNoSW26Iyu2rtTYbrL3onnbqby mi5EcVZqbJQDizocfmSTZUeQNyjTWSYM7+cZE=
- Domainkey-signature: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=message-id:date:from:to:subject:in-reply-to:mime-version :content-type:content-transfer-encoding:content-disposition :references; b=DZMUtN/l5XFe7SjuJUVz4vG7PPfHp8tmDyei87Ph9PxeyG1ga8B773gmuFN7GnGu1+ +AKF44EtMPmSIhB/PEcc+lygDahu4RnYXqDuvHcs8dx+VY8eVJoKiBMRnGeRhnicdlkE 1+0tCh/zO/OWLWFaNNTElqTzboxUrUPGDAfi8=
I agree with what you say.
I agree that it would be useful to have a FFT module that has both a
builtin implementation (Tim's code is the beginning of it) and that
allows to use FFTW as optional backend.
As Keir notes, FFTW however is GPL so at least GPL-incompatible users
would always be using our builtin implementation.
We need extension to dynamic-size; of course that can't be done with
template metaprogramming but would involve a dynamic loop; we could
still try to "peel" it partly into a fixed-size version. Gael's idea
is exactly that and of course it applies not only to dynamic size but
also to too large fixed sizes.
Regarding complex numbers vectorization: this means that we'd need to
provide our own ei_complex type... and that we'd need to become much
more careful than we are now as our current code, at many places, is
based around the assumption that PacketSize==1 implies no
vectorization; also some of our vectorized paths (like in Dot.h)
assume real numbers as for now we don't vectorize complex numbers.
2008/11/28 Gael Guennebaud <gael.guennebaud@xxxxxxxxx>:
> here are my thoughts.
> First of all, I'm 100% OK to provide a Fourier module in Eigen. My issues
> are more about how to do it. Looking at FFTW benchmarks, it is pretty clear
> to me that we *must* provide support for it. So, in the same vein of the
> other modules, I would provide a default implementation based on Tim's
> lightweight source code such that we can make the dependency to FFTW
> optional for most use cases of this module. Of course this requires to
> extend a bit Tim's code to support at least: real data, dynamic sizes, 1D
> and 2D data (I guess all are trivial to add). IMO non pot sizes is not
> About current Tim's code, what about adding a specialization of
> Radix2_Decimation<> for sizes greater than, let's say 32, which would use
> normal recursive calls untill the size is 32. This would automatically add
> support for dynamic sizes, and significantly reduce the size of the
> generated code.
> About the vectorization of complex types, this is indeed perfectly doable.
> However note that for std::complex<double> the size of the packet would be 1
> ! So, somehow, this means that for double it is rather the std::complex's
> operators that should be updated to use the SSE instruction set....
> On Fri, Nov 28, 2008 at 3:12 AM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
>> 2008/11/27 Tim Molteno <tim@xxxxxxxxxxxxxxxxxxx>:
>> > I'd love to vectorize this better. I'll look into this, but ran into the
>> > issue of addressing sub-vectors in-place. I could not work out how to do
>> > it
>> > in Eigen. How do I carry out operations like
>> > data[1:N] += data[N+1:2N]
>> assuming that data is a vector of length 2N :
>> data.start<N>() += data.end<N>();
>> data.start<N>() += data.segment<N>(N);
>> (I think that segment only got that name after beta1 so you should use
>> trunk and "make doc" instead of using the online documentation).
>> > In other words, can I create a vector object that refers into part of
>> > another vector object?
>> > I tried using MatrixBlock objects, but could not get it to work?
>> Yes Eigen::Block is exactly what you need here and the above-mentioned
>> start/end/segment methods return objects of type Block.
>> > I think that FFT's are inherently non-local operations.
>> OK, so be it, then :)
>> > There might be
>> > better ways of doing this, but as it stands this code is pretty fast for
>> > an
>> > FFT and extremely small and lightweight, so I'd rather vectorize the
>> > apply
>> > function first as this is where most of the computation is happening.
>> Good to know.
>> As I said for now Eigen doesn't vectorize operations involve complex
>> numbers but I don't see anything making it impossible in theory, I'll
>> think about that when I can.
>> So, from Matthias's and your answers, I gather that there is no
>> generalization to non-complex numbers or to non-power-of-two sizes.
>> What about dynamic-size: are there, in your experience, use cases
>> where the size N isn't known at compile-time and it's be useful to
>> have a FFT able to handle any size N as a runtime parameter? If yes,
>> that's also something that we want.
>> More nitpicking:
>> Why are fft() and ifft() methods of a class CFFT ? I'd suggest either
>> making them static methods, or global functions. In Eigen, we tend to
>> do the latter as it makes for shorter calling syntax.
>> Let's give us a bit more time to determine whether we want FFT inside
>> Eigen -- don't want to sound un-enthusiastic, it's rather that I'm
>> currently busy and I don't want to make the wrong decision.
>> Anybody reading this, opinions?
>> > Cheers
>> > Tim
>> > ---