Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

BUG: Force `npymath to respect npy_longdouble` #20465

Merged
merged 1 commit into from Nov 26, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
29 changes: 23 additions & 6 deletions numpy/core/include/numpy/npy_common.h
Expand Up @@ -356,14 +356,31 @@ typedef unsigned long npy_ulonglong;
typedef unsigned char npy_bool;
#define NPY_FALSE 0
#define NPY_TRUE 1


/*
* `NPY_SIZEOF_LONGDOUBLE` isn't usually equal to sizeof(long double).
* In some certain cases, it may forced to be equal to sizeof(double)
* even against the compiler implementation and the same goes for
* `complex long double`.
*
* Therefore, avoid `long double`, use `npy_longdouble` instead,
* and when it comes to standard math functions make sure of using
* the double version when `NPY_SIZEOF_LONGDOUBLE` == `NPY_SIZEOF_DOUBLE`.
* For example:
* npy_longdouble *ptr, x;
* #if NPY_SIZEOF_LONGDOUBLE == NPY_SIZEOF_DOUBLE
* npy_longdouble r = modf(x, ptr);
* #else
* npy_longdouble r = modfl(x, ptr);
* #endif
*
* See https://github.com/numpy/numpy/issues/20348
*/
#if NPY_SIZEOF_LONGDOUBLE == NPY_SIZEOF_DOUBLE
typedef double npy_longdouble;
#define NPY_LONGDOUBLE_FMT "g"
#define NPY_LONGDOUBLE_FMT "g"
typedef double npy_longdouble;
#else
typedef long double npy_longdouble;
#define NPY_LONGDOUBLE_FMT "Lg"
#define NPY_LONGDOUBLE_FMT "Lg"
typedef long double npy_longdouble;
#endif

#ifndef Py_USING_UNICODE
Expand Down
2 changes: 1 addition & 1 deletion numpy/core/src/multiarray/_multiarray_tests.c.src
Expand Up @@ -2093,7 +2093,7 @@ PrintFloat_Printf_g(PyObject *obj, int precision)
}
else if (PyArray_IsScalar(obj, LongDouble)) {
npy_longdouble x = PyArrayScalar_VAL(obj, LongDouble);
PyOS_snprintf(str, sizeof(str), "%.*Lg", precision, x);
PyOS_snprintf(str, sizeof(str), "%.*" NPY_LONGDOUBLE_FMT, precision, x);
}
else{
double val = PyFloat_AsDouble(obj);
Expand Down
40 changes: 26 additions & 14 deletions numpy/core/src/npymath/npy_math_internal.h.src
Expand Up @@ -480,10 +480,16 @@ NPY_INPLACE @type@ npy_frexp@c@(@type@ x, int* exp)

/**begin repeat
* #type = npy_longdouble, npy_double, npy_float#
* #TYPE = LONGDOUBLE, DOUBLE, FLOAT#
* #c = l,,f#
* #C = L,,F#
*/

#undef NPY__FP_SFX
#if NPY_SIZEOF_@TYPE@ == NPY_SIZEOF_DOUBLE
#define NPY__FP_SFX(X) X
#else
#define NPY__FP_SFX(X) NPY_CAT(X, @c@)
#endif
/*
* On arm64 macOS, there's a bug with sin, cos, and tan where they don't
* raise "invalid" when given INFINITY as input.
Expand All @@ -509,7 +515,7 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x)
return (x - x);
}
#endif
return @kind@@c@(x);
return NPY__FP_SFX(@kind@)(x);
}
#endif

Expand All @@ -524,7 +530,7 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x)
#ifdef HAVE_@KIND@@C@
NPY_INPLACE @type@ npy_@kind@@c@(@type@ x, @type@ y)
{
return @kind@@c@(x, y);
return NPY__FP_SFX(@kind@)(x, y);
}
#endif
/**end repeat1**/
Expand Down Expand Up @@ -555,21 +561,21 @@ npy_@kind@@c@(@type@ x, @type@ y)
#ifdef HAVE_MODF@C@
NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr)
{
return modf@c@(x, iptr);
return NPY__FP_SFX(modf)(x, iptr);
}
#endif

#ifdef HAVE_LDEXP@C@
NPY_INPLACE @type@ npy_ldexp@c@(@type@ x, int exp)
{
return ldexp@c@(x, exp);
return NPY__FP_SFX(ldexp)(x, exp);
}
#endif

#ifdef HAVE_FREXP@C@
NPY_INPLACE @type@ npy_frexp@c@(@type@ x, int* exp)
{
return frexp@c@(x, exp);
return NPY__FP_SFX(frexp)(x, exp);
}
#endif

Expand All @@ -592,10 +598,10 @@ NPY_INPLACE @type@ npy_cbrt@c@(@type@ x)
#else
NPY_INPLACE @type@ npy_cbrt@c@(@type@ x)
{
return cbrt@c@(x);
return NPY__FP_SFX(cbrt)(x);
}
#endif

#undef NPY__FP_SFX
/**end repeat**/


Expand All @@ -605,10 +611,16 @@ NPY_INPLACE @type@ npy_cbrt@c@(@type@ x)

/**begin repeat
* #type = npy_float, npy_double, npy_longdouble#
* #TYPE = FLOAT, DOUBLE, LONGDOUBLE#
* #c = f, ,l#
* #C = F, ,L#
*/

#undef NPY__FP_SFX
#if NPY_SIZEOF_@TYPE@ == NPY_SIZEOF_DOUBLE
#define NPY__FP_SFX(X) X
#else
#define NPY__FP_SFX(X) NPY_CAT(X, @c@)
#endif
@type@ npy_heaviside@c@(@type@ x, @type@ h0)
{
if (npy_isnan(x)) {
Expand All @@ -625,10 +637,10 @@ NPY_INPLACE @type@ npy_cbrt@c@(@type@ x)
}
}

#define LOGE2 NPY_LOGE2@c@
#define LOG2E NPY_LOG2E@c@
#define RAD2DEG (180.0@c@/NPY_PI@c@)
#define DEG2RAD (NPY_PI@c@/180.0@c@)
#define LOGE2 NPY__FP_SFX(NPY_LOGE2)
#define LOG2E NPY__FP_SFX(NPY_LOG2E)
#define RAD2DEG (NPY__FP_SFX(180.0)/NPY__FP_SFX(NPY_PI))
#define DEG2RAD (NPY__FP_SFX(NPY_PI)/NPY__FP_SFX(180.0))

NPY_INPLACE @type@ npy_rad2deg@c@(@type@ x)
{
Expand Down Expand Up @@ -783,7 +795,7 @@ npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus)
#undef LOG2E
#undef RAD2DEG
#undef DEG2RAD

#undef NPY__FP_SFX
/**end repeat**/

/**begin repeat
Expand Down