[SciPy-Dev] Accuracy of single-precision FFT

[Anne Archibald]

I'm not sure exactly what you're listing here, but I wrote a test
program to compare single- and double-precision FFTs using FFTW; the
root-mean-squared fractional errors are:
1 1.10991e-09
17 9.26845e-08
37 1.16509e-07
97 1.20773e-07
313 1.33648e-07
701 1.55432e-07
1447 1.75132e-07
2011 1.80602e-07
3469 1.93587e-07

We can't use FFTW, unfortunately, but it's clear that prime-size FFTs
can do much better than we are. (Though what exactly is the definition
of the errors, below?) FFTW does this, IIRC, by not actually doing
prime-size FFTs; they have some scheme by which the prime-size FFT is
viewed as a convolution, then padded and computed with a
nicely-divisible-size FFT.

Anne

On 24 June 2010 12:22, Pauli Virtanen <pav at iki.fi> wrote:
> Thu, 24 Jun 2010 22:31:52 +0800, Ralf Gommers wrote:
> [clip: significant errors in float32 fft in Scipy]
>>>    1 0.0
>>>   17 4.76837e-07
>>>   37 2.98023e-06
>>>   97 0.000104427
>>>  313 0.000443935
>>>  701 0.00112867
>>> 1447 0.00620008
>>> 2011 0.0138307
>>> 3469 0.16958
>>>
>>> So even decimal=4 would fail for 97 already. For larger primes the FFT
>>> should be slower but not less accurate, right?
>>
>> Any opinion on this? Is it easily fixable? This is the last thing
>> holding up 0.8.0 I think, can we mark it knownfail for that or does
>> anyone think it's important enough to delay the release for?
>
> IIRC, single precision (float32) FFT is a new feature in Scipy 0.8, and
> was not present in earlier releases. I think Numpy and previous versions
> of Scipy were doing the FFT all the time in double precision (check
> this!).
>
> There are now two possibilities:
>
> 1) the single precision FFT in Scipy works incorrectly,
>
> 2) the single precision FFT in Scipy works correctly, but the precision
>   unavoidably sucks for large arrays.
>
> I guess (2) is more likely here.
>
> Moreover, I would guess the error is the largest in the high-frequency
> components, but remains low in at low frequencies.
>
>    ***
>
> The main question for the release, IMO, is:
>
>        Should the single-precision FFT implementation be used
>        by default if input data is single precision?
>
> I'm not completely convinced it should be, since if (2) is true, then
> this may surprise some people who expect more accuracy.
>
> Perhaps it should be put as a keyword option instead, or just split the
> functionality to different functions. I think we should maybe delay the
> release until we are sure which way is best to take here.
>
>    ***
>
> David, any comments?
>
> --
> Pauli Virtanen
>
> _______________________________________________
> SciPy-Dev mailing list
> SciPy-Dev at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-dev
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: fftwaccuracy.c
Type: text/x-csrc
Size: 1632 bytes
Desc: not available
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20100624/2f9da6c8/attachment.c>