From eb2dadb0f6684e361962d374e6211f1333ad656a Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Sun, 16 Aug 2026 11:34:42 -0700 Subject: [PATCH 1/9] Remove `_Py_ADJUST_ERANGE*()` functions. --- Include/internal/pycore_pymath.h | 49 -------------------------------- Objects/complexobject.c | 15 +++++++--- Objects/floatobject.c | 12 +++++--- 3 files changed, 19 insertions(+), 57 deletions(-) diff --git a/Include/internal/pycore_pymath.h b/Include/internal/pycore_pymath.h index 532c5ceafb56395..f7c244a70cf9d62 100644 --- a/Include/internal/pycore_pymath.h +++ b/Include/internal/pycore_pymath.h @@ -8,55 +8,6 @@ extern "C" { # error "this header requires Py_BUILD_CORE define" #endif - -/* _Py_ADJUST_ERANGE1(x) - * _Py_ADJUST_ERANGE2(x, y) - * Set errno to 0 before calling a libm function, and invoke one of these - * macros after, passing the function result(s) (_Py_ADJUST_ERANGE2 is useful - * for functions returning complex results). This makes two kinds of - * adjustments to errno: (A) If it looks like the platform libm set - * errno=ERANGE due to underflow, clear errno. (B) If it looks like the - * platform libm overflowed but didn't set errno, force errno to ERANGE. In - * effect, we're trying to force a useful implementation of C89 errno - * behavior. - * Caution: - * This isn't reliable. C99 no longer requires libm to set errno under - * any exceptional condition, but does require +- HUGE_VAL return - * values on overflow. A 754 box *probably* maps HUGE_VAL to a - * double infinity, and we're cool if that's so, unless the input - * was an infinity and an infinity is the expected result. A C89 - * system sets errno to ERANGE, so we check for that too. We're - * out of luck if a C99 754 box doesn't map HUGE_VAL to +Inf, or - * if the returned result is a NaN, or if a C89 box returns HUGE_VAL - * in non-overflow cases. - */ -static inline void _Py_ADJUST_ERANGE1(double x) -{ - if (errno == 0) { - if (x == INFINITY || x == -INFINITY) { - errno = ERANGE; - } - } - else if (errno == ERANGE && x == 0.0) { - errno = 0; - } -} - -static inline void _Py_ADJUST_ERANGE2(double x, double y) -{ - if (x == INFINITY || x == -INFINITY || - y == INFINITY || y == -INFINITY) - { - if (errno == 0) { - errno = ERANGE; - } - } - else if (errno == ERANGE) { - errno = 0; - } -} - - //--- HAVE_PY_SET_53BIT_PRECISION macro ------------------------------------ // // The functions _Py_dg_strtod() and _Py_dg_dtoa() in Python/dtoa.c (which are diff --git a/Objects/complexobject.c b/Objects/complexobject.c index 3612c2699a557db..085bdc51abaef31 100644 --- a/Objects/complexobject.c +++ b/Objects/complexobject.c @@ -11,7 +11,6 @@ #include "pycore_freelist.h" // _Py_FREELIST_FREE(), _Py_FREELIST_POP() #include "pycore_long.h" // _PyLong_GetZero() #include "pycore_object.h" // _PyObject_Init() -#include "pycore_pymath.h" // _Py_ADJUST_ERANGE2() #define _PyComplexObject_CAST(op) ((PyComplexObject *)(op)) @@ -333,7 +332,16 @@ _Py_c_pow(Py_complex a, Py_complex b) r.real = len*cos(phase); r.imag = len*sin(phase); - _Py_ADJUST_ERANGE2(r.real, r.imag); + if (isinf(r.real) || isinf(r.imag)) { + if (errno == 0) { + /* It looks like overflow, but libm didn't set errno. */ + errno = ERANGE; + } + } + else if (errno == ERANGE) { + /* It looks like libm set errno because of underflow. */ + errno = 0; + } } return r; } @@ -753,7 +761,6 @@ complex_pow(PyObject *v, PyObject *w, PyObject *z) // a faster and more accurate algorithm. if (b.imag == 0.0 && b.real == floor(b.real) && fabs(b.real) <= 100.0) { p = c_powi(a, (long)b.real); - _Py_ADJUST_ERANGE2(p.real, p.imag); } else { p = _Py_c_pow(a, b); @@ -764,7 +771,7 @@ complex_pow(PyObject *v, PyObject *w, PyObject *z) "zero to a negative or complex power"); return NULL; } - else if (errno == ERANGE) { + else if (isinf(p.real) || isinf(p.imag)) { PyErr_SetString(PyExc_OverflowError, "complex exponentiation"); return NULL; diff --git a/Objects/floatobject.c b/Objects/floatobject.c index 17e6a729dcd83fc..cb72007205f64e4 100644 --- a/Objects/floatobject.c +++ b/Objects/floatobject.c @@ -787,13 +787,17 @@ float_pow(PyObject *v, PyObject *w, PyObject *z) */ errno = 0; ix = pow(iv, iw); - _Py_ADJUST_ERANGE1(ix); if (negate_result) ix = -ix; - if (errno != 0) { - /* We don't expect any errno value other than ERANGE, but - * the range of libm bugs appears unbounded. + if (isinf(ix)) { + PyErr_SetString(PyExc_OverflowError, + "Numerical result out of range"); + return NULL; + } + if (ix != 0.0 && errno != 0) { + /* We don't expect any errno value other than ERANGE (for underflow, + * which should mean ix == 0.0), but libm may have bugs. */ PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError : PyExc_ValueError); From 087a86ea4fdd74405f419b0878b1c1e97baa68dd Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Tue, 18 Aug 2026 17:04:39 -0700 Subject: [PATCH 2/9] Remove `errno` from floatobject and complexobject more aggressively --- Objects/complexobject.c | 13 +------------ Objects/floatobject.c | 23 ++++------------------- 2 files changed, 5 insertions(+), 31 deletions(-) diff --git a/Objects/complexobject.c b/Objects/complexobject.c index 085bdc51abaef31..f51d13937fbdd30 100644 --- a/Objects/complexobject.c +++ b/Objects/complexobject.c @@ -331,17 +331,6 @@ _Py_c_pow(Py_complex a, Py_complex b) } r.real = len*cos(phase); r.imag = len*sin(phase); - - if (isinf(r.real) || isinf(r.imag)) { - if (errno == 0) { - /* It looks like overflow, but libm didn't set errno. */ - errno = ERANGE; - } - } - else if (errno == ERANGE) { - /* It looks like libm set errno because of underflow. */ - errno = 0; - } } return r; } @@ -773,7 +762,7 @@ complex_pow(PyObject *v, PyObject *w, PyObject *z) } else if (isinf(p.real) || isinf(p.imag)) { PyErr_SetString(PyExc_OverflowError, - "complex exponentiation"); + "complex exponentiation result out of range"); return NULL; } return PyComplex_FromCComplex(p); diff --git a/Objects/floatobject.c b/Objects/floatobject.c index cb72007205f64e4..59c549bf38d97fb 100644 --- a/Objects/floatobject.c +++ b/Objects/floatobject.c @@ -785,22 +785,13 @@ float_pow(PyObject *v, PyObject *w, PyObject *z) * positive and not equal to 1.0. We finally allow * the platform pow to step in and do the rest. */ - errno = 0; ix = pow(iv, iw); if (negate_result) ix = -ix; if (isinf(ix)) { PyErr_SetString(PyExc_OverflowError, - "Numerical result out of range"); - return NULL; - } - if (ix != 0.0 && errno != 0) { - /* We don't expect any errno value other than ERANGE (for underflow, - * which should mean ix == 0.0), but libm may have bugs. - */ - PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError : - PyExc_ValueError); + "float exponentiation result out of range"); return NULL; } return PyFloat_FromDouble(ix); @@ -840,20 +831,14 @@ float_is_integer_impl(PyObject *self) /*[clinic end generated code: output=7112acf95a4d31ea input=311810d3f777e10d]*/ { double x = PyFloat_AsDouble(self); - PyObject *o; if (x == -1.0 && PyErr_Occurred()) return NULL; if (!isfinite(x)) Py_RETURN_FALSE; - errno = 0; - o = (floor(x) == x) ? Py_True : Py_False; - if (errno != 0) { - PyErr_SetFromErrno(errno == ERANGE ? PyExc_OverflowError : - PyExc_ValueError); - return NULL; - } - return Py_NewRef(o); + if (floor(x) == x) + Py_RETURN_TRUE; + Py_RETURN_FALSE; } /*[clinic input] From de5f38230d8781b3fcea1e657c3f6cc9b5fe8cb7 Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Tue, 18 Aug 2026 22:38:19 -0700 Subject: [PATCH 3/9] Ignore `errno` in C API complex pow overflow test --- Lib/test/test_capi/test_complex.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/Lib/test/test_capi/test_complex.py b/Lib/test/test_capi/test_complex.py index c3189a67cc7e2d3..964fa8a872b5d83 100644 --- a/Lib/test/test_capi/test_complex.py +++ b/Lib/test/test_capi/test_complex.py @@ -271,10 +271,8 @@ def test_py_c_pow(self): self.assertEqual(_py_c_pow(0j, -1)[1], errno.EDOM) self.assertEqual(_py_c_pow(0j, 1j)[1], errno.EDOM) max_num = DBL_MAX+1j - self.assertEqual(_py_c_pow(max_num, max_num), - (complex(INF, INF), errno.ERANGE)) - self.assertEqual(_py_c_pow(max_num, 2), - (complex(INF, INF), errno.ERANGE)) + self.assertEqual(_py_c_pow(max_num, max_num)[0], complex(INF, INF)) + self.assertEqual(_py_c_pow(max_num, 2)[0], complex(INF, INF)) def test_py_c_abs(self): From ed0cd08023ae0131bcee1efa3843e1f82eb687ef Mon Sep 17 00:00:00 2001 From: "hpkfft.com" Date: Tue, 18 Aug 2026 22:44:27 -0700 Subject: [PATCH 4/9] Update `float_is_integer_impl()` Co-authored-by: Sergey B Kirpichev --- Objects/floatobject.c | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/Objects/floatobject.c b/Objects/floatobject.c index 59c549bf38d97fb..e1ef3b937421d69 100644 --- a/Objects/floatobject.c +++ b/Objects/floatobject.c @@ -836,9 +836,7 @@ float_is_integer_impl(PyObject *self) return NULL; if (!isfinite(x)) Py_RETURN_FALSE; - if (floor(x) == x) - Py_RETURN_TRUE; - Py_RETURN_FALSE; + return PyBool_FromLong(floor(x) == x); } /*[clinic input] From fbb466e50e5b8286e58354decdd5d61756540977 Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Wed, 19 Aug 2026 14:51:43 -0700 Subject: [PATCH 5/9] Revert "Ignore `errno` in C API complex pow overflow test" This reverts commit de5f38230d8781b3fcea1e657c3f6cc9b5fe8cb7. --- Lib/test/test_capi/test_complex.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/Lib/test/test_capi/test_complex.py b/Lib/test/test_capi/test_complex.py index 964fa8a872b5d83..c3189a67cc7e2d3 100644 --- a/Lib/test/test_capi/test_complex.py +++ b/Lib/test/test_capi/test_complex.py @@ -271,8 +271,10 @@ def test_py_c_pow(self): self.assertEqual(_py_c_pow(0j, -1)[1], errno.EDOM) self.assertEqual(_py_c_pow(0j, 1j)[1], errno.EDOM) max_num = DBL_MAX+1j - self.assertEqual(_py_c_pow(max_num, max_num)[0], complex(INF, INF)) - self.assertEqual(_py_c_pow(max_num, 2)[0], complex(INF, INF)) + self.assertEqual(_py_c_pow(max_num, max_num), + (complex(INF, INF), errno.ERANGE)) + self.assertEqual(_py_c_pow(max_num, 2), + (complex(INF, INF), errno.ERANGE)) def test_py_c_abs(self): From 980cc961d5c93f2583a656f2771a54c8e8445452 Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Wed, 19 Aug 2026 14:59:34 -0700 Subject: [PATCH 6/9] Another idea to avoid errno in complex_pow --- Objects/complexobject.c | 54 +++++++++++++++++++++++++++++++---------- 1 file changed, 41 insertions(+), 13 deletions(-) diff --git a/Objects/complexobject.c b/Objects/complexobject.c index f51d13937fbdd30..db7715322f04ad4 100644 --- a/Objects/complexobject.c +++ b/Objects/complexobject.c @@ -305,8 +305,8 @@ _Py_rc_quot(double a, Py_complex b) #pragma optimize("", on) #endif -Py_complex -_Py_c_pow(Py_complex a, Py_complex b) +static Py_complex +c_pow(Py_complex a, Py_complex b, int *e) { Py_complex r; double vabs,len,at,phase; @@ -316,7 +316,7 @@ _Py_c_pow(Py_complex a, Py_complex b) } else if (a.real == 0. && a.imag == 0.) { if (b.imag != 0. || b.real < 0.) - errno = EDOM; + *e = EDOM; r.real = 0.; r.imag = 0.; } @@ -331,10 +331,26 @@ _Py_c_pow(Py_complex a, Py_complex b) } r.real = len*cos(phase); r.imag = len*sin(phase); + if ((isinf(r.real) || isinf(r.imag)) + && isfinite(a.real) && isfinite(a.imag) + && isfinite(b.real) && isfinite(b.imag)) { + *e = ERANGE; + } } return r; } +Py_complex +_Py_c_pow(Py_complex a, Py_complex b) +{ + Py_complex r; + int e = 0; + int saved_errno = errno; + r = c_pow(a, b, &e); + errno = (e) ? e : saved_errno; + return r; +} + static Py_complex c_powu(Py_complex x, long n) { @@ -352,13 +368,25 @@ c_powu(Py_complex x, long n) } static Py_complex -c_powi(Py_complex x, long n) +c_powi(Py_complex x, long n, int *e) { - if (n > 0) - return c_powu(x,n); - else - return _Py_c_quot(c_1, c_powu(x,-n)); + Py_complex r; + if (x.real == 0. && x.imag == 0. && n != 0) { + if (n < 0) + *e = EDOM; + r.real = 0.; + r.imag = 0.; + } + else { + r = (n < 0) ? _Py_rc_quot(1.0, c_powu(x, -n)) : c_powu(x, n); + } + + if ((isinf(r.real) || isinf(r.imag)) + && isfinite(x.real) && isfinite(x.imag)) { + *e = ERANGE; + } + return r; } double @@ -745,22 +773,22 @@ complex_pow(PyObject *v, PyObject *w, PyObject *z) PyErr_SetString(PyExc_ValueError, "complex modulo"); return NULL; } - errno = 0; + int e = 0; // Check whether the exponent has a small integer value, and if so use // a faster and more accurate algorithm. if (b.imag == 0.0 && b.real == floor(b.real) && fabs(b.real) <= 100.0) { - p = c_powi(a, (long)b.real); + p = c_powi(a, (long)b.real, &e); } else { - p = _Py_c_pow(a, b); + p = c_pow(a, b, &e); } - if (errno == EDOM) { + if (e == EDOM) { PyErr_SetString(PyExc_ZeroDivisionError, "zero to a negative or complex power"); return NULL; } - else if (isinf(p.real) || isinf(p.imag)) { + else if (e == ERANGE) { PyErr_SetString(PyExc_OverflowError, "complex exponentiation result out of range"); return NULL; From eb6f38e89a87fc6e5f98311183a65a4a62265e5c Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Thu, 20 Aug 2026 23:53:43 -0700 Subject: [PATCH 7/9] Code as benchmarked for comment in Issue #156145 --- Objects/complexobject.c | 192 ++++++++++++++++++++-------------------- 1 file changed, 94 insertions(+), 98 deletions(-) diff --git a/Objects/complexobject.c b/Objects/complexobject.c index db7715322f04ad4..c26d87a42182060 100644 --- a/Objects/complexobject.c +++ b/Objects/complexobject.c @@ -165,8 +165,8 @@ _Py_rc_prod(double a, Py_complex b) #ifdef _M_ARM64 #pragma optimize("", off) #endif -Py_complex -_Py_c_quot(Py_complex a, Py_complex b) +static Py_complex +c_quot(Py_complex a, Py_complex b) { /****************************************************************** This was the original algorithm. It's grossly prone to spurious @@ -196,16 +196,11 @@ _Py_c_quot(Py_complex a, Py_complex b) if (abs_breal >= abs_bimag) { /* divide tops and bottom by b.real */ - if (abs_breal == 0.0) { - errno = EDOM; - r.real = r.imag = 0.0; - } - else { - const double ratio = b.imag / b.real; - const double denom = b.real + b.imag * ratio; - r.real = (a.real + a.imag * ratio) / denom; - r.imag = (a.imag - a.real * ratio) / denom; - } + /* If b.real == b.imag == 0.0, r.real and r.imag will be NaN. */ + const double ratio = b.imag / b.real; + const double denom = b.real + b.imag * ratio; + r.real = (a.real + a.imag * ratio) / denom; + r.imag = (a.imag - a.real * ratio) / denom; } else if (abs_bimag >= abs_breal) { /* divide tops and bottom by b.imag */ @@ -245,21 +240,30 @@ _Py_c_quot(Py_complex a, Py_complex b) } Py_complex -_Py_cr_quot(Py_complex a, double b) +_Py_c_quot(Py_complex a, Py_complex b) { - Py_complex r = a; - if (b) { - r.real /= b; - r.imag /= b; + Py_complex r; + if (b.real == 0.0 && b.imag == 0.0) { + errno = EDOM; + r.real = r.imag = 0.0; // As documented, so shall it be done. } else { - errno = EDOM; - r.real = r.imag = 0.0; + r = c_quot(a, b); } return r; } -/* an equivalent of _Py_c_quot() function, when 1st argument is real */ +Py_complex +_Py_cr_quot(Py_complex a, double b) +{ + Py_complex r = a; + // May overflow or may result in NaN. + r.real /= b; + r.imag /= b; + return r; +} + +/* an equivalent of c_quot() function, when 1st argument is real */ Py_complex _Py_rc_quot(double a, Py_complex b) { @@ -268,16 +272,10 @@ _Py_rc_quot(double a, Py_complex b) const double abs_bimag = b.imag < 0 ? -b.imag : b.imag; if (abs_breal >= abs_bimag) { - if (abs_breal == 0.0) { - errno = EDOM; - r.real = r.imag = 0.0; - } - else { - const double ratio = b.imag / b.real; - const double denom = b.real + b.imag * ratio; - r.real = a / denom; - r.imag = (-a * ratio) / denom; - } + const double ratio = b.imag / b.real; + const double denom = b.real + b.imag * ratio; + r.real = a / denom; + r.imag = (-a * ratio) / denom; } else if (abs_bimag >= abs_breal) { const double ratio = b.real / b.imag; @@ -305,6 +303,22 @@ _Py_rc_quot(double a, Py_complex b) #pragma optimize("", on) #endif +static Py_complex +c_powu(Py_complex x, long n) +{ + Py_complex r, p; + long mask = 1; + r = c_1; + p = x; + while (mask > 0 && n >= mask) { + if (n & mask) + r = _Py_c_prod(r,p); + mask <<= 1; + p = _Py_c_prod(p,p); + } + return r; +} + static Py_complex c_pow(Py_complex a, Py_complex b, int *e) { @@ -320,6 +334,17 @@ c_pow(Py_complex a, Py_complex b, int *e) r.real = 0.; r.imag = 0.; } + else if (b.imag == 0.0 && b.real == floor(b.real) + && fabs(b.real) <= 100.0) { + // The exponent is a small integer value, so use a faster and + // more accurate algorithm. + long n = (long)b.real; + r = (n < 0) ? _Py_rc_quot(1.0, c_powu(a, -n)) : c_powu(a, n); + if ((isinf(r.real) || isinf(r.imag)) + && isfinite(a.real) && isfinite(a.imag)) { + *e = ERANGE; + } + } else { vabs = hypot(a.real,a.imag); len = pow(vabs,b.real); @@ -347,52 +372,14 @@ _Py_c_pow(Py_complex a, Py_complex b) int e = 0; int saved_errno = errno; r = c_pow(a, b, &e); + // Sets errno = EDOM on invalid, ERANGE on overflow, else unchanged. errno = (e) ? e : saved_errno; return r; } -static Py_complex -c_powu(Py_complex x, long n) +static double +c_abs(Py_complex z, int *e) { - Py_complex r, p; - long mask = 1; - r = c_1; - p = x; - while (mask > 0 && n >= mask) { - if (n & mask) - r = _Py_c_prod(r,p); - mask <<= 1; - p = _Py_c_prod(p,p); - } - return r; -} - -static Py_complex -c_powi(Py_complex x, long n, int *e) -{ - Py_complex r; - - if (x.real == 0. && x.imag == 0. && n != 0) { - if (n < 0) - *e = EDOM; - r.real = 0.; - r.imag = 0.; - } - else { - r = (n < 0) ? _Py_rc_quot(1.0, c_powu(x, -n)) : c_powu(x, n); - } - - if ((isinf(r.real) || isinf(r.imag)) - && isfinite(x.real) && isfinite(x.imag)) { - *e = ERANGE; - } - return r; -} - -double -_Py_c_abs(Py_complex z) -{ - /* sets errno = ERANGE on overflow; otherwise errno = 0 */ double result; if (!isfinite(z.real) || !isfinite(z.imag)) { @@ -401,12 +388,10 @@ _Py_c_abs(Py_complex z) NaN. */ if (isinf(z.real)) { result = fabs(z.real); - errno = 0; return result; } if (isinf(z.imag)) { result = fabs(z.imag); - errno = 0; return result; } /* either the real or imaginary part is a NaN, @@ -415,12 +400,22 @@ _Py_c_abs(Py_complex z) } result = hypot(z.real, z.imag); if (!isfinite(result)) - errno = ERANGE; - else - errno = 0; + *e = ERANGE; return result; } +double +_Py_c_abs(Py_complex z) +{ + double r; + int e = 0; + int saved_errno = errno; + r = c_abs(z, &e); + // Sets errno = ERANGE on overflow, else errno is unchanged. + errno = (e) ? e : saved_errno; + return r; +} + static PyObject * complex_subtype_from_c_complex(PyTypeObject *type, Py_complex cval) { @@ -718,15 +713,17 @@ real_to_complex(PyObject **pobj, Py_complex *pc) See C11's Annex G, sections G.5.1 and G.5.2. */ -#define COMPLEX_BINOP(NAME, FUNC) \ +#define COMPLEX_BINOP(NAME, FUNC, IS_DIVISION) \ static PyObject * \ complex_##NAME(PyObject *v, PyObject *w) \ { \ Py_complex a; \ - errno = 0; \ if (PyComplex_Check(w)) { \ Py_complex b = ((PyComplexObject *)w)->cval; \ if (PyComplex_Check(v)) { \ + if (IS_DIVISION && b.real == 0.0 \ + && b.imag == 0.0) \ + goto DivByZero; \ a = ((PyComplexObject *)v)->cval; \ a = _Py_c_##FUNC(a, b); \ } \ @@ -734,6 +731,9 @@ real_to_complex(PyObject **pobj, Py_complex *pc) return v; \ } \ else { \ + if (IS_DIVISION && b.real == 0.0 \ + && b.imag == 0.0) \ + goto DivByZero; \ a = _Py_rc_##FUNC(a.real, b); \ } \ } \ @@ -746,20 +746,22 @@ real_to_complex(PyObject **pobj, Py_complex *pc) if (real_to_double(&w, &b) < 0) { \ return w; \ } \ + if (IS_DIVISION && b == 0.0) \ + goto DivByZero; \ a = _Py_cr_##FUNC(a, b); \ } \ - if (errno == EDOM) { \ - PyErr_SetString(PyExc_ZeroDivisionError, \ - "division by zero"); \ - return NULL; \ - } \ return PyComplex_FromCComplex(a); \ - } + \ + DivByZero: \ + PyErr_SetString(PyExc_ZeroDivisionError, \ + "division by zero"); \ + return NULL; \ + } -COMPLEX_BINOP(add, sum) -COMPLEX_BINOP(mul, prod) -COMPLEX_BINOP(sub, diff) -COMPLEX_BINOP(div, quot) +COMPLEX_BINOP(add, sum, 0) +COMPLEX_BINOP(mul, prod, 0) +COMPLEX_BINOP(sub, diff, 0) +COMPLEX_BINOP(div, quot, 1) static PyObject * complex_pow(PyObject *v, PyObject *w, PyObject *z) @@ -774,14 +776,7 @@ complex_pow(PyObject *v, PyObject *w, PyObject *z) return NULL; } int e = 0; - // Check whether the exponent has a small integer value, and if so use - // a faster and more accurate algorithm. - if (b.imag == 0.0 && b.real == floor(b.real) && fabs(b.real) <= 100.0) { - p = c_powi(a, (long)b.real, &e); - } - else { - p = c_pow(a, b, &e); - } + p = c_pow(a, b, &e); if (e == EDOM) { PyErr_SetString(PyExc_ZeroDivisionError, @@ -820,8 +815,9 @@ static PyObject * complex_abs(PyObject *op) { PyComplexObject *v = _PyComplexObject_CAST(op); - double result = _Py_c_abs(v->cval); - if (errno == ERANGE) { + int e = 0; + double result = c_abs(v->cval, &e); + if (e == ERANGE) { PyErr_SetString(PyExc_OverflowError, "absolute value too large"); return NULL; From eb4aee6c1c52c1a5ec59d4bab5b60a5eb8682fd4 Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Wed, 26 Aug 2026 10:31:41 -0700 Subject: [PATCH 8/9] Add complex tests from PR 118000 by Sergey. Co-authored-by: Sergey B Kirpichev --- Lib/test/test_complex.py | 61 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 61 insertions(+) diff --git a/Lib/test/test_complex.py b/Lib/test/test_complex.py index bb307191dffcc14..4ac9ba5def4544a 100644 --- a/Lib/test/test_complex.py +++ b/Lib/test/test_complex.py @@ -84,6 +84,10 @@ def assertClose(self, x, y, eps=1e-9): # check that relative difference < eps self.assertTrue(abs(x-y)/abs(y) < eps) + def assertSameSign(self, x, y): + if copysign(1., x) != copysign(1., y): + self.fail(f'{x!r} and {y!r} have different signs') + def check_div(self, x, y): """Compute complex z=x*y, and check that z/x==y and z/y==x.""" z = x * y @@ -446,6 +450,63 @@ def test_pow_with_small_integer_exponents(self): self.assertEqual(str(float_pow), str(int_pow)) self.assertEqual(str(complex_pow), str(int_pow)) + # Check that complex numbers with special components + # are correctly handled. + values = [complex(x, y) + for x in [5, -5, +0.0, -0.0, INF, -INF, NAN] + for y in [12, -12, +0.0, -0.0, INF, -INF, NAN]] + for c in values: + with self.subTest(value=c): + self.assertComplexesAreIdentical(c**0, complex(1, +0.0)) + self.assertComplexesAreIdentical(c**1, c) + self.assertComplexesAreIdentical(c**2, c*c) + self.assertComplexesAreIdentical(c**3, c*(c*c)) + self.assertComplexesAreIdentical(c**3, (c*c)*c) + if not c: + continue + for n in range(1, 9): + with self.subTest(exponent=-n): + self.assertComplexesAreIdentical(c**-n, 1/(c**n)) + + # Special cases for complex division. + for x in [+2, -2]: + for y in [+0.0, -0.0]: + c = complex(x, y) + with self.subTest(value=c): + self.assertComplexesAreIdentical(c**-1, complex(1/x, -y)) + c = complex(y, x) + with self.subTest(value=c): + self.assertComplexesAreIdentical(c**-1, complex(y, -1/x)) + for x in [+INF, -INF]: + for y in [+1, -1]: + c = complex(x, y) + with self.subTest(value=c): + self.assertComplexesAreIdentical(c**-1, complex(1/x, -0.0*y)) + self.assertComplexesAreIdentical(c**-2, complex(0.0, -y/x)) + c = complex(y, x) + with self.subTest(value=c): + self.assertComplexesAreIdentical(c**-1, complex(+0.0*y, -1/x)) + self.assertComplexesAreIdentical(c**-2, complex(-0.0, -y/x)) + + # Test that zeros have the same sign as small non-zero values. + eps = 1e-11 + pairs = [(complex(x, y), complex(x, copysign(0.0, y))) + for x in [+1, -1] for y in [+eps, -eps]] + pairs += [(complex(y, x), complex(copysign(0.0, y), x)) + for x in [+1, -1] for y in [+eps, -eps]] + for c1, c2 in pairs: + for n in exponents: + with self.subTest(value=c1, exponent=n): + r1 = c1**n + r2 = c2**n + self.assertClose(r1, r2) + self.assertSameSign(r1.real, r2.real) + self.assertSameSign(r1.imag, r2.imag) + self.assertNotEqual(r1.real, 0.0) + if n != 0: + self.assertNotEqual(r1.imag, 0.0) + self.assertTrue(r2.real == 0.0 or r2.imag == 0.0) + def test_boolcontext(self): for i in range(100): self.assertTrue(complex(random() + 1e-6, random() + 1e-6)) From d697320d1e7b4cf4438792a8c679cbf0bc192afd Mon Sep 17 00:00:00 2001 From: Paul Caprioli Date: Fri, 28 Aug 2026 21:33:21 -0700 Subject: [PATCH 9/9] Rewrite c_powi() --- Objects/complexobject.c | 34 +++++++++++++++++----------------- 1 file changed, 17 insertions(+), 17 deletions(-) diff --git a/Objects/complexobject.c b/Objects/complexobject.c index c26d87a42182060..2e649f6d7849772 100644 --- a/Objects/complexobject.c +++ b/Objects/complexobject.c @@ -25,8 +25,6 @@ class complex "PyComplexObject *" "&PyComplex_Type" /* elementary operations on complex numbers */ -static Py_complex c_1 = {1., 0.}; - Py_complex _Py_c_sum(Py_complex a, Py_complex b) { @@ -304,18 +302,22 @@ _Py_rc_quot(double a, Py_complex b) #endif static Py_complex -c_powu(Py_complex x, long n) +c_powi(Py_complex x, long n) { - Py_complex r, p; - long mask = 1; - r = c_1; - p = x; - while (mask > 0 && n >= mask) { - if (n & mask) - r = _Py_c_prod(r,p); - mask <<= 1; - p = _Py_c_prod(p,p); + Py_complex r = x; + long absn = (n < 0) ? -n : n; + + assert(absn > 0); + if (--absn & 0x1) + r = _Py_c_prod(r, x); + while (absn >>= 1) { + x = _Py_c_prod(x, x); + if (absn & 0x1) + r = _Py_c_prod(r, x); } + + if (n < 0) + r = _Py_rc_quot(1.0, r); return r; } @@ -328,18 +330,16 @@ c_pow(Py_complex a, Py_complex b, int *e) r.real = 1.; r.imag = 0.; } - else if (a.real == 0. && a.imag == 0.) { - if (b.imag != 0. || b.real < 0.) - *e = EDOM; + else if (a.real == 0. && a.imag == 0. && (b.imag != 0. || b.real < 0.)) { r.real = 0.; r.imag = 0.; + *e = EDOM; } else if (b.imag == 0.0 && b.real == floor(b.real) && fabs(b.real) <= 100.0) { // The exponent is a small integer value, so use a faster and // more accurate algorithm. - long n = (long)b.real; - r = (n < 0) ? _Py_rc_quot(1.0, c_powu(a, -n)) : c_powu(a, n); + r = c_powi(a, (long)b.real); if ((isinf(r.real) || isinf(r.imag)) && isfinite(a.real) && isfinite(a.imag)) { *e = ERANGE;