0

I have been translating IDL code into Python and noticed differences in the end products. The differences are not negligible either (e.g. via IDL I get 19.03 while Python yields 19.16). I confirmed that the two scripts yield the same values (they typically differ in the fourth, fifth, or sixth decimal) up to the point where I start doing multiple array operations. I suspected the precision may be the cause (both Python and IDL arrays are type=FLOAT). I did a simple experiment and I see significant differences here too.

IDL

a = 0.01
b = 0.0

for r = 1,1000 do begin
   b += a
endfor
c = a * 1000

print,b
    >>> 10.0001
print,c
    >>> 10.0000

Python

a = 0.01
b = 0.00

for r in range(1000):
    b += a
c = a * 1000

print(b)
    >>> 9.999999999999831
print(c)
    >>> 10.0

Granted, the difference is still small in this example, but Python is clearly far closer to truth than is IDL. I expected the results to be identical since both languages are using FLOAT precision. The bottom line is the errors are propagating in both languages in different ways, which yields differing results. My questions are:

  1. is there a difference in precision between the two languages (i.e. does FLOAT mean the same thing in IDL as it does in Python, I THINK it does)?

  2. is there any way to reconcile the precision differences?

  3. I do not have much experience with IDL; am I missing something obvious?

PS:

As I was writing this post popped up. I am seeing the same issue as that OP (Python is correct).

IDL

>>> 3015/0.0002529821940697729
>>> 11917835.

Python

>>> 3015/0.0002529821940697729
>>> 11917834.814763514
DavidH
  • 415
  • 4
  • 21
tnknepp
  • 5,888
  • 6
  • 43
  • 57

1 Answers1

6

FLOAT in IDL is single precision (32 bit), and DOUBLE is double precision (64 bit).

In python (and numpy), the default floating point type is double precision (64 bit).

You can recreate the single precision IDL calculation using numpy as follows:

In [9]: import numpy as np

In [10]: a = np.array(0.01, dtype=np.float32)

In [11]: b = np.array(0.0, dtype=np.float32)

In [12]: for r in range(1000):
    ...:     b += a
    ...:     

In [13]: c = a*1000

In [14]: print(b)
10.000133514404297

In [15]: print(c)
9.99999977648

You can use the FORMAT argument of the IDL print statement to print more digits of a and b, and compare them to the values shown above.

Warren Weckesser
  • 110,654
  • 19
  • 194
  • 214
  • That is very helpful. I always assumed that "FLOAT" meant the same thing across languages. It is good to know that that is NOT the case. – tnknepp May 05 '17 at 17:50
  • @tnkneep not only is it not consistent accross languages, but primitives aren't even consistent across *platforms* in the same language. Double might not even be 64bits or might have higher precision in calculation depending on your hardware (some intel processesors for example have 80bit floating point units that perform floating point calculations even if it is stored in 64bits internally) – Krupip May 05 '17 at 17:59