[Numpy-discussion] Numpy float precision vs Python list float issue

[David Cournapeau]

On Tue, Apr 21, 2009 at 12:49 AM, Rob Clewley <rob.clewley at gmail.com> wrote:
> On Mon, Apr 20, 2009 at 10:48 AM, David Cournapeau
> <david at ar.media.kyoto-u.ac.jp> wrote:
>> Rob Clewley wrote:
>>> David,
>>>
>>> I'm confused about your reply. I don't think Ruben was only asking why
>>> you'd ever get non-zero error after the forward and inverse transform,
>>> but why his implementation using lists gives zero error but using
>>> arrays he gets something of order 1e-15.
>>>
>>
>> That's more likely just an accident. Forward + inverse = id is the
>> surprising thing, actually. In any numerical package, if you do
>> ifft(fft(a)), you will not recover a exactly for any non trivial size.
>> For example, with floating point numbers, the order in which you do
>> operations matters, so:
> <SNIP ARITHMETIC>
>> Will give you different values for d and c, even if you "on paper",
>> those are exactly the same. For those reasons, it is virtually
>> impossible to have exactly the same values for two different
>> implementations of the same algorithm. As long as the difference is
>> small (if the reconstruction error falls in the 1e-15 range, it is
>> mostly likely the case), it should not matter,
>
> I understand the numerical mathematics behind this very well but my
> point is that his two algorithms appear to be identical (same
> operations, same order), he simply uses lists in one and arrays in the
> other. It's not like he used vectorization or other array-related
> operations - he uses for loops in both cases. Of course I agree that
> 1e-15 error should be acceptable, but that's not the point. I think
> there is legitimate curiosity in wondering why there is any difference
> between using the two data types in exactly the same algorithm.

Yes, it is legitimate and healthy to worry about the difference - but
the surprising thing really is the list behavior when you are used to
numerical computation :) And I maintain that the algorithms are not
the same in both operations. For once, the operation of using arrays
on the data do not give the same data in both cases, you can see right
away that m and ml are not the same, e.g.

print ml - morig

shows that the internal representation is not exactly the same.

David