On Tue, Oct 11, 2011 at 12:17 PM, Matthew Brett
<matthew.brett@gmail.com> wrote:
Hi,
While struggling with floating point precision, I ran into this:
In [52]: a = 2**54+3
In [53]: a
Out[53]: 18014398509481987L
In [54]: np.float128(a)
Out[54]: 18014398509481988.0
In [55]: np.float128(a)-1
Out[55]: 18014398509481987.0
The line above tells us that float128 can exactly represent 2**54+3,
but the line above that says that np.float128(2**54+3) rounds upwards
as if it were a float64:
In [59]: np.float64(a)
Out[59]: 18014398509481988.0
In [60]: np.float64(a)-1
Out[60]: 18014398509481988.0
Similarly:
In [66]: np.float128('1e308')
Out[66]: 1.000000000000000011e+308
In [67]: np.float128('1e309')
Out[67]: inf
Is it possible that float64 is being used somewhere in float128 casting?
The problem is probably in specifying the values. Python doesn't support long double and I expect python integers to be converted to doubles, then cast to long double. The only way to get around this is probably using string representations of the numbers, and I don't know how well/consistently numpy does that at the moment. If it calls python to do the job, then double is probably what is returned. It doesn't help on my system:
In [1]: float128("18014398509481987.0")
Out[1]: 18014398509481988.0
Chuck