From 0c3a5617ff4ac4a2a3fb32cbbcea3c8f30e1e05b Mon Sep 17 00:00:00 2001 From: Lysandros Nikolaou Date: Wed, 7 Feb 2024 12:18:39 +0100 Subject: [PATCH] Rename include folder & rename symbols to npymath_* - Rename include folder so that `#include ` has to be used. - Rename npymath_* symbols, so that there are no conflicts with old symbols in NumPy --- meson.options | 2 +- src/arm64_exports.c | 8 +- src/halffloat.c | 160 +++-- src/ieee754.c.src | 42 +- .../npymath => libnpymath}/_config.h.in | 0 .../{numpy/npymath => libnpymath}/block.h | 12 +- .../{numpy/npymath => libnpymath}/common.h | 51 +- .../{numpy/npymath => libnpymath}/config.h | 8 +- .../{numpy/npymath => libnpymath}/cpu.h | 8 +- .../{numpy/npymath => libnpymath}/endian.h | 8 +- .../{numpy/npymath => libnpymath}/fpmath.h | 12 +- src/include/libnpymath/halffloat.h | 96 +++ .../npymath => libnpymath}/internal.h.src | 240 +++---- .../{numpy/npymath => libnpymath}/meson.build | 2 +- .../npymath => libnpymath}/meson_config.h.in | 4 +- src/include/libnpymath/npy_math.h | 502 +++++++++++++++ .../libnpymath}/npy_math_common.h | 4 +- .../{numpy/npymath => libnpymath}/private.h | 14 +- .../{numpy/npymath => libnpymath}/utils.h | 6 +- src/include/numpy/halffloat.h | 70 --- src/include/numpy/npy_math.h | 517 --------------- src/meson.build | 33 +- src/npy_math.c | 2 +- src/npy_math_complex.c.src | 594 +++++++++--------- 24 files changed, 1197 insertions(+), 1198 deletions(-) rename src/include/{numpy/npymath => libnpymath}/_config.h.in (100%) rename src/include/{numpy/npymath => libnpymath}/block.h (95%) rename src/include/{numpy/npymath => libnpymath}/common.h (89%) rename src/include/{numpy/npymath => libnpymath}/config.h (93%) rename src/include/{numpy/npymath => libnpymath}/cpu.h (97%) rename src/include/{numpy/npymath => libnpymath}/endian.h (95%) rename src/include/{numpy/npymath => libnpymath}/fpmath.h (82%) create mode 100644 src/include/libnpymath/halffloat.h rename src/include/{numpy/npymath => libnpymath}/internal.h.src (71%) rename src/include/{numpy/npymath => libnpymath}/meson.build (90%) rename src/include/{numpy/npymath => libnpymath}/meson_config.h.in (95%) create mode 100644 src/include/libnpymath/npy_math.h rename src/{ => include/libnpymath}/npy_math_common.h (66%) rename src/include/{numpy/npymath => libnpymath}/private.h (98%) rename src/include/{numpy/npymath => libnpymath}/utils.h (89%) delete mode 100644 src/include/numpy/halffloat.h delete mode 100644 src/include/numpy/npy_math.h diff --git a/meson.options b/meson.options index 5fb8b46..d04ea56 100644 --- a/meson.options +++ b/meson.options @@ -2,7 +2,7 @@ option('numpy_build', type: 'boolean', value: false, description: 'Whether libnpymath is getting built as part of numpy') option('install_lib', type: 'boolean', value: true, description: 'Whether to install the libnpymath library') -option('headers_dir', type: 'string', value: 'libnpymath', +option('headers_dir', type: 'string', value: '', description: 'Installation directory for the libnpymath headers') option('lib_dir', type: 'string', value: 'libnpymath', description: 'Installation directory for the libnpymath library') diff --git a/src/arm64_exports.c b/src/arm64_exports.c index dfa8054..1e24326 100644 --- a/src/arm64_exports.c +++ b/src/arm64_exports.c @@ -11,19 +11,19 @@ * This file is actually compiled as part of the main module. */ -double npy_asinh(double x) { +double npymath_asinh(double x) { return asinh(x); } -double npy_copysign(double y, double x) { +double npymath_copysign(double y, double x) { return copysign(y, x); } -double npy_log1p(double x) { +double npymath_log1p(double x) { return log1p(x); } -double npy_nextafter(double x, double y) { +double npymath_nextafter(double x, double y) { return nextafter(x, y); } diff --git a/src/halffloat.c b/src/halffloat.c index c71ed4f..97bac7c 100644 --- a/src/halffloat.c +++ b/src/halffloat.c @@ -1,18 +1,16 @@ -#define NPY_NO_DEPRECATED_API NPY_API_VERSION - -#include "numpy/halffloat.h" +#include "libnpymath/halffloat.h" /* * This chooses between 'ties to even' and 'ties away from zero'. */ -#define NPY_HALF_ROUND_TIES_TO_EVEN 1 +#define NPYMATH_HALF_ROUND_TIES_TO_EVEN 1 /* * If these are 1, the conversions try to trigger underflow, * overflow, and invalid exceptions in the FP system when needed. */ -#define NPY_HALF_GENERATE_OVERFLOW 1 -#define NPY_HALF_GENERATE_UNDERFLOW 1 -#define NPY_HALF_GENERATE_INVALID 1 +#define NPYMATH_HALF_GENERATE_OVERFLOW 1 +#define NPYMATH_HALF_GENERATE_UNDERFLOW 1 +#define NPYMATH_HALF_GENERATE_INVALID 1 /* ******************************************************************** @@ -20,7 +18,7 @@ ******************************************************************** */ -float npy_half_to_float(npymath_half h) +float npymath_half_to_float(npymath_half h) { #if defined(NPYMATH_HAVE_FP16) float ret; @@ -36,12 +34,12 @@ float npy_half_to_float(npymath_half h) return vec_extract(vf32, 0); #else union { float ret; npymath_uint32 retbits; } conv; - conv.retbits = npy_halfbits_to_floatbits(h); + conv.retbits = npymath_halfbits_to_floatbits(h); return conv.ret; #endif } -double npy_half_to_double(npymath_half h) +double npymath_half_to_double(npymath_half h) { #if defined(NPYMATH_HAVE_AVX512FP16) double ret; @@ -55,12 +53,12 @@ double npy_half_to_double(npymath_half h) return vec_extract(vf64, 0); #else union { double ret; npymath_uint64 retbits; } conv; - conv.retbits = npy_halfbits_to_doublebits(h); + conv.retbits = npymath_halfbits_to_doublebits(h); return conv.ret; #endif } -npymath_half npy_float_to_half(float f) +npymath_half npymath_float_to_half(float f) { #if defined(NPYMATH_HAVE_FP16) __m128 mf = _mm_load_ss(&f); @@ -73,12 +71,12 @@ npymath_half npy_float_to_half(float f) #else union { float f; npymath_uint32 fbits; } conv; conv.f = f; - return npy_floatbits_to_halfbits(conv.fbits); + return npymath_floatbits_to_halfbits(conv.fbits); #endif } -npymath_half npy_double_to_half(double d) +npymath_half npymath_double_to_half(double d) { #if defined(NPYMATH_HAVE_AVX512FP16) __m128d md = _mm_load_sd(&f); @@ -91,50 +89,50 @@ npymath_half npy_double_to_half(double d) #else union { double d; npymath_uint64 dbits; } conv; conv.d = d; - return npy_doublebits_to_halfbits(conv.dbits); + return npymath_doublebits_to_halfbits(conv.dbits); #endif } -int npy_half_iszero(npymath_half h) +int npymath_half_iszero(npymath_half h) { return (h&0x7fff) == 0; } -int npy_half_isnan(npymath_half h) +int npymath_half_isnan(npymath_half h) { return ((h&0x7c00u) == 0x7c00u) && ((h&0x03ffu) != 0x0000u); } -int npy_half_isinf(npymath_half h) +int npymath_half_isinf(npymath_half h) { return ((h&0x7fffu) == 0x7c00u); } -int npy_half_isfinite(npymath_half h) +int npymath_half_isfinite(npymath_half h) { return ((h&0x7c00u) != 0x7c00u); } -int npy_half_signbit(npymath_half h) +int npymath_half_signbit(npymath_half h) { return (h&0x8000u) != 0; } -npymath_half npy_half_spacing(npymath_half h) +npymath_half npymath_half_spacing(npymath_half h) { npymath_half ret; npymath_uint16 h_exp = h&0x7c00u; npymath_uint16 h_sig = h&0x03ffu; if (h_exp == 0x7c00u) { -#if NPY_HALF_GENERATE_INVALID - npy_set_floatstatus_invalid(); +#if NPYMATH_HALF_GENERATE_INVALID + npymath_set_floatstatus_invalid(); #endif - ret = NPY_HALF_NAN; + ret = NPYMATH_HALF_NAN; } else if (h == 0x7bffu) { -#if NPY_HALF_GENERATE_OVERFLOW - npy_set_floatstatus_overflow(); +#if NPYMATH_HALF_GENERATE_OVERFLOW + npymath_set_floatstatus_overflow(); #endif - ret = NPY_HALF_PINF; + ret = NPYMATH_HALF_PINF; } else if ((h&0x8000u) && h_sig == 0) { /* Negative boundary case */ if (h_exp > 0x2c00u) { /* If result is normalized */ ret = h_exp - 0x2c00u; @@ -154,20 +152,20 @@ npymath_half npy_half_spacing(npymath_half h) return ret; } -npymath_half npy_half_copysign(npymath_half x, npymath_half y) +npymath_half npymath_half_copysign(npymath_half x, npymath_half y) { return (x&0x7fffu) | (y&0x8000u); } -npymath_half npy_half_nextafter(npymath_half x, npymath_half y) +npymath_half npymath_half_nextafter(npymath_half x, npymath_half y) { npymath_half ret; - if (npy_half_isnan(x) || npy_half_isnan(y)) { - ret = NPY_HALF_NAN; - } else if (npy_half_eq_nonan(x, y)) { + if (npymath_half_isnan(x) || npymath_half_isnan(y)) { + ret = NPYMATH_HALF_NAN; + } else if (npymath_half_eq_nonan(x, y)) { ret = x; - } else if (npy_half_iszero(x)) { + } else if (npymath_half_iszero(x)) { ret = (y&0x8000u) + 1; /* Smallest subnormal half */ } else if (!(x&0x8000u)) { /* x > 0 */ if ((npymath_int16)x > (npymath_int16)y) { /* x > y */ @@ -182,21 +180,21 @@ npymath_half npy_half_nextafter(npymath_half x, npymath_half y) ret = x+1; } } -#if NPY_HALF_GENERATE_OVERFLOW - if (npy_half_isinf(ret) && npy_half_isfinite(x)) { - npy_set_floatstatus_overflow(); +#if NPYMATH_HALF_GENERATE_OVERFLOW + if (npymath_half_isinf(ret) && npymath_half_isfinite(x)) { + npymath_set_floatstatus_overflow(); } #endif return ret; } -int npy_half_eq_nonan(npymath_half h1, npymath_half h2) +int npymath_half_eq_nonan(npymath_half h1, npymath_half h2) { return (h1 == h2 || ((h1 | h2) & 0x7fff) == 0); } -int npy_half_eq(npymath_half h1, npymath_half h2) +int npymath_half_eq(npymath_half h1, npymath_half h2) { /* * The equality cases are as follows: @@ -204,16 +202,16 @@ int npy_half_eq(npymath_half h1, npymath_half h2) * - If the values are equal, equal. * - If the values are both signed zeros, equal. */ - return (!npy_half_isnan(h1) && !npy_half_isnan(h2)) && + return (!npymath_half_isnan(h1) && !npymath_half_isnan(h2)) && (h1 == h2 || ((h1 | h2) & 0x7fff) == 0); } -int npy_half_ne(npymath_half h1, npymath_half h2) +int npymath_half_ne(npymath_half h1, npymath_half h2) { - return !npy_half_eq(h1, h2); + return !npymath_half_eq(h1, h2); } -int npy_half_lt_nonan(npymath_half h1, npymath_half h2) +int npymath_half_lt_nonan(npymath_half h1, npymath_half h2) { int sign_h1 = (h1 & 0x8000u) == 0x8000u; int sign_h2 = (h2 & 0x8000u) == 0x8000u; @@ -229,17 +227,17 @@ int npy_half_lt_nonan(npymath_half h1, npymath_half h2) : sign_h1 && ((h1 | h2) != 0x8000u); } -int npy_half_lt(npymath_half h1, npymath_half h2) +int npymath_half_lt(npymath_half h1, npymath_half h2) { - return (!npy_half_isnan(h1) && !npy_half_isnan(h2)) && npy_half_lt_nonan(h1, h2); + return (!npymath_half_isnan(h1) && !npymath_half_isnan(h2)) && npymath_half_lt_nonan(h1, h2); } -int npy_half_gt(npymath_half h1, npymath_half h2) +int npymath_half_gt(npymath_half h1, npymath_half h2) { - return npy_half_lt(h2, h1); + return npymath_half_lt(h2, h1); } -int npy_half_le_nonan(npymath_half h1, npymath_half h2) +int npymath_half_le_nonan(npymath_half h1, npymath_half h2) { int sign_h1 = (h1 & 0x8000u) == 0x8000u; int sign_h2 = (h2 & 0x8000u) == 0x8000u; @@ -255,25 +253,25 @@ int npy_half_le_nonan(npymath_half h1, npymath_half h2) : sign_h1 || ((h1 | h2) == 0x8000u); } -int npy_half_le(npymath_half h1, npymath_half h2) +int npymath_half_le(npymath_half h1, npymath_half h2) { - return (!npy_half_isnan(h1) && !npy_half_isnan(h2)) && npy_half_le_nonan(h1, h2); + return (!npymath_half_isnan(h1) && !npymath_half_isnan(h2)) && npymath_half_le_nonan(h1, h2); } -int npy_half_ge(npymath_half h1, npymath_half h2) +int npymath_half_ge(npymath_half h1, npymath_half h2) { - return npy_half_le(h2, h1); + return npymath_half_le(h2, h1); } -npymath_half npy_half_divmod(npymath_half h1, npymath_half h2, npymath_half *modulus) +npymath_half npymath_half_divmod(npymath_half h1, npymath_half h2, npymath_half *modulus) { - float fh1 = npy_half_to_float(h1); - float fh2 = npy_half_to_float(h2); + float fh1 = npymath_half_to_float(h1); + float fh2 = npymath_half_to_float(h2); float div, mod; - div = npy_divmodf(fh1, fh2, &mod); - *modulus = npy_float_to_half(mod); - return npy_float_to_half(div); + div = npymath_divmodf(fh1, fh2, &mod); + *modulus = npymath_float_to_half(mod); + return npymath_float_to_half(div); } @@ -284,7 +282,7 @@ npymath_half npy_half_divmod(npymath_half h1, npymath_half h2, npymath_half *mod ******************************************************************** */ -npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) +npymath_uint16 npymath_floatbits_to_halfbits(npymath_uint32 f) { npymath_uint32 f_exp, f_sig; npymath_uint16 h_sgn, h_exp, h_sig; @@ -311,8 +309,8 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) } } else { /* overflow to signed inf */ -#if NPY_HALF_GENERATE_OVERFLOW - npy_set_floatstatus_overflow(); +#if NPYMATH_HALF_GENERATE_OVERFLOW + npymath_set_floatstatus_overflow(); #endif return (npymath_uint16) (h_sgn + 0x7c00u); } @@ -325,10 +323,10 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) * exponents all convert to signed zero half-floats. */ if (f_exp < 0x33000000u) { -#if NPY_HALF_GENERATE_UNDERFLOW +#if NPYMATH_HALF_GENERATE_UNDERFLOW /* If f != 0, it underflowed to 0 */ if ((f&0x7fffffff) != 0) { - npy_set_floatstatus_underflow(); + npymath_set_floatstatus_underflow(); } #endif return h_sgn; @@ -336,10 +334,10 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) /* Make the subnormal significand */ f_exp >>= 23; f_sig = (0x00800000u + (f&0x007fffffu)); -#if NPY_HALF_GENERATE_UNDERFLOW +#if NPYMATH_HALF_GENERATE_UNDERFLOW /* If it's not exactly represented, it underflowed */ if ((f_sig&(((npymath_uint32)1 << (126 - f_exp)) - 1)) != 0) { - npy_set_floatstatus_underflow(); + npymath_set_floatstatus_underflow(); } #endif /* @@ -350,7 +348,7 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) */ f_sig >>= (113 - f_exp); /* Handle rounding by adding 1 to the bit beyond half precision */ -#if NPY_HALF_ROUND_TIES_TO_EVEN +#if NPYMATH_HALF_ROUND_TIES_TO_EVEN /* * If the last bit in the half significand is 0 (already even), and * the remaining bit pattern is 1000...0, then we do not add one @@ -377,7 +375,7 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) h_exp = (npymath_uint16) ((f_exp - 0x38000000u) >> 13); /* Handle rounding by adding 1 to the bit beyond half precision */ f_sig = (f&0x007fffffu); -#if NPY_HALF_ROUND_TIES_TO_EVEN +#if NPYMATH_HALF_ROUND_TIES_TO_EVEN /* * If the last bit in the half significand is 0 (already even), and * the remaining bit pattern is 1000...0, then we do not add one @@ -396,10 +394,10 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) * correct result. h_exp may increment to 15, at greatest, in * which case the result overflows to a signed inf. */ -#if NPY_HALF_GENERATE_OVERFLOW +#if NPYMATH_HALF_GENERATE_OVERFLOW h_sig += h_exp; if (h_sig == 0x7c00u) { - npy_set_floatstatus_overflow(); + npymath_set_floatstatus_overflow(); } return h_sgn + h_sig; #else @@ -407,7 +405,7 @@ npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f) #endif } -npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) +npymath_uint16 npymath_doublebits_to_halfbits(npymath_uint64 d) { npymath_uint64 d_exp, d_sig; npymath_uint16 h_sgn, h_exp, h_sig; @@ -434,8 +432,8 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) } } else { /* overflow to signed inf */ -#if NPY_HALF_GENERATE_OVERFLOW - npy_set_floatstatus_overflow(); +#if NPYMATH_HALF_GENERATE_OVERFLOW + npymath_set_floatstatus_overflow(); #endif return h_sgn + 0x7c00u; } @@ -448,10 +446,10 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) * exponents all convert to signed zero half-floats. */ if (d_exp < 0x3e60000000000000ULL) { -#if NPY_HALF_GENERATE_UNDERFLOW +#if NPYMATH_HALF_GENERATE_UNDERFLOW /* If d != 0, it underflowed to 0 */ if ((d&0x7fffffffffffffffULL) != 0) { - npy_set_floatstatus_underflow(); + npymath_set_floatstatus_underflow(); } #endif return h_sgn; @@ -459,10 +457,10 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) /* Make the subnormal significand */ d_exp >>= 52; d_sig = (0x0010000000000000ULL + (d&0x000fffffffffffffULL)); -#if NPY_HALF_GENERATE_UNDERFLOW +#if NPYMATH_HALF_GENERATE_UNDERFLOW /* If it's not exactly represented, it underflowed */ if ((d_sig&(((npymath_uint64)1 << (1051 - d_exp)) - 1)) != 0) { - npy_set_floatstatus_underflow(); + npymath_set_floatstatus_underflow(); } #endif /* @@ -476,7 +474,7 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) assert(d_exp - 998 >= 0); d_sig <<= (d_exp - 998); /* Handle rounding by adding 1 to the bit beyond half precision */ -#if NPY_HALF_ROUND_TIES_TO_EVEN +#if NPYMATH_HALF_ROUND_TIES_TO_EVEN /* * If the last bit in the half significand is 0 (already even), and * the remaining bit pattern is 1000...0, then we do not add one @@ -501,7 +499,7 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) h_exp = (npymath_uint16) ((d_exp - 0x3f00000000000000ULL) >> 42); /* Handle rounding by adding 1 to the bit beyond half precision */ d_sig = (d&0x000fffffffffffffULL); -#if NPY_HALF_ROUND_TIES_TO_EVEN +#if NPYMATH_HALF_ROUND_TIES_TO_EVEN /* * If the last bit in the half significand is 0 (already even), and * the remaining bit pattern is 1000...0, then we do not add one @@ -521,10 +519,10 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) * correct result. h_exp may increment to 15, at greatest, in * which case the result overflows to a signed inf. */ -#if NPY_HALF_GENERATE_OVERFLOW +#if NPYMATH_HALF_GENERATE_OVERFLOW h_sig += h_exp; if (h_sig == 0x7c00u) { - npy_set_floatstatus_overflow(); + npymath_set_floatstatus_overflow(); } return h_sgn + h_sig; #else @@ -532,7 +530,7 @@ npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d) #endif } -npymath_uint32 npy_halfbits_to_floatbits(npymath_uint16 h) +npymath_uint32 npymath_halfbits_to_floatbits(npymath_uint16 h) { npymath_uint16 h_exp, h_sig; npymath_uint32 f_sgn, f_exp, f_sig; @@ -564,7 +562,7 @@ npymath_uint32 npy_halfbits_to_floatbits(npymath_uint16 h) } } -npymath_uint64 npy_halfbits_to_doublebits(npymath_uint16 h) +npymath_uint64 npymath_halfbits_to_doublebits(npymath_uint16 h) { npymath_uint16 h_exp, h_sig; npymath_uint64 d_sgn, d_exp, d_sig; diff --git a/src/ieee754.c.src b/src/ieee754.c.src index ad6ede5..fc6627d 100644 --- a/src/ieee754.c.src +++ b/src/ieee754.c.src @@ -4,9 +4,9 @@ * * Low-level routines related to IEEE-754 format */ -#include "npy_math_common.h" -#include "numpy/npymath/private.h" -#include "numpy/npymath/utils.h" +#include "libnpymath/npy_math_common.h" +#include "libnpymath/private.h" +#include "libnpymath/utils.h" /* The below code is provided for compilers which do not yet provide C11 compatibility (gcc 4.5 and older) @@ -301,24 +301,24 @@ static npymath_longdouble _nextl(npymath_longdouble x, int p) * #SUFF = F,,L# * #type = npymath_float, npymath_double, npymath_longdouble# */ -@type@ npy_spacing@suff@(@type@ x) +@type@ npymath_spacing@suff@(@type@ x) { /* XXX: npy isnan/isinf may be optimized by bit twiddling */ - if (npy_isinf(x)) { - return NPY_NAN@SUFF@; + if (npymath_isinf(x)) { + return NPYMATH_NAN@SUFF@; } return _next@suff@(x, 1) - x; } /**end repeat**/ -int npy_clear_floatstatus() { +int npymath_clear_floatstatus() { char x=0; - return npy_clear_floatstatus_barrier(&x); + return npymath_clear_floatstatus_barrier(&x); } -int npy_get_floatstatus() { +int npymath_get_floatstatus() { char x=0; - return npy_get_floatstatus_barrier(&x); + return npymath_get_floatstatus_barrier(&x); } @@ -354,7 +354,7 @@ int npy_get_floatstatus() { #endif -int npy_get_floatstatus_barrier(char* param) +int npymath_get_floatstatus_barrier(char* param) { int fpstatus = fetestexcept(FE_DIVBYZERO | FE_OVERFLOW | FE_UNDERFLOW | FE_INVALID); @@ -365,16 +365,16 @@ int npy_get_floatstatus_barrier(char* param) volatile char NPYMATH_UNUSED(c) = *(char*)param; } - return ((FE_DIVBYZERO & fpstatus) ? NPY_FPE_DIVIDEBYZERO : 0) | - ((FE_OVERFLOW & fpstatus) ? NPY_FPE_OVERFLOW : 0) | - ((FE_UNDERFLOW & fpstatus) ? NPY_FPE_UNDERFLOW : 0) | - ((FE_INVALID & fpstatus) ? NPY_FPE_INVALID : 0); + return ((FE_DIVBYZERO & fpstatus) ? NPYMATH_FPE_DIVIDEBYZERO : 0) | + ((FE_OVERFLOW & fpstatus) ? NPYMATH_FPE_OVERFLOW : 0) | + ((FE_UNDERFLOW & fpstatus) ? NPYMATH_FPE_UNDERFLOW : 0) | + ((FE_INVALID & fpstatus) ? NPYMATH_FPE_INVALID : 0); } -int npy_clear_floatstatus_barrier(char * param) +int npymath_clear_floatstatus_barrier(char * param) { /* testing float status is 50-100 times faster than clearing on x86 */ - int fpstatus = npy_get_floatstatus_barrier(param); + int fpstatus = npymath_get_floatstatus_barrier(param); if (fpstatus != 0) { feclearexcept(FE_DIVBYZERO | FE_OVERFLOW | FE_UNDERFLOW | FE_INVALID); @@ -384,22 +384,22 @@ int npy_clear_floatstatus_barrier(char * param) } -void npy_set_floatstatus_divbyzero(void) +void npymath_set_floatstatus_divbyzero(void) { feraiseexcept(FE_DIVBYZERO); } -void npy_set_floatstatus_overflow(void) +void npymath_set_floatstatus_overflow(void) { feraiseexcept(FE_OVERFLOW); } -void npy_set_floatstatus_underflow(void) +void npymath_set_floatstatus_underflow(void) { feraiseexcept(FE_UNDERFLOW); } -void npy_set_floatstatus_invalid(void) +void npymath_set_floatstatus_invalid(void) { feraiseexcept(FE_INVALID); } diff --git a/src/include/numpy/npymath/_config.h.in b/src/include/libnpymath/_config.h.in similarity index 100% rename from src/include/numpy/npymath/_config.h.in rename to src/include/libnpymath/_config.h.in diff --git a/src/include/numpy/npymath/block.h b/src/include/libnpymath/block.h similarity index 95% rename from src/include/numpy/npymath/block.h rename to src/include/libnpymath/block.h index 8a70fc8..2f5f96b 100644 --- a/src/include/numpy/npymath/block.h +++ b/src/include/libnpymath/block.h @@ -1,9 +1,9 @@ -#ifndef NPYMATH_CONFIG_H_ -#define NPYMATH_CONFIG_H_ +#ifndef LIBNPYMATH_BLOCK_H_ +#define LIBNPYMATH_BLOCK_H_ -#include "numpy/npymath/meson_config.h" -#include "numpy/npymath/config.h" -#include "numpy/npymath/utils.h" +#include "libnpymath/meson_config.h" +#include "libnpymath/config.h" +#include "libnpymath/utils.h" /* blocklist */ @@ -183,4 +183,4 @@ #endif /* defined(__GLIBC) */ #endif /* defined(NPYMATH_HAVE_FEATURES_H) */ -#endif /* NPYMATH_CONFIG_H_ */ +#endif /* LIBNPYMATH_BLOCK_H_ */ diff --git a/src/include/numpy/npymath/common.h b/src/include/libnpymath/common.h similarity index 89% rename from src/include/numpy/npymath/common.h rename to src/include/libnpymath/common.h index 124bd24..9a94760 100644 --- a/src/include/numpy/npymath/common.h +++ b/src/include/libnpymath/common.h @@ -1,13 +1,13 @@ -#ifndef NPYMATH_COMMON_H_ -#define NPYMATH_COMMON_H_ +#ifndef LIBNPYMATH_COMMON_H_ +#define LIBNPYMATH_COMMON_H_ /* need Python.h for Py_intptr_t */ #include /* numpconfig.h is auto-generated */ -#include "numpy/npymath/config.h" +#include "libnpymath/config.h" #if defined(NPYMATH_HAVE_CONFIG_H) -#include "numpy/npymath/block.h" +#include "libnpymath/block.h" #endif /* @@ -166,42 +166,42 @@ typedef longdouble_t _Complex npymath_clongdouble; #define NPYMATH_BITSOF_LONGDOUBLE (NPYMATH_SIZEOF_LONGDOUBLE * CHAR_BIT) #if NPYMATH_BITSOF_LONG == 8 -#define NPYMATH_INT8 NPY_LONG +#define NPYMATH_INT8 typedef unsigned long npymath_uint8; #elif NPYMATH_BITSOF_LONG == 16 -#define NPYMATH_INT16 NPY_LONG +#define NPYMATH_INT16 typedef long npymath_int16; typedef unsigned long npymath_uint16; #elif NPYMATH_BITSOF_LONG == 32 -#define NPYMATH_INT32 NPY_LONG +#define NPYMATH_INT32 typedef long npymath_int32; typedef unsigned long npymath_uint32; #elif NPYMATH_BITSOF_LONG == 64 -#define NPYMATH_INT64 NPY_LONG +#define NPYMATH_INT64 typedef long npymath_int64; typedef unsigned long npymath_uint64; #endif #if NPYMATH_BITSOF_LONGLONG == 8 # ifndef NPYMATH_INT8 -# define NPYMATH_INT8 NPY_LONGLONG +# define NPYMATH_INT8 typedef npymath_ulonglong npymath_uint8; # endif #elif NPYMATH_BITSOF_LONGLONG == 16 # ifndef NPYMATH_INT16 -# define NPYMATH_INT16 NPY_LONGLONG +# define NPYMATH_INT16 typedef npymath_longlong npymath_int16; typedef npymath_ulonglong npymath_uint16; # endif #elif NPYMATH_BITSOF_LONGLONG == 32 # ifndef NPYMATH_INT32 -# define NPYMATH_INT32 NPY_LONGLONG +# define NPYMATH_INT32 typedef npymath_longlong npymath_int32; typedef npymath_ulonglong npymath_uint32; # endif #elif NPYMATH_BITSOF_LONGLONG == 64 # ifndef NPYMATH_INT64 -# define NPYMATH_INT64 NPY_LONGLONG +# define NPYMATH_INT64 typedef npymath_longlong npymath_int64; typedef npymath_ulonglong npymath_uint64; # endif @@ -209,24 +209,24 @@ typedef longdouble_t _Complex npymath_clongdouble; #if NPYMATH_BITSOF_INT == 8 #ifndef NPYMATH_INT8 -#define NPYMATH_INT8 NPY_INT +#define NPYMATH_INT8 typedef unsigned int npymath_uint8; #endif #elif NPYMATH_BITSOF_INT == 16 #ifndef NPYMATH_INT16 -#define NPYMATH_INT16 NPY_INT +#define NPYMATH_INT16 typedef int npymath_int16; typedef unsigned int npymath_uint16; #endif #elif NPYMATH_BITSOF_INT == 32 #ifndef NPYMATH_INT32 -#define NPYMATH_INT32 NPY_INT +#define NPYMATH_INT32 typedef int npymath_int32; typedef unsigned int npymath_uint32; #endif #elif NPYMATH_BITSOF_INT == 64 #ifndef NPYMATH_INT64 -#define NPYMATH_INT64 NPY_INT +#define NPYMATH_INT64 typedef int npymath_int64; typedef unsigned int npymath_uint64; #endif @@ -234,24 +234,24 @@ typedef longdouble_t _Complex npymath_clongdouble; #if NPYMATH_BITSOF_SHORT == 8 #ifndef NPYMATH_INT8 -#define NPYMATH_INT8 NPY_SHORT +#define NPYMATH_INT8 typedef unsigned short npymath_uint8; #endif #elif NPYMATH_BITSOF_SHORT == 16 #ifndef NPYMATH_INT16 -#define NPYMATH_INT16 NPY_SHORT +#define NPYMATH_INT16 typedef short npymath_int16; typedef unsigned short npymath_uint16; #endif #elif NPYMATH_BITSOF_SHORT == 32 #ifndef NPYMATH_INT32 -#define NPYMATH_INT32 NPY_SHORT +#define NPYMATH_INT32 typedef short npymath_int32; typedef unsigned short npymath_uint32; #endif #elif NPYMATH_BITSOF_SHORT == 64 #ifndef NPYMATH_INT64 -#define NPYMATH_INT64 NPY_SHORT +#define NPYMATH_INT64 typedef short npymath_int64; typedef unsigned short npymath_uint64; #endif @@ -260,33 +260,32 @@ typedef longdouble_t _Complex npymath_clongdouble; #if NPYMATH_BITSOF_CHAR == 8 #ifndef NPYMATH_INT8 -#define NPYMATH_INT8 NPY_BYTE +#define NPYMATH_INT8 typedef unsigned char npymath_uint8; #endif #elif NPYMATH_BITSOF_CHAR == 16 #ifndef NPYMATH_INT16 -#define NPYMATH_INT16 NPY_BYTE +#define NPYMATH_INT16 typedef signed char npymath_int16; typedef unsigned char npymath_uint16; #endif #elif NPYMATH_BITSOF_CHAR == 32 #ifndef NPYMATH_INT32 -#define NPYMATH_INT32 NPY_BYTE +#define NPYMATH_INT32 typedef signed char npymath_int32; typedef unsigned char npymath_uint32; #endif #elif NPYMATH_BITSOF_CHAR == 64 #ifndef NPYMATH_INT64 -#define NPYMATH_INT64 NPY_BYTE +#define NPYMATH_INT64 typedef signed char npymath_int64; typedef unsigned char npymath_uint64; #endif #endif /* half/float16 isn't a floating-point type in C */ -#define NPY_FLOAT16 NPY_HALF typedef npymath_uint16 npymath_half; /* End of typedefs for numarray style bit-width names */ -#endif /* NPYMATH_COMMON_H_ */ +#endif /* LIBNPYMATH_COMMON_H_ */ diff --git a/src/include/numpy/npymath/config.h b/src/include/libnpymath/config.h similarity index 93% rename from src/include/numpy/npymath/config.h rename to src/include/libnpymath/config.h index 9e2c16b..f09247e 100644 --- a/src/include/numpy/npymath/config.h +++ b/src/include/libnpymath/config.h @@ -1,7 +1,7 @@ -#ifndef NPYMATHCONFIG_H_ -#define NPYMATHCONFIG_H_ +#ifndef LIBNPYMATH_CONFIG_H_ +#define LIBNPYMATH_CONFIG_H_ -#include "numpy/npymath/_config.h" +#include "libnpymath/_config.h" /* * On Mac OS X, because there is only one configuration stage for all the archs @@ -51,4 +51,4 @@ #endif #endif -#endif /* NPYMATHCONFIG_H_ */ +#endif /* LIBNPYMATH_CONFIG_H_ */ diff --git a/src/include/numpy/npymath/cpu.h b/src/include/libnpymath/cpu.h similarity index 97% rename from src/include/numpy/npymath/cpu.h rename to src/include/libnpymath/cpu.h index aa6999d..5c6fec4 100644 --- a/src/include/numpy/npymath/cpu.h +++ b/src/include/libnpymath/cpu.h @@ -21,10 +21,10 @@ * NPYMATH_CPU_LOONGARCH * NPYMATH_CPU_WASM */ -#ifndef NPYMATH_CPU_H_ -#define NPYMATH_CPU_H_ +#ifndef LIBNPYMATH_CPU_H_ +#define LIBNPYMATH_CPU_H_ -#include "numpy/npymath/config.h" +#include "libnpymath/config.h" #if defined( __i386__ ) || defined(i386) || defined(_M_IX86) /* @@ -114,4 +114,4 @@ information about your platform (OS, CPU and compiler) #endif -#endif /* NPYMATH_CPU_H_ */ +#endif /* LIBNPYMATH_CPU_H_ */ diff --git a/src/include/numpy/npymath/endian.h b/src/include/libnpymath/endian.h similarity index 95% rename from src/include/numpy/npymath/endian.h rename to src/include/libnpymath/endian.h index 64c7771..912341f 100644 --- a/src/include/numpy/npymath/endian.h +++ b/src/include/libnpymath/endian.h @@ -1,5 +1,5 @@ -#ifndef NPYMATH_ENDIAN_H_ -#define NPYMATH_ENDIAN_H_ +#ifndef LIBNPYMATH_ENDIAN_H_ +#define LIBNPYMATH_ENDIAN_H_ /* * NPYMATH_BYTE_ORDER is set to the same value as BYTE_ORDER set by glibc in @@ -32,7 +32,7 @@ #ifndef NPYMATH_BYTE_ORDER /* Set endianness info using target CPU */ - #include "numpy/npymath/cpu.h" + #include "libnpymath/cpu.h" #define NPYMATH_LITTLE_ENDIAN 1234 #define NPYMATH_BIG_ENDIAN 4321 @@ -74,4 +74,4 @@ #endif -#endif /* NPYMATH_ENDIAN_H_ */ +#endif /* LIBNPYMATH_ENDIAN_H_ */ diff --git a/src/include/numpy/npymath/fpmath.h b/src/include/libnpymath/fpmath.h similarity index 82% rename from src/include/numpy/npymath/fpmath.h rename to src/include/libnpymath/fpmath.h index 9c3add4..5fb4978 100644 --- a/src/include/numpy/npymath/fpmath.h +++ b/src/include/libnpymath/fpmath.h @@ -1,10 +1,10 @@ -#ifndef NPYMATH_FPMATH_H_ -#define NPYMATH_FPMATH_H_ +#ifndef LIBNPYMATH_FPMATH_H_ +#define LIBNPYMATH_FPMATH_H_ -#include "numpy/npymath/block.h" +#include "libnpymath/block.h" -#include "numpy/npymath/cpu.h" -#include "numpy/npymath/common.h" +#include "libnpymath/cpu.h" +#include "libnpymath/common.h" #if !(defined(NPYMATH_HAVE_LDOUBLE_IEEE_QUAD_BE) || \ defined(NPYMATH_HAVE_LDOUBLE_IEEE_QUAD_LE) || \ @@ -26,4 +26,4 @@ #define NPYMATH_HAVE_LDOUBLE_DOUBLE_DOUBLE_BE #endif -#endif /* NPYMATH_FPMATH_H_ */ +#endif /* LIBNPYMATH_FPMATH_H_ */ diff --git a/src/include/libnpymath/halffloat.h b/src/include/libnpymath/halffloat.h new file mode 100644 index 0000000..f5348f3 --- /dev/null +++ b/src/include/libnpymath/halffloat.h @@ -0,0 +1,96 @@ +#ifndef LIBNPYMATH_HALFFLOAT_H_ +#define LIBNPYMATH_HALFFLOAT_H_ + +#include +#include "libnpymath/npy_math.h" + +#ifdef __cplusplus +extern "C" { +#endif + +/* + * Half-precision routines + */ + +/* Conversions */ +float npymath_half_to_float(npymath_half h); +double npymath_half_to_double(npymath_half h); +npymath_half npymath_float_to_half(float f); +npymath_half npymath_double_to_half(double d); +#define NPYMATH_HALF_TO_FLOAT(x) npymath_half_to_float(x) +#define NPYMATH_HALF_TO_DOUBLE(x) npymath_half_to_double(x) +#define NPYMATH_FLOAT_TO_HALF(x) npymath_float_to_half(x) +#define NPYMATH_DOUBLE_TO_HALF(x) npymath_double_to_half(x) +/* Comparisons */ +int npymath_half_eq(npymath_half h1, npymath_half h2); +int npymath_half_ne(npymath_half h1, npymath_half h2); +int npymath_half_le(npymath_half h1, npymath_half h2); +int npymath_half_lt(npymath_half h1, npymath_half h2); +int npymath_half_ge(npymath_half h1, npymath_half h2); +int npymath_half_gt(npymath_half h1, npymath_half h2); +#define NPYMATH_HALF_EQ(h1, h2) npymath_half_eq(h1, h2) +#define NPYMATH_HALF_NE(h1, h2) npymath_half_ne(h1, h2) +#define NPYMATH_HALF_LE(h1, h2) npymath_half_le(h1, h2) +#define NPYMATH_HALF_LT(h1, h2) npymath_half_lt(h1, h2) +#define NPYMATH_HALF_GE(h1, h2) npymath_half_ge(h1, h2) +#define NPYMATH_HALF_GT(h1, h2) npymath_half_gt(h1, h2) +/* faster *_nonan variants for when you know h1 and h2 are not NaN */ +int npymath_half_eq_nonan(npymath_half h1, npymath_half h2); +int npymath_half_lt_nonan(npymath_half h1, npymath_half h2); +int npymath_half_le_nonan(npymath_half h1, npymath_half h2); +#define NPYMATH_HALF_EQ_NONAN(h1, h2) npymath_half_eq_nonan(h1, h2) +#define NPYMATH_HALF_LT_NONAN(h1, h2) npymath_half_lt_nonan(h1, h2) +#define NPYMATH_HALF_LE_NONAN(h1, h2) npymath_half_le_nonan(h1, h2) +/* Miscellaneous functions */ +int npymath_half_iszero(npymath_half h); +int npymath_half_isnan(npymath_half h); +int npymath_half_isinf(npymath_half h); +int npymath_half_isfinite(npymath_half h); +int npymath_half_signbit(npymath_half h); +npymath_half npymath_half_copysign(npymath_half x, npymath_half y); +npymath_half npymath_half_spacing(npymath_half h); +npymath_half npymath_half_nextafter(npymath_half x, npymath_half y); +npymath_half npymath_half_divmod(npymath_half x, npymath_half y, npymath_half *modulus); +#define NPYMATH_HALF_ISZERO(h) npymath_half_iszero(h) +#define NPYMATH_HALF_ISNAN(h) npymath_half_isnan(h) +#define NPYMATH_HALF_ISINF(h) npymath_half_isinf(h) +#define NPYMATH_HALF_ISFINITE(h) npymath_half_isfinite(h) +#define NPYMATH_HALF_SIGNBIT(h) npymath_half_signbit(h) +#define NPYMATH_HALF_COPYSIGN(x, y) npymath_half_copysign(x, y) +#define NPYMATH_HALF_SPACING(h) npymath_half_spacing(h) +#define NPYMATH_HALF_NEXTAFTER(x, y) npymath_half_nextafter(x, y) +#define NPYMATH_HALF_DIVMOD(x, y, modulus) npymath_half_divmod(x, y, modulus) + +/* + * Half-precision constants + */ + +#define NPYMATH_HALF_ZERO (0x0000u) +#define NPYMATH_HALF_PZERO (0x0000u) +#define NPYMATH_HALF_NZERO (0x8000u) +#define NPYMATH_HALF_ONE (0x3c00u) +#define NPYMATH_HALF_NEGONE (0xbc00u) +#define NPYMATH_HALF_PINF (0x7c00u) +#define NPYMATH_HALF_NINF (0xfc00u) +#define NPYMATH_HALF_NAN (0x7e00u) + +#define NPYMATH_MAX_HALF (0x7bffu) + +/* + * Bit-level conversions + */ + +npymath_uint16 npymath_floatbits_to_halfbits(npymath_uint32 f); +npymath_uint16 npymath_doublebits_to_halfbits(npymath_uint64 d); +npymath_uint32 npymath_halfbits_to_floatbits(npymath_uint16 h); +npymath_uint64 npymath_halfbits_to_doublebits(npymath_uint16 h); +#define NPYMATH_FLOATBITS_TO_HALFBITS(x) npymath_floatbits_to_halfbits(x) +#define NPYMATH_DOUBLEBITS_TO_HALFBITS(x) npymath_doublebits_to_halfbits(x) +#define NPYMATH_HALFBITS_TO_FLOATBITS(x) npymath_halfbits_to_floatbits(x) +#define NPYMATH_HALFBITS_TO_DOUBLEBITS(x) npymath_halfbits_to_doublebits(x) + +#ifdef __cplusplus +} +#endif + +#endif /* LIBNPYMATH_HALFFLOAT_H_ */ diff --git a/src/include/numpy/npymath/internal.h.src b/src/include/libnpymath/internal.h.src similarity index 71% rename from src/include/numpy/npymath/internal.h.src rename to src/include/libnpymath/internal.h.src index e0757fa..74b07be 100644 --- a/src/include/numpy/npymath/internal.h.src +++ b/src/include/libnpymath/internal.h.src @@ -53,10 +53,10 @@ * is preserved. * ==================================================== */ -#ifndef NPYMATH_INTERNAL_H_ -#define NPYMATH_INTERNAL_H_ +#ifndef LIBNPYMATH_INTERNAL_H_ +#define LIBNPYMATH_INTERNAL_H_ -#include "numpy/npymath/private.h" +#include "libnpymath/private.h" #ifdef _MSC_VER # include // for __popcnt #endif @@ -96,17 +96,17 @@ static const npymath_uint64 MAGIC64[] = {0x5555555555555555ull, 0x33333333333333 /* Taken from FreeBSD mlib, adapted for numpy * * XXX: we could be a bit faster by reusing high/low words for inf/nan - * classification instead of calling npy_isinf/npy_isnan: we should have some + * classification instead of calling npymath_isinf/npymath_isnan: we should have some * macros for this, though, instead of doing it manually */ -NPY_INPLACE double npy_log2(double x) +NPYMATH_INPLACE double npymath_log2(double x) { #ifndef NPYMATH_BLOCK_LOG2 return log2(x); #else - if (!npy_isfinite(x) || x <= 0.) { + if (!npymath_isfinite(x) || x <= 0.) { /* special value result */ - return npy_log(x); + return npymath_log(x); } else { /* @@ -134,12 +134,12 @@ NPY_INPLACE double npy_log2(double x) /* Taken from FreeBSD mlib, adapted for numpy * * XXX: we could be a bit faster by reusing high/low words for inf/nan - * classification instead of calling npy_isinf/npy_isnan: we should have some + * classification instead of calling npymath_isinf/npymath_isnan: we should have some * macros for this, though, instead of doing it manually */ /* XXX: we should have this in npy_math.h */ -#define NPY_DBL_EPSILON 1.2246467991473531772E-16 -NPY_INPLACE double npy_atan2(double y, double x) +#define NPYMATH_DBL_EPSILON 1.2246467991473531772E-16 +NPYMATH_INPLACE double npymath_atan2(double y, double x) { #ifndef NPYMATH_BLOCK_ATAN2 return atan2(y, x); @@ -154,67 +154,67 @@ NPY_INPLACE double npy_atan2(double y, double x) iy = hy & 0x7fffffff; /* if x or y is nan, return nan */ - if (npy_isnan(x * y)) { + if (npymath_isnan(x * y)) { return x + y; } if (x == 1.0) { - return npy_atan(y); + return npymath_atan(y); } - m = 2 * (npy_signbit((x)) != 0) + (npy_signbit((y)) != 0); + m = 2 * (npymath_signbit((x)) != 0) + (npymath_signbit((y)) != 0); if (y == 0.0) { switch(m) { case 0: case 1: return y; /* atan(+-0,+anything)=+-0 */ - case 2: return NPY_PI;/* atan(+0,-anything) = pi */ - case 3: return -NPY_PI;/* atan(-0,-anything) =-pi */ + case 2: return NPYMATH_PI;/* atan(+0,-anything) = pi */ + case 3: return -NPYMATH_PI;/* atan(-0,-anything) =-pi */ } } if (x == 0.0) { - return y > 0 ? NPY_PI_2 : -NPY_PI_2; + return y > 0 ? NPYMATH_PI_2 : -NPYMATH_PI_2; } - if (npy_isinf(x)) { - if (npy_isinf(y)) { + if (npymath_isinf(x)) { + if (npymath_isinf(y)) { switch(m) { - case 0: return NPY_PI_4;/* atan(+INF,+INF) */ - case 1: return -NPY_PI_4;/* atan(-INF,+INF) */ - case 2: return 3.0*NPY_PI_4;/*atan(+INF,-INF)*/ - case 3: return -3.0*NPY_PI_4;/*atan(-INF,-INF)*/ + case 0: return NPYMATH_PI_4;/* atan(+INF,+INF) */ + case 1: return -NPYMATH_PI_4;/* atan(-INF,+INF) */ + case 2: return 3.0*NPYMATH_PI_4;/*atan(+INF,-INF)*/ + case 3: return -3.0*NPYMATH_PI_4;/*atan(-INF,-INF)*/ } } else { switch(m) { - case 0: return NPY_PZERO; /* atan(+...,+INF) */ - case 1: return NPY_NZERO; /* atan(-...,+INF) */ - case 2: return NPY_PI; /* atan(+...,-INF) */ - case 3: return -NPY_PI; /* atan(-...,-INF) */ + case 0: return NPYMATH_PZERO; /* atan(+...,+INF) */ + case 1: return NPYMATH_NZERO; /* atan(-...,+INF) */ + case 2: return NPYMATH_PI; /* atan(+...,-INF) */ + case 3: return -NPYMATH_PI; /* atan(-...,-INF) */ } } } - if (npy_isinf(y)) { - return y > 0 ? NPY_PI_2 : -NPY_PI_2; + if (npymath_isinf(y)) { + return y > 0 ? NPYMATH_PI_2 : -NPYMATH_PI_2; } /* compute y/x */ k = (iy - ix) >> 20; if (k > 60) { /* |y/x| > 2**60 */ - z = NPY_PI_2 + 0.5 * NPY_DBL_EPSILON; + z = NPYMATH_PI_2 + 0.5 * NPYMATH_DBL_EPSILON; m &= 1; } else if (hx < 0 && k < -60) { z = 0.0; /* 0 > |y|/x > -2**-60 */ } else { - z = npy_atan(npy_fabs(y/x)); /* safe to do y/x */ + z = npymath_atan(npymath_fabs(y/x)); /* safe to do y/x */ } switch (m) { case 0: return z ; /* atan(+,+) */ case 1: return -z ; /* atan(-,+) */ - case 2: return NPY_PI - (z - NPY_DBL_EPSILON);/* atan(+,-) */ + case 2: return NPYMATH_PI - (z - NPYMATH_DBL_EPSILON);/* atan(+,-) */ default: /* case 3 */ - return (z - NPY_DBL_EPSILON) - NPY_PI;/* atan(-,-) */ + return (z - NPYMATH_DBL_EPSILON) - NPYMATH_PI;/* atan(-,-) */ } #endif } @@ -223,23 +223,23 @@ NPY_INPLACE double npy_atan2(double y, double x) -NPY_INPLACE double npy_hypot(double x, double y) +NPYMATH_INPLACE double npymath_hypot(double x, double y) { #ifndef NPYMATH_BLOCK_HYPOT return hypot(x, y); #else double yx; - if (npy_isinf(x) || npy_isinf(y)) { - return NPY_INFINITY; + if (npymath_isinf(x) || npymath_isinf(y)) { + return NPYMATH_INFINITY; } - if (npy_isnan(x) || npy_isnan(y)) { - return NPY_NAN; + if (npymath_isnan(x) || npymath_isnan(y)) { + return NPYMATH_NAN; } - x = npy_fabs(x); - y = npy_fabs(y); + x = npymath_fabs(x); + y = npymath_fabs(y); if (x < y) { double temp = x; x = y; @@ -250,7 +250,7 @@ NPY_INPLACE double npy_hypot(double x, double y) } else { yx = y/x; - return x*npy_sqrt(1.+yx*yx); + return x*npymath_sqrt(1.+yx*yx); } #endif } @@ -285,11 +285,11 @@ NPY_INPLACE double npy_hypot(double x, double y) * #c = l,,f# * #C = L,,F# */ -#undef NPY__FP_SFX +#undef NPYMATH__FP_SFX #if NPYMATH_SIZEOF_@TYPE@ == NPYMATH_SIZEOF_DOUBLE - #define NPY__FP_SFX(X) X + #define NPYMATH__FP_SFX(X) X #else - #define NPY__FP_SFX(X) NPYMATH_CAT(X, @c@) + #define NPYMATH__FP_SFX(X) NPYMATH_CAT(X, @c@) #endif /* * On arm64 macOS, there's a bug with sin, cos, and tan where they don't @@ -305,14 +305,14 @@ NPY_INPLACE double npy_hypot(double x, double y) * #kind = sin,cos,tan# * #TRIG_WORKAROUND = WORKAROUND_APPLE_TRIG_BUG*3# */ -NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_@kind@@c@(@type@ x) { #if @TRIG_WORKAROUND@ - if (!npy_isfinite(x)) { + if (!npymath_isfinite(x)) { return (x - x); } #endif - return NPY__FP_SFX(@kind@)(x); + return NPYMATH__FP_SFX(@kind@)(x); } /**end repeat1**/ @@ -329,11 +329,11 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) * #c = f,l# * #C = F,L# */ -#undef NPY__FP_SFX +#undef NPYMATH__FP_SFX #if NPYMATH_SIZEOF_@TYPE@ == NPYMATH_SIZEOF_DOUBLE - #define NPY__FP_SFX(X) X + #define NPYMATH__FP_SFX(X) X #else - #define NPY__FP_SFX(X) NPYMATH_CAT(X, @c@) + #define NPYMATH__FP_SFX(X) NPYMATH_CAT(X, @c@) #endif /**begin repeat1 @@ -345,16 +345,16 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) #undef @kind@@c@ #endif #ifdef NPYMATH_BLOCK_@KIND@@C@ -NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_@kind@@c@(@type@ x) { - return (@type@) npy_@kind@((double)x); + return (@type@) npymath_@kind@((double)x); } #endif #ifndef NPYMATH_BLOCK_@KIND@@C@ -NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_@kind@@c@(@type@ x) { - return NPY__FP_SFX(@kind@)(x); + return NPYMATH__FP_SFX(@kind@)(x); } #endif @@ -369,16 +369,16 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x) #undef @kind@@c@ #endif #ifdef NPYMATH_BLOCK_@KIND@@C@ -NPY_INPLACE @type@ npy_@kind@@c@(@type@ x, @type@ y) +NPYMATH_INPLACE @type@ npymath_@kind@@c@(@type@ x, @type@ y) { - return (@type@) npy_@kind@((double)x, (double) y); + return (@type@) npymath_@kind@((double)x, (double) y); } #endif #ifndef NPYMATH_BLOCK_@KIND@@C@ -NPY_INPLACE @type@ npy_@kind@@c@(@type@ x, @type@ y) +NPYMATH_INPLACE @type@ npymath_@kind@@c@(@type@ x, @type@ y) { - return NPY__FP_SFX(@kind@)(x, y); + return NPYMATH__FP_SFX(@kind@)(x, y); } #endif /**end repeat1**/ @@ -387,19 +387,19 @@ NPY_INPLACE @type@ npy_@kind@@c@(@type@ x, @type@ y) #undef modf@c@ #endif #ifdef NPYMATH_BLOCK_MODF@C@ -NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr) +NPYMATH_INPLACE @type@ npymath_modf@c@(@type@ x, @type@ *iptr) { double niptr; - double y = npy_modf((double)x, &niptr); + double y = npymath_modf((double)x, &niptr); *iptr = (@type@) niptr; return (@type@) y; } #endif #ifndef NPYMATH_BLOCK_MODF@C@ -NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr) +NPYMATH_INPLACE @type@ npymath_modf@c@(@type@ x, @type@ *iptr) { - return NPY__FP_SFX(modf)(x, iptr); + return NPYMATH__FP_SFX(modf)(x, iptr); } #endif @@ -407,7 +407,7 @@ NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr) /**end repeat**/ -#undef NPY__FP_SFX +#undef NPYMATH__FP_SFX /* @@ -420,16 +420,16 @@ NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr) * #c = f, ,l# * #C = F, ,L# */ -#undef NPY__FP_SFX +#undef NPYMATH__FP_SFX #if NPYMATH_SIZEOF_@TYPE@ == NPYMATH_SIZEOF_DOUBLE - #define NPY__FP_SFX(X) X + #define NPYMATH__FP_SFX(X) X #else - #define NPY__FP_SFX(X) NPYMATH_CAT(X, @c@) + #define NPYMATH__FP_SFX(X) NPYMATH_CAT(X, @c@) #endif -@type@ npy_heaviside@c@(@type@ x, @type@ h0) +@type@ npymath_heaviside@c@(@type@ x, @type@ h0) { - if (npy_isnan(x)) { - return (@type@) NPY_NAN; + if (npymath_isnan(x)) { + return (@type@) NPYMATH_NAN; } else if (x == 0) { return h0; @@ -442,32 +442,32 @@ NPY_INPLACE @type@ npy_modf@c@(@type@ x, @type@ *iptr) } } -#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)) +#define LOGE2 NPYMATH__FP_SFX(NPYMATH_LOGE2) +#define LOG2E NPYMATH__FP_SFX(NPYMATH_LOG2E) +#define RAD2DEG (NPYMATH__FP_SFX(180.0)/NPYMATH__FP_SFX(NPYMATH_PI)) +#define DEG2RAD (NPYMATH__FP_SFX(NPYMATH_PI)/NPYMATH__FP_SFX(180.0)) -NPY_INPLACE @type@ npy_rad2deg@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_rad2deg@c@(@type@ x) { return x*RAD2DEG; } -NPY_INPLACE @type@ npy_deg2rad@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_deg2rad@c@(@type@ x) { return x*DEG2RAD; } -NPY_INPLACE @type@ npy_log2_1p@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_log2_1p@c@(@type@ x) { - return LOG2E*npy_log1p@c@(x); + return LOG2E*npymath_log1p@c@(x); } -NPY_INPLACE @type@ npy_exp2_m1@c@(@type@ x) +NPYMATH_INPLACE @type@ npymath_exp2_m1@c@(@type@ x) { - return npy_expm1@c@(LOGE2*x); + return npymath_expm1@c@(LOGE2*x); } -NPY_INPLACE @type@ npy_logaddexp@c@(@type@ x, @type@ y) +NPYMATH_INPLACE @type@ npymath_logaddexp@c@(@type@ x, @type@ y) { if (x == y) { /* Handles infinities of the same sign without warnings */ @@ -476,10 +476,10 @@ NPY_INPLACE @type@ npy_logaddexp@c@(@type@ x, @type@ y) else { const @type@ tmp = x - y; if (tmp > 0) { - return x + npy_log1p@c@(npy_exp@c@(-tmp)); + return x + npymath_log1p@c@(npymath_exp@c@(-tmp)); } else if (tmp <= 0) { - return y + npy_log1p@c@(npy_exp@c@(tmp)); + return y + npymath_log1p@c@(npymath_exp@c@(tmp)); } else { /* NaNs */ @@ -488,7 +488,7 @@ NPY_INPLACE @type@ npy_logaddexp@c@(@type@ x, @type@ y) } } -NPY_INPLACE @type@ npy_logaddexp2@c@(@type@ x, @type@ y) +NPYMATH_INPLACE @type@ npymath_logaddexp2@c@(@type@ x, @type@ y) { if (x == y) { /* Handles infinities of the same sign without warnings */ @@ -497,10 +497,10 @@ NPY_INPLACE @type@ npy_logaddexp2@c@(@type@ x, @type@ y) else { const @type@ tmp = x - y; if (tmp > 0) { - return x + npy_log2_1p@c@(npy_exp2@c@(-tmp)); + return x + npymath_log2_1p@c@(npymath_exp2@c@(-tmp)); } else if (tmp <= 0) { - return y + npy_log2_1p@c@(npy_exp2@c@(tmp)); + return y + npymath_log2_1p@c@(npymath_exp2@c@(tmp)); } else { /* NaNs */ @@ -511,10 +511,10 @@ NPY_INPLACE @type@ npy_logaddexp2@c@(@type@ x, @type@ y) /* * Wrapper function for remainder edge cases - * Internally calls npy_divmod* + * Internally calls npymath_divmod* */ -NPY_INPLACE @type@ -npy_remainder@c@(@type@ a, @type@ b) +NPYMATH_INPLACE @type@ +npymath_remainder@c@(@type@ a, @type@ b) { @type@ mod; if (NPYMATH_UNLIKELY(!b)) { @@ -523,16 +523,16 @@ npy_remainder@c@(@type@ a, @type@ b) * result (always NaN). `divmod` may set additional FPE for the * division by zero creating an inf. */ - mod = npy_fmod@c@(a, b); + mod = npymath_fmod@c@(a, b); } else { - npy_divmod@c@(a, b, &mod); + npymath_divmod@c@(a, b, &mod); } return mod; } -NPY_INPLACE @type@ -npy_floor_divide@c@(@type@ a, @type@ b) { +NPYMATH_INPLACE @type@ +npymath_floor_divide@c@(@type@ a, @type@ b) { @type@ div, mod; if (NPYMATH_UNLIKELY(!b)) { /* @@ -543,7 +543,7 @@ npy_floor_divide@c@(@type@ a, @type@ b) { div = a / b; } else { - div = npy_divmod@c@(a, b, &mod); + div = npymath_divmod@c@(a, b, &mod); } return div; } @@ -553,12 +553,12 @@ npy_floor_divide@c@(@type@ a, @type@ b) { * * The implementation is mostly copied from cpython 3.5. */ -NPY_INPLACE @type@ -npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus) +NPYMATH_INPLACE @type@ +npymath_divmod@c@(@type@ a, @type@ b, @type@ *modulus) { @type@ div, mod, floordiv; - mod = npy_fmod@c@(a, b); + mod = npymath_fmod@c@(a, b); if (NPYMATH_UNLIKELY(!b)) { /* b == 0 (not NaN): return result of fmod. For IEEE is nan */ *modulus = mod; @@ -577,18 +577,18 @@ npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus) } else { /* if mod is zero ensure correct sign */ - mod = npy_copysign@c@(0, b); + mod = npymath_copysign@c@(0, b); } /* snap quotient to nearest integral value */ if (div) { - floordiv = npy_floor@c@(div); + floordiv = npymath_floor@c@(div); if (isgreater(div - floordiv, 0.5@c@)) floordiv += 1.0@c@; } else { /* if div is zero ensure correct sign */ - floordiv = npy_copysign@c@(0, a/b); + floordiv = npymath_copysign@c@(0, a/b); } *modulus = mod; @@ -599,7 +599,7 @@ npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus) #undef LOG2E #undef RAD2DEG #undef DEG2RAD -#undef NPY__FP_SFX +#undef NPYMATH__FP_SFX /**end repeat**/ /**begin repeat @@ -607,8 +607,8 @@ npy_divmod@c@(@type@ a, @type@ b, @type@ *modulus) * #type = npymath_uint, npymath_ulong, npymath_ulonglong# * #c = u,ul,ull# */ -NPY_INPLACE @type@ -npy_gcd@c@(@type@ a, @type@ b) +NPYMATH_INPLACE @type@ +npymath_gcd@c@(@type@ a, @type@ b) { @type@ c; while (a != 0) { @@ -619,10 +619,10 @@ npy_gcd@c@(@type@ a, @type@ b) return b; } -NPY_INPLACE @type@ -npy_lcm@c@(@type@ a, @type@ b) +NPYMATH_INPLACE @type@ +npymath_lcm@c@(@type@ a, @type@ b) { - @type@ gcd = npy_gcd@c@(a, b); + @type@ gcd = npymath_gcd@c@(a, b); return gcd == 0 ? 0 : a / gcd * b; } /**end repeat**/ @@ -633,10 +633,10 @@ npy_lcm@c@(@type@ a, @type@ b) * #c = (,l,ll)*2# * #func=gcd*3,lcm*3# */ -NPY_INPLACE @type@ -npy_@func@@c@(@type@ a, @type@ b) +NPYMATH_INPLACE @type@ +npymath_@func@@c@(@type@ a, @type@ b) { - return npy_@func@u@c@(a < 0 ? -a : a, b < 0 ? -b : b); + return npymath_@func@u@c@(a < 0 ? -a : a, b < 0 ? -b : b); } /**end repeat**/ @@ -653,8 +653,8 @@ npy_@func@@c@(@type@ a, @type@ b) * #u = u,# * #is_signed = 0,1# */ -NPY_INPLACE npymath_@u@@type@ -npy_lshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) +NPYMATH_INPLACE npymath_@u@@type@ +npymath_lshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) { if (NPYMATH_LIKELY((size_t)b < sizeof(a) * CHAR_BIT)) { return a << b; @@ -663,8 +663,8 @@ npy_lshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) return 0; } } -NPY_INPLACE npymath_@u@@type@ -npy_rshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) +NPYMATH_INPLACE npymath_@u@@type@ +npymath_rshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) { if (NPYMATH_LIKELY((size_t)b < sizeof(a) * CHAR_BIT)) { return a >> b; @@ -700,8 +700,8 @@ npy_rshift@u@@c@(npymath_@u@@type@ a, npymath_@u@@type@ b) #endif -NPY_INPLACE uint8_t -npy_popcount_parallel@c@(npymath_@type@ a) +NPYMATH_INPLACE uint8_t +npymath_popcount_parallel@c@(npymath_@type@ a) { a = a - ((a >> 1) & (npymath_@type@) TO_BITS_LEN(MAGIC)[0]); a = ((a & (npymath_@type@) TO_BITS_LEN(MAGIC)[1])) + ((a >> 2) & (npymath_@type@) TO_BITS_LEN(MAGIC)[1]); @@ -709,8 +709,8 @@ npy_popcount_parallel@c@(npymath_@type@ a) return (npymath_@type@) (a * (npymath_@type@) TO_BITS_LEN(MAGIC)[3]) >> ((NPYMATH_SIZEOF_@STYPE@ - 1) * CHAR_BIT); } -NPY_INPLACE uint8_t -npy_popcountu@c@(npymath_@type@ a) +NPYMATH_INPLACE uint8_t +npymath_popcountu@c@(npymath_@type@ a) { /* use built-in popcount if present, else use our implementation */ #if (defined(__clang__) || defined(__GNUC__)) && NPYMATH_BITSOF_@STYPE@ >= 32 @@ -726,7 +726,7 @@ npy_popcountu@c@(npymath_@type@ a) return __popcnt32(left) + __popcnt32(right); #endif #else - return npy_popcount_parallel@c@(a); + return npymath_popcount_parallel@c@(a); #endif } /**end repeat**/ @@ -736,12 +736,12 @@ npy_popcountu@c@(npymath_@type@ a) * #type = byte, short, int, long, longlong# * #c = hh, h, , l, ll# */ -NPY_INPLACE uint8_t -npy_popcount@c@(npymath_@type@ a) +NPYMATH_INPLACE uint8_t +npymath_popcount@c@(npymath_@type@ a) { /* Return popcount of abs(a) */ - return npy_popcountu@c@(a < 0 ? -a : a); + return npymath_popcountu@c@(a < 0 ? -a : a); } /**end repeat**/ -#endif +#endif /* LIBNPYMATH_INTERNAL_H_ */ diff --git a/src/include/numpy/npymath/meson.build b/src/include/libnpymath/meson.build similarity index 90% rename from src/include/numpy/npymath/meson.build rename to src/include/libnpymath/meson.build index c3a8351..4a5f0f0 100644 --- a/src/include/numpy/npymath/meson.build +++ b/src/include/libnpymath/meson.build @@ -21,7 +21,7 @@ if get_option('install_lib') command: [src_file_cli, '@INPUT@', '-o', '@OUTPUT@'], # This installs the header in the numpy include folder install: true, - install_dir: install_dir / get_option('headers_dir') / 'npymath', + install_dir: install_dir / get_option('headers_dir') / 'libnpymath', ) else internal_h = custom_target( diff --git a/src/include/numpy/npymath/meson_config.h.in b/src/include/libnpymath/meson_config.h.in similarity index 95% rename from src/include/numpy/npymath/meson_config.h.in rename to src/include/libnpymath/meson_config.h.in index f309b85..267a49d 100644 --- a/src/include/numpy/npymath/meson_config.h.in +++ b/src/include/libnpymath/meson_config.h.in @@ -78,6 +78,6 @@ /* #undef inline */ #endif -#ifndef NPYMATH_CONFIG_H_ -#error npymath_config.h should never be included directly, include npymath_config.h instead +#ifndef LIBNPYMATH_BLOCK_H_ +#error meson_config.h should never be included directly, include block.h instead #endif diff --git a/src/include/libnpymath/npy_math.h b/src/include/libnpymath/npy_math.h new file mode 100644 index 0000000..be9f261 --- /dev/null +++ b/src/include/libnpymath/npy_math.h @@ -0,0 +1,502 @@ +#ifndef LIBNPYMATH_NPY_MATH_H_ +#define LIBNPYMATH_NPY_MATH_H_ + +#include "libnpymath/common.h" + +#include + +#if defined(NPYMATH_INLINE_MATH) +#define NPYMATH_INPLACE static inline +#else +#define NPYMATH_INPLACE +#endif + + +#ifdef __cplusplus +extern "C" { +#endif + +#define PyArray_MAX(a,b) (((a)>(b))?(a):(b)) +#define PyArray_MIN(a,b) (((a)<(b))?(a):(b)) + +/* + * NAN and INFINITY like macros (same behavior as glibc for NAN, same as C99 + * for INFINITY) + * + * XXX: I should test whether INFINITY and NAN are available on the platform + */ +static inline float __npymath_inff(void) +{ + const union { npymath_uint32 __i; float __f;} __bint = {0x7f800000UL}; + return __bint.__f; +} + +static inline float __npymath_nanf(void) +{ + const union { npymath_uint32 __i; float __f;} __bint = {0x7fc00000UL}; + return __bint.__f; +} + +static inline float __npymath_pzerof(void) +{ + const union { npymath_uint32 __i; float __f;} __bint = {0x00000000UL}; + return __bint.__f; +} + +static inline float __npymath_nzerof(void) +{ + const union { npymath_uint32 __i; float __f;} __bint = {0x80000000UL}; + return __bint.__f; +} + +#define NPYMATH_INFINITYF __npymath_inff() +#define NPYMATH_NANF __npymath_nanf() +#define NPYMATH_PZEROF __npymath_pzerof() +#define NPYMATH_NZEROF __npymath_nzerof() + +#define NPYMATH_INFINITY ((npymath_double)NPYMATH_INFINITYF) +#define NPYMATH_NAN ((npymath_double)NPYMATH_NANF) +#define NPYMATH_PZERO ((npymath_double)NPYMATH_PZEROF) +#define NPYMATH_NZERO ((npymath_double)NPYMATH_NZEROF) + +#define NPYMATH_INFINITYL ((npymath_longdouble)NPYMATH_INFINITYF) +#define NPYMATH_NANL ((npymath_longdouble)NPYMATH_NANF) +#define NPYMATH_PZEROL ((npymath_longdouble)NPYMATH_PZEROF) +#define NPYMATH_NZEROL ((npymath_longdouble)NPYMATH_NZEROF) + +/* + * Useful constants + */ +#define NPYMATH_E 2.718281828459045235360287471352662498 /* e */ +#define NPYMATH_LOG2E 1.442695040888963407359924681001892137 /* log_2 e */ +#define NPYMATH_LOG10E 0.434294481903251827651128918916605082 /* log_10 e */ +#define NPYMATH_LOGE2 0.693147180559945309417232121458176568 /* log_e 2 */ +#define NPYMATH_LOGE10 2.302585092994045684017991454684364208 /* log_e 10 */ +#define NPYMATH_PI 3.141592653589793238462643383279502884 /* pi */ +#define NPYMATH_PI_2 1.570796326794896619231321691639751442 /* pi/2 */ +#define NPYMATH_PI_4 0.785398163397448309615660845819875721 /* pi/4 */ +#define NPYMATH_1_PI 0.318309886183790671537767526745028724 /* 1/pi */ +#define NPYMATH_2_PI 0.636619772367581343075535053490057448 /* 2/pi */ +#define NPYMATH_EULER 0.577215664901532860606512090082402431 /* Euler constant */ +#define NPYMATH_SQRT2 1.414213562373095048801688724209698079 /* sqrt(2) */ +#define NPYMATH_SQRT1_2 0.707106781186547524400844362104849039 /* 1/sqrt(2) */ + +#define NPYMATH_Ef 2.718281828459045235360287471352662498F /* e */ +#define NPYMATH_LOG2Ef 1.442695040888963407359924681001892137F /* log_2 e */ +#define NPYMATH_LOG10Ef 0.434294481903251827651128918916605082F /* log_10 e */ +#define NPYMATH_LOGE2f 0.693147180559945309417232121458176568F /* log_e 2 */ +#define NPYMATH_LOGE10f 2.302585092994045684017991454684364208F /* log_e 10 */ +#define NPYMATH_PIf 3.141592653589793238462643383279502884F /* pi */ +#define NPYMATH_PI_2f 1.570796326794896619231321691639751442F /* pi/2 */ +#define NPYMATH_PI_4f 0.785398163397448309615660845819875721F /* pi/4 */ +#define NPYMATH_1_PIf 0.318309886183790671537767526745028724F /* 1/pi */ +#define NPYMATH_2_PIf 0.636619772367581343075535053490057448F /* 2/pi */ +#define NPYMATH_EULERf 0.577215664901532860606512090082402431F /* Euler constant */ +#define NPYMATH_SQRT2f 1.414213562373095048801688724209698079F /* sqrt(2) */ +#define NPYMATH_SQRT1_2f 0.707106781186547524400844362104849039F /* 1/sqrt(2) */ + +#define NPYMATH_El 2.718281828459045235360287471352662498L /* e */ +#define NPYMATH_LOG2El 1.442695040888963407359924681001892137L /* log_2 e */ +#define NPYMATH_LOG10El 0.434294481903251827651128918916605082L /* log_10 e */ +#define NPYMATH_LOGE2l 0.693147180559945309417232121458176568L /* log_e 2 */ +#define NPYMATH_LOGE10l 2.302585092994045684017991454684364208L /* log_e 10 */ +#define NPYMATH_PIl 3.141592653589793238462643383279502884L /* pi */ +#define NPYMATH_PI_2l 1.570796326794896619231321691639751442L /* pi/2 */ +#define NPYMATH_PI_4l 0.785398163397448309615660845819875721L /* pi/4 */ +#define NPYMATH_1_PIl 0.318309886183790671537767526745028724L /* 1/pi */ +#define NPYMATH_2_PIl 0.636619772367581343075535053490057448L /* 2/pi */ +#define NPYMATH_EULERl 0.577215664901532860606512090082402431L /* Euler constant */ +#define NPYMATH_SQRT2l 1.414213562373095048801688724209698079L /* sqrt(2) */ +#define NPYMATH_SQRT1_2l 0.707106781186547524400844362104849039L /* 1/sqrt(2) */ + +/* + * Integer functions. + */ +NPYMATH_INPLACE npymath_uint npymath_gcdu(npymath_uint a, npymath_uint b); +NPYMATH_INPLACE npymath_uint npymath_lcmu(npymath_uint a, npymath_uint b); +NPYMATH_INPLACE npymath_ulong npymath_gcdul(npymath_ulong a, npymath_ulong b); +NPYMATH_INPLACE npymath_ulong npymath_lcmul(npymath_ulong a, npymath_ulong b); +NPYMATH_INPLACE npymath_ulonglong npymath_gcdull(npymath_ulonglong a, npymath_ulonglong b); +NPYMATH_INPLACE npymath_ulonglong npymath_lcmull(npymath_ulonglong a, npymath_ulonglong b); + +NPYMATH_INPLACE npymath_int npymath_gcd(npymath_int a, npymath_int b); +NPYMATH_INPLACE npymath_int npymath_lcm(npymath_int a, npymath_int b); +NPYMATH_INPLACE npymath_long npymath_gcdl(npymath_long a, npymath_long b); +NPYMATH_INPLACE npymath_long npymath_lcml(npymath_long a, npymath_long b); +NPYMATH_INPLACE npymath_longlong npymath_gcdll(npymath_longlong a, npymath_longlong b); +NPYMATH_INPLACE npymath_longlong npymath_lcmll(npymath_longlong a, npymath_longlong b); + +NPYMATH_INPLACE npymath_ubyte npymath_rshiftuhh(npymath_ubyte a, npymath_ubyte b); +NPYMATH_INPLACE npymath_ubyte npymath_lshiftuhh(npymath_ubyte a, npymath_ubyte b); +NPYMATH_INPLACE npymath_ushort npymath_rshiftuh(npymath_ushort a, npymath_ushort b); +NPYMATH_INPLACE npymath_ushort npymath_lshiftuh(npymath_ushort a, npymath_ushort b); +NPYMATH_INPLACE npymath_uint npymath_rshiftu(npymath_uint a, npymath_uint b); +NPYMATH_INPLACE npymath_uint npymath_lshiftu(npymath_uint a, npymath_uint b); +NPYMATH_INPLACE npymath_ulong npymath_rshiftul(npymath_ulong a, npymath_ulong b); +NPYMATH_INPLACE npymath_ulong npymath_lshiftul(npymath_ulong a, npymath_ulong b); +NPYMATH_INPLACE npymath_ulonglong npymath_rshiftull(npymath_ulonglong a, npymath_ulonglong b); +NPYMATH_INPLACE npymath_ulonglong npymath_lshiftull(npymath_ulonglong a, npymath_ulonglong b); + +NPYMATH_INPLACE npymath_byte npymath_rshifthh(npymath_byte a, npymath_byte b); +NPYMATH_INPLACE npymath_byte npymath_lshifthh(npymath_byte a, npymath_byte b); +NPYMATH_INPLACE npymath_short npymath_rshifth(npymath_short a, npymath_short b); +NPYMATH_INPLACE npymath_short npymath_lshifth(npymath_short a, npymath_short b); +NPYMATH_INPLACE npymath_int npymath_rshift(npymath_int a, npymath_int b); +NPYMATH_INPLACE npymath_int npymath_lshift(npymath_int a, npymath_int b); +NPYMATH_INPLACE npymath_long npymath_rshiftl(npymath_long a, npymath_long b); +NPYMATH_INPLACE npymath_long npymath_lshiftl(npymath_long a, npymath_long b); +NPYMATH_INPLACE npymath_longlong npymath_rshiftll(npymath_longlong a, npymath_longlong b); +NPYMATH_INPLACE npymath_longlong npymath_lshiftll(npymath_longlong a, npymath_longlong b); + +NPYMATH_INPLACE uint8_t npymath_popcountuhh(npymath_ubyte a); +NPYMATH_INPLACE uint8_t npymath_popcountuh(npymath_ushort a); +NPYMATH_INPLACE uint8_t npymath_popcountu(npymath_uint a); +NPYMATH_INPLACE uint8_t npymath_popcountul(npymath_ulong a); +NPYMATH_INPLACE uint8_t npymath_popcountull(npymath_ulonglong a); +NPYMATH_INPLACE uint8_t npymath_popcounthh(npymath_byte a); +NPYMATH_INPLACE uint8_t npymath_popcounth(npymath_short a); +NPYMATH_INPLACE uint8_t npymath_popcount(npymath_int a); +NPYMATH_INPLACE uint8_t npymath_popcountl(npymath_long a); +NPYMATH_INPLACE uint8_t npymath_popcountll(npymath_longlong a); + +/* + * C99 double math funcs that need fixups or are blocklist-able + */ +NPYMATH_INPLACE double npymath_sin(double x); +NPYMATH_INPLACE double npymath_cos(double x); +NPYMATH_INPLACE double npymath_tan(double x); +NPYMATH_INPLACE double npymath_hypot(double x, double y); +NPYMATH_INPLACE double npymath_log2(double x); +NPYMATH_INPLACE double npymath_atan2(double x, double y); + +/* Mandatory C99 double math funcs, no blocklisting or fixups */ +/* defined for legacy reasons, should be deprecated at some point */ +#define npymath_sinh sinh +#define npymath_cosh cosh +#define npymath_tanh tanh +#define npymath_asin asin +#define npymath_acos acos +#define npymath_atan atan +#define npymath_log log +#define npymath_log10 log10 +#define npymath_cbrt cbrt +#define npymath_fabs fabs +#define npymath_ceil ceil +#define npymath_fmod fmod +#define npymath_floor floor +#define npymath_expm1 expm1 +#define npymath_log1p log1p +#define npymath_acosh acosh +#define npymath_asinh asinh +#define npymath_atanh atanh +#define npymath_rint rint +#define npymath_trunc trunc +#define npymath_exp2 exp2 +#define npymath_frexp frexp +#define npymath_ldexp ldexp +#define npymath_copysign copysign +#define npymath_exp exp +#define npymath_sqrt sqrt +#define npymath_pow pow +#define npymath_modf modf +#define npymath_nextafter nextafter + +double npymath_spacing(double x); + +/* + * IEEE 754 fpu handling + */ + +/* use builtins to avoid function calls in tight loops + * only available if npymath_config.h is available (= numpys own build) */ +#ifdef NPYMATH_HAVE___BUILTIN_ISNAN + #define npymath_isnan(x) __builtin_isnan(x) +#else + #define npymath_isnan(x) isnan(x) +#endif + + +/* only available if npymath_config.h is available (= numpys own build) */ +#ifdef NPYMATH_HAVE___BUILTIN_ISFINITE + #define npymath_isfinite(x) __builtin_isfinite(x) +#else + #define npymath_isfinite(x) isfinite((x)) +#endif + +/* only available if npymath_config.h is available (= numpys own build) */ +#ifdef NPYMATH_HAVE___BUILTIN_ISINF + #define npymath_isinf(x) __builtin_isinf(x) +#else + #define npymath_isinf(x) isinf((x)) +#endif + +#define npymath_signbit(x) signbit((x)) + +/* + * float C99 math funcs that need fixups or are blocklist-able + */ +NPYMATH_INPLACE float npymath_sinf(float x); +NPYMATH_INPLACE float npymath_cosf(float x); +NPYMATH_INPLACE float npymath_tanf(float x); +NPYMATH_INPLACE float npymath_expf(float x); +NPYMATH_INPLACE float npymath_sqrtf(float x); +NPYMATH_INPLACE float npymath_hypotf(float x, float y); +NPYMATH_INPLACE float npymath_log2f(float x); +NPYMATH_INPLACE float npymath_atan2f(float x, float y); +NPYMATH_INPLACE float npymath_powf(float x, float y); +NPYMATH_INPLACE float npymath_modff(float x, float* y); + +/* Mandatory C99 float math funcs, no blocklisting or fixups */ +/* defined for legacy reasons, should be deprecated at some point */ + +#define npymath_sinhf sinhf +#define npymath_coshf coshf +#define npymath_tanhf tanhf +#define npymath_asinf asinf +#define npymath_acosf acosf +#define npymath_atanf atanf +#define npymath_logf logf +#define npymath_log10f log10f +#define npymath_cbrtf cbrtf +#define npymath_fabsf fabsf +#define npymath_ceilf ceilf +#define npymath_fmodf fmodf +#define npymath_floorf floorf +#define npymath_expm1f expm1f +#define npymath_log1pf log1pf +#define npymath_asinhf asinhf +#define npymath_acoshf acoshf +#define npymath_atanhf atanhf +#define npymath_rintf rintf +#define npymath_truncf truncf +#define npymath_exp2f exp2f +#define npymath_frexpf frexpf +#define npymath_ldexpf ldexpf +#define npymath_copysignf copysignf +#define npymath_nextafterf nextafterf + +float npymath_spacingf(float x); + +/* + * long double C99 double math funcs that need fixups or are blocklist-able + */ +NPYMATH_INPLACE npymath_longdouble npymath_sinl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_cosl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_tanl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_expl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_sqrtl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_hypotl(npymath_longdouble x, npymath_longdouble y); +NPYMATH_INPLACE npymath_longdouble npymath_log2l(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_atan2l(npymath_longdouble x, npymath_longdouble y); +NPYMATH_INPLACE npymath_longdouble npymath_powl(npymath_longdouble x, npymath_longdouble y); +NPYMATH_INPLACE npymath_longdouble npymath_modfl(npymath_longdouble x, npymath_longdouble* y); + +/* Mandatory C99 double math funcs, no blocklisting or fixups */ +/* defined for legacy reasons, should be deprecated at some point */ +#define npymath_sinhl sinhl +#define npymath_coshl coshl +#define npymath_tanhl tanhl +#define npymath_fabsl fabsl +#define npymath_floorl floorl +#define npymath_ceill ceill +#define npymath_rintl rintl +#define npymath_truncl truncl +#define npymath_cbrtl cbrtl +#define npymath_log10l log10l +#define npymath_logl logl +#define npymath_expm1l expm1l +#define npymath_asinl asinl +#define npymath_acosl acosl +#define npymath_atanl atanl +#define npymath_asinhl asinhl +#define npymath_acoshl acoshl +#define npymath_atanhl atanhl +#define npymath_log1pl log1pl +#define npymath_exp2l exp2l +#define npymath_fmodl fmodl +#define npymath_frexpl frexpl +#define npymath_ldexpl ldexpl +#define npymath_copysignl copysignl +#define npymath_nextafterl nextafterl + +npymath_longdouble npymath_spacingl(npymath_longdouble x); + +/* + * Non standard functions + */ +NPYMATH_INPLACE double npymath_deg2rad(double x); +NPYMATH_INPLACE double npymath_rad2deg(double x); +NPYMATH_INPLACE double npymath_logaddexp(double x, double y); +NPYMATH_INPLACE double npymath_logaddexp2(double x, double y); +NPYMATH_INPLACE double npymath_divmod(double x, double y, double *modulus); +NPYMATH_INPLACE double npymath_heaviside(double x, double h0); + +NPYMATH_INPLACE float npymath_deg2radf(float x); +NPYMATH_INPLACE float npymath_rad2degf(float x); +NPYMATH_INPLACE float npymath_logaddexpf(float x, float y); +NPYMATH_INPLACE float npymath_logaddexp2f(float x, float y); +NPYMATH_INPLACE float npymath_divmodf(float x, float y, float *modulus); +NPYMATH_INPLACE float npymath_heavisidef(float x, float h0); + +NPYMATH_INPLACE npymath_longdouble npymath_deg2radl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_rad2degl(npymath_longdouble x); +NPYMATH_INPLACE npymath_longdouble npymath_logaddexpl(npymath_longdouble x, npymath_longdouble y); +NPYMATH_INPLACE npymath_longdouble npymath_logaddexp2l(npymath_longdouble x, npymath_longdouble y); +NPYMATH_INPLACE npymath_longdouble npymath_divmodl(npymath_longdouble x, npymath_longdouble y, + npymath_longdouble *modulus); +NPYMATH_INPLACE npymath_longdouble npymath_heavisidel(npymath_longdouble x, npymath_longdouble h0); + +#define npymath_degrees npymath_rad2deg +#define npymath_degreesf npymath_rad2degf +#define npymath_degreesl npymath_rad2degl + +#define npymath_radians npymath_deg2rad +#define npymath_radiansf npymath_deg2radf +#define npymath_radiansl npymath_deg2radl + +/* + * Complex declarations + */ + +npymath_double npymath_creal(const npymath_cdouble *z); +void npymath_csetreal(npymath_cdouble *z, npymath_double r); +npymath_double npymath_cimag(const npymath_cdouble *z); +void npymath_csetimag(npymath_cdouble *z, npymath_double i); +npymath_float npymath_crealf(const npymath_cfloat *z); +void npymath_csetrealf(npymath_cfloat *z, npymath_float r); +npymath_float npymath_cimagf(const npymath_cfloat *z); +void npymath_csetimagf(npymath_cfloat *z, npymath_float i); +npymath_longdouble npymath_creall(const npymath_clongdouble *z); +void npymath_csetreall(npymath_clongdouble *z, npymath_longdouble r); +npymath_longdouble npymath_cimagl(const npymath_clongdouble *z); +void npymath_csetimagl(npymath_clongdouble *z, npymath_longdouble i); + +npymath_cdouble npymath_cpack(npymath_double x, npymath_double y); +npymath_cfloat npymath_cpackf(npymath_float x, npymath_float y); +npymath_clongdouble npymath_cpackl(npymath_longdouble x, npymath_longdouble y); + +/* + * Double precision complex functions + */ +double npymath_cabs(npymath_cdouble z); +double npymath_carg(npymath_cdouble z); + +npymath_cdouble npymath_cexp(npymath_cdouble z); +npymath_cdouble npymath_clog(npymath_cdouble z); +npymath_cdouble npymath_cpow(npymath_cdouble x, npymath_cdouble y); + +npymath_cdouble npymath_csqrt(npymath_cdouble z); + +npymath_cdouble npymath_ccos(npymath_cdouble z); +npymath_cdouble npymath_csin(npymath_cdouble z); +npymath_cdouble npymath_ctan(npymath_cdouble z); + +npymath_cdouble npymath_ccosh(npymath_cdouble z); +npymath_cdouble npymath_csinh(npymath_cdouble z); +npymath_cdouble npymath_ctanh(npymath_cdouble z); + +npymath_cdouble npymath_cacos(npymath_cdouble z); +npymath_cdouble npymath_casin(npymath_cdouble z); +npymath_cdouble npymath_catan(npymath_cdouble z); + +npymath_cdouble npymath_cacosh(npymath_cdouble z); +npymath_cdouble npymath_casinh(npymath_cdouble z); +npymath_cdouble npymath_catanh(npymath_cdouble z); + +/* + * Single precision complex functions + */ +float npymath_cabsf(npymath_cfloat z); +float npymath_cargf(npymath_cfloat z); + +npymath_cfloat npymath_cexpf(npymath_cfloat z); +npymath_cfloat npymath_clogf(npymath_cfloat z); +npymath_cfloat npymath_cpowf(npymath_cfloat x, npymath_cfloat y); + +npymath_cfloat npymath_csqrtf(npymath_cfloat z); + +npymath_cfloat npymath_ccosf(npymath_cfloat z); +npymath_cfloat npymath_csinf(npymath_cfloat z); +npymath_cfloat npymath_ctanf(npymath_cfloat z); + +npymath_cfloat npymath_ccoshf(npymath_cfloat z); +npymath_cfloat npymath_csinhf(npymath_cfloat z); +npymath_cfloat npymath_ctanhf(npymath_cfloat z); + +npymath_cfloat npymath_cacosf(npymath_cfloat z); +npymath_cfloat npymath_casinf(npymath_cfloat z); +npymath_cfloat npymath_catanf(npymath_cfloat z); + +npymath_cfloat npymath_cacoshf(npymath_cfloat z); +npymath_cfloat npymath_casinhf(npymath_cfloat z); +npymath_cfloat npymath_catanhf(npymath_cfloat z); + + +/* + * Extended precision complex functions + */ +npymath_longdouble npymath_cabsl(npymath_clongdouble z); +npymath_longdouble npymath_cargl(npymath_clongdouble z); + +npymath_clongdouble npymath_cexpl(npymath_clongdouble z); +npymath_clongdouble npymath_clogl(npymath_clongdouble z); +npymath_clongdouble npymath_cpowl(npymath_clongdouble x, npymath_clongdouble y); + +npymath_clongdouble npymath_csqrtl(npymath_clongdouble z); + +npymath_clongdouble npymath_ccosl(npymath_clongdouble z); +npymath_clongdouble npymath_csinl(npymath_clongdouble z); +npymath_clongdouble npymath_ctanl(npymath_clongdouble z); + +npymath_clongdouble npymath_ccoshl(npymath_clongdouble z); +npymath_clongdouble npymath_csinhl(npymath_clongdouble z); +npymath_clongdouble npymath_ctanhl(npymath_clongdouble z); + +npymath_clongdouble npymath_cacosl(npymath_clongdouble z); +npymath_clongdouble npymath_casinl(npymath_clongdouble z); +npymath_clongdouble npymath_catanl(npymath_clongdouble z); + +npymath_clongdouble npymath_cacoshl(npymath_clongdouble z); +npymath_clongdouble npymath_casinhl(npymath_clongdouble z); +npymath_clongdouble npymath_catanhl(npymath_clongdouble z); + + +/* + * Functions that set the floating point error + * status word. + */ + +/* + * platform-dependent code translates floating point + * status to an integer sum of these values + */ +#define NPYMATH_FPE_DIVIDEBYZERO 1 +#define NPYMATH_FPE_OVERFLOW 2 +#define NPYMATH_FPE_UNDERFLOW 4 +#define NPYMATH_FPE_INVALID 8 + +int npymath_clear_floatstatus_barrier(char*); +int npymath_get_floatstatus_barrier(char*); +/* + * use caution with these - clang and gcc8.1 are known to reorder calls + * to this form of the function which can defeat the check. The _barrier + * form of the call is preferable, where the argument is + * (char*)&local_variable + */ +int npymath_clear_floatstatus(void); +int npymath_get_floatstatus(void); + +void npymath_set_floatstatus_divbyzero(void); +void npymath_set_floatstatus_overflow(void); +void npymath_set_floatstatus_underflow(void); +void npymath_set_floatstatus_invalid(void); + +#ifdef __cplusplus +} +#endif + +#ifdef NPYMATH_INLINE_MATH +#include "libnpymath/internal.h" +#endif + +#endif /* LIBNPYMATH_NPY_MATH_H_ */ diff --git a/src/npy_math_common.h b/src/include/libnpymath/npy_math_common.h similarity index 66% rename from src/npy_math_common.h rename to src/include/libnpymath/npy_math_common.h index da38e4b..c8f35c3 100644 --- a/src/npy_math_common.h +++ b/src/include/libnpymath/npy_math_common.h @@ -5,5 +5,5 @@ #include #include -#include "numpy/npymath/block.h" -#include "numpy/npy_math.h" +#include "libnpymath/block.h" +#include "libnpymath/npy_math.h" diff --git a/src/include/numpy/npymath/private.h b/src/include/libnpymath/private.h similarity index 98% rename from src/include/numpy/npymath/private.h rename to src/include/libnpymath/private.h index 66e791c..eb16159 100644 --- a/src/include/numpy/npymath/private.h +++ b/src/include/libnpymath/private.h @@ -15,8 +15,8 @@ * $FreeBSD$ */ -#ifndef _NPYMATH_PRIVATE_H_ -#define _NPYMATH_PRIVATE_H_ +#ifndef LIBNPYMATH_PRIVATE_H_ +#define LIBNPYMATH_PRIVATE_H_ #include #ifdef __cplusplus @@ -27,11 +27,11 @@ using std::isless; #include #endif -#include "numpy/npymath/fpmath.h" +#include "libnpymath/fpmath.h" -#include "numpy/npy_math.h" -#include "numpy/npymath/endian.h" -#include "numpy/npymath/common.h" +#include "libnpymath/npy_math.h" +#include "libnpymath/endian.h" +#include "libnpymath/common.h" /* * The original fdlibm code used statements like: @@ -478,4 +478,4 @@ do { \ #endif /* !NPYMATH_HAVE_LDOUBLE_DOUBLE_DOUBLE_* */ -#endif /* !_NPYMATH_PRIVATE_H_ */ +#endif /* LIBNPYMATH_PRIVATE_H_ */ diff --git a/src/include/numpy/npymath/utils.h b/src/include/libnpymath/utils.h similarity index 89% rename from src/include/numpy/npymath/utils.h rename to src/include/libnpymath/utils.h index 98a1054..65e6e70 100644 --- a/src/include/numpy/npymath/utils.h +++ b/src/include/libnpymath/utils.h @@ -1,5 +1,5 @@ -#ifndef NPYMATH_UTILS_H_ -#define NPYMATH_UTILS_H_ +#ifndef LIBNPYMATH_UTILS_H_ +#define LIBNPYMATH_UTILS_H_ #ifndef __COMP_NPYMATH_UNUSED #if defined(__GNUC__) @@ -22,4 +22,4 @@ #define NPYMATH_CAT_(a, b) NPYMATH_CAT__(a, b) #define NPYMATH_CAT(a, b) NPYMATH_CAT_(a, b) -#endif /* NPYMATH_UTILS_H_ */ +#endif /* LIBNPYMATH_UTILS_H_ */ diff --git a/src/include/numpy/halffloat.h b/src/include/numpy/halffloat.h deleted file mode 100644 index 99ac35d..0000000 --- a/src/include/numpy/halffloat.h +++ /dev/null @@ -1,70 +0,0 @@ -#ifndef NPYMATH_HALFFLOAT_H_ -#define NPYMATH_HALFFLOAT_H_ - -#include -#include "numpy/npy_math.h" - -#ifdef __cplusplus -extern "C" { -#endif - -/* - * Half-precision routines - */ - -/* Conversions */ -float npy_half_to_float(npymath_half h); -double npy_half_to_double(npymath_half h); -npymath_half npy_float_to_half(float f); -npymath_half npy_double_to_half(double d); -/* Comparisons */ -int npy_half_eq(npymath_half h1, npymath_half h2); -int npy_half_ne(npymath_half h1, npymath_half h2); -int npy_half_le(npymath_half h1, npymath_half h2); -int npy_half_lt(npymath_half h1, npymath_half h2); -int npy_half_ge(npymath_half h1, npymath_half h2); -int npy_half_gt(npymath_half h1, npymath_half h2); -/* faster *_nonan variants for when you know h1 and h2 are not NaN */ -int npy_half_eq_nonan(npymath_half h1, npymath_half h2); -int npy_half_lt_nonan(npymath_half h1, npymath_half h2); -int npy_half_le_nonan(npymath_half h1, npymath_half h2); -/* Miscellaneous functions */ -int npy_half_iszero(npymath_half h); -int npy_half_isnan(npymath_half h); -int npy_half_isinf(npymath_half h); -int npy_half_isfinite(npymath_half h); -int npy_half_signbit(npymath_half h); -npymath_half npy_half_copysign(npymath_half x, npymath_half y); -npymath_half npy_half_spacing(npymath_half h); -npymath_half npy_half_nextafter(npymath_half x, npymath_half y); -npymath_half npy_half_divmod(npymath_half x, npymath_half y, npymath_half *modulus); - -/* - * Half-precision constants - */ - -#define NPY_HALF_ZERO (0x0000u) -#define NPY_HALF_PZERO (0x0000u) -#define NPY_HALF_NZERO (0x8000u) -#define NPY_HALF_ONE (0x3c00u) -#define NPY_HALF_NEGONE (0xbc00u) -#define NPY_HALF_PINF (0x7c00u) -#define NPY_HALF_NINF (0xfc00u) -#define NPY_HALF_NAN (0x7e00u) - -#define NPY_MAX_HALF (0x7bffu) - -/* - * Bit-level conversions - */ - -npymath_uint16 npy_floatbits_to_halfbits(npymath_uint32 f); -npymath_uint16 npy_doublebits_to_halfbits(npymath_uint64 d); -npymath_uint32 npy_halfbits_to_floatbits(npymath_uint16 h); -npymath_uint64 npy_halfbits_to_doublebits(npymath_uint16 h); - -#ifdef __cplusplus -} -#endif - -#endif /* NPYMATH_HALFFLOAT_H_ */ diff --git a/src/include/numpy/npy_math.h b/src/include/numpy/npy_math.h deleted file mode 100644 index f3b53d3..0000000 --- a/src/include/numpy/npy_math.h +++ /dev/null @@ -1,517 +0,0 @@ -#ifndef NPY_MATH_H_ -#define NPY_MATH_H_ - -#include "numpy/npymath/common.h" - -#include - -/* By adding static inline specifiers to npy_math function definitions when - appropriate, compiler is given the opportunity to optimize */ -#if defined(NPYMATH_INLINE_MATH) -#define NPY_INPLACE static inline -#else -#define NPY_INPLACE -#endif - - -#ifdef __cplusplus -extern "C" { -#endif - -#define PyArray_MAX(a,b) (((a)>(b))?(a):(b)) -#define PyArray_MIN(a,b) (((a)<(b))?(a):(b)) - -/* - * NAN and INFINITY like macros (same behavior as glibc for NAN, same as C99 - * for INFINITY) - * - * XXX: I should test whether INFINITY and NAN are available on the platform - */ -static inline float __npy_inff(void) -{ - const union { npymath_uint32 __i; float __f;} __bint = {0x7f800000UL}; - return __bint.__f; -} - -static inline float __npy_nanf(void) -{ - const union { npymath_uint32 __i; float __f;} __bint = {0x7fc00000UL}; - return __bint.__f; -} - -static inline float __npy_pzerof(void) -{ - const union { npymath_uint32 __i; float __f;} __bint = {0x00000000UL}; - return __bint.__f; -} - -static inline float __npy_nzerof(void) -{ - const union { npymath_uint32 __i; float __f;} __bint = {0x80000000UL}; - return __bint.__f; -} - -#define NPY_INFINITYF __npy_inff() -#define NPY_NANF __npy_nanf() -#define NPY_PZEROF __npy_pzerof() -#define NPY_NZEROF __npy_nzerof() - -#define NPY_INFINITY ((npymath_double)NPY_INFINITYF) -#define NPY_NAN ((npymath_double)NPY_NANF) -#define NPY_PZERO ((npymath_double)NPY_PZEROF) -#define NPY_NZERO ((npymath_double)NPY_NZEROF) - -#define NPY_INFINITYL ((npymath_longdouble)NPY_INFINITYF) -#define NPY_NANL ((npymath_longdouble)NPY_NANF) -#define NPY_PZEROL ((npymath_longdouble)NPY_PZEROF) -#define NPY_NZEROL ((npymath_longdouble)NPY_NZEROF) - -/* - * Useful constants - */ -#define NPY_E 2.718281828459045235360287471352662498 /* e */ -#define NPY_LOG2E 1.442695040888963407359924681001892137 /* log_2 e */ -#define NPY_LOG10E 0.434294481903251827651128918916605082 /* log_10 e */ -#define NPY_LOGE2 0.693147180559945309417232121458176568 /* log_e 2 */ -#define NPY_LOGE10 2.302585092994045684017991454684364208 /* log_e 10 */ -#define NPY_PI 3.141592653589793238462643383279502884 /* pi */ -#define NPY_PI_2 1.570796326794896619231321691639751442 /* pi/2 */ -#define NPY_PI_4 0.785398163397448309615660845819875721 /* pi/4 */ -#define NPY_1_PI 0.318309886183790671537767526745028724 /* 1/pi */ -#define NPY_2_PI 0.636619772367581343075535053490057448 /* 2/pi */ -#define NPY_EULER 0.577215664901532860606512090082402431 /* Euler constant */ -#define NPY_SQRT2 1.414213562373095048801688724209698079 /* sqrt(2) */ -#define NPY_SQRT1_2 0.707106781186547524400844362104849039 /* 1/sqrt(2) */ - -#define NPY_Ef 2.718281828459045235360287471352662498F /* e */ -#define NPY_LOG2Ef 1.442695040888963407359924681001892137F /* log_2 e */ -#define NPY_LOG10Ef 0.434294481903251827651128918916605082F /* log_10 e */ -#define NPY_LOGE2f 0.693147180559945309417232121458176568F /* log_e 2 */ -#define NPY_LOGE10f 2.302585092994045684017991454684364208F /* log_e 10 */ -#define NPY_PIf 3.141592653589793238462643383279502884F /* pi */ -#define NPY_PI_2f 1.570796326794896619231321691639751442F /* pi/2 */ -#define NPY_PI_4f 0.785398163397448309615660845819875721F /* pi/4 */ -#define NPY_1_PIf 0.318309886183790671537767526745028724F /* 1/pi */ -#define NPY_2_PIf 0.636619772367581343075535053490057448F /* 2/pi */ -#define NPY_EULERf 0.577215664901532860606512090082402431F /* Euler constant */ -#define NPY_SQRT2f 1.414213562373095048801688724209698079F /* sqrt(2) */ -#define NPY_SQRT1_2f 0.707106781186547524400844362104849039F /* 1/sqrt(2) */ - -#define NPY_El 2.718281828459045235360287471352662498L /* e */ -#define NPY_LOG2El 1.442695040888963407359924681001892137L /* log_2 e */ -#define NPY_LOG10El 0.434294481903251827651128918916605082L /* log_10 e */ -#define NPY_LOGE2l 0.693147180559945309417232121458176568L /* log_e 2 */ -#define NPY_LOGE10l 2.302585092994045684017991454684364208L /* log_e 10 */ -#define NPY_PIl 3.141592653589793238462643383279502884L /* pi */ -#define NPY_PI_2l 1.570796326794896619231321691639751442L /* pi/2 */ -#define NPY_PI_4l 0.785398163397448309615660845819875721L /* pi/4 */ -#define NPY_1_PIl 0.318309886183790671537767526745028724L /* 1/pi */ -#define NPY_2_PIl 0.636619772367581343075535053490057448L /* 2/pi */ -#define NPY_EULERl 0.577215664901532860606512090082402431L /* Euler constant */ -#define NPY_SQRT2l 1.414213562373095048801688724209698079L /* sqrt(2) */ -#define NPY_SQRT1_2l 0.707106781186547524400844362104849039L /* 1/sqrt(2) */ - -/* - * Integer functions. - */ -NPY_INPLACE npymath_uint npy_gcdu(npymath_uint a, npymath_uint b); -NPY_INPLACE npymath_uint npy_lcmu(npymath_uint a, npymath_uint b); -NPY_INPLACE npymath_ulong npy_gcdul(npymath_ulong a, npymath_ulong b); -NPY_INPLACE npymath_ulong npy_lcmul(npymath_ulong a, npymath_ulong b); -NPY_INPLACE npymath_ulonglong npy_gcdull(npymath_ulonglong a, npymath_ulonglong b); -NPY_INPLACE npymath_ulonglong npy_lcmull(npymath_ulonglong a, npymath_ulonglong b); - -NPY_INPLACE npymath_int npy_gcd(npymath_int a, npymath_int b); -NPY_INPLACE npymath_int npy_lcm(npymath_int a, npymath_int b); -NPY_INPLACE npymath_long npy_gcdl(npymath_long a, npymath_long b); -NPY_INPLACE npymath_long npy_lcml(npymath_long a, npymath_long b); -NPY_INPLACE npymath_longlong npy_gcdll(npymath_longlong a, npymath_longlong b); -NPY_INPLACE npymath_longlong npy_lcmll(npymath_longlong a, npymath_longlong b); - -NPY_INPLACE npymath_ubyte npy_rshiftuhh(npymath_ubyte a, npymath_ubyte b); -NPY_INPLACE npymath_ubyte npy_lshiftuhh(npymath_ubyte a, npymath_ubyte b); -NPY_INPLACE npymath_ushort npy_rshiftuh(npymath_ushort a, npymath_ushort b); -NPY_INPLACE npymath_ushort npy_lshiftuh(npymath_ushort a, npymath_ushort b); -NPY_INPLACE npymath_uint npy_rshiftu(npymath_uint a, npymath_uint b); -NPY_INPLACE npymath_uint npy_lshiftu(npymath_uint a, npymath_uint b); -NPY_INPLACE npymath_ulong npy_rshiftul(npymath_ulong a, npymath_ulong b); -NPY_INPLACE npymath_ulong npy_lshiftul(npymath_ulong a, npymath_ulong b); -NPY_INPLACE npymath_ulonglong npy_rshiftull(npymath_ulonglong a, npymath_ulonglong b); -NPY_INPLACE npymath_ulonglong npy_lshiftull(npymath_ulonglong a, npymath_ulonglong b); - -NPY_INPLACE npymath_byte npy_rshifthh(npymath_byte a, npymath_byte b); -NPY_INPLACE npymath_byte npy_lshifthh(npymath_byte a, npymath_byte b); -NPY_INPLACE npymath_short npy_rshifth(npymath_short a, npymath_short b); -NPY_INPLACE npymath_short npy_lshifth(npymath_short a, npymath_short b); -NPY_INPLACE npymath_int npy_rshift(npymath_int a, npymath_int b); -NPY_INPLACE npymath_int npy_lshift(npymath_int a, npymath_int b); -NPY_INPLACE npymath_long npy_rshiftl(npymath_long a, npymath_long b); -NPY_INPLACE npymath_long npy_lshiftl(npymath_long a, npymath_long b); -NPY_INPLACE npymath_longlong npy_rshiftll(npymath_longlong a, npymath_longlong b); -NPY_INPLACE npymath_longlong npy_lshiftll(npymath_longlong a, npymath_longlong b); - -NPY_INPLACE uint8_t npy_popcountuhh(npymath_ubyte a); -NPY_INPLACE uint8_t npy_popcountuh(npymath_ushort a); -NPY_INPLACE uint8_t npy_popcountu(npymath_uint a); -NPY_INPLACE uint8_t npy_popcountul(npymath_ulong a); -NPY_INPLACE uint8_t npy_popcountull(npymath_ulonglong a); -NPY_INPLACE uint8_t npy_popcounthh(npymath_byte a); -NPY_INPLACE uint8_t npy_popcounth(npymath_short a); -NPY_INPLACE uint8_t npy_popcount(npymath_int a); -NPY_INPLACE uint8_t npy_popcountl(npymath_long a); -NPY_INPLACE uint8_t npy_popcountll(npymath_longlong a); - -/* - * C99 double math funcs that need fixups or are blocklist-able - */ -NPY_INPLACE double npy_sin(double x); -NPY_INPLACE double npy_cos(double x); -NPY_INPLACE double npy_tan(double x); -NPY_INPLACE double npy_hypot(double x, double y); -NPY_INPLACE double npy_log2(double x); -NPY_INPLACE double npy_atan2(double x, double y); - -/* Mandatory C99 double math funcs, no blocklisting or fixups */ -/* defined for legacy reasons, should be deprecated at some point */ -#define npy_sinh sinh -#define npy_cosh cosh -#define npy_tanh tanh -#define npy_asin asin -#define npy_acos acos -#define npy_atan atan -#define npy_log log -#define npy_log10 log10 -#define npy_cbrt cbrt -#define npy_fabs fabs -#define npy_ceil ceil -#define npy_fmod fmod -#define npy_floor floor -#define npy_expm1 expm1 -#define npy_log1p log1p -#define npy_acosh acosh -#define npy_asinh asinh -#define npy_atanh atanh -#define npy_rint rint -#define npy_trunc trunc -#define npy_exp2 exp2 -#define npy_frexp frexp -#define npy_ldexp ldexp -#define npy_copysign copysign -#define npy_exp exp -#define npy_sqrt sqrt -#define npy_pow pow -#define npy_modf modf -#define npy_nextafter nextafter - -double npy_spacing(double x); - -/* - * IEEE 754 fpu handling - */ - -/* use builtins to avoid function calls in tight loops - * only available if npymath_config.h is available (= numpys own build) */ -#ifdef NPYMATH_HAVE___BUILTIN_ISNAN - #define npy_isnan(x) __builtin_isnan(x) -#else - #define npy_isnan(x) isnan(x) -#endif - - -/* only available if npymath_config.h is available (= numpys own build) */ -#ifdef NPYMATH_HAVE___BUILTIN_ISFINITE - #define npy_isfinite(x) __builtin_isfinite(x) -#else - #define npy_isfinite(x) isfinite((x)) -#endif - -/* only available if npymath_config.h is available (= numpys own build) */ -#ifdef NPYMATH_HAVE___BUILTIN_ISINF - #define npy_isinf(x) __builtin_isinf(x) -#else - #define npy_isinf(x) isinf((x)) -#endif - -#define npy_signbit(x) signbit((x)) - -/* - * float C99 math funcs that need fixups or are blocklist-able - */ -NPY_INPLACE float npy_sinf(float x); -NPY_INPLACE float npy_cosf(float x); -NPY_INPLACE float npy_tanf(float x); -NPY_INPLACE float npy_expf(float x); -NPY_INPLACE float npy_sqrtf(float x); -NPY_INPLACE float npy_hypotf(float x, float y); -NPY_INPLACE float npy_log2f(float x); -NPY_INPLACE float npy_atan2f(float x, float y); -NPY_INPLACE float npy_powf(float x, float y); -NPY_INPLACE float npy_modff(float x, float* y); - -/* Mandatory C99 float math funcs, no blocklisting or fixups */ -/* defined for legacy reasons, should be deprecated at some point */ - -#define npy_sinhf sinhf -#define npy_coshf coshf -#define npy_tanhf tanhf -#define npy_asinf asinf -#define npy_acosf acosf -#define npy_atanf atanf -#define npy_logf logf -#define npy_log10f log10f -#define npy_cbrtf cbrtf -#define npy_fabsf fabsf -#define npy_ceilf ceilf -#define npy_fmodf fmodf -#define npy_floorf floorf -#define npy_expm1f expm1f -#define npy_log1pf log1pf -#define npy_asinhf asinhf -#define npy_acoshf acoshf -#define npy_atanhf atanhf -#define npy_rintf rintf -#define npy_truncf truncf -#define npy_exp2f exp2f -#define npy_frexpf frexpf -#define npy_ldexpf ldexpf -#define npy_copysignf copysignf -#define npy_nextafterf nextafterf - -float npy_spacingf(float x); - -/* - * long double C99 double math funcs that need fixups or are blocklist-able - */ -NPY_INPLACE npymath_longdouble npy_sinl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_cosl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_tanl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_expl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_sqrtl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_hypotl(npymath_longdouble x, npymath_longdouble y); -NPY_INPLACE npymath_longdouble npy_log2l(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_atan2l(npymath_longdouble x, npymath_longdouble y); -NPY_INPLACE npymath_longdouble npy_powl(npymath_longdouble x, npymath_longdouble y); -NPY_INPLACE npymath_longdouble npy_modfl(npymath_longdouble x, npymath_longdouble* y); - -/* Mandatory C99 double math funcs, no blocklisting or fixups */ -/* defined for legacy reasons, should be deprecated at some point */ -#define npy_sinhl sinhl -#define npy_coshl coshl -#define npy_tanhl tanhl -#define npy_fabsl fabsl -#define npy_floorl floorl -#define npy_ceill ceill -#define npy_rintl rintl -#define npy_truncl truncl -#define npy_cbrtl cbrtl -#define npy_log10l log10l -#define npy_logl logl -#define npy_expm1l expm1l -#define npy_asinl asinl -#define npy_acosl acosl -#define npy_atanl atanl -#define npy_asinhl asinhl -#define npy_acoshl acoshl -#define npy_atanhl atanhl -#define npy_log1pl log1pl -#define npy_exp2l exp2l -#define npy_fmodl fmodl -#define npy_frexpl frexpl -#define npy_ldexpl ldexpl -#define npy_copysignl copysignl -#define npy_nextafterl nextafterl - -npymath_longdouble npy_spacingl(npymath_longdouble x); - -/* - * Non standard functions - */ -NPY_INPLACE double npy_deg2rad(double x); -NPY_INPLACE double npy_rad2deg(double x); -NPY_INPLACE double npy_logaddexp(double x, double y); -NPY_INPLACE double npy_logaddexp2(double x, double y); -NPY_INPLACE double npy_divmod(double x, double y, double *modulus); -NPY_INPLACE double npy_heaviside(double x, double h0); - -NPY_INPLACE float npy_deg2radf(float x); -NPY_INPLACE float npy_rad2degf(float x); -NPY_INPLACE float npy_logaddexpf(float x, float y); -NPY_INPLACE float npy_logaddexp2f(float x, float y); -NPY_INPLACE float npy_divmodf(float x, float y, float *modulus); -NPY_INPLACE float npy_heavisidef(float x, float h0); - -NPY_INPLACE npymath_longdouble npy_deg2radl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_rad2degl(npymath_longdouble x); -NPY_INPLACE npymath_longdouble npy_logaddexpl(npymath_longdouble x, npymath_longdouble y); -NPY_INPLACE npymath_longdouble npy_logaddexp2l(npymath_longdouble x, npymath_longdouble y); -NPY_INPLACE npymath_longdouble npy_divmodl(npymath_longdouble x, npymath_longdouble y, - npymath_longdouble *modulus); -NPY_INPLACE npymath_longdouble npy_heavisidel(npymath_longdouble x, npymath_longdouble h0); - -#define npy_degrees npy_rad2deg -#define npy_degreesf npy_rad2degf -#define npy_degreesl npy_rad2degl - -#define npy_radians npy_deg2rad -#define npy_radiansf npy_deg2radf -#define npy_radiansl npy_deg2radl - -/* - * Complex declarations - */ - -npymath_double npy_creal(const npymath_cdouble *z); -void npy_csetreal(npymath_cdouble *z, npymath_double r); -npymath_double npy_cimag(const npymath_cdouble *z); -void npy_csetimag(npymath_cdouble *z, npymath_double i); -npymath_float npy_crealf(const npymath_cfloat *z); -void npy_csetrealf(npymath_cfloat *z, npymath_float r); -npymath_float npy_cimagf(const npymath_cfloat *z); -void npy_csetimagf(npymath_cfloat *z, npymath_float i); -npymath_longdouble npy_creall(const npymath_clongdouble *z); -void npy_csetreall(npymath_clongdouble *z, npymath_longdouble r); -npymath_longdouble npy_cimagl(const npymath_clongdouble *z); -void npy_csetimagl(npymath_clongdouble *z, npymath_longdouble i); - -#define NPY_CREAL(z) npy_creal(z) -#define NPY_CREALF(z) npy_crealf(z) -#define NPY_CREALL(z) npy_creall(z) -#define NPY_CIMAG(z) npy_cimag(z) -#define NPY_CIMAGF(z) npy_cimagf(z) -#define NPY_CIMAGL(z) npy_cimagl(z) -#define NPY_CSETREAL(z, r) npy_csetreal(z, r) -#define NPY_CSETIMAG(z, i) npy_csetimag(z, i) -#define NPY_CSETREALF(z, r) npy_csetrealf(z, r) -#define NPY_CSETIMAGF(z, i) npy_csetimagf(z, i) -#define NPY_CSETREALL(z, r) npy_csetreall(z, r) -#define NPY_CSETIMAGL(z, i) npy_csetimagl(z, i) - -npymath_cdouble npy_cpack(npymath_double x, npymath_double y); -npymath_cfloat npy_cpackf(npymath_float x, npymath_float y); -npymath_clongdouble npy_cpackl(npymath_longdouble x, npymath_longdouble y); - -/* - * Double precision complex functions - */ -double npy_cabs(npymath_cdouble z); -double npy_carg(npymath_cdouble z); - -npymath_cdouble npy_cexp(npymath_cdouble z); -npymath_cdouble npy_clog(npymath_cdouble z); -npymath_cdouble npy_cpow(npymath_cdouble x, npymath_cdouble y); - -npymath_cdouble npy_csqrt(npymath_cdouble z); - -npymath_cdouble npy_ccos(npymath_cdouble z); -npymath_cdouble npy_csin(npymath_cdouble z); -npymath_cdouble npy_ctan(npymath_cdouble z); - -npymath_cdouble npy_ccosh(npymath_cdouble z); -npymath_cdouble npy_csinh(npymath_cdouble z); -npymath_cdouble npy_ctanh(npymath_cdouble z); - -npymath_cdouble npy_cacos(npymath_cdouble z); -npymath_cdouble npy_casin(npymath_cdouble z); -npymath_cdouble npy_catan(npymath_cdouble z); - -npymath_cdouble npy_cacosh(npymath_cdouble z); -npymath_cdouble npy_casinh(npymath_cdouble z); -npymath_cdouble npy_catanh(npymath_cdouble z); - -/* - * Single precision complex functions - */ -float npy_cabsf(npymath_cfloat z); -float npy_cargf(npymath_cfloat z); - -npymath_cfloat npy_cexpf(npymath_cfloat z); -npymath_cfloat npy_clogf(npymath_cfloat z); -npymath_cfloat npy_cpowf(npymath_cfloat x, npymath_cfloat y); - -npymath_cfloat npy_csqrtf(npymath_cfloat z); - -npymath_cfloat npy_ccosf(npymath_cfloat z); -npymath_cfloat npy_csinf(npymath_cfloat z); -npymath_cfloat npy_ctanf(npymath_cfloat z); - -npymath_cfloat npy_ccoshf(npymath_cfloat z); -npymath_cfloat npy_csinhf(npymath_cfloat z); -npymath_cfloat npy_ctanhf(npymath_cfloat z); - -npymath_cfloat npy_cacosf(npymath_cfloat z); -npymath_cfloat npy_casinf(npymath_cfloat z); -npymath_cfloat npy_catanf(npymath_cfloat z); - -npymath_cfloat npy_cacoshf(npymath_cfloat z); -npymath_cfloat npy_casinhf(npymath_cfloat z); -npymath_cfloat npy_catanhf(npymath_cfloat z); - - -/* - * Extended precision complex functions - */ -npymath_longdouble npy_cabsl(npymath_clongdouble z); -npymath_longdouble npy_cargl(npymath_clongdouble z); - -npymath_clongdouble npy_cexpl(npymath_clongdouble z); -npymath_clongdouble npy_clogl(npymath_clongdouble z); -npymath_clongdouble npy_cpowl(npymath_clongdouble x, npymath_clongdouble y); - -npymath_clongdouble npy_csqrtl(npymath_clongdouble z); - -npymath_clongdouble npy_ccosl(npymath_clongdouble z); -npymath_clongdouble npy_csinl(npymath_clongdouble z); -npymath_clongdouble npy_ctanl(npymath_clongdouble z); - -npymath_clongdouble npy_ccoshl(npymath_clongdouble z); -npymath_clongdouble npy_csinhl(npymath_clongdouble z); -npymath_clongdouble npy_ctanhl(npymath_clongdouble z); - -npymath_clongdouble npy_cacosl(npymath_clongdouble z); -npymath_clongdouble npy_casinl(npymath_clongdouble z); -npymath_clongdouble npy_catanl(npymath_clongdouble z); - -npymath_clongdouble npy_cacoshl(npymath_clongdouble z); -npymath_clongdouble npy_casinhl(npymath_clongdouble z); -npymath_clongdouble npy_catanhl(npymath_clongdouble z); - - -/* - * Functions that set the floating point error - * status word. - */ - -/* - * platform-dependent code translates floating point - * status to an integer sum of these values - */ -#define NPY_FPE_DIVIDEBYZERO 1 -#define NPY_FPE_OVERFLOW 2 -#define NPY_FPE_UNDERFLOW 4 -#define NPY_FPE_INVALID 8 - -int npy_clear_floatstatus_barrier(char*); -int npy_get_floatstatus_barrier(char*); -/* - * use caution with these - clang and gcc8.1 are known to reorder calls - * to this form of the function which can defeat the check. The _barrier - * form of the call is preferable, where the argument is - * (char*)&local_variable - */ -int npy_clear_floatstatus(void); -int npy_get_floatstatus(void); - -void npy_set_floatstatus_divbyzero(void); -void npy_set_floatstatus_overflow(void); -void npy_set_floatstatus_underflow(void); -void npy_set_floatstatus_invalid(void); - -#ifdef __cplusplus -} -#endif - -#ifdef NPYMATH_INLINE_MATH -#include "numpy/npymath/internal.h" -#endif - -#endif /* NPY_MATH_H_ */ diff --git a/src/meson.build b/src/meson.build index 53fa560..ef828ce 100644 --- a/src/meson.build +++ b/src/meson.build @@ -73,9 +73,6 @@ foreach symbol_type: complex_types_to_check cdata.set(symbol_type[0], cc.sizeof(symbol_type[1], prefix: '#include ')) endforeach -# Mandatory functions: if not found, fail the build -# Some of these can still be blocklisted if the C99 implementation -# is buggy, see numpy/_core/src/common/npy_config.h mandatory_math_funcs = [ 'sin', 'cos', 'tan', 'sinh', 'cosh', 'tanh', 'fabs', 'floor', 'ceil', 'sqrt', 'log10', 'log', 'exp', 'asin', @@ -251,7 +248,7 @@ if cc.has_header('sys/endian.h') cdata.set10('NPYMATH_HAVE_SYS_ENDIAN_H', true) endif -subdir('include/numpy/npymath') +subdir('include/libnpymath') # Build npymath static library # ---------------------------- @@ -300,26 +297,20 @@ if get_option('install_lib') py.install_sources( [ - 'include/numpy/npy_math.h', - 'include/numpy/halffloat.h', - ], - subdir: get_option('headers_dir'), - ) - - py.install_sources( - [ - 'include/numpy/npymath/common.h', - 'include/numpy/npymath/config.h', + 'include/libnpymath/npy_math.h', + 'include/libnpymath/halffloat.h', + 'include/libnpymath/common.h', + 'include/libnpymath/config.h', _config_h, meson_config_h, - 'include/numpy/npymath/block.h', - 'include/numpy/npymath/cpu.h', - 'include/numpy/npymath/endian.h', - 'include/numpy/npymath/fpmath.h', - 'include/numpy/npymath/utils.h', - 'include/numpy/npymath/private.h', + 'include/libnpymath/block.h', + 'include/libnpymath/cpu.h', + 'include/libnpymath/endian.h', + 'include/libnpymath/fpmath.h', + 'include/libnpymath/utils.h', + 'include/libnpymath/private.h', ], - subdir: get_option('headers_dir') / 'npymath', + subdir: get_option('headers_dir') / 'libnpymath', ) else libnpymath = library( diff --git a/src/npy_math.c b/src/npy_math.c index 8583972..b6c5c96 100644 --- a/src/npy_math.c +++ b/src/npy_math.c @@ -6,4 +6,4 @@ */ #define NPYMATH_INLINE_MATH 0 -#include "numpy/npymath/internal.h" +#include "libnpymath/internal.h" diff --git a/src/npy_math_complex.c.src b/src/npy_math_complex.c.src index 8241477..d346e42 100644 --- a/src/npy_math_complex.c.src +++ b/src/npy_math_complex.c.src @@ -31,9 +31,9 @@ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF * SUCH DAMAGE. */ -#include "npy_math_common.h" -#include "numpy/npymath/private.h" -#include "numpy/npymath/utils.h" +#include "libnpymath/npy_math_common.h" +#include "libnpymath/private.h" +#include "libnpymath/utils.h" /* * Hack inherited from BSD, the intent is to set the FPU inexact @@ -82,11 +82,11 @@ static inline cmul@c@(@ctype@ a, @ctype@ b) { @type@ ar, ai, br, bi; - ar = npy_creal@c@(&a); - ai = npy_cimag@c@(&a); - br = npy_creal@c@(&b); - bi = npy_cimag@c@(&b); - return npy_cpack@c@(ar*br - ai*bi, ar*bi + ai*br); + ar = npymath_creal@c@(&a); + ai = npymath_cimag@c@(&a); + br = npymath_creal@c@(&b); + bi = npymath_cimag@c@(&b); + return npymath_cpack@c@(ar*br - ai*bi, ar*bi + ai*br); } static inline @@ -94,28 +94,28 @@ static inline cdiv@c@(@ctype@ a, @ctype@ b) { @type@ ar, ai, br, bi, abs_br, abs_bi; - ar = npy_creal@c@(&a); - ai = npy_cimag@c@(&a); - br = npy_creal@c@(&b); - bi = npy_cimag@c@(&b); - abs_br = npy_fabs@c@(br); - abs_bi = npy_fabs@c@(bi); + ar = npymath_creal@c@(&a); + ai = npymath_cimag@c@(&a); + br = npymath_creal@c@(&b); + bi = npymath_cimag@c@(&b); + abs_br = npymath_fabs@c@(br); + abs_bi = npymath_fabs@c@(bi); if (abs_br >= abs_bi) { if (abs_br == 0 && abs_bi == 0) { /* divide by zeros should yield a complex inf or nan */ - return npy_cpack@c@(ar/abs_br, ai/abs_bi); + return npymath_cpack@c@(ar/abs_br, ai/abs_bi); } else { @type@ rat = bi/br; @type@ scl = 1.0@C@/(br+bi*rat); - return npy_cpack@c@((ar + ai*rat)*scl, (ai - ar*rat)*scl); + return npymath_cpack@c@((ar + ai*rat)*scl, (ai - ar*rat)*scl); } } else { @type@ rat = br/bi; @type@ scl = 1.0@C@/(bi + br*rat); - return npy_cpack@c@((ar*rat + ai)*scl, (ai*rat - ar)*scl); + return npymath_cpack@c@((ar*rat + ai)*scl, (ai*rat - ar)*scl); } } @@ -124,51 +124,51 @@ cdiv@c@(@ctype@ a, @ctype@ b) *=========================================================*/ @type@ -npy_creal@c@(const @ctype@ *z) +npymath_creal@c@(const @ctype@ *z) { return ((@type@ *) z)[0]; } void -npy_csetreal@c@(@ctype@ *z, const @type@ r) +npymath_csetreal@c@(@ctype@ *z, const @type@ r) { ((@type@ *) z)[0] = r; } @type@ -npy_cimag@c@(const @ctype@ *z) +npymath_cimag@c@(const @ctype@ *z) { return ((@type@ *) z)[1]; } void -npy_csetimag@c@(@ctype@ *z, const @type@ i) +npymath_csetimag@c@(@ctype@ *z, const @type@ i) { ((@type@ *) z)[1] = i; } @ctype@ -npy_cpack@c@(@type@ x, @type@ y) +npymath_cpack@c@(@type@ x, @type@ y) { @ctype@ z; - npy_csetreal@c@(&z, x); - npy_csetimag@c@(&z, y); + npymath_csetreal@c@(&z, x); + npymath_csetimag@c@(&z, y); return z; } #ifndef NPYMATH_HAVE_CABS@C@ @type@ -npy_cabs@c@(@ctype@ z) +npymath_cabs@c@(@ctype@ z) { - return npy_hypot@c@(npy_creal@c@(&z), npy_cimag@c@(&z)); + return npymath_hypot@c@(npymath_creal@c@(&z), npymath_cimag@c@(&z)); } #endif #ifndef NPYMATH_HAVE_CARG@C@ @type@ -npy_carg@c@(@ctype@ z) +npymath_carg@c@(@ctype@ z) { - return npy_atan2@c@(npy_cimag@c@(&z), npy_creal@c@(&z)); + return npymath_atan2@c@(npymath_cimag@c@(&z), npymath_creal@c@(&z)); } #endif @@ -202,7 +202,7 @@ npy_carg@c@(@ctype@ z) static @ctype@ -_npy_scaled_cexp@c@(@type@ x, @type@ y, npymath_int expt) +_npymath_scaled_cexp@c@(@type@ x, @type@ y, npymath_int expt) { #if @precision@ == 1 const npymath_int k = 235; @@ -213,17 +213,17 @@ _npy_scaled_cexp@c@(@type@ x, @type@ y, npymath_int expt) #if @precision@ == 3 const npymath_int k = 19547; #endif - const @type@ kln2 = k * NPY_LOGE2@c@; + const @type@ kln2 = k * NPYMATH_LOGE2@c@; @type@ mant, mantcos, mantsin; npymath_int ex, excos, exsin; - mant = npy_frexp@c@(npy_exp@c@(x - kln2), &ex); - mantcos = npy_frexp@c@(npy_cos@c@(y), &excos); - mantsin = npy_frexp@c@(npy_sin@c@(y), &exsin); + mant = npymath_frexp@c@(npymath_exp@c@(x - kln2), &ex); + mantcos = npymath_frexp@c@(npymath_cos@c@(y), &excos); + mantsin = npymath_frexp@c@(npymath_sin@c@(y), &exsin); expt += ex + k; - return npy_cpack@c@( npy_ldexp@c@(mant * mantcos, expt + excos), - npy_ldexp@c@(mant * mantsin, expt + exsin)); + return npymath_cpack@c@( npymath_ldexp@c@(mant * mantcos, expt + excos), + npymath_ldexp@c@(mant * mantsin, expt + exsin)); } #endif @@ -231,72 +231,72 @@ _npy_scaled_cexp@c@(@type@ x, @type@ y, npymath_int expt) #ifndef NPYMATH_HAVE_CEXP@C@ @ctype@ -npy_cexp@c@(@ctype@ z) +npymath_cexp@c@(@ctype@ z) { @type@ x, c, s; @type@ r, i; @ctype@ ret; - r = npy_creal@c@(&z); - i = npy_cimag@c@(&z); + r = npymath_creal@c@(&z); + i = npymath_cimag@c@(&z); - if (npy_isfinite(r)) { + if (npymath_isfinite(r)) { if (r >= SCALED_CEXP_LOWER@C@ && r <= SCALED_CEXP_UPPER@C@) { - ret = _npy_scaled_cexp@c@(r, i, 0); + ret = _npymath_scaled_cexp@c@(r, i, 0); } else { - x = npy_exp@c@(r); + x = npymath_exp@c@(r); - c = npy_cos@c@(i); - s = npy_sin@c@(i); + c = npymath_cos@c@(i); + s = npymath_sin@c@(i); - if (npy_isfinite(i)) { - ret = npy_cpack@c@(x * c, x * s); + if (npymath_isfinite(i)) { + ret = npymath_cpack@c@(x * c, x * s); } else { - ret = npy_cpack@c@(NPY_NAN@C@, npy_copysign@c@(NPY_NAN@C@, i)); + ret = npymath_cpack@c@(NPYMATH_NAN@C@, npymath_copysign@c@(NPYMATH_NAN@C@, i)); } } } - else if (npy_isnan(r)) { + else if (npymath_isnan(r)) { /* r is nan */ if (i == 0) { ret = z; } else { - ret = npy_cpack@c@(r, npy_copysign@c@(NPY_NAN@C@, i)); + ret = npymath_cpack@c@(r, npymath_copysign@c@(NPYMATH_NAN@C@, i)); } } else { /* r is +- inf */ if (r > 0) { if (i == 0) { - ret = npy_cpack@c@(r, i); + ret = npymath_cpack@c@(r, i); } - else if (npy_isfinite(i)) { - c = npy_cos@c@(i); - s = npy_sin@c@(i); + else if (npymath_isfinite(i)) { + c = npymath_cos@c@(i); + s = npymath_sin@c@(i); - ret = npy_cpack@c@(r * c, r * s); + ret = npymath_cpack@c@(r * c, r * s); } else { /* x = +inf, y = +-inf | nan */ - npy_set_floatstatus_invalid(); - ret = npy_cpack@c@(r, NPY_NAN@C@); + npymath_set_floatstatus_invalid(); + ret = npymath_cpack@c@(r, NPYMATH_NAN@C@); } } else { - if (npy_isfinite(i)) { - x = npy_exp@c@(r); - c = npy_cos@c@(i); - s = npy_sin@c@(i); + if (npymath_isfinite(i)) { + x = npymath_exp@c@(r); + c = npymath_cos@c@(i); + s = npymath_sin@c@(i); - ret = npy_cpack@c@(x * c, x * s); + ret = npymath_cpack@c@(x * c, x * s); } else { /* x = -inf, y = nan | +i inf */ - ret = npy_cpack@c@(0, 0); + ret = npymath_cpack@c@(0, 0); } } } @@ -334,86 +334,86 @@ npy_cexp@c@(@ctype@ z) * is fine here. */ @ctype@ -npy_clog@c@(@ctype@ z) +npymath_clog@c@(@ctype@ z) { - @type@ ax = npy_fabs@c@(npy_creal@c@(&z)); - @type@ ay = npy_fabs@c@(npy_cimag@c@(&z)); + @type@ ax = npymath_fabs@c@(npymath_creal@c@(&z)); + @type@ ay = npymath_fabs@c@(npymath_cimag@c@(&z)); @type@ rr, ri; if (ax > @TMAX@/4 || ay > @TMAX@/4) { - rr = npy_log@c@(npy_hypot@c@(ax/2, ay/2)) + NPY_LOGE2@c@; + rr = npymath_log@c@(npymath_hypot@c@(ax/2, ay/2)) + NPYMATH_LOGE2@c@; } else if (ax < @TMIN@ && ay < @TMIN@) { if (ax > 0 || ay > 0) { /* catch cases where hypot(ax, ay) is subnormal */ - rr = npy_log@c@(npy_hypot@c@(npy_ldexp@c@(ax, @TMANT_DIG@), - npy_ldexp@c@(ay, @TMANT_DIG@))) - @TMANT_DIG@*NPY_LOGE2@c@; + rr = npymath_log@c@(npymath_hypot@c@(npymath_ldexp@c@(ax, @TMANT_DIG@), + npymath_ldexp@c@(ay, @TMANT_DIG@))) - @TMANT_DIG@*NPYMATH_LOGE2@c@; } else { /* log(+/-0 +/- 0i) */ /* raise divide-by-zero floating point exception */ - rr = -1.0@c@ / npy_creal@c@(&z); - rr = npy_copysign@c@(rr, -1); - ri = npy_carg@c@(z); - return npy_cpack@c@(rr, ri); + rr = -1.0@c@ / npymath_creal@c@(&z); + rr = npymath_copysign@c@(rr, -1); + ri = npymath_carg@c@(z); + return npymath_cpack@c@(rr, ri); } } else { - @type@ h = npy_hypot@c@(ax, ay); + @type@ h = npymath_hypot@c@(ax, ay); if (0.71 <= h && h <= 1.73) { @type@ am = ax > ay ? ax : ay; /* max(ax, ay) */ @type@ an = ax > ay ? ay : ax; /* min(ax, ay) */ - rr = npy_log1p@c@((am-1)*(am+1)+an*an)/2; + rr = npymath_log1p@c@((am-1)*(am+1)+an*an)/2; } else { - rr = npy_log@c@(h); + rr = npymath_log@c@(h); } } - ri = npy_carg@c@(z); + ri = npymath_carg@c@(z); - return npy_cpack@c@(rr, ri); + return npymath_cpack@c@(rr, ri); } #endif #ifndef NPYMATH_HAVE_CSQRT@C@ /* We risk spurious overflow for components >= DBL_MAX / (1 + sqrt(2)). */ -#define THRESH (@TMAX@ / (1 + NPY_SQRT2@c@)) +#define THRESH (@TMAX@ / (1 + NPYMATH_SQRT2@c@)) @ctype@ -npy_csqrt@c@(@ctype@ z) +npymath_csqrt@c@(@ctype@ z) { @ctype@ result; @type@ a, b; @type@ t; int scale; - a = npy_creal@c@(&z); - b = npy_cimag@c@(&z); + a = npymath_creal@c@(&z); + b = npymath_cimag@c@(&z); /* Handle special cases. */ if (a == 0 && b == 0) { - return (npy_cpack@c@(0, b)); + return (npymath_cpack@c@(0, b)); } - if (npy_isinf(b)) { - return (npy_cpack@c@(NPY_INFINITY@C@, b)); + if (npymath_isinf(b)) { + return (npymath_cpack@c@(NPYMATH_INFINITY@C@, b)); } - if (npy_isnan(a)) { + if (npymath_isnan(a)) { t = (b - b) / (b - b); /* raise invalid if b is not a NaN */ - return (npy_cpack@c@(a, t)); /* return NaN + NaN i */ + return (npymath_cpack@c@(a, t)); /* return NaN + NaN i */ } - if (npy_isinf(a)) { + if (npymath_isinf(a)) { /* * csqrt(inf + NaN i) = inf + NaN i * csqrt(inf + y i) = inf + 0 i * csqrt(-inf + NaN i) = NaN +- inf i * csqrt(-inf + y i) = 0 + inf i */ - if (npy_signbit(a)) { - return (npy_cpack@c@(npy_fabs@c@(b - b), npy_copysign@c@(a, b))); + if (npymath_signbit(a)) { + return (npymath_cpack@c@(npymath_fabs@c@(b - b), npymath_copysign@c@(a, b))); } else { - return (npy_cpack@c@(a, npy_copysign@c@(b - b, b))); + return (npymath_cpack@c@(a, npymath_copysign@c@(b - b, b))); } } /* @@ -422,7 +422,7 @@ npy_csqrt@c@(@ctype@ z) */ /* Scale to avoid overflow. */ - if (npy_fabs@c@(a) >= THRESH || npy_fabs@c@(b) >= THRESH) { + if (npymath_fabs@c@(a) >= THRESH || npymath_fabs@c@(b) >= THRESH) { a *= 0.25; b *= 0.25; scale = 1; @@ -433,17 +433,17 @@ npy_csqrt@c@(@ctype@ z) /* Algorithm 312, CACM vol 10, Oct 1967. */ if (a >= 0) { - t = npy_sqrt@c@((a + npy_hypot@c@(a, b)) * 0.5@c@); - result = npy_cpack@c@(t, b / (2 * t)); + t = npymath_sqrt@c@((a + npymath_hypot@c@(a, b)) * 0.5@c@); + result = npymath_cpack@c@(t, b / (2 * t)); } else { - t = npy_sqrt@c@((-a + npy_hypot@c@(a, b)) * 0.5@c@); - result = npy_cpack@c@(npy_fabs@c@(b) / (2 * t), npy_copysign@c@(t, b)); + t = npymath_sqrt@c@((-a + npymath_hypot@c@(a, b)) * 0.5@c@); + result = npymath_cpack@c@(npymath_fabs@c@(b) / (2 * t), npymath_copysign@c@(t, b)); } /* Rescale. */ if (scale) { - return (npy_cpack@c@(npy_creal@c@(&result) * 2, npy_cimag@c@(&result))); + return (npymath_cpack@c@(npymath_creal@c@(&result) * 2, npymath_cimag@c@(&result))); } else { return (result); @@ -457,7 +457,7 @@ npy_csqrt@c@(@ctype@ z) * integer powers, but in the body use cpow if it is available. */ -/* private function for use in npy_pow{f, ,l} */ +/* private function for use in npymath_pow{f, ,l} */ #ifdef NPYMATH_HAVE_CPOW@C@ static @ctype@ sys_cpow@c@(@ctype@ x, @ctype@ y) @@ -468,13 +468,13 @@ sys_cpow@c@(@ctype@ x, @ctype@ y) @ctype@ -npy_cpow@c@ (@ctype@ a, @ctype@ b) +npymath_cpow@c@ (@ctype@ a, @ctype@ b) { Py_intptr_t n; - @type@ ar = npy_creal@c@(&a); - @type@ br = npy_creal@c@(&b); - @type@ ai = npy_cimag@c@(&a); - @type@ bi = npy_cimag@c@(&b); + @type@ ar = npymath_creal@c@(&a); + @type@ br = npymath_creal@c@(&b); + @type@ ai = npymath_cimag@c@(&a); + @type@ bi = npymath_cimag@c@(&b); @ctype@ r; /* @@ -483,7 +483,7 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) * If a is also zero then 0^0 is best defined as 1. */ if (br == 0. && bi == 0.) { - return npy_cpack@c@(1., 0.); + return npymath_cpack@c@(1., 0.); } /* case 0^b * If a is a complex zero (ai=ar=0), then the result depends @@ -499,7 +499,7 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) * real and imaginary part. */ if (br > 0) { - return npy_cpack@c@(0., 0.); + return npymath_cpack@c@(0., 0.); } /* else we are in the case where the * real part of b is negative (br<0). @@ -508,17 +508,17 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) */ /* Raise invalid value by calling inf - inf*/ - volatile @type@ tmp = NPY_INFINITY@C@; - tmp -= NPY_INFINITY@C@; + volatile @type@ tmp = NPYMATH_INFINITY@C@; + tmp -= NPYMATH_INFINITY@C@; ar = tmp; - r = npy_cpack@c@(NPY_NAN@C@, NPY_NAN@C@); + r = npymath_cpack@c@(NPYMATH_NAN@C@, NPYMATH_NAN@C@); return r; } if (bi == 0 && br > -100 && br < 100 && (n=(Py_intptr_t)br) == br) { if (n == 1) { /* unroll: handle inf better */ - return npy_cpack@c@(ar, ai); + return npymath_cpack@c@(ar, ai); } else if (n == 2) { /* unroll: handle inf better */ @@ -535,7 +535,7 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) n = -n; } aa = c_1@c@; - p = npy_cpack@c@(ar, ai); + p = npymath_cpack@c@(ar, ai); while (1) { if (n & mask) { aa = cmul@c@(aa,p); @@ -546,7 +546,7 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) } p = cmul@c@(p,p); } - r = npy_cpack@c@(npy_creal@c@(&aa), npy_cimag@c@(&aa)); + r = npymath_cpack@c@(npymath_creal@c@(&aa), npymath_cimag@c@(&aa)); if (br < 0) { r = cdiv@c@(c_1@c@, r); } @@ -559,11 +559,11 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) #else { - @ctype@ loga = npy_clog@c@(a); + @ctype@ loga = npymath_clog@c@(a); - ar = npy_creal@c@(&loga); - ai = npy_cimag@c@(&loga); - return npy_cexp@c@(npy_cpack@c@(ar*br - ai*bi, ar*bi + ai*br)); + ar = npymath_creal@c@(&loga); + ai = npymath_cimag@c@(&loga); + return npymath_cexp@c@(npymath_cpack@c@(ar*br - ai*bi, ar*bi + ai*br)); } #endif @@ -572,30 +572,30 @@ npy_cpow@c@ (@ctype@ a, @ctype@ b) #ifndef NPYMATH_HAVE_CCOS@C@ @ctype@ -npy_ccos@c@(@ctype@ z) +npymath_ccos@c@(@ctype@ z) { /* ccos(z) = ccosh(I * z) */ - return npy_ccosh@c@(npy_cpack@c@(-npy_cimag@c@(&z), npy_creal@c@(&z))); + return npymath_ccosh@c@(npymath_cpack@c@(-npymath_cimag@c@(&z), npymath_creal@c@(&z))); } #endif #ifndef NPYMATH_HAVE_CSIN@C@ @ctype@ -npy_csin@c@(@ctype@ z) +npymath_csin@c@(@ctype@ z) { /* csin(z) = -I * csinh(I * z) */ - z = npy_csinh@c@(npy_cpack@c@(-npy_cimag@c@(&z), npy_creal@c@(&z))); - return npy_cpack@c@(npy_cimag@c@(&z), -npy_creal@c@(&z)); + z = npymath_csinh@c@(npymath_cpack@c@(-npymath_cimag@c@(&z), npymath_creal@c@(&z))); + return npymath_cpack@c@(npymath_cimag@c@(&z), -npymath_creal@c@(&z)); } #endif #ifndef NPYMATH_HAVE_CTAN@C@ @ctype@ -npy_ctan@c@(@ctype@ z) +npymath_ctan@c@(@ctype@ z) { /* ctan(z) = -I * ctanh(I * z) */ - z = npy_ctanh@c@(npy_cpack@c@(-npy_cimag@c@(&z), npy_creal@c@(&z))); - return (npy_cpack@c@(npy_cimag@c@(&z), -npy_creal@c@(&z))); + z = npymath_ctanh@c@(npymath_cpack@c@(-npymath_cimag@c@(&z), npymath_creal@c@(&z))); + return (npymath_cpack@c@(npymath_cimag@c@(&z), -npymath_creal@c@(&z))); } #endif @@ -616,7 +616,7 @@ npy_ctan@c@(@ctype@ z) * although the exact value assigned to CCOSH_BIG is not so important */ @ctype@ -npy_ccosh@c@(@ctype@ z) +npymath_ccosh@c@(@ctype@ z) { #if @precision@ == 1 const npymath_float CCOSH_BIG = 9.0f; @@ -639,41 +639,41 @@ npy_ccosh@c@(@ctype@ z) @type@ x, y, h, absx; npymath_int xfinite, yfinite; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); - xfinite = npy_isfinite(x); - yfinite = npy_isfinite(y); + xfinite = npymath_isfinite(x); + yfinite = npymath_isfinite(y); /* Handle the nearly-non-exceptional cases where x and y are finite. */ if (xfinite && yfinite) { if (y == 0) { - return npy_cpack@c@(npy_cosh@c@(x), x * y); + return npymath_cpack@c@(npymath_cosh@c@(x), x * y); } - absx = npy_fabs@c@(x); + absx = npymath_fabs@c@(x); if (absx < CCOSH_BIG) { /* small x: normal case */ - return npy_cpack@c@(npy_cosh@c@(x) * npy_cos@c@(y), - npy_sinh@c@(x) * npy_sin@c@(y)); + return npymath_cpack@c@(npymath_cosh@c@(x) * npymath_cos@c@(y), + npymath_sinh@c@(x) * npymath_sin@c@(y)); } /* |x| >= 22, so cosh(x) ~= exp(|x|) */ if (absx < SCALED_CEXP_LOWER@C@) { /* x < 710: exp(|x|) won't overflow */ - h = npy_exp@c@(absx) * 0.5@c@; - return npy_cpack@c@(h * npy_cos@c@(y), - npy_copysign@c@(h, x) * npy_sin@c@(y)); + h = npymath_exp@c@(absx) * 0.5@c@; + return npymath_cpack@c@(h * npymath_cos@c@(y), + npymath_copysign@c@(h, x) * npymath_sin@c@(y)); } else if (absx < SCALED_CEXP_UPPER@C@) { /* x < 1455: scale to avoid overflow */ - z = _npy_scaled_cexp@c@(absx, y, -1); - return npy_cpack@c@(npy_creal@c@(&z), - npy_cimag@c@(&z) * npy_copysign@c@(1, x)); + z = _npymath_scaled_cexp@c@(absx, y, -1); + return npymath_cpack@c@(npymath_creal@c@(&z), + npymath_cimag@c@(&z) * npymath_copysign@c@(1, x)); } else { /* x >= 1455: the result always overflows */ h = CCOSH_HUGE * x; - return npy_cpack@c@(h * h * npy_cos@c@(y), h * npy_sin@c@(y)); + return npymath_cpack@c@(h * h * npymath_cos@c@(y), h * npymath_sin@c@(y)); } } @@ -687,7 +687,7 @@ npy_ccosh@c@(@ctype@ z) * the same as d(NaN). */ if (x == 0 && !yfinite) { - return npy_cpack@c@(y - y, npy_copysign@c@(0, x * (y - y))); + return npymath_cpack@c@(y - y, npymath_copysign@c@(0, x * (y - y))); } /* @@ -697,7 +697,7 @@ npy_ccosh@c@(@ctype@ z) * The sign of 0 in the result is unspecified. */ if (y == 0 && !xfinite) { - return npy_cpack@c@(x * x, npy_copysign@c@(0, x) * y); + return npymath_cpack@c@(x * x, npymath_copysign@c@(0, x) * y); } /* @@ -709,7 +709,7 @@ npy_ccosh@c@(@ctype@ z) * nonzero x. Choice = don't raise (except for signaling NaNs). */ if (xfinite && !yfinite) { - return npy_cpack@c@(y - y, x * (y - y)); + return npymath_cpack@c@(y - y, x * (y - y)); } /* @@ -721,11 +721,11 @@ npy_ccosh@c@(@ctype@ z) * * cosh(+-Inf + I y) = +Inf cos(y) +- I Inf sin(y) */ - if (npy_isinf(x)) { + if (npymath_isinf(x)) { if (!yfinite) { - return npy_cpack@c@(x * x, x * (y - y)); + return npymath_cpack@c@(x * x, x * (y - y)); } - return npy_cpack@c@((x * x) * npy_cos@c@(y), x * npy_sin@c@(y)); + return npymath_cpack@c@((x * x) * npymath_cos@c@(y), x * npymath_sin@c@(y)); } /* @@ -739,7 +739,7 @@ npy_ccosh@c@(@ctype@ z) * Optionally raises the invalid floating-point exception for finite * nonzero y. Choice = don't raise (except for signaling NaNs). */ - return npy_cpack@c@((x * x) * (y - y), (x + x) * (y - y)); + return npymath_cpack@c@((x * x) * (y - y), (x + x) * (y - y)); } #undef CCOSH_BIG #undef CCOSH_HUGE @@ -758,7 +758,7 @@ npy_ccosh@c@(@ctype@ z) * These values and the return value were taken from n1124.pdf. */ @ctype@ -npy_csinh@c@(@ctype@ z) +npymath_csinh@c@(@ctype@ z) { #if @precision@ == 1 const npymath_float CSINH_BIG = 9.0f; @@ -781,41 +781,41 @@ npy_csinh@c@(@ctype@ z) @type@ x, y, h, absx; npymath_int xfinite, yfinite; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); - xfinite = npy_isfinite(x); - yfinite = npy_isfinite(y); + xfinite = npymath_isfinite(x); + yfinite = npymath_isfinite(y); /* Handle the nearly-non-exceptional cases where x and y are finite. */ if (xfinite && yfinite) { if (y == 0) { - return npy_cpack@c@(npy_sinh@c@(x), y); + return npymath_cpack@c@(npymath_sinh@c@(x), y); } - absx = npy_fabs@c@(x); + absx = npymath_fabs@c@(x); if (absx < CSINH_BIG) { /* small x: normal case */ - return npy_cpack@c@(npy_sinh@c@(x) * npy_cos@c@(y), - npy_cosh@c@(x) * npy_sin@c@(y)); + return npymath_cpack@c@(npymath_sinh@c@(x) * npymath_cos@c@(y), + npymath_cosh@c@(x) * npymath_sin@c@(y)); } /* |x| >= 22, so cosh(x) ~= exp(|x|) */ if (absx < SCALED_CEXP_LOWER@C@) { /* x < 710: exp(|x|) won't overflow */ - h = npy_exp@c@(npy_fabs@c@(x)) * 0.5@c@; - return npy_cpack@c@(npy_copysign@c@(h, x) * npy_cos@c@(y), - h * npy_sin@c@(y)); + h = npymath_exp@c@(npymath_fabs@c@(x)) * 0.5@c@; + return npymath_cpack@c@(npymath_copysign@c@(h, x) * npymath_cos@c@(y), + h * npymath_sin@c@(y)); } else if (x < SCALED_CEXP_UPPER@C@) { /* x < 1455: scale to avoid overflow */ - z = _npy_scaled_cexp@c@(absx, y, -1); - return npy_cpack@c@(npy_creal@c@(&z) * npy_copysign@c@(1, x), - npy_cimag@c@(&z)); + z = _npymath_scaled_cexp@c@(absx, y, -1); + return npymath_cpack@c@(npymath_creal@c@(&z) * npymath_copysign@c@(1, x), + npymath_cimag@c@(&z)); } else { /* x >= 1455: the result always overflows */ h = CSINH_HUGE * x; - return npy_cpack@c@(h * npy_cos@c@(y), h * h * npy_sin@c@(y)); + return npymath_cpack@c@(h * npymath_cos@c@(y), h * h * npymath_sin@c@(y)); } } @@ -829,7 +829,7 @@ npy_csinh@c@(@ctype@ z) * the same as d(NaN). */ if (x == 0 && !yfinite) { - return npy_cpack@c@(npy_copysign@c@(0, x * (y - y)), y - y); + return npymath_cpack@c@(npymath_copysign@c@(0, x * (y - y)), y - y); } /* @@ -838,10 +838,10 @@ npy_csinh@c@(@ctype@ z) * sinh(NaN +- I 0) = d(NaN) + I +-0. */ if (y == 0 && !xfinite) { - if (npy_isnan(x)) { + if (npymath_isnan(x)) { return z; } - return npy_cpack@c@(x, npy_copysign@c@(0, y)); + return npymath_cpack@c@(x, npymath_copysign@c@(0, y)); } /* @@ -853,7 +853,7 @@ npy_csinh@c@(@ctype@ z) * nonzero x. Choice = don't raise (except for signaling NaNs). */ if (xfinite && !yfinite) { - return npy_cpack@c@(y - y, x * (y - y)); + return npymath_cpack@c@(y - y, x * (y - y)); } /* @@ -867,12 +867,12 @@ npy_csinh@c@(@ctype@ z) * * sinh(+-Inf + I y) = +-Inf cos(y) + I Inf sin(y) */ - if (!xfinite && !npy_isnan(x)) { + if (!xfinite && !npymath_isnan(x)) { if (!yfinite) { - return npy_cpack@c@(x * x, x * (y - y)); + return npymath_cpack@c@(x * x, x * (y - y)); } - return npy_cpack@c@(x * npy_cos@c@(y), - NPY_INFINITY@C@ * npy_sin@c@(y)); + return npymath_cpack@c@(x * npymath_cos@c@(y), + NPYMATH_INFINITY@C@ * npymath_sin@c@(y)); } /* @@ -886,7 +886,7 @@ npy_csinh@c@(@ctype@ z) * Optionally raises the invalid floating-point exception for finite * nonzero y. Choice = don't raise (except for signaling NaNs). */ - return npy_cpack@c@((x * x) * (y - y), (x + x) * (y - y)); + return npymath_cpack@c@((x * x) * (y - y), (x + x) * (y - y)); } #undef CSINH_BIG #undef CSINH_HUGE @@ -939,13 +939,13 @@ npy_csinh@c@(@ctype@ z) #define TANHL_HUGE 42.0L @ctype@ -npy_ctanh@c@(@ctype@ z) +npymath_ctanh@c@(@ctype@ z) { @type@ x, y; @type@ t, beta, s, rho, denom; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); /* * ctanh(NaN + i 0) = NaN + i 0 @@ -963,22 +963,22 @@ npy_ctanh@c@(@ctype@ z) * case is only needed to avoid a spurious invalid exception when * y is infinite. */ - if (!npy_isfinite(x)) { - if (npy_isnan(x)) { - return npy_cpack@c@(x, (y == 0 ? y : x * y)); + if (!npymath_isfinite(x)) { + if (npymath_isnan(x)) { + return npymath_cpack@c@(x, (y == 0 ? y : x * y)); } - return npy_cpack@c@(npy_copysign@c@(1,x), - npy_copysign@c@(0, - npy_isinf(y) ? - y : npy_sin@c@(y) * npy_cos@c@(y))); + return npymath_cpack@c@(npymath_copysign@c@(1,x), + npymath_copysign@c@(0, + npymath_isinf(y) ? + y : npymath_sin@c@(y) * npymath_cos@c@(y))); } /* * ctanh(x + i NAN) = NaN + i NaN * ctanh(x +- i Inf) = NaN + i NaN */ - if (!npy_isfinite(y)) { - return (npy_cpack@c@(y - y, y - y)); + if (!npymath_isfinite(y)) { + return (npymath_cpack@c@(y - y, y - y)); } /* @@ -986,20 +986,20 @@ npy_ctanh@c@(@ctype@ z) * approximation sinh^2(huge) ~= exp(2*huge) / 4. * We use a modified formula to avoid spurious overflow. */ - if (npy_fabs@c@(x) >= TANH@C@_HUGE) { - @type@ exp_mx = npy_exp@c@(-npy_fabs@c@(x)); - return npy_cpack@c@(npy_copysign@c@(1, x), - 4 * npy_sin@c@(y) * npy_cos@c@(y) * + if (npymath_fabs@c@(x) >= TANH@C@_HUGE) { + @type@ exp_mx = npymath_exp@c@(-npymath_fabs@c@(x)); + return npymath_cpack@c@(npymath_copysign@c@(1, x), + 4 * npymath_sin@c@(y) * npymath_cos@c@(y) * exp_mx * exp_mx); } /* Kahan's algorithm */ - t = npy_tan@c@(y); + t = npymath_tan@c@(y); beta = 1 + t * t; /* = 1 / cos^2(y) */ - s = npy_sinh@c@(x); - rho = npy_sqrt@c@(1 + s * s); /* = cosh(x) */ + s = npymath_sinh@c@(x); + rho = npymath_sqrt@c@(1 + s * s); /* = cosh(x) */ denom = 1 + beta * s * s; - return (npy_cpack@c@((beta * rho * s) / denom, t / denom)); + return (npymath_cpack@c@((beta * rho * s) / denom, t / denom)); } #undef TANH_HUGE #undef TANHF_HUGE @@ -1106,8 +1106,8 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, @type@ R, S, A; /* A, B, R, and S are as in Hull et al. */ @type@ Am1, Amy; /* A-1, A-y. */ - R = npy_hypot@c@(x, y + 1); /* |z+I| */ - S = npy_hypot@c@(x, y - 1); /* |z-I| */ + R = npymath_hypot@c@(x, y + 1); /* |z+I| */ + S = npymath_hypot@c@(x, y - 1); /* |z-I| */ /* A = (|z+I| + |z-I|) / 2 */ A = (R + S) / 2; @@ -1130,32 +1130,32 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, * fp is of order x^2, and fm = x/2. * A = 1 (inexactly). */ - *rx = npy_sqrt@c@(x); + *rx = npymath_sqrt@c@(x); } - else if (x >= @TEPS@ * npy_fabs@c@(y - 1)) { + else if (x >= @TEPS@ * npymath_fabs@c@(y - 1)) { /* * Underflow will not occur because * x >= DBL_EPSILON^2/128 >= FOUR_SQRT_MIN */ Am1 = _f@c@(x, 1 + y, R) + _f@c@(x, 1 - y, S); - *rx = npy_log1p@c@(Am1 + npy_sqrt@c@(Am1 * (A + 1))); + *rx = npymath_log1p@c@(Am1 + npymath_sqrt@c@(Am1 * (A + 1))); } else if (y < 1) { /* * fp = x*x/(1+y)/4, fm = x*x/(1-y)/4, and * A = 1 (inexactly). */ - *rx = x / npy_sqrt@c@((1 - y) * (1 + y)); + *rx = x / npymath_sqrt@c@((1 - y) * (1 + y)); } else { /* if (y > 1) */ /* * A-1 = y-1 (inexactly). */ - *rx = npy_log1p@c@((y - 1) + npy_sqrt@c@((y - 1) * (y + 1))); + *rx = npymath_log1p@c@((y - 1) + npymath_sqrt@c@((y - 1) * (y + 1))); } } else { - *rx = npy_log@c@(A + npy_sqrt@c@(A * A - 1)); + *rx = npymath_log@c@(A + npymath_sqrt@c@(A * A - 1)); } *new_y = y; @@ -1187,9 +1187,9 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, * fp is of order x^2, and fm = x/2. * A = 1 (inexactly). */ - *sqrt_A2my2 = npy_sqrt@c@(x) * npy_sqrt@c@((A + y) / 2); + *sqrt_A2my2 = npymath_sqrt@c@(x) * npymath_sqrt@c@((A + y) / 2); } - else if (x >= @TEPS@ * npy_fabs@c@(y - 1)) { + else if (x >= @TEPS@ * npymath_fabs@c@(y - 1)) { /* * Underflow will not occur because * x >= DBL_EPSILON/128 >= FOUR_SQRT_MIN @@ -1197,7 +1197,7 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, * x >= DBL_EPSILON^2 >= FOUR_SQRT_MIN */ Amy = _f@c@(x, y + 1, R) + _f@c@(x, y - 1, S); - *sqrt_A2my2 = npy_sqrt@c@(Amy * (A + y)); + *sqrt_A2my2 = npymath_sqrt@c@(Amy * (A + y)); } else if (y > 1) { /* @@ -1208,7 +1208,7 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, * scaling should avoid any underflow problems. */ *sqrt_A2my2 = x * (4 / @TEPS@ / @TEPS@) * y / - npy_sqrt@c@((y + 1) * (y - 1)); + npymath_sqrt@c@((y + 1) * (y - 1)); *new_y = y * (4 / @TEPS@ / @TEPS@); } else { /* if (y < 1) */ @@ -1216,7 +1216,7 @@ _do_hard_work@c@(@type@ x, @type@ y, @type@ *rx, * fm = 1-y >= DBL_EPSILON, fp is of order x^2, and * A = 1 (inexactly). */ - *sqrt_A2my2 = npy_sqrt@c@((1 - y) * (1 + y)); + *sqrt_A2my2 = npymath_sqrt@c@((1 - y) * (1 + y)); } } } @@ -1247,8 +1247,8 @@ _clog_for_large_values@c@(@type@ x, @type@ y, #endif @type@ ax, ay, t; - ax = npy_fabs@c@(x); - ay = npy_fabs@c@(y); + ax = npymath_fabs@c@(x); + ay = npymath_fabs@c@(y); if (ax < ay) { t = ax; ax = ay; @@ -1263,25 +1263,25 @@ _clog_for_large_values@c@(@type@ x, @type@ y, * this method is still poor since it is unnecessarily slow. */ if (ax > @TMAX@ / 2) { - *rr = npy_log@c@(npy_hypot@c@(x / NPY_E@c@, y / NPY_E@c@)) + 1; + *rr = npymath_log@c@(npymath_hypot@c@(x / NPYMATH_E@c@, y / NPYMATH_E@c@)) + 1; } /* * Avoid overflow when x or y is large. Avoid underflow when x or * y is small. */ else if (ax > QUARTER_SQRT_MAX || ay < SQRT_MIN) { - *rr = npy_log@c@(npy_hypot@c@(x, y)); + *rr = npymath_log@c@(npymath_hypot@c@(x, y)); } else { - *rr = npy_log@c@(ax * ax + ay * ay) / 2; + *rr = npymath_log@c@(ax * ax + ay * ay) / 2; } - *ri = npy_atan2@c@(y, x); + *ri = npymath_atan2@c@(y, x); } #endif #ifndef NPYMATH_HAVE_CACOS@C@ @ctype@ -npy_cacos@c@(@ctype@ z) +npymath_cacos@c@(@ctype@ z) { #if @precision@ == 1 /* this is sqrt(6*EPS) */ @@ -1298,109 +1298,109 @@ npy_cacos@c@(@ctype@ z) const volatile npymath_longdouble pio2_lo = 2.710505431213761085e-20l; #endif const @type@ RECIP_EPSILON = 1.0@c@ / @TEPS@; - const @type@ pio2_hi = NPY_PI_2@c@; + const @type@ pio2_hi = NPYMATH_PI_2@c@; @type@ x, y, ax, ay, wx, wy, rx, ry, B, sqrt_A2mx2, new_x; npymath_int sx, sy; npymath_int B_is_usable; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); - sx = npy_signbit(x); - sy = npy_signbit(y); - ax = npy_fabs@c@(x); - ay = npy_fabs@c@(y); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); + sx = npymath_signbit(x); + sy = npymath_signbit(y); + ax = npymath_fabs@c@(x); + ay = npymath_fabs@c@(y); - if (npy_isnan(x) || npy_isnan(y)) { + if (npymath_isnan(x) || npymath_isnan(y)) { /* cacos(+-Inf + I*NaN) = NaN + I*opt(-)Inf */ - if (npy_isinf(x)) { - return npy_cpack@c@(y + y, -NPY_INFINITY@C@); + if (npymath_isinf(x)) { + return npymath_cpack@c@(y + y, -NPYMATH_INFINITY@C@); } /* cacos(NaN + I*+-Inf) = NaN + I*-+Inf */ - if (npy_isinf(y)) { - return npy_cpack@c@(x + x, -y); + if (npymath_isinf(y)) { + return npymath_cpack@c@(x + x, -y); } /* cacos(0 + I*NaN) = PI/2 + I*NaN with inexact */ if (x == 0) { - return npy_cpack@c@(pio2_hi + pio2_lo, y + y); + return npymath_cpack@c@(pio2_hi + pio2_lo, y + y); } /* * All other cases involving NaN return NaN + I*NaN. * C99 leaves it optional whether to raise invalid if one of * the arguments is not NaN, so we opt not to raise it. */ - return npy_cpack@c@(NPY_NAN@C@, NPY_NAN@C@); + return npymath_cpack@c@(NPYMATH_NAN@C@, NPYMATH_NAN@C@); } if (ax > RECIP_EPSILON || ay > RECIP_EPSILON) { /* clog...() will raise inexact unless x or y is infinite. */ _clog_for_large_values@c@(x, y, &wx, &wy); - rx = npy_fabs@c@(wy); - ry = wx + NPY_LOGE2@c@; + rx = npymath_fabs@c@(wy); + ry = wx + NPYMATH_LOGE2@c@; if (sy == 0) { ry = -ry; } - return npy_cpack@c@(rx, ry); + return npymath_cpack@c@(rx, ry); } /* Avoid spuriously raising inexact for z = 1. */ if (x == 1 && y == 0) { - return npy_cpack@c@(0, -y); + return npymath_cpack@c@(0, -y); } /* All remaining cases are inexact. */ raise_inexact(); if (ax < SQRT_6_EPSILON / 4 && ay < SQRT_6_EPSILON / 4) { - return npy_cpack@c@(pio2_hi - (x - pio2_lo), -y); + return npymath_cpack@c@(pio2_hi - (x - pio2_lo), -y); } _do_hard_work@c@(ay, ax, &ry, &B_is_usable, &B, &sqrt_A2mx2, &new_x); if (B_is_usable) { if (sx == 0) { - rx = npy_acos@c@(B); + rx = npymath_acos@c@(B); } else { - rx = npy_acos@c@(-B); + rx = npymath_acos@c@(-B); } } else { if (sx == 0) { - rx = npy_atan2@c@(sqrt_A2mx2, new_x); + rx = npymath_atan2@c@(sqrt_A2mx2, new_x); } else { - rx = npy_atan2@c@(sqrt_A2mx2, -new_x); + rx = npymath_atan2@c@(sqrt_A2mx2, -new_x); } } if (sy == 0) { ry = -ry; } - return npy_cpack@c@(rx, ry); + return npymath_cpack@c@(rx, ry); } #endif #ifndef NPYMATH_HAVE_CASIN@C@ @ctype@ -npy_casin@c@(@ctype@ z) +npymath_casin@c@(@ctype@ z) { /* casin(z) = I * conj( casinh(I * conj(z)) ) */ - z = npy_casinh@c@(npy_cpack@c@(npy_cimag@c@(&z), npy_creal@c@(&z))); - return npy_cpack@c@(npy_cimag@c@(&z), npy_creal@c@(&z)); + z = npymath_casinh@c@(npymath_cpack@c@(npymath_cimag@c@(&z), npymath_creal@c@(&z))); + return npymath_cpack@c@(npymath_cimag@c@(&z), npymath_creal@c@(&z)); } #endif #ifndef NPYMATH_HAVE_CATAN@C@ @ctype@ -npy_catan@c@(@ctype@ z) +npymath_catan@c@(@ctype@ z) { /* catan(z) = I * conj( catanh(I * conj(z)) ) */ - z = npy_catanh@c@(npy_cpack@c@(npy_cimag@c@(&z), npy_creal@c@(&z))); - return npy_cpack@c@(npy_cimag@c@(&z), npy_creal@c@(&z)); + z = npymath_catanh@c@(npymath_cpack@c@(npymath_cimag@c@(&z), npymath_creal@c@(&z))); + return npymath_cpack@c@(npymath_cimag@c@(&z), npymath_creal@c@(&z)); } #endif #ifndef NPYMATH_HAVE_CACOSH@C@ @ctype@ -npy_cacosh@c@(@ctype@ z) +npymath_cacosh@c@(@ctype@ z) { /* * cacosh(z) = I*cacos(z) or -I*cacos(z) @@ -1409,23 +1409,23 @@ npy_cacosh@c@(@ctype@ z) @ctype@ w; @type@ rx, ry; - w = npy_cacos@c@(z); - rx = npy_creal@c@(&w); - ry = npy_cimag@c@(&w); + w = npymath_cacos@c@(z); + rx = npymath_creal@c@(&w); + ry = npymath_cimag@c@(&w); /* cacosh(NaN + I*NaN) = NaN + I*NaN */ - if (npy_isnan(rx) && npy_isnan(ry)) { - return npy_cpack@c@(ry, rx); + if (npymath_isnan(rx) && npymath_isnan(ry)) { + return npymath_cpack@c@(ry, rx); } /* cacosh(NaN + I*+-Inf) = +Inf + I*NaN */ /* cacosh(+-Inf + I*NaN) = +Inf + I*NaN */ - if (npy_isnan(rx)) { - return npy_cpack@c@(npy_fabs@c@(ry), rx); + if (npymath_isnan(rx)) { + return npymath_cpack@c@(npymath_fabs@c@(ry), rx); } /* cacosh(0 + I*NaN) = NaN + I*NaN */ - if (npy_isnan(ry)) { - return npy_cpack@c@(ry, ry); + if (npymath_isnan(ry)) { + return npymath_cpack@c@(ry, ry); } - return npy_cpack@c@(npy_fabs@c@(ry), npy_copysign@c@(rx, npy_cimag@c@(&z))); + return npymath_cpack@c@(npymath_fabs@c@(ry), npymath_copysign@c@(rx, npymath_cimag@c@(&z))); } #endif @@ -1439,7 +1439,7 @@ npy_cacosh@c@(@ctype@ z) * as z -> infinity, uniformly in y */ @ctype@ -npy_casinh@c@(@ctype@ z) +npymath_casinh@c@(@ctype@ z) { #if @precision@ == 1 /* this is sqrt(6*EPS) */ @@ -1455,43 +1455,43 @@ npy_casinh@c@(@ctype@ z) @type@ x, y, ax, ay, wx, wy, rx, ry, B, sqrt_A2my2, new_y; npymath_int B_is_usable; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); - ax = npy_fabs@c@(x); - ay = npy_fabs@c@(y); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); + ax = npymath_fabs@c@(x); + ay = npymath_fabs@c@(y); - if (npy_isnan(x) || npy_isnan(y)) { + if (npymath_isnan(x) || npymath_isnan(y)) { /* casinh(+-Inf + I*NaN) = +-Inf + I*NaN */ - if (npy_isinf(x)) { - return npy_cpack@c@(x, y + y); + if (npymath_isinf(x)) { + return npymath_cpack@c@(x, y + y); } /* casinh(NaN + I*+-Inf) = opt(+-)Inf + I*NaN */ - if (npy_isinf(y)) { - return npy_cpack@c@(y, x + x); + if (npymath_isinf(y)) { + return npymath_cpack@c@(y, x + x); } /* casinh(NaN + I*0) = NaN + I*0 */ if (y == 0) { - return npy_cpack@c@(x + x, y); + return npymath_cpack@c@(x + x, y); } /* * All other cases involving NaN return NaN + I*NaN. * C99 leaves it optional whether to raise invalid if one of * the arguments is not NaN, so we opt not to raise it. */ - return npy_cpack@c@(NPY_NAN@C@, NPY_NAN@C@); + return npymath_cpack@c@(NPYMATH_NAN@C@, NPYMATH_NAN@C@); } if (ax > RECIP_EPSILON || ay > RECIP_EPSILON) { /* clog...() will raise inexact unless x or y is infinite. */ - if (npy_signbit(x) == 0) { + if (npymath_signbit(x) == 0) { _clog_for_large_values@c@(x, y, &wx, &wy); - wx += NPY_LOGE2@c@; + wx += NPYMATH_LOGE2@c@; } else { _clog_for_large_values@c@(-x, -y, &wx, &wy); - wx += NPY_LOGE2@c@; + wx += NPYMATH_LOGE2@c@; } - return npy_cpack@c@(npy_copysign@c@(wx, x), npy_copysign@c@(wy, y)); + return npymath_cpack@c@(npymath_copysign@c@(wx, x), npymath_copysign@c@(wy, y)); } /* Avoid spuriously raising inexact for z = 0. */ @@ -1508,12 +1508,12 @@ npy_casinh@c@(@ctype@ z) _do_hard_work@c@(ax, ay, &rx, &B_is_usable, &B, &sqrt_A2my2, &new_y); if (B_is_usable) { - ry = npy_asin@c@(B); + ry = npymath_asin@c@(B); } else { - ry = npy_atan2@c@(new_y, sqrt_A2my2); + ry = npymath_atan2@c@(new_y, sqrt_A2my2); } - return npy_cpack@c@(npy_copysign@c@(rx, x), npy_copysign@c@(ry, y)); + return npymath_cpack@c@(npymath_copysign@c@(rx, x), npymath_copysign@c@(ry, y)); } #endif @@ -1573,7 +1573,7 @@ _real_part_reciprocalf(npymath_float x, npymath_float y) ix = hx & 0x7f800000; GET_FLOAT_WORD(hy, y); iy = hy & 0x7f800000; - if (ix - iy >= CUTOFF << 23 || npy_isinf(x)) { + if (ix - iy >= CUTOFF << 23 || npymath_isinf(x)) { return (1 / x); } if (iy - ix >= CUTOFF << 23) { @@ -1610,7 +1610,7 @@ _real_part_reciprocal(npymath_double x, npymath_double y) ix = hx & 0x7ff00000; GET_HIGH_WORD(hy, y); iy = hy & 0x7ff00000; - if (ix - iy >= CUTOFF << 20 || npy_isinf(x)) { + if (ix - iy >= CUTOFF << 20 || npymath_isinf(x)) { /* +-Inf -> +-0 is special */ return (1 / x); } @@ -1649,7 +1649,7 @@ _real_part_reciprocall(npymath_longdouble x, ix = GET_LDOUBLE_EXP(ux); uy.e = y; iy = GET_LDOUBLE_EXP(uy); - if (ix - iy >= CUTOFF || npy_isinf(x)) { + if (ix - iy >= CUTOFF || npymath_isinf(x)) { return (1/x); } if (iy - ix >= CUTOFF) { @@ -1681,7 +1681,7 @@ _real_part_reciprocall(npymath_longdouble x, #endif @ctype@ -npy_catanh@c@(@ctype@ z) +npymath_catanh@c@(@ctype@ z) { #if @precision@ == 1 /* this is sqrt(3*EPS) */ @@ -1698,45 +1698,45 @@ npy_catanh@c@(@ctype@ z) const volatile npymath_longdouble pio2_lo = 2.710505431213761085e-20l; #endif const @type@ RECIP_EPSILON = 1.0@c@ / @TEPS@; - const @type@ pio2_hi = NPY_PI_2@c@; + const @type@ pio2_hi = NPYMATH_PI_2@c@; @type@ x, y, ax, ay, rx, ry; - x = npy_creal@c@(&z); - y = npy_cimag@c@(&z); - ax = npy_fabs@c@(x); - ay = npy_fabs@c@(y); + x = npymath_creal@c@(&z); + y = npymath_cimag@c@(&z); + ax = npymath_fabs@c@(x); + ay = npymath_fabs@c@(y); /* This helps handle many cases. */ if (y == 0 && ax <= 1) { - return npy_cpack@c@(npy_atanh@c@(x), y); + return npymath_cpack@c@(npymath_atanh@c@(x), y); } /* To ensure the same accuracy as atan(), and to filter out z = 0. */ if (x == 0) { - return npy_cpack@c@(x, npy_atan@c@(y)); + return npymath_cpack@c@(x, npymath_atan@c@(y)); } - if (npy_isnan(x) || npy_isnan(y)) { + if (npymath_isnan(x) || npymath_isnan(y)) { /* catanh(+-Inf + I*NaN) = +-0 + I*NaN */ - if (npy_isinf(x)) { - return npy_cpack@c@(npy_copysign@c@(0, x), y + y); + if (npymath_isinf(x)) { + return npymath_cpack@c@(npymath_copysign@c@(0, x), y + y); } /* catanh(NaN + I*+-Inf) = sign(NaN)0 + I*+-PI/2 */ - if (npy_isinf(y)) { - return npy_cpack@c@(npy_copysign@c@(0, x), - npy_copysign@c@(pio2_hi + pio2_lo, y)); + if (npymath_isinf(y)) { + return npymath_cpack@c@(npymath_copysign@c@(0, x), + npymath_copysign@c@(pio2_hi + pio2_lo, y)); } /* * All other cases involving NaN return NaN + I*NaN. * C99 leaves it optional whether to raise invalid if one of * the arguments is not NaN, so we opt not to raise it. */ - return npy_cpack@c@(NPY_NAN@C@, NPY_NAN@C@); + return npymath_cpack@c@(NPYMATH_NAN@C@, NPYMATH_NAN@C@); } if (ax > RECIP_EPSILON || ay > RECIP_EPSILON) { - return npy_cpack@c@(_real_part_reciprocal@c@(x, y), - npy_copysign@c@(pio2_hi + pio2_lo, y)); + return npymath_cpack@c@(_real_part_reciprocal@c@(x, y), + npymath_copysign@c@(pio2_hi + pio2_lo, y)); } if (ax < SQRT_3_EPSILON / 2 && ay < SQRT_3_EPSILON / 2) { @@ -1750,23 +1750,23 @@ npy_catanh@c@(@ctype@ z) } if (ax == 1 && ay < @TEPS@) { - rx = (NPY_LOGE2@c@ - npy_log@c@(ay)) / 2; + rx = (NPYMATH_LOGE2@c@ - npymath_log@c@(ay)) / 2; } else { - rx = npy_log1p@c@(4 * ax / _sum_squares@c@(ax - 1, ay)) / 4; + rx = npymath_log1p@c@(4 * ax / _sum_squares@c@(ax - 1, ay)) / 4; } if (ax == 1) { - ry = npy_atan2@c@(2, -ay) / 2; + ry = npymath_atan2@c@(2, -ay) / 2; } else if (ay < @TEPS@) { - ry = npy_atan2@c@(2 * ay, (1 - ax) * (1 + ax)) / 2; + ry = npymath_atan2@c@(2 * ay, (1 - ax) * (1 + ax)) / 2; } else { - ry = npy_atan2@c@(2 * ay, (1 - ax) * (1 + ax) - ay * ay) / 2; + ry = npymath_atan2@c@(2 * ay, (1 - ax) * (1 + ax) - ay * ay) / 2; } - return npy_cpack@c@(npy_copysign@c@(rx, x), npy_copysign@c@(ry, y)); + return npymath_cpack@c@(npymath_copysign@c@(rx, x), npymath_copysign@c@(ry, y)); } #endif /**end repeat**/ @@ -1788,7 +1788,7 @@ npy_catanh@c@(@ctype@ z) */ #ifdef NPYMATH_HAVE_@KIND@@C@ @type@ -npy_@kind@@c@(@ctype@ z) +npymath_@kind@@c@(@ctype@ z) { return @kind@@c@(z); } @@ -1803,7 +1803,7 @@ npy_@kind@@c@(@ctype@ z) */ #ifdef NPYMATH_HAVE_@KIND@@C@ @ctype@ -npy_@kind@@c@(@ctype@ z) +npymath_@kind@@c@(@ctype@ z) { return @kind@@c@(z); }