Re: FFT module thoughts, was Re: [eigen] Eigen 3.0.5 Could NOT find FFTW

[Mark Borgerding <mark@xxxxxxxxxxxxxx>]

I tried to keep this in mind when writing originally.  It does not 
handle it, but it might not be too hard.
There are at least two Achille's heels,  both in 
kiss_cpx_fft::make_twiddles  ( file==ei_kissfft_impl.h )
1. The exp() function, as others have mentioned.
2. The calculation of pi, as done via acos(-1)

Stepping back, what make_twiddles does is make one full cycle of a 
complex sinusoid with nfft points and unity magnitude.
This problem easily decomposes into a primitive root-finding problem.  
This can be done very efficiently in most cases for the FFT, since 
people tend to use a length that is a multiple of small primes (2,3,5)
(Let me know if this is not clear and I can provide an example.)
In other words, if we had a root-finding function

  template<typename Scalar>
  Scalar  nth_root( const RealScalar & x, int n); // maybe implemented 
via Newton-Raphson ?

then make_twiddles could be specialized to use it for types other than 
complex<float|double|long double>
The real_twiddles function also contains some trig calls -- these could 
derived in a similar way.




On 04/05/2012 07:04 AM, Pavel Holoborodko wrote:
> Yes, I have my own multiple precision complex exp() in std space, as
> well as other math. functions. Seems to work fine in all other
> contexts (including Eigen itself), but doesn't help with FFT :(
>
> On Apr 5, 2012, at 19:52, "Björn Piltz"<bjornpiltz@xxxxxxxxxxxxxx>  wrote:
>
> > I think std::exp() might be a candidate.
> >
> > Try adding
> >      namespace std {  std::complex<Real>  exp(const std::complex<Real>&); }
> > to your file to see if it's the culprit.
> >
> > Björn