Loss of precision with np.str_ input to np.longdouble

While trying to debug issues in a pulsar timing package, a very odd bug appeared:

import numpy as np
mjd = '50000.012345678901234'
np.longdouble(mjd)
# 50000.012345678901234
np.longdouble(np.str_(mjd))
# 50000.0123456789006
np.str_(mjd)
# '50000.012345678901234'

So, apparently python string and numpy string get converted differently...

p.s. Extends as expected to array/astype:

np.array([mjd, mjd], dtype='f16')[0]
# 50000.012345678901234
np.array([mjd, mjd]).astype('f16')[0]
# 50000.0123456789006

Work around np.longdouble bug.

I'm likely either responsible for a regression here, or responsible for forgetting to implement the precise behavior in all cases.

Is this present in master?

I think this may have existed for a while: the same problem is present in 1.17.4, 1.18.1, as well as current master. So, likely not your responsibility!

(We had not seen this before because in other code we used to return python string scalars instead of numpy string scalars - and we never tested with arrays.)

Also, should probably have noted on top: it is definitely a case of the numpy string first being transformed to float64 before being turned into float128, while the python string is correctly parsed at high precision:

np.longdouble(float(mjd))
# 50000.0123456789006

IIRC, the high precision parsing is fairly new. It was a long standing issue...

My first gut feeling would be that we have diverging implementation: One for SETITEM that knows about strings, and one for the casting machinery that does not. But, I have not dug into it, seems like a pretty annoying issue in practice.
Yeah, as Chuck said it is "fairly" new, although fairly by now is also a few years probably :).

The tests that started to fail were only put in half a year ago, so it is quite possible they happened to pass because of the "fairly new" ability to parse python string correctly!

OK, it arrived in #6199 (2015) and the improvements were made, IIRC, for pulsar data. So not so recent.

@aarchiba

Ah, then let me ping @aarchiba here as well.

Work around np.longdouble bug.

I went through numpy's parsing code and fixed what I could - there are situations where I couldn't figure out how to avoid going through a python float (f-strings, for example, IIRC), but this seems like one I just missed - I almost never use numpy strings. Fixing it should just be a matter of tracing the dispatch through to see where conversion occurs.

Its a bit confusing, but as far as I can tell the issue is that we are simply have too much logic. Deleting a bunch of lines (that I do not understand why they are needed to begin with) seems to fix the issue:

The @TYPE@_setitem functions really have to be able (and for all I can see do) handle all Python objects. But these cast functions seem to assume that they have to help out the _setitem function. I think deleting it is fine and we just have to write a test to ensure the correctness? (Should open a PR soon).

Details

diff --git a/numpy/core/src/multiarray/arraytypes.c.src b/numpy/core/src/multiarray/arraytypes.c.src
index c16e0f311..5e07f0df4 100644
--- a/numpy/core/src/multiarray/arraytypes.c.src
+++ b/numpy/core/src/multiarray/arraytypes.c.src
@@ -1506,70 +1506,8 @@ OBJECT_to_@TOTYPE@(void *input, void *output, npy_intp n,
  * #oskip = 1*18,(PyArray_DESCR(aop)->elsize)*3,1*2,
  *          1*18,(PyArray_DESCR(aop)->elsize)*3,1*2,
  *          1*18,(PyArray_DESCR(aop)->elsize)*3,1*2#
- * #convert = 1*18, 0*3, 1*2,
- *            1*18, 0*3, 1*2,
- *            0*23#
- * #convstr = (Long*9, Long*2, Float*4, Complex*3, Tuple*3, Long*2)*3#
  */
 
-#if @convert@
-
-#define IS_@from@
-
-static void
-@from@_to_@to@(void *input, void *output, npy_intp n,
-        void *vaip, void *aop)
-{
-    @fromtyp@ *ip = input;
-    @totyp@ *op = output;
-    PyArrayObject *aip = vaip;
-
-    npy_intp i;
-    int skip = PyArray_DESCR(aip)->elsize;
-    int oskip = @oskip@;
-
-    for (i = 0; i < n; i++, ip+=skip, op+=oskip) {
-        PyObject *new;
-        PyObject *temp = PyArray_Scalar(ip, PyArray_DESCR(aip), (PyObject *)aip);
-        if (temp == NULL) {
-            return;
-        }
-
-#if defined(IS_STRING)
-        /* Work around some Python 3K */
-        new = PyUnicode_FromEncodedObject(temp, "ascii", "strict");
-        Py_DECREF(temp);
-        temp = new;
-        if (temp == NULL) {
-            return;
-        }
-#endif
-        /* convert from Python object to needed one */
-        {
-            PyObject *args;
-
-            /* call out to the Python builtin given by convstr */
-            args = Py_BuildValue("(N)", temp);
-            new = Py@convstr@_Type.tp_new(&Py@convstr@_Type, args, NULL);
-            Py_DECREF(args);
-            temp = new;
-            if (temp == NULL) {
-                return;
-            }
-        }
-
-        if (@to@_setitem(temp, op, aop)) {
-            Py_DECREF(temp);
-            return;
-        }
-        Py_DECREF(temp);
-    }
-}
-
-#undef IS_@from@
-
-#else
-
 static void
 @from@_to_@to@(void *input, void *output, npy_intp n,
         void *vaip, void *aop)
@@ -1595,7 +1533,6 @@ static void
     }
 }
 
-#endif
 
 /**end repeat**/

Just as a side note, I guess clongdouble does not have correct parsing code currently.

BUG,MAINT: Remove incorrect special case in string to number casts

I tried to fix this in ~~gh-15608~~ gh-15766, one related issue I noticed that casting a longdouble to string, currently uses a string length of 32, which is OK for 80bit extended precision integers, but not in general. I am thinking of simply upping it to 48 to at least make this less of a problem.

Just curious what you guys think:

Is increasing to 48 correct (maybe assuming IEEE float128)
Can we do that change in a minor release? (i.e. a big enough bug with little actual change)
What happens with double-double, I assume that the double-double range of representable numbers is so weird, that it is not actually possible to print them reasonably for all numbers.

BUG,MAINT: Remove incorrect special case in string to number casts

Navigation Menu

Repository navigation

Loss of precision with np.str_ input to np.longdouble #15608

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Loss of precision with np.str_ input to np.longdouble #15608

Description

Activity

eric-wieser commented on Feb 19, 2020

mhvk commented on Feb 19, 2020

mhvk commented on Feb 19, 2020

charris commented on Feb 19, 2020

seberg commented on Feb 19, 2020

mhvk commented on Feb 19, 2020

charris commented on Feb 19, 2020

mhvk commented on Feb 19, 2020

aarchiba commented on Feb 25, 2020

seberg commented on Mar 16, 2020

seberg commented on Mar 16, 2020

seberg commented on Mar 17, 2020

Metadata

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Issue actions