Floating point arithmetic and Fourier collocation

[pasmith]

>>> from numpy import exp, pi
>>> exp(1j*pi)
(-1+1.2246467991473532e-16j)

The fact that the imaginary part of this is not zero is wrecking a fourier collocation scheme for a nonlinear PDE with periodic boundary conditions: the coefficient corresponding to the Nyquist frequency, which should be purely real, does not remain so. This had no bearing on the physical quantity, because scipy.fft.rfft apparently ignores the imaginary part of this coefficient, but it was disconcerting to look at the spectrum of what was apparently a uniform final steady state and find that the magnitude of this coefficient was on the order of $10^{100}$.

 

[pbuk]

It's that pesky roundoff error again. A simpler illustration is

>>> from math import sin, pi
>>> sin(pi)
1.2246467991473532e-16

Don't go out into the world of floating point arithmetic without your big boy trousers on.