summaryrefslogtreecommitdiff
path: root/numpy
diff options
context:
space:
mode:
Diffstat (limited to 'numpy')
-rw-r--r--numpy/core/setup_common.py12
-rw-r--r--numpy/core/src/multiarray/dragon4.c251
-rw-r--r--numpy/core/src/multiarray/dragon4.h4
-rw-r--r--numpy/core/src/npymath/ieee754.c.src2
-rw-r--r--numpy/core/src/npymath/npy_math_private.h4
-rw-r--r--numpy/core/src/private/npy_fpmath.h12
-rw-r--r--numpy/core/src/umath/_umath_tests.c.src46
-rw-r--r--numpy/core/src/umath/override.c12
-rw-r--r--numpy/core/src/umath/ufunc_object.c144
-rw-r--r--numpy/core/tests/test_scalarprint.py66
-rw-r--r--numpy/core/tests/test_ufunc.py63
-rw-r--r--numpy/core/tests/test_umath.py1
-rw-r--r--numpy/ctypeslib.py133
-rw-r--r--numpy/lib/tests/test_ufunclike.py4
-rw-r--r--numpy/testing/_private/utils.py76
-rw-r--r--numpy/testing/tests/test_utils.py18
-rw-r--r--numpy/tests/test_ctypeslib.py77
17 files changed, 735 insertions, 190 deletions
diff --git a/numpy/core/setup_common.py b/numpy/core/setup_common.py
index 450fb1b9d..70a43046c 100644
--- a/numpy/core/setup_common.py
+++ b/numpy/core/setup_common.py
@@ -335,9 +335,9 @@ _MOTOROLA_EXTENDED_12B = ['300', '031', '000', '000', '353', '171',
_IEEE_QUAD_PREC_BE = ['300', '031', '326', '363', '105', '100', '000', '000',
'000', '000', '000', '000', '000', '000', '000', '000']
_IEEE_QUAD_PREC_LE = _IEEE_QUAD_PREC_BE[::-1]
-_DOUBLE_DOUBLE_BE = (['301', '235', '157', '064', '124', '000', '000', '000'] +
+_IBM_DOUBLE_DOUBLE_BE = (['301', '235', '157', '064', '124', '000', '000', '000'] +
['000'] * 8)
-_DOUBLE_DOUBLE_LE = (['000', '000', '000', '124', '064', '157', '235', '301'] +
+_IBM_DOUBLE_DOUBLE_LE = (['000', '000', '000', '124', '064', '157', '235', '301'] +
['000'] * 8)
def long_double_representation(lines):
@@ -381,10 +381,10 @@ def long_double_representation(lines):
return 'IEEE_QUAD_BE'
elif read[8:-8] == _IEEE_QUAD_PREC_LE:
return 'IEEE_QUAD_LE'
- elif read[8:-8] == _DOUBLE_DOUBLE_BE:
- return 'DOUBLE_DOUBLE_BE'
- elif read[8:-8] == _DOUBLE_DOUBLE_LE:
- return 'DOUBLE_DOUBLE_LE'
+ elif read[8:-8] == _IBM_DOUBLE_DOUBLE_LE:
+ return 'IBM_DOUBLE_DOUBLE_LE'
+ elif read[8:-8] == _IBM_DOUBLE_DOUBLE_BE:
+ return 'IBM_DOUBLE_DOUBLE_BE'
# if the content was 8 bytes, left with 32-8-8 = 16 bytes
elif read[:16] == _BEFORE_SEQ:
if read[16:-8] == _IEEE_DOUBLE_LE:
diff --git a/numpy/core/src/multiarray/dragon4.c b/numpy/core/src/multiarray/dragon4.c
index 8f7dc0e2a..c14653ac5 100644
--- a/numpy/core/src/multiarray/dragon4.c
+++ b/numpy/core/src/multiarray/dragon4.c
@@ -215,7 +215,9 @@ BigInt_Set_uint64(BigInt *i, npy_uint64 val)
}
}
-#if defined(HAVE_LDOUBLE_IEEE_QUAD_LE)
+#if (defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE) || \
+ defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE) || \
+ defined(HAVE_LDOUBLE_IEEE_QUAD_LE))
static void
BigInt_Set_2x_uint64(BigInt *i, npy_uint64 hi, npy_uint64 lo)
{
@@ -247,7 +249,7 @@ BigInt_Set_2x_uint64(BigInt *i, npy_uint64 hi, npy_uint64 lo)
i->blocks[0] = lo & bitmask_u64(32);
}
}
-#endif /* QUAD */
+#endif /* DOUBLE_DOUBLE and QUAD */
static void
BigInt_Set_uint32(BigInt *i, npy_uint32 val)
@@ -2810,6 +2812,251 @@ Dragon4_PrintFloat_IEEE_binary128(
}
#endif /* HAVE_LDOUBLE_IEEE_QUAD_LE */
+#if (defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE) || \
+ defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE))
+/*
+ * IBM extended precision 128-bit floating-point format, aka IBM double-dobule
+ *
+ * IBM's double-double type is a pair of IEEE binary64 values, which you add
+ * together to get a total value. The exponents are arranged so that the lower
+ * double is about 2^52 times smaller than the high one, and the nearest
+ * float64 value is simply the upper double, in which case the pair is
+ * considered "normalized" (not to confuse with "normal" and "subnormal"
+ * binary64 values). We assume normalized values. You can see the glibc's
+ * printf on ppc does so too by constructing un-normalized values to get
+ * strange behavior from the OS printf:
+ *
+ * >>> from numpy.core._multiarray_tests import format_float_OSprintf_g
+ * >>> x = np.array([0.3,0.3], dtype='f8').view('f16')[0]
+ * >>> format_float_OSprintf_g(x, 2)
+ * 0.30
+ * >>> format_float_OSprintf_g(2*x, 2)
+ * 1.20
+ *
+ * If we don't assume normalization, x should really print as 0.6.
+ *
+ * For normalized values gcc assumes that the total mantissa is no
+ * more than 106 bits (53+53), so we can drop bits from the second double which
+ * would be pushed past 106 when left-shifting by its exponent, as happens
+ * sometimes. (There has been debate about this, see
+ * https://gcc.gnu.org/bugzilla/show_bug.cgi?format=multiple&id=70117,
+ * https://sourceware.org/bugzilla/show_bug.cgi?id=22752 )
+ *
+ * Note: This function is for the IBM-double-double which is a pair of IEEE
+ * binary64 floats, like on ppc64 systems. This is *not* the hexadecimal
+ * IBM-double-double type, which is a pair of IBM hexadecimal64 floats.
+ *
+ * See also:
+ * https://gcc.gnu.org/wiki/Ieee128PowerPCA
+ * https://www.ibm.com/support/knowledgecenter/en/ssw_aix_71/com.ibm.aix.genprogc/128bit_long_double_floating-point_datatype.htm
+ */
+static npy_uint32
+Dragon4_PrintFloat_IBM_double_double(
+ Dragon4_Scratch *scratch, FloatVal128 val128, Dragon4_Options *opt)
+{
+ char *buffer = scratch->repr;
+ npy_uint32 bufferSize = sizeof(scratch->repr);
+ BigInt *bigints = scratch->bigints;
+
+ npy_uint32 floatExponent1, floatExponent2;
+ npy_uint64 floatMantissa1, floatMantissa2;
+ npy_uint32 floatSign1, floatSign2;
+
+ npy_uint64 mantissa1, mantissa2;
+ npy_int32 exponent1, exponent2;
+ int shift;
+ npy_uint32 mantissaBit;
+ npy_bool hasUnequalMargins;
+ char signbit = '\0';
+
+ if (bufferSize == 0) {
+ return 0;
+ }
+
+ if (bufferSize == 1) {
+ buffer[0] = '\0';
+ return 0;
+ }
+
+ /* deconstruct the floating point values */
+ floatMantissa1 = val128.hi & bitmask_u64(52);
+ floatExponent1 = (val128.hi >> 52) & bitmask_u32(11);
+ floatSign1 = (val128.hi >> 63) != 0;
+
+ floatMantissa2 = val128.lo & bitmask_u64(52);
+ floatExponent2 = (val128.lo >> 52) & bitmask_u32(11);
+ floatSign2 = (val128.lo >> 63) != 0;
+
+ /* output the sign using 1st float's sign */
+ if (floatSign1) {
+ signbit = '-';
+ }
+ else if (opt->sign) {
+ signbit = '+';
+ }
+
+ /* we only need to look at the first float for inf/nan */
+ if (floatExponent1 == bitmask_u32(11)) {
+ return PrintInfNan(buffer, bufferSize, floatMantissa1, 13, signbit);
+ }
+
+ /* else this is a number */
+
+ /* Factor the 1st value into its parts, see binary64 for comments. */
+ if (floatExponent1 == 0) {
+ /*
+ * If the first number is a subnormal value, the 2nd has to be 0 for
+ * the float128 to be normalized, so we can ignore it. In this case
+ * the float128 only has the precision of a single binary64 value.
+ */
+ mantissa1 = floatMantissa1;
+ exponent1 = 1 - 1023 - 52;
+ mantissaBit = LogBase2_64(mantissa1);
+ hasUnequalMargins = NPY_FALSE;
+
+ BigInt_Set_uint64(&bigints[0], mantissa1);
+ }
+ else {
+ mantissa1 = (1ull << 52) | floatMantissa1;
+ exponent1 = floatExponent1 - 1023 - 52;
+ mantissaBit = 52 + 53;
+
+ /*
+ * Computing hasUnequalMargins and mantissaBit:
+ * This is a little trickier than for IEEE formats.
+ *
+ * When both doubles are "normal" it is clearer since we can think of
+ * it as an IEEE type with a 106 bit mantissa. This value can never
+ * have "unequal" margins because of the implied 1 bit in the 2nd
+ * value. (unequal margins only happen when the mantissa has a value
+ * like "10000000000...", all zeros except the implied bit at the
+ * start, since the next lowest number has a different exponent).
+ * mantissaBits will always be 52+53 in this case.
+ *
+ * If the 1st number is a very small normal, and the 2nd is subnormal
+ * and not 2^52 times smaller, the number behaves like a subnormal
+ * overall, where the upper number just adds some bits on the left.
+ * Like usual subnormals, it has "equal" margins. The slightly tricky
+ * thing is that the number of mantissaBits varies. It will be 52
+ * (from lower double) plus a variable number depending on the upper
+ * number's exponent. We recompute the number of bits in the shift
+ * calculation below, because the shift will be equal to the number of
+ * lost bits.
+ *
+ * We can get unequal margins only if the first value has all-0
+ * mantissa (except implied bit), and the second value is exactly 0. As
+ * a special exception the smallest normal value (smallest exponent, 0
+ * mantissa) should have equal margins, since it is "next to" a
+ * subnormal value.
+ */
+
+ /* factor the 2nd value into its parts */
+ if (floatExponent2 != 0) {
+ mantissa2 = (1ull << 52) | floatMantissa2;
+ exponent2 = floatExponent2 - 1023 - 52;
+ hasUnequalMargins = NPY_FALSE;
+ }
+ else {
+ /* shift exp by one so that leading mantissa bit is still bit 53 */
+ mantissa2 = floatMantissa2 << 1;
+ exponent2 = - 1023 - 52;
+ hasUnequalMargins = (floatExponent1 != 1) && (floatMantissa1 == 0)
+ && (floatMantissa2 == 0);
+ }
+
+ /*
+ * The 2nd val's exponent might not be exactly 52 smaller than the 1st,
+ * it can vary a little bit. So do some shifting of the low mantissa,
+ * so that the total mantissa is equivalent to bits 53 to 0 of the
+ * first double immediately followed by bits 53 to 0 of the second.
+ */
+ shift = exponent1 - exponent2 - 53;
+ if (shift > 0) {
+ /* shift more than 64 is undefined behavior */
+ mantissa2 = shift < 64 ? mantissa2 >> shift : 0;
+ }
+ else if (shift < 0) {
+ /*
+ * This only happens if the 2nd value is subnormal.
+ * We expect that shift > -64, but check it anyway
+ */
+ mantissa2 = -shift < 64 ? mantissa2 << -shift : 0;
+ }
+
+ /*
+ * If the low double is a different sign from the high double,
+ * rearrange so that the total mantissa is the sum of the two
+ * mantissas, instead of a subtraction.
+ * hi - lo -> (hi-1) + (1-lo), where lo < 1
+ */
+ if (floatSign1 != floatSign2 && mantissa2 != 0) {
+ mantissa1--;
+ mantissa2 = (1ull << 53) - mantissa2;
+ }
+
+ /*
+ * Compute the number of bits if we are in the subnormal range.
+ * The value "shift" happens to be exactly the number of lost bits.
+ * Also, shift the bits so that the least significant bit is at
+ * bit position 0, like a typical subnormal. After this exponent1
+ * should always be 2^-1022
+ */
+ if (shift < 0) {
+ mantissa2 = (mantissa2 >> -shift) | (mantissa1 << (53 + shift));
+ mantissa1 = mantissa1 >> -shift;
+ mantissaBit = mantissaBit -(-shift);
+ exponent1 -= shift;
+ DEBUG_ASSERT(exponent1 == -1022);
+ }
+
+ /*
+ * set up the BigInt mantissa, by shifting the parts as needed
+ * We can use | instead of + since the mantissas should not overlap
+ */
+ BigInt_Set_2x_uint64(&bigints[0], mantissa1 >> 11,
+ (mantissa1 << 53) | (mantissa2));
+ exponent1 = exponent1 - 53;
+ }
+
+ return Format_floatbits(buffer, bufferSize, bigints, exponent1,
+ signbit, mantissaBit, hasUnequalMargins, opt);
+}
+
+#if defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE)
+static npy_uint32
+Dragon4_PrintFloat_IBM_double_double_le(
+ Dragon4_Scratch *scratch, npy_float128 *value, Dragon4_Options *opt)
+{
+ FloatVal128 val128;
+ FloatUnion128 buf128;
+
+ buf128.floatingPoint = *value;
+ val128.lo = buf128.integer.a;
+ val128.hi = buf128.integer.b;
+
+ return Dragon4_PrintFloat_IBM_double_double(scratch, val128, opt);
+}
+#endif /* HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE */
+
+#if defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE)
+static npy_uint32
+Dragon4_PrintFloat_IBM_double_double_be(
+ Dragon4_Scratch *scratch, npy_float128 *value, Dragon4_Options *opt)
+{
+ FloatVal128 val128;
+ FloatUnion128 buf128;
+
+ buf128.floatingPoint = *value;
+ val128.hi = buf128.integer.a;
+ val128.lo = buf128.integer.b;
+
+ return Dragon4_PrintFloat_IBM_double_double(scratch, val128, opt);
+}
+
+#endif /* HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE */
+
+#endif /* HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE | HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE */
+
#endif /* NPY_FLOAT128 */
diff --git a/numpy/core/src/multiarray/dragon4.h b/numpy/core/src/multiarray/dragon4.h
index 176b79f69..383a0949d 100644
--- a/numpy/core/src/multiarray/dragon4.h
+++ b/numpy/core/src/multiarray/dragon4.h
@@ -75,6 +75,10 @@
#define NPY_LONGDOUBLE_BINFMT_NAME Intel_extended128
#elif defined(HAVE_LDOUBLE_MOTOROLA_EXTENDED_12_BYTES_BE)
#define NPY_LONGDOUBLE_BINFMT_NAME Motorola_extended96
+#elif defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE)
+ #define NPY_LONGDOUBLE_BINFMT_NAME IBM_double_double_le
+#elif defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE)
+ #define NPY_LONGDOUBLE_BINFMT_NAME IBM_double_double_be
#else
#error No long double representation defined
#endif
diff --git a/numpy/core/src/npymath/ieee754.c.src b/numpy/core/src/npymath/ieee754.c.src
index fcfd15ab4..ea1792887 100644
--- a/numpy/core/src/npymath/ieee754.c.src
+++ b/numpy/core/src/npymath/ieee754.c.src
@@ -798,7 +798,7 @@ int npy_clear_floatstatus_barrier(char *param)
#else
-int npy_get_floatstatus_barrier(char NPY_UNUSED(*param))
+int npy_get_floatstatus_barrier(char *NPY_UNUSED(param))
{
return 0;
}
diff --git a/numpy/core/src/npymath/npy_math_private.h b/numpy/core/src/npymath/npy_math_private.h
index 7ee564239..e4a919db6 100644
--- a/numpy/core/src/npymath/npy_math_private.h
+++ b/numpy/core/src/npymath/npy_math_private.h
@@ -434,8 +434,8 @@ do { \
typedef npy_uint32 ldouble_sign_t;
#endif
-#if !defined(HAVE_LDOUBLE_DOUBLE_DOUBLE_BE) && \
- !defined(HAVE_LDOUBLE_DOUBLE_DOUBLE_LE)
+#if !defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE) && \
+ !defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE)
/* Get the sign bit of x. x should be of type IEEEl2bitsrep */
#define GET_LDOUBLE_SIGN(x) \
(((x).a[LDBL_SIGN_INDEX] & LDBL_SIGN_MASK) >> LDBL_SIGN_SHIFT)
diff --git a/numpy/core/src/private/npy_fpmath.h b/numpy/core/src/private/npy_fpmath.h
index e1521de3b..dbb3fb23d 100644
--- a/numpy/core/src/private/npy_fpmath.h
+++ b/numpy/core/src/private/npy_fpmath.h
@@ -14,9 +14,17 @@
defined(HAVE_LDOUBLE_INTEL_EXTENDED_16_BYTES_LE) || \
defined(HAVE_LDOUBLE_INTEL_EXTENDED_12_BYTES_LE) || \
defined(HAVE_LDOUBLE_MOTOROLA_EXTENDED_12_BYTES_BE) || \
- defined(HAVE_LDOUBLE_DOUBLE_DOUBLE_BE) || \
- defined(HAVE_LDOUBLE_DOUBLE_DOUBLE_LE))
+ defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE) || \
+ defined(HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE))
#error No long double representation defined
#endif
+/* for back-compat, also keep old name for double-double */
+#ifdef HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_LE
+ #define HAVE_LDOUBLE_DOUBLE_DOUBLE_LE
+#endif
+#ifdef HAVE_LDOUBLE_IBM_DOUBLE_DOUBLE_BE
+ #define HAVE_LDOUBLE_DOUBLE_DOUBLE_BE
+#endif
+
#endif
diff --git a/numpy/core/src/umath/_umath_tests.c.src b/numpy/core/src/umath/_umath_tests.c.src
index 2a74c1aaa..fcbdbe330 100644
--- a/numpy/core/src/umath/_umath_tests.c.src
+++ b/numpy/core/src/umath/_umath_tests.c.src
@@ -253,6 +253,38 @@ static void
/**end repeat**/
+char *cumsum_signature = "(i)->(i)";
+
+/*
+ * This implements the function
+ * out[n] = sum_i^n in[i]
+ */
+
+/**begin repeat
+
+ #TYPE=LONG,DOUBLE#
+ #typ=npy_long,npy_double#
+*/
+
+static void
+@TYPE@_cumsum(char **args, npy_intp *dimensions, npy_intp *steps, void *NPY_UNUSED(func))
+{
+ INIT_OUTER_LOOP_2
+ npy_intp di = dimensions[0];
+ npy_intp i;
+ npy_intp is=steps[0], os=steps[1];
+ BEGIN_OUTER_LOOP_2
+ char *ip=args[0], *op=args[1];
+ @typ@ cumsum = 0;
+ for (i = 0; i < di; i++, ip += is, op += os) {
+ cumsum += (*(@typ@ *)ip);
+ *(@typ@ *)op = cumsum;
+ }
+ END_OUTER_LOOP
+}
+
+/**end repeat**/
+
static PyUFuncGenericFunction inner1d_functions[] = { LONG_inner1d, DOUBLE_inner1d };
static void * inner1d_data[] = { (void *)NULL, (void *)NULL };
@@ -270,6 +302,10 @@ static void *eucldiean_pdist_data[] = { (void *)NULL, (void *)NULL };
static char euclidean_pdist_signatures[] = { NPY_FLOAT, NPY_FLOAT,
NPY_DOUBLE, NPY_DOUBLE };
+static PyUFuncGenericFunction cumsum_functions[] = { LONG_cumsum, DOUBLE_cumsum };
+static void * cumsum_data[] = { (void *)NULL, (void *)NULL };
+static char cumsum_signatures[] = { NPY_LONG, NPY_LONG, NPY_DOUBLE, NPY_DOUBLE };
+
static int
addUfuncs(PyObject *dictionary) {
@@ -321,6 +357,16 @@ addUfuncs(PyObject *dictionary) {
}
PyDict_SetItemString(dictionary, "euclidean_pdist", f);
Py_DECREF(f);
+ f = PyUFunc_FromFuncAndDataAndSignature(cumsum_functions,
+ cumsum_data, cumsum_signatures,
+ 2, 1, 1, PyUFunc_None, "cumsum",
+ "Cumulative sum of the input (n)->(n)\n",
+ 0, cumsum_signature);
+ if (f == NULL) {
+ return -1;
+ }
+ PyDict_SetItemString(dictionary, "cumsum", f);
+ Py_DECREF(f);
f = PyUFunc_FromFuncAndDataAndSignature(inner1d_functions, inner1d_data,
inner1d_signatures, 2, 2, 1, PyUFunc_None, "inner1d_no_doc",
NULL,
diff --git a/numpy/core/src/umath/override.c b/numpy/core/src/umath/override.c
index c298fe315..c0bc47b7b 100644
--- a/numpy/core/src/umath/override.c
+++ b/numpy/core/src/umath/override.c
@@ -51,6 +51,7 @@ normalize___call___args(PyUFuncObject *ufunc, PyObject *args,
npy_intp nin = ufunc->nin;
npy_intp nout = ufunc->nout;
npy_intp nargs = PyTuple_GET_SIZE(args);
+ npy_intp nkwds = PyDict_Size(*normal_kwds);
PyObject *obj;
if (nargs < nin) {
@@ -74,7 +75,7 @@ normalize___call___args(PyUFuncObject *ufunc, PyObject *args,
/* If we have more args than nin, they must be the output variables.*/
if (nargs > nin) {
- if(PyDict_GetItemString(*normal_kwds, "out")) {
+ if(nkwds > 0 && PyDict_GetItemString(*normal_kwds, "out")) {
PyErr_Format(PyExc_TypeError,
"argument given by name ('out') and position "
"(%"NPY_INTP_FMT")", nin);
@@ -112,8 +113,15 @@ normalize___call___args(PyUFuncObject *ufunc, PyObject *args,
Py_DECREF(obj);
}
}
+ /* gufuncs accept either 'axes' or 'axis', but not both */
+ if (nkwds >= 2 && (PyDict_GetItemString(*normal_kwds, "axis") &&
+ PyDict_GetItemString(*normal_kwds, "axes"))) {
+ PyErr_SetString(PyExc_TypeError,
+ "cannot specify both 'axis' and 'axes'");
+ return -1;
+ }
/* finally, ufuncs accept 'sig' or 'signature' normalize to 'signature' */
- return normalize_signature_keyword(*normal_kwds);
+ return nkwds == 0 ? 0 : normalize_signature_keyword(*normal_kwds);
}
static int
diff --git a/numpy/core/src/umath/ufunc_object.c b/numpy/core/src/umath/ufunc_object.c
index 631999bc1..9b03a7916 100644
--- a/numpy/core/src/umath/ufunc_object.c
+++ b/numpy/core/src/umath/ufunc_object.c
@@ -566,11 +566,12 @@ get_ufunc_arguments(PyUFuncObject *ufunc,
NPY_ORDER *out_order,
NPY_CASTING *out_casting,
PyObject **out_extobj,
- PyObject **out_typetup,
- int *out_subok,
- PyArrayObject **out_wheremask,
- PyObject **out_axes,
- int *out_keepdims)
+ PyObject **out_typetup, /* type: Tuple[np.dtype] */
+ int *out_subok, /* bool */
+ PyArrayObject **out_wheremask, /* PyArray of bool */
+ PyObject **out_axes, /* type: List[Tuple[T]] */
+ PyObject **out_axis, /* type: T */
+ int *out_keepdims) /* bool */
{
int i, nargs;
int nin = ufunc->nin;
@@ -588,6 +589,9 @@ get_ufunc_arguments(PyUFuncObject *ufunc,
if (out_axes != NULL) {
*out_axes = NULL;
}
+ if (out_axis != NULL) {
+ *out_axis = NULL;
+ }
if (out_wheremask != NULL) {
*out_wheremask = NULL;
}
@@ -827,9 +831,23 @@ get_ufunc_arguments(PyUFuncObject *ufunc,
case 'a':
/* possible axes argument for generalized ufunc */
if (out_axes != NULL && strcmp(str, "axes") == 0) {
+ if (out_axis != NULL && *out_axis != NULL) {
+ PyErr_SetString(PyExc_TypeError,
+ "cannot specify both 'axis' and 'axes'");
+ goto fail;
+ }
Py_INCREF(value);
*out_axes = value;
-
+ bad_arg = 0;
+ }
+ else if (out_axis != NULL && strcmp(str, "axis") == 0) {
+ if (out_axes != NULL && *out_axes != NULL) {
+ PyErr_SetString(PyExc_TypeError,
+ "cannot specify both 'axis' and 'axes'");
+ goto fail;
+ }
+ Py_INCREF(value);
+ *out_axis = value;
bad_arg = 0;
}
break;
@@ -1045,6 +1063,10 @@ fail:
Py_XDECREF(*out_axes);
*out_axes = NULL;
}
+ if (out_axis != NULL) {
+ Py_XDECREF(*out_axis);
+ *out_axis = NULL;
+ }
return -1;
}
@@ -1891,6 +1913,27 @@ _has_output_coredims(PyUFuncObject *ufunc) {
}
/*
+ * Check whether the gufunc can be used with axis, i.e., that there is only
+ * a single, shared core dimension (which means that operands either have
+ * that dimension, or have no core dimensions). Returns 0 if all is fine,
+ * and sets an error and returns -1 if not.
+ */
+static int
+_check_axis_support(PyUFuncObject *ufunc) {
+ if (ufunc->core_num_dim_ix != 1) {
+ PyErr_Format(PyExc_TypeError,
+ "%s: axis can only be used with a single shared core "
+ "dimension, not with the %d distinct ones implied by "
+ "signature %s.",
+ ufunc_get_name_cstr(ufunc),
+ ufunc->core_num_dim_ix,
+ ufunc->core_signature);
+ return -1;
+ }
+ return 0;
+}
+
+/*
* Check whether the gufunc can be used with keepdims, i.e., that all its
* input arguments have the same number of core dimension, and all output
* arguments have no core dimensions. Returns 0 if all is fine, and sets
@@ -1905,7 +1948,7 @@ _check_keepdims_support(PyUFuncObject *ufunc) {
if (ufunc->core_num_dims[i] != (i < nin ? input_core_dims : 0)) {
PyErr_Format(PyExc_TypeError,
"%s does not support keepdims: its signature %s requires "
- "that %s %d has %d core dimensions, but keepdims can only "
+ "%s %d to have %d core dimensions, but keepdims can only "
"be used when all inputs have the same number of core "
"dimensions and all outputs have no core dimensions.",
ufunc_get_name_cstr(ufunc),
@@ -1931,8 +1974,7 @@ static int
_parse_axes_arg(PyUFuncObject *ufunc, int core_num_dims[], PyObject *axes,
PyArrayObject **op, int broadcast_ndim, int **remap_axis) {
int nin = ufunc->nin;
- int nout = ufunc->nout;
- int nop = nin + nout;
+ int nop = ufunc->nargs;
int iop, list_size;
if (!PyList_Check(axes)) {
@@ -2050,6 +2092,59 @@ _parse_axes_arg(PyUFuncObject *ufunc, int core_num_dims[], PyObject *axes,
return 0;
}
+/*
+ * Simplified version of the above, using axis to fill the remap_axis
+ * array, which maps default to actual axes for each operand, indexed as
+ * as remap_axis[iop][iaxis]. The default axis order has first all broadcast
+ * axes and then the core axes the gufunc operates on.
+ *
+ * Returns 0 on success, and -1 on failure
+ */
+static int
+_parse_axis_arg(PyUFuncObject *ufunc, int core_num_dims[], PyObject *axis,
+ PyArrayObject **op, int broadcast_ndim, int **remap_axis) {
+ int nop = ufunc->nargs;
+ int iop, axis_int;
+
+ axis_int = PyArray_PyIntAsInt(axis);
+ if (error_converting(axis_int)) {
+ return -1;
+ }
+
+ for (iop = 0; iop < nop; ++iop) {
+ int axis, op_ndim, op_axis;
+
+ /* _check_axis_support ensures core_num_dims is 0 or 1 */
+ if (core_num_dims[iop] == 0) {
+ remap_axis[iop] = NULL;
+ continue;
+ }
+ if (op[iop]) {
+ op_ndim = PyArray_NDIM(op[iop]);
+ }
+ else {
+ op_ndim = broadcast_ndim + 1;
+ }
+ op_axis = axis_int; /* ensure we don't modify axis_int */
+ if (check_and_adjust_axis(&op_axis, op_ndim) < 0) {
+ return -1;
+ }
+ /* Are we actually remapping away from last axis? */
+ if (op_axis == op_ndim - 1) {
+ remap_axis[iop] = NULL;
+ continue;
+ }
+ remap_axis[iop][op_ndim - 1] = op_axis;
+ for (axis = 0; axis < op_axis; axis++) {
+ remap_axis[iop][axis] = axis;
+ }
+ for (axis = op_axis; axis < op_ndim - 1; axis++) {
+ remap_axis[iop][axis] = axis + 1;
+ }
+ } /* end of for(iop) loop over operands */
+ return 0;
+}
+
#define REMAP_AXIS(iop, axis) ((remap_axis != NULL && \
remap_axis[iop] != NULL)? \
remap_axis[iop][axis] : axis)
@@ -2245,7 +2340,8 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
NPY_ORDER order = NPY_KEEPORDER;
/* Use the default assignment casting rule */
NPY_CASTING casting = NPY_DEFAULT_ASSIGN_CASTING;
- PyObject *extobj = NULL, *type_tup = NULL, *axes = NULL;
+ /* other possible keyword arguments */
+ PyObject *extobj = NULL, *type_tup = NULL, *axes = NULL, *axis = NULL;
int keepdims = -1;
if (ufunc == NULL) {
@@ -2270,10 +2366,12 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
NPY_UF_DBG_PRINT("Getting arguments\n");
- /* Get all the arguments */
+ /*
+ * Get all the arguments.
+ */
retval = get_ufunc_arguments(ufunc, args, kwds,
op, &order, &casting, &extobj,
- &type_tup, &subok, NULL, &axes, &keepdims);
+ &type_tup, &subok, NULL, &axes, &axis, &keepdims);
if (retval < 0) {
goto fail;
}
@@ -2288,6 +2386,12 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
goto fail;
}
}
+ if (axis != NULL) {
+ retval = _check_axis_support(ufunc);
+ if (retval < 0) {
+ goto fail;
+ }
+ }
/*
* If keepdims is set and true, signal all dimensions will be the same.
*/
@@ -2356,7 +2460,7 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
}
/* Possibly remap axes. */
- if (axes) {
+ if (axes != NULL || axis != NULL) {
remap_axis = PyArray_malloc(sizeof(remap_axis[0]) * nop);
remap_axis_memory = PyArray_malloc(sizeof(remap_axis_memory[0]) *
nop * NPY_MAXDIMS);
@@ -2367,8 +2471,14 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
for (i=0; i < nop; i++) {
remap_axis[i] = remap_axis_memory + i * NPY_MAXDIMS;
}
- retval = _parse_axes_arg(ufunc, core_num_dims, axes, op, broadcast_ndim,
- remap_axis);
+ if (axis) {
+ retval = _parse_axis_arg(ufunc, core_num_dims, axis, op,
+ broadcast_ndim, remap_axis);
+ }
+ else {
+ retval = _parse_axes_arg(ufunc, core_num_dims, axes, op,
+ broadcast_ndim, remap_axis);
+ }
if(retval < 0) {
goto fail;
}
@@ -2714,6 +2824,7 @@ PyUFunc_GeneralizedFunction(PyUFuncObject *ufunc,
Py_XDECREF(type_tup);
Py_XDECREF(extobj);
Py_XDECREF(axes);
+ Py_XDECREF(axis);
Py_XDECREF(full_args.in);
Py_XDECREF(full_args.out);
@@ -2735,6 +2846,7 @@ fail:
Py_XDECREF(type_tup);
Py_XDECREF(extobj);
Py_XDECREF(axes);
+ Py_XDECREF(axis);
Py_XDECREF(full_args.in);
Py_XDECREF(full_args.out);
PyArray_free(remap_axis_memory);
@@ -2812,7 +2924,7 @@ PyUFunc_GenericFunction(PyUFuncObject *ufunc,
/* Get all the arguments */
retval = get_ufunc_arguments(ufunc, args, kwds,
op, &order, &casting, &extobj,
- &type_tup, &subok, &wheremask, NULL, NULL);
+ &type_tup, &subok, &wheremask, NULL, NULL, NULL);
if (retval < 0) {
goto fail;
}
diff --git a/numpy/core/tests/test_scalarprint.py b/numpy/core/tests/test_scalarprint.py
index a20ec9f74..472ff691d 100644
--- a/numpy/core/tests/test_scalarprint.py
+++ b/numpy/core/tests/test_scalarprint.py
@@ -5,10 +5,12 @@
from __future__ import division, absolute_import, print_function
import code, sys
+import platform
+import pytest
+
from tempfile import TemporaryFile
import numpy as np
-from numpy.testing import assert_, assert_equal, suppress_warnings
-
+from numpy.testing import assert_, assert_equal, suppress_warnings, dec
class TestRealScalars(object):
def test_str(self):
@@ -250,6 +252,66 @@ class TestRealScalars(object):
"1.2" if tp != np.float16 else "1.2002")
assert_equal(fpos(tp('1.'), trim='-'), "1")
+ @pytest.mark.skipif(not platform.machine().startswith("ppc64"),
+ reason="only applies to ppc float128 values")
+ def test_ppc64_ibm_double_double128(self):
+ # check that the precision decreases once we get into the subnormal
+ # range. Unlike float64, this starts around 1e-292 instead of 1e-308,
+ # which happens when the first double is normal and the second is
+ # subnormal.
+ x = np.float128('2.123123123123123123123123123123123e-286')
+ got = [str(x/np.float128('2e' + str(i))) for i in range(0,40)]
+ expected = [
+ "1.06156156156156156156156156156157e-286",
+ "1.06156156156156156156156156156158e-287",
+ "1.06156156156156156156156156156159e-288",
+ "1.0615615615615615615615615615616e-289",
+ "1.06156156156156156156156156156157e-290",
+ "1.06156156156156156156156156156156e-291",
+ "1.0615615615615615615615615615616e-292",
+ "1.0615615615615615615615615615615e-293",
+ "1.061561561561561561561561561562e-294",
+ "1.06156156156156156156156156155e-295",
+ "1.0615615615615615615615615616e-296",
+ "1.06156156156156156156156156e-297",
+ "1.06156156156156156156156157e-298",
+ "1.0615615615615615615615616e-299",
+ "1.06156156156156156156156e-300",
+ "1.06156156156156156156155e-301",
+ "1.0615615615615615615616e-302",
+ "1.061561561561561561562e-303",
+ "1.06156156156156156156e-304",
+ "1.0615615615615615618e-305",
+ "1.06156156156156156e-306",
+ "1.06156156156156157e-307",
+ "1.0615615615615616e-308",
+ "1.06156156156156e-309",
+ "1.06156156156157e-310",
+ "1.0615615615616e-311",
+ "1.06156156156e-312",
+ "1.06156156154e-313",
+ "1.0615615616e-314",
+ "1.06156156e-315",
+ "1.06156155e-316",
+ "1.061562e-317",
+ "1.06156e-318",
+ "1.06155e-319",
+ "1.0617e-320",
+ "1.06e-321",
+ "1.04e-322",
+ "1e-323",
+ "0.0",
+ "0.0"]
+ assert_equal(got, expected)
+
+ # Note: we follow glibc behavior, but it (or gcc) might not be right.
+ # In particular we can get two values that print the same but are not
+ # equal:
+ a = np.float128('2')/np.float128('3')
+ b = np.float128(str(a))
+ assert_equal(str(a), str(b))
+ assert_(a != b)
+
def float32_roundtrip(self):
# gh-9360
x = np.float32(1024 - 2**-14)
diff --git a/numpy/core/tests/test_ufunc.py b/numpy/core/tests/test_ufunc.py
index ac4346461..ef9ced354 100644
--- a/numpy/core/tests/test_ufunc.py
+++ b/numpy/core/tests/test_ufunc.py
@@ -726,6 +726,50 @@ class TestUfunc(object):
# should be able to deal with bad unrelated kwargs.
assert_raises(TypeError, mm, z, z, axes=[0, 1], parrot=True)
+ def test_axis_argument(self):
+ # inner1d signature: '(i),(i)->()'
+ inner1d = umt.inner1d
+ a = np.arange(27.).reshape((3, 3, 3))
+ b = np.arange(10., 19.).reshape((3, 1, 3))
+ c = inner1d(a, b)
+ assert_array_equal(c, (a * b).sum(-1))
+ c = inner1d(a, b, axis=-1)
+ assert_array_equal(c, (a * b).sum(-1))
+ out = np.zeros_like(c)
+ d = inner1d(a, b, axis=-1, out=out)
+ assert_(d is out)
+ assert_array_equal(d, c)
+ c = inner1d(a, b, axis=0)
+ assert_array_equal(c, (a * b).sum(0))
+ # Sanity checks on innerwt and cumsum.
+ a = np.arange(6).reshape((2, 3))
+ b = np.arange(10, 16).reshape((2, 3))
+ w = np.arange(20, 26).reshape((2, 3))
+ assert_array_equal(umt.innerwt(a, b, w, axis=0),
+ np.sum(a * b * w, axis=0))
+ assert_array_equal(umt.cumsum(a, axis=0), np.cumsum(a, axis=0))
+ assert_array_equal(umt.cumsum(a, axis=-1), np.cumsum(a, axis=-1))
+ out = np.empty_like(a)
+ b = umt.cumsum(a, out=out, axis=0)
+ assert_(out is b)
+ assert_array_equal(b, np.cumsum(a, axis=0))
+ b = umt.cumsum(a, out=out, axis=1)
+ assert_(out is b)
+ assert_array_equal(b, np.cumsum(a, axis=-1))
+ # Check errors.
+ # Cannot pass in both axis and axes.
+ assert_raises(TypeError, inner1d, a, b, axis=0, axes=[0, 0])
+ # Not an integer.
+ assert_raises(TypeError, inner1d, a, b, axis=[0])
+ # more than 1 core dimensions.
+ mm = umt.matrix_multiply
+ assert_raises(TypeError, mm, a, b, axis=1)
+ # Output wrong size in axis.
+ out = np.empty((1, 2, 3), dtype=a.dtype)
+ assert_raises(ValueError, umt.cumsum, a, out=out, axis=0)
+ # Regular ufuncs should not accept axis.
+ assert_raises(TypeError, np.add, 1., 1., axis=0)
+
def test_keepdims_argument(self):
# inner1d signature: '(i),(i)->()'
inner1d = umt.inner1d
@@ -741,7 +785,15 @@ class TestUfunc(object):
d = inner1d(a, b, keepdims=True, out=out)
assert_(d is out)
assert_array_equal(d, c)
- # Now combined with axes.
+ # Now combined with axis and axes.
+ c = inner1d(a, b, axis=-1, keepdims=False)
+ assert_array_equal(c, (a * b).sum(-1, keepdims=False))
+ c = inner1d(a, b, axis=-1, keepdims=True)
+ assert_array_equal(c, (a * b).sum(-1, keepdims=True))
+ c = inner1d(a, b, axis=0, keepdims=False)
+ assert_array_equal(c, (a * b).sum(0, keepdims=False))
+ c = inner1d(a, b, axis=0, keepdims=True)
+ assert_array_equal(c, (a * b).sum(0, keepdims=True))
c = inner1d(a, b, axes=[(-1,), (-1,), ()], keepdims=False)
assert_array_equal(c, (a * b).sum(-1))
c = inner1d(a, b, axes=[(-1,), (-1,), (-1,)], keepdims=True)
@@ -785,10 +837,12 @@ class TestUfunc(object):
w = np.arange(20, 26).reshape((2, 3))
assert_array_equal(umt.innerwt(a, b, w, keepdims=True),
np.sum(a * b * w, axis=-1, keepdims=True))
+ assert_array_equal(umt.innerwt(a, b, w, axis=0, keepdims=True),
+ np.sum(a * b * w, axis=0, keepdims=True))
# Check errors.
# Not a boolean
assert_raises(TypeError, inner1d, a, b, keepdims='true')
- # 1 core dimension only.
+ # More than 1 core dimension, and core output dimensions.
mm = umt.matrix_multiply
assert_raises(TypeError, mm, a, b, keepdims=True)
assert_raises(TypeError, mm, a, b, keepdims=False)
@@ -886,6 +940,11 @@ class TestUfunc(object):
# An output array is required to determine p with signature (n,d)->(p)
assert_raises(ValueError, umt.euclidean_pdist, a)
+ def test_cumsum(self):
+ a = np.arange(10)
+ result = umt.cumsum(a)
+ assert_array_equal(result, a.cumsum())
+
def test_object_logical(self):
a = np.array([3, None, True, False, "test", ""], dtype=object)
assert_equal(np.logical_or(a, None),
diff --git a/numpy/core/tests/test_umath.py b/numpy/core/tests/test_umath.py
index 93ec73094..3c0d1759a 100644
--- a/numpy/core/tests/test_umath.py
+++ b/numpy/core/tests/test_umath.py
@@ -1810,6 +1810,7 @@ class TestSpecialMethods(object):
assert_raises(TypeError, np.multiply, a)
assert_raises(TypeError, np.multiply, a, a, a, a)
assert_raises(TypeError, np.multiply, a, a, sig='a', signature='a')
+ assert_raises(TypeError, ncu_tests.inner1d, a, a, axis=0, axes=[0, 0])
# reduce, positional args
res = np.multiply.reduce(a, 'axis0', 'dtype0', 'out0', 'keep0')
diff --git a/numpy/ctypeslib.py b/numpy/ctypeslib.py
index b8457c78b..9d71adbdb 100644
--- a/numpy/ctypeslib.py
+++ b/numpy/ctypeslib.py
@@ -319,120 +319,47 @@ def ndpointer(dtype=None, ndim=None, shape=None, flags=None):
_pointer_type_cache[(dtype, shape, ndim, num)] = klass
return klass
-if ctypes is not None:
- ct = ctypes
- ################################################################
- # simple types
-
- # maps the numpy typecodes like '<f8' to simple ctypes types like
- # c_double. Filled in by prep_simple.
- _typecodes = {}
-
- def prep_simple(simple_type, dtype):
- """Given a ctypes simple type, construct and attach an
- __array_interface__ property to it if it does not yet have one.
- """
- try: simple_type.__array_interface__
- except AttributeError: pass
- else: return
-
- typestr = _dtype(dtype).str
- _typecodes[typestr] = simple_type
-
- def __array_interface__(self):
- return {'descr': [('', typestr)],
- '__ref': self,
- 'strides': None,
- 'shape': (),
- 'version': 3,
- 'typestr': typestr,
- 'data': (ct.addressof(self), False),
- }
-
- simple_type.__array_interface__ = property(__array_interface__)
+def _get_typecodes():
+ """ Return a dictionary mapping __array_interface__ formats to ctypes types """
+ ct = ctypes
simple_types = [
- ((ct.c_byte, ct.c_short, ct.c_int, ct.c_long, ct.c_longlong), "i"),
- ((ct.c_ubyte, ct.c_ushort, ct.c_uint, ct.c_ulong, ct.c_ulonglong), "u"),
- ((ct.c_float, ct.c_double), "f"),
+ ct.c_byte, ct.c_short, ct.c_int, ct.c_long, ct.c_longlong,
+ ct.c_ubyte, ct.c_ushort, ct.c_uint, ct.c_ulong, ct.c_ulonglong,
+ ct.c_float, ct.c_double,
]
- # Prep that numerical ctypes types:
- for types, code in simple_types:
- for tp in types:
- prep_simple(tp, "%c%d" % (code, ct.sizeof(tp)))
+ return {_dtype(ctype).str: ctype for ctype in simple_types}
- ################################################################
- # array types
- _ARRAY_TYPE = type(ct.c_int * 1)
+def _ctype_ndarray(element_type, shape):
+ """ Create an ndarray of the given element type and shape """
+ for dim in shape[::-1]:
+ element_type = element_type * dim
+ return element_type
- def prep_array(array_type):
- """Given a ctypes array type, construct and attach an
- __array_interface__ property to it if it does not yet have one.
- """
- try: array_type.__array_interface__
- except AttributeError: pass
- else: return
-
- shape = []
- ob = array_type
- while type(ob) is _ARRAY_TYPE:
- shape.append(ob._length_)
- ob = ob._type_
- shape = tuple(shape)
- ai = ob().__array_interface__
- descr = ai['descr']
- typestr = ai['typestr']
-
- def __array_interface__(self):
- return {'descr': descr,
- '__ref': self,
- 'strides': None,
- 'shape': shape,
- 'version': 3,
- 'typestr': typestr,
- 'data': (ct.addressof(self), False),
- }
-
- array_type.__array_interface__ = property(__array_interface__)
-
- def prep_pointer(pointer_obj, shape):
- """Given a ctypes pointer object, construct and
- attach an __array_interface__ property to it if it does not
- yet have one.
- """
- try: pointer_obj.__array_interface__
- except AttributeError: pass
- else: return
-
- contents = pointer_obj.contents
- dtype = _dtype(type(contents))
-
- inter = {'version': 3,
- 'typestr': dtype.str,
- 'data': (ct.addressof(contents), False),
- 'shape': shape}
-
- pointer_obj.__array_interface__ = inter
- ################################################################
- # public functions
+if ctypes is not None:
+ _typecodes = _get_typecodes()
def as_array(obj, shape=None):
- """Create a numpy array from a ctypes array or a ctypes POINTER.
+ """
+ Create a numpy array from a ctypes array or POINTER.
+
The numpy array shares the memory with the ctypes object.
- The size parameter must be given if converting from a ctypes POINTER.
- The size parameter is ignored if converting from a ctypes array
+ The shape parameter must be given if converting from a ctypes POINTER.
+ The shape parameter is ignored if converting from a ctypes array
"""
- tp = type(obj)
- try: tp.__array_interface__
- except AttributeError:
- if hasattr(obj, 'contents'):
- prep_pointer(obj, shape)
- else:
- prep_array(tp)
+ if isinstance(obj, ctypes._Pointer):
+ # convert pointers to an array of the desired shape
+ if shape is None:
+ raise TypeError(
+ 'as_array() requires a shape argument when called on a '
+ 'pointer')
+ p_arr_type = ctypes.POINTER(_ctype_ndarray(obj._type_, shape))
+ obj = ctypes.cast(obj, p_arr_type).contents
+
return array(obj, copy=False)
def as_ctypes(obj):
@@ -446,9 +373,7 @@ if ctypes is not None:
addr, readonly = ai["data"]
if readonly:
raise TypeError("readonly arrays unsupported")
- tp = _typecodes[ai["typestr"]]
- for dim in ai["shape"][::-1]:
- tp = tp * dim
+ tp = _ctype_ndarray(_typecodes[ai["typestr"]], ai["shape"])
result = tp.from_address(addr)
result.__keep = ai
return result
diff --git a/numpy/lib/tests/test_ufunclike.py b/numpy/lib/tests/test_ufunclike.py
index ad006fe17..5604b3744 100644
--- a/numpy/lib/tests/test_ufunclike.py
+++ b/numpy/lib/tests/test_ufunclike.py
@@ -55,6 +55,10 @@ class TestUfunclike(object):
obj.metadata = self.metadata
return obj
+ def __array_finalize__(self, obj):
+ self.metadata = getattr(obj, 'metadata', None)
+ return self
+
a = nx.array([1.1, -1.1])
m = MyArray(a, metadata='foo')
f = ufl.fix(m)
diff --git a/numpy/testing/_private/utils.py b/numpy/testing/_private/utils.py
index 528d28b06..032c4a116 100644
--- a/numpy/testing/_private/utils.py
+++ b/numpy/testing/_private/utils.py
@@ -686,7 +686,7 @@ def assert_array_compare(comparison, x, y, err_msg='', verbose=True,
header='', precision=6, equal_nan=True,
equal_inf=True):
__tracebackhide__ = True # Hide traceback for py.test
- from numpy.core import array, isnan, isinf, any, inf
+ from numpy.core import array, isnan, inf, bool_
x = array(x, copy=False, subok=True)
y = array(y, copy=False, subok=True)
@@ -696,17 +696,28 @@ def assert_array_compare(comparison, x, y, err_msg='', verbose=True,
def istime(x):
return x.dtype.char in "Mm"
- def chk_same_position(x_id, y_id, hasval='nan'):
- """Handling nan/inf: check that x and y have the nan/inf at the same
- locations."""
- try:
- assert_array_equal(x_id, y_id)
- except AssertionError:
+ def func_assert_same_pos(x, y, func=isnan, hasval='nan'):
+ """Handling nan/inf: combine results of running func on x and y,
+ checking that they are True at the same locations."""
+ # Both the != True comparison here and the cast to bool_ at
+ # the end are done to deal with `masked`, which cannot be
+ # compared usefully, and for which .all() yields masked.
+ x_id = func(x)
+ y_id = func(y)
+ if (x_id == y_id).all() != True:
msg = build_err_msg([x, y],
err_msg + '\nx and y %s location mismatch:'
% (hasval), verbose=verbose, header=header,
names=('x', 'y'), precision=precision)
raise AssertionError(msg)
+ # If there is a scalar, then here we know the array has the same
+ # flag as it everywhere, so we should return the scalar flag.
+ if x_id.ndim == 0:
+ return bool_(x_id)
+ elif y_id.ndim == 0:
+ return bool_(y_id)
+ else:
+ return y_id
try:
cond = (x.shape == () or y.shape == ()) or x.shape == y.shape
@@ -719,49 +730,32 @@ def assert_array_compare(comparison, x, y, err_msg='', verbose=True,
names=('x', 'y'), precision=precision)
raise AssertionError(msg)
+ flagged = bool_(False)
if isnumber(x) and isnumber(y):
- has_nan = has_inf = False
if equal_nan:
- x_isnan, y_isnan = isnan(x), isnan(y)
- # Validate that NaNs are in the same place
- has_nan = any(x_isnan) or any(y_isnan)
- if has_nan:
- chk_same_position(x_isnan, y_isnan, hasval='nan')
+ flagged = func_assert_same_pos(x, y, func=isnan, hasval='nan')
if equal_inf:
- x_isinf, y_isinf = isinf(x), isinf(y)
- # Validate that infinite values are in the same place
- has_inf = any(x_isinf) or any(y_isinf)
- if has_inf:
- # Check +inf and -inf separately, since they are different
- chk_same_position(x == +inf, y == +inf, hasval='+inf')
- chk_same_position(x == -inf, y == -inf, hasval='-inf')
-
- if has_nan and has_inf:
- x = x[~(x_isnan | x_isinf)]
- y = y[~(y_isnan | y_isinf)]
- elif has_nan:
- x = x[~x_isnan]
- y = y[~y_isnan]
- elif has_inf:
- x = x[~x_isinf]
- y = y[~y_isinf]
-
- # Only do the comparison if actual values are left
- if x.size == 0:
- return
+ flagged |= func_assert_same_pos(x, y,
+ func=lambda xy: xy == +inf,
+ hasval='+inf')
+ flagged |= func_assert_same_pos(x, y,
+ func=lambda xy: xy == -inf,
+ hasval='-inf')
elif istime(x) and istime(y):
# If one is datetime64 and the other timedelta64 there is no point
if equal_nan and x.dtype.type == y.dtype.type:
- x_isnat, y_isnat = isnat(x), isnat(y)
-
- if any(x_isnat) or any(y_isnat):
- chk_same_position(x_isnat, y_isnat, hasval="NaT")
+ flagged = func_assert_same_pos(x, y, func=isnat, hasval="NaT")
- if any(x_isnat) or any(y_isnat):
- x = x[~x_isnat]
- y = y[~y_isnat]
+ if flagged.ndim > 0:
+ x, y = x[~flagged], y[~flagged]
+ # Only do the comparison if actual values are left
+ if x.size == 0:
+ return
+ elif flagged:
+ # no sense doing comparison if everything is flagged.
+ return
val = comparison(x, y)
diff --git a/numpy/testing/tests/test_utils.py b/numpy/testing/tests/test_utils.py
index 602cdf5f2..465c217d4 100644
--- a/numpy/testing/tests/test_utils.py
+++ b/numpy/testing/tests/test_utils.py
@@ -151,6 +151,17 @@ class TestArrayEqual(_GenericTest):
self._test_not_equal(c, b)
assert_equal(len(l), 1)
+ def test_masked_nan_inf(self):
+ # Regression test for gh-11121
+ a = np.ma.MaskedArray([3., 4., 6.5], mask=[False, True, False])
+ b = np.array([3., np.nan, 6.5])
+ self._test_equal(a, b)
+ self._test_equal(b, a)
+ a = np.ma.MaskedArray([3., 4., 6.5], mask=[True, False, False])
+ b = np.array([np.inf, 4., 6.5])
+ self._test_equal(a, b)
+ self._test_equal(b, a)
+
class TestBuildErrorMessage(object):
@@ -390,6 +401,9 @@ class TestArrayAlmostEqual(_GenericTest):
# comparison operators, not on them being able to store booleans
# (which, e.g., astropy Quantity cannot usefully do). See gh-8452.
class MyArray(np.ndarray):
+ def __eq__(self, other):
+ return super(MyArray, self).__eq__(other).view(np.ndarray)
+
def __lt__(self, other):
return super(MyArray, self).__lt__(other).view(np.ndarray)
@@ -489,6 +503,9 @@ class TestAlmostEqual(_GenericTest):
# comparison operators, not on them being able to store booleans
# (which, e.g., astropy Quantity cannot usefully do). See gh-8452.
class MyArray(np.ndarray):
+ def __eq__(self, other):
+ return super(MyArray, self).__eq__(other).view(np.ndarray)
+
def __lt__(self, other):
return super(MyArray, self).__lt__(other).view(np.ndarray)
@@ -650,6 +667,7 @@ class TestArrayAssertLess(object):
assert_raises(AssertionError, lambda: self._assert_func(-ainf, -x))
self._assert_func(-ainf, x)
+
@pytest.mark.skip(reason="The raises decorator depends on Nose")
class TestRaises(object):
diff --git a/numpy/tests/test_ctypeslib.py b/numpy/tests/test_ctypeslib.py
index 0f0d6dbc4..75ce9c8ca 100644
--- a/numpy/tests/test_ctypeslib.py
+++ b/numpy/tests/test_ctypeslib.py
@@ -4,9 +4,9 @@ import sys
import pytest
import numpy as np
-from numpy.ctypeslib import ndpointer, load_library
+from numpy.ctypeslib import ndpointer, load_library, as_array
from numpy.distutils.misc_util import get_shared_lib_extension
-from numpy.testing import assert_, assert_raises
+from numpy.testing import assert_, assert_array_equal, assert_raises, assert_equal
try:
cdll = None
@@ -21,11 +21,12 @@ try:
except ImportError:
_HAS_CTYPE = False
+
+@pytest.mark.skipif(not _HAS_CTYPE,
+ reason="ctypes not available in this python")
+@pytest.mark.skipif(sys.platform == 'cygwin',
+ reason="Known to fail on cygwin")
class TestLoadLibrary(object):
- @pytest.mark.skipif(not _HAS_CTYPE,
- reason="ctypes not available in this python")
- @pytest.mark.skipif(sys.platform == 'cygwin',
- reason="Known to fail on cygwin")
def test_basic(self):
try:
# Should succeed
@@ -35,10 +36,6 @@ class TestLoadLibrary(object):
" (import error was: %s)" % str(e))
print(msg)
- @pytest.mark.skipif(not _HAS_CTYPE,
- reason="ctypes not available in this python")
- @pytest.mark.skipif(sys.platform == 'cygwin',
- reason="Known to fail on cygwin")
def test_basic2(self):
# Regression for #801: load_library with a full library name
# (including extension) does not work.
@@ -54,6 +51,7 @@ class TestLoadLibrary(object):
" (import error was: %s)" % str(e))
print(msg)
+
class TestNdpointer(object):
def test_dtype(self):
dt = np.intc
@@ -113,3 +111,62 @@ class TestNdpointer(object):
a1 = ndpointer(dtype=np.float64)
a2 = ndpointer(dtype=np.float64)
assert_(a1 == a2)
+
+
+@pytest.mark.skipif(not _HAS_CTYPE,
+ reason="ctypes not available on this python installation")
+class TestAsArray(object):
+ def test_array(self):
+ from ctypes import c_int
+
+ pair_t = c_int * 2
+ a = as_array(pair_t(1, 2))
+ assert_equal(a.shape, (2,))
+ assert_array_equal(a, np.array([1, 2]))
+ a = as_array((pair_t * 3)(pair_t(1, 2), pair_t(3, 4), pair_t(5, 6)))
+ assert_equal(a.shape, (3, 2))
+ assert_array_equal(a, np.array([[1, 2], [3, 4], [5, 6]]))
+
+ def test_pointer(self):
+ from ctypes import c_int, cast, POINTER
+
+ p = cast((c_int * 10)(*range(10)), POINTER(c_int))
+
+ a = as_array(p, shape=(10,))
+ assert_equal(a.shape, (10,))
+ assert_array_equal(a, np.arange(10))
+
+ a = as_array(p, shape=(2, 5))
+ assert_equal(a.shape, (2, 5))
+ assert_array_equal(a, np.arange(10).reshape((2, 5)))
+
+ # shape argument is required
+ assert_raises(TypeError, as_array, p)
+
+ def test_struct_array_pointer(self):
+ from ctypes import c_int16, Structure, pointer
+
+ class Struct(Structure):
+ _fields_ = [('a', c_int16)]
+
+ Struct3 = 3 * Struct
+
+ c_array = (2 * Struct3)(
+ Struct3(Struct(a=1), Struct(a=2), Struct(a=3)),
+ Struct3(Struct(a=4), Struct(a=5), Struct(a=6))
+ )
+
+ expected = np.array([
+ [(1,), (2,), (3,)],
+ [(4,), (5,), (6,)],
+ ], dtype=[('a', np.int16)])
+
+ def check(x):
+ assert_equal(x.dtype, expected.dtype)
+ assert_equal(x, expected)
+
+ # all of these should be equivalent
+ check(as_array(c_array))
+ check(as_array(pointer(c_array), shape=()))
+ check(as_array(pointer(c_array[0]), shape=(2,)))
+ check(as_array(pointer(c_array[0][0]), shape=(2, 3)))