[SciPy-User] Single precision FFT insufficiently accurate.

Pauli Virtanen pav@iki...
Mon Jun 28 06:45:07 CDT 2010


Mon, 28 Jun 2010 13:21:25 +0200, Sebastian Haase wrote:
[clip]
> What size of error are talking about anyway ..? 

We are talking about 0.1% ... 5% relative error,

	import numpy as np
	from scipy.fftpack import fft, ifft

	x = np.random.rand(2011).astype(np.float32)
	np.linalg.norm(x - ifft(fft(x))) / np.linalg.norm(x)
	# -> 0.001  norm-2 relative error

	x = np.random.rand(2012).astype(np.float32)
	np.linalg.norm(x - ifft(fft(x))) / np.linalg.norm(x)
	# -> 6e-5  norm-2 relative error

	x = np.random.rand(8923).astype(np.float32)
	np.linalg.norm(x - ifft(fft(x))) / np.linalg.norm(x)
	# -> 0.03  norm-2 relative error

	x = np.random.rand(8925).astype(np.float32)
	np.linalg.norm(x - ifft(fft(x))) / np.linalg.norm(x)
	# -> 2.4902545e-07 norm-2 relative error

So for "difficult" cases the error is up to several orders of magnitude 
larger than for the "easy" cases.

> Personally I would
> leave it in, and make a note in the doc-string about expected precision
> error for non multiple-2 dimensions for single precision float.
> Maybe one could (for now) even append an option for
> "workaroundFloat32PrecisionLoss"

Several percent errors are not something I'd like to leave for the users 
to sort out by themselves, even if mentioned in the documentation.

I would perhaps rather drop the feature in 0.8 and wait for a proper fix 
in 0.9 (hopefully later this year), than add keyword arguments that we 
have to deprecate later on.

-- 
Pauli Virtanen



More information about the SciPy-User mailing list