Major fixes and new features

venv/lib/python3.12/site-packages/mypyc/lib-rt/float_ops.c

@@ -0,0 +1,192 @@
+// Float primitive operations
+//
+// These are registered in mypyc.primitives.float_ops.
+
+#include <Python.h>
+#include "CPy.h"
+
+
+static double CPy_DomainError(void) {
+    PyErr_SetString(PyExc_ValueError, "math domain error");
+    return CPY_FLOAT_ERROR;
+}
+
+static double CPy_MathRangeError(void) {
+    PyErr_SetString(PyExc_OverflowError, "math range error");
+    return CPY_FLOAT_ERROR;
+}
+
+double CPyFloat_FromTagged(CPyTagged x) {
+    if (CPyTagged_CheckShort(x)) {
+        return CPyTagged_ShortAsSsize_t(x);
+    }
+    double result = PyFloat_AsDouble(CPyTagged_LongAsObject(x));
+    if (unlikely(result == -1.0) && PyErr_Occurred()) {
+        return CPY_FLOAT_ERROR;
+    }
+    return result;
+}
+
+double CPyFloat_Sin(double x) {
+    double v = sin(x);
+    if (unlikely(isnan(v)) && !isnan(x)) {
+        return CPy_DomainError();
+    }
+    return v;
+}
+
+double CPyFloat_Cos(double x) {
+    double v = cos(x);
+    if (unlikely(isnan(v)) && !isnan(x)) {
+        return CPy_DomainError();
+    }
+    return v;
+}
+
+double CPyFloat_Tan(double x) {
+    if (unlikely(isinf(x))) {
+        return CPy_DomainError();
+    }
+    return tan(x);
+}
+
+double CPyFloat_Sqrt(double x) {
+    if (x < 0.0) {
+        return CPy_DomainError();
+    }
+    return sqrt(x);
+}
+
+double CPyFloat_Exp(double x) {
+    double v = exp(x);
+    if (unlikely(v == INFINITY) && x != INFINITY) {
+        return CPy_MathRangeError();
+    }
+    return v;
+}
+
+double CPyFloat_Log(double x) {
+    if (x <= 0.0) {
+        return CPy_DomainError();
+    }
+    return log(x);
+}
+
+CPyTagged CPyFloat_Floor(double x) {
+    double v = floor(x);
+    return CPyTagged_FromFloat(v);
+}
+
+CPyTagged CPyFloat_Ceil(double x) {
+    double v = ceil(x);
+    return CPyTagged_FromFloat(v);
+}
+
+bool CPyFloat_IsInf(double x) {
+    return isinf(x) != 0;
+}
+
+bool CPyFloat_IsNaN(double x) {
+    return isnan(x) != 0;
+}
+
+// From CPython 3.10.0, Objects/floatobject.c
+static void
+_float_div_mod(double vx, double wx, double *floordiv, double *mod)
+{
+    double div;
+    *mod = fmod(vx, wx);
+    /* fmod is typically exact, so vx-mod is *mathematically* an
+       exact multiple of wx.  But this is fp arithmetic, and fp
+       vx - mod is an approximation; the result is that div may
+       not be an exact integral value after the division, although
+       it will always be very close to one.
+    */
+    div = (vx - *mod) / wx;
+    if (*mod) {
+        /* ensure the remainder has the same sign as the denominator */
+        if ((wx < 0) != (*mod < 0)) {
+            *mod += wx;
+            div -= 1.0;
+        }
+    }
+    else {
+        /* the remainder is zero, and in the presence of signed zeroes
+           fmod returns different results across platforms; ensure
+           it has the same sign as the denominator. */
+        *mod = copysign(0.0, wx);
+    }
+    /* snap quotient to nearest integral value */
+    if (div) {
+        *floordiv = floor(div);
+        if (div - *floordiv > 0.5) {
+            *floordiv += 1.0;
+        }
+    }
+    else {
+        /* div is zero - get the same sign as the true quotient */
+        *floordiv = copysign(0.0, vx / wx); /* zero w/ sign of vx/wx */
+    }
+}
+
+double CPyFloat_FloorDivide(double x, double y) {
+    double mod, floordiv;
+    if (y == 0) {
+        PyErr_SetString(PyExc_ZeroDivisionError, "float floor division by zero");
+        return CPY_FLOAT_ERROR;
+    }
+    _float_div_mod(x, y, &floordiv, &mod);
+    return floordiv;
+}
+
+// Adapted from CPython 3.10.7
+double CPyFloat_Pow(double x, double y) {
+    if (!isfinite(x) || !isfinite(y)) {
+        if (isnan(x))
+            return y == 0.0 ? 1.0 : x; /* NaN**0 = 1 */
+        else if (isnan(y))
+            return x == 1.0 ? 1.0 : y; /* 1**NaN = 1 */
+        else if (isinf(x)) {
+            int odd_y = isfinite(y) && fmod(fabs(y), 2.0) == 1.0;
+            if (y > 0.0)
+                return odd_y ? x : fabs(x);
+            else if (y == 0.0)
+                return 1.0;
+            else /* y < 0. */
+                return odd_y ? copysign(0.0, x) : 0.0;
+        }
+        else if (isinf(y)) {
+            if (fabs(x) == 1.0)
+                return 1.0;
+            else if (y > 0.0 && fabs(x) > 1.0)
+                return y;
+            else if (y < 0.0 && fabs(x) < 1.0) {
+                #if PY_VERSION_HEX < 0x030B0000
+                if (x == 0.0) { /* 0**-inf: divide-by-zero */
+                    return CPy_DomainError();
+                }
+                #endif
+                return -y; /* result is +inf */
+            } else
+                return 0.0;
+        }
+    }
+    double r = pow(x, y);
+    if (!isfinite(r)) {
+        if (isnan(r)) {
+            return CPy_DomainError();
+        }
+        /*
+           an infinite result here arises either from:
+           (A) (+/-0.)**negative (-> divide-by-zero)
+           (B) overflow of x**y with x and y finite
+        */
+        else if (isinf(r)) {
+            if (x == 0.0)
+                return CPy_DomainError();
+            else
+                return CPy_MathRangeError();
+        }
+    }
+    return r;
+}