From 6e39fa19555043588910d10f1fe677cf6b04d77e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Matthias=20G=C3=B6rgens?= Date: Sat, 20 May 2023 04:03:49 +0800 Subject: [PATCH] gh-94906: Support multiple steps in math.nextafter (#103881) This PR updates `math.nextafter` to add a new `steps` argument. The behaviour is as though `math.nextafter` had been called `steps` times in succession. --------- Co-authored-by: Mark Dickinson --- Doc/library/math.rst | 9 +- .../pycore_global_objects_fini_generated.h | 1 + Include/internal/pycore_global_strings.h | 1 + .../internal/pycore_runtime_init_generated.h | 1 + .../internal/pycore_unicodeobject_generated.h | 3 + Lib/test/test_math.py | 20 +++- Lib/test/test_math_property.py | 41 +++++++ ...2-07-16-17-15-29.gh-issue-94906.C4G8DG.rst | 1 + Modules/clinic/mathmodule.c.h | 55 +++++++-- Modules/mathmodule.c | 109 +++++++++++++++++- 10 files changed, 223 insertions(+), 18 deletions(-) create mode 100644 Lib/test/test_math_property.py create mode 100644 Misc/NEWS.d/next/Library/2022-07-16-17-15-29.gh-issue-94906.C4G8DG.rst diff --git a/Doc/library/math.rst b/Doc/library/math.rst index 797f32408ea..9e58b552576 100644 --- a/Doc/library/math.rst +++ b/Doc/library/math.rst @@ -224,11 +224,11 @@ Number-theoretic and representation functions of *x* and are floats. -.. function:: nextafter(x, y) +.. function:: nextafter(x, y, steps=1) - Return the next floating-point value after *x* towards *y*. + Return the floating-point value *steps* steps after *x* towards *y*. - If *x* is equal to *y*, return *y*. + If *x* is equal to *y*, return *y*, unless *steps* is zero. Examples: @@ -239,6 +239,9 @@ Number-theoretic and representation functions See also :func:`math.ulp`. + .. versionchanged:: 3.12 + Added the *steps* argument. + .. versionadded:: 3.9 .. function:: perm(n, k=None) diff --git a/Include/internal/pycore_global_objects_fini_generated.h b/Include/internal/pycore_global_objects_fini_generated.h index eeaa2ad96c9..8ca3545d8b3 100644 --- a/Include/internal/pycore_global_objects_fini_generated.h +++ b/Include/internal/pycore_global_objects_fini_generated.h @@ -1193,6 +1193,7 @@ _PyStaticObjects_CheckRefcnt(PyInterpreterState *interp) { _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(stdin)); _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(stdout)); _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(step)); + _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(steps)); _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(store_name)); _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(strategy)); _PyStaticObject_CheckRefcnt((PyObject *)&_Py_ID(strftime)); diff --git a/Include/internal/pycore_global_strings.h b/Include/internal/pycore_global_strings.h index 5cc790d126b..8e429bbfa26 100644 --- a/Include/internal/pycore_global_strings.h +++ b/Include/internal/pycore_global_strings.h @@ -681,6 +681,7 @@ struct _Py_global_strings { STRUCT_FOR_ID(stdin) STRUCT_FOR_ID(stdout) STRUCT_FOR_ID(step) + STRUCT_FOR_ID(steps) STRUCT_FOR_ID(store_name) STRUCT_FOR_ID(strategy) STRUCT_FOR_ID(strftime) diff --git a/Include/internal/pycore_runtime_init_generated.h b/Include/internal/pycore_runtime_init_generated.h index 0cb24a92dff..3edf076696d 100644 --- a/Include/internal/pycore_runtime_init_generated.h +++ b/Include/internal/pycore_runtime_init_generated.h @@ -1187,6 +1187,7 @@ extern "C" { INIT_ID(stdin), \ INIT_ID(stdout), \ INIT_ID(step), \ + INIT_ID(steps), \ INIT_ID(store_name), \ INIT_ID(strategy), \ INIT_ID(strftime), \ diff --git a/Include/internal/pycore_unicodeobject_generated.h b/Include/internal/pycore_unicodeobject_generated.h index fe2479beb8a..0e1f71798a6 100644 --- a/Include/internal/pycore_unicodeobject_generated.h +++ b/Include/internal/pycore_unicodeobject_generated.h @@ -1884,6 +1884,9 @@ _PyUnicode_InitStaticStrings(PyInterpreterState *interp) { string = &_Py_ID(step); assert(_PyUnicode_CheckConsistency(string, 1)); _PyUnicode_InternInPlace(interp, &string); + string = &_Py_ID(steps); + assert(_PyUnicode_CheckConsistency(string, 1)); + _PyUnicode_InternInPlace(interp, &string); string = &_Py_ID(store_name); assert(_PyUnicode_CheckConsistency(string, 1)); _PyUnicode_InternInPlace(interp, &string); diff --git a/Lib/test/test_math.py b/Lib/test/test_math.py index f282434c9a3..2bda6101216 100644 --- a/Lib/test/test_math.py +++ b/Lib/test/test_math.py @@ -2296,11 +2296,20 @@ class MathTests(unittest.TestCase): float.fromhex('0x1.fffffffffffffp-1')) self.assertEqual(math.nextafter(1.0, INF), float.fromhex('0x1.0000000000001p+0')) + self.assertEqual(math.nextafter(1.0, -INF, steps=1), + float.fromhex('0x1.fffffffffffffp-1')) + self.assertEqual(math.nextafter(1.0, INF, steps=1), + float.fromhex('0x1.0000000000001p+0')) + self.assertEqual(math.nextafter(1.0, -INF, steps=3), + float.fromhex('0x1.ffffffffffffdp-1')) + self.assertEqual(math.nextafter(1.0, INF, steps=3), + float.fromhex('0x1.0000000000003p+0')) # x == y: y is returned - self.assertEqual(math.nextafter(2.0, 2.0), 2.0) - self.assertEqualSign(math.nextafter(-0.0, +0.0), +0.0) - self.assertEqualSign(math.nextafter(+0.0, -0.0), -0.0) + for steps in range(1, 5): + self.assertEqual(math.nextafter(2.0, 2.0, steps=steps), 2.0) + self.assertEqualSign(math.nextafter(-0.0, +0.0, steps=steps), +0.0) + self.assertEqualSign(math.nextafter(+0.0, -0.0, steps=steps), -0.0) # around 0.0 smallest_subnormal = sys.float_info.min * sys.float_info.epsilon @@ -2325,6 +2334,11 @@ class MathTests(unittest.TestCase): self.assertIsNaN(math.nextafter(1.0, NAN)) self.assertIsNaN(math.nextafter(NAN, NAN)) + self.assertEqual(1.0, math.nextafter(1.0, INF, steps=0)) + with self.assertRaises(ValueError): + math.nextafter(1.0, INF, steps=-1) + + @requires_IEEE_754 def test_ulp(self): self.assertEqual(math.ulp(1.0), sys.float_info.epsilon) diff --git a/Lib/test/test_math_property.py b/Lib/test/test_math_property.py new file mode 100644 index 00000000000..7d51aa17b4c --- /dev/null +++ b/Lib/test/test_math_property.py @@ -0,0 +1,41 @@ +import functools +import unittest +from math import isnan, nextafter +from test.support import requires_IEEE_754 +from test.support.hypothesis_helper import hypothesis + +floats = hypothesis.strategies.floats +integers = hypothesis.strategies.integers + + +def assert_equal_float(x, y): + assert isnan(x) and isnan(y) or x == y + + +def via_reduce(x, y, steps): + return functools.reduce(nextafter, [y] * steps, x) + + +class NextafterTests(unittest.TestCase): + @requires_IEEE_754 + @hypothesis.given( + x=floats(), + y=floats(), + steps=integers(min_value=0, max_value=2**16)) + def test_count(self, x, y, steps): + assert_equal_float(via_reduce(x, y, steps), + nextafter(x, y, steps=steps)) + + @requires_IEEE_754 + @hypothesis.given( + x=floats(), + y=floats(), + a=integers(min_value=0), + b=integers(min_value=0)) + def test_addition_commutes(self, x, y, a, b): + first = nextafter(x, y, steps=a) + second = nextafter(first, y, steps=b) + combined = nextafter(x, y, steps=a+b) + hypothesis.note(f"{first} -> {second} == {combined}") + + assert_equal_float(second, combined) diff --git a/Misc/NEWS.d/next/Library/2022-07-16-17-15-29.gh-issue-94906.C4G8DG.rst b/Misc/NEWS.d/next/Library/2022-07-16-17-15-29.gh-issue-94906.C4G8DG.rst new file mode 100644 index 00000000000..663343371d1 --- /dev/null +++ b/Misc/NEWS.d/next/Library/2022-07-16-17-15-29.gh-issue-94906.C4G8DG.rst @@ -0,0 +1 @@ +Support multiple steps in :func:`math.nextafter`. Patch by Shantanu Jain and Matthias Gorgens. diff --git a/Modules/clinic/mathmodule.c.h b/Modules/clinic/mathmodule.c.h index bc5bbceb4c9..c16c1b08398 100644 --- a/Modules/clinic/mathmodule.c.h +++ b/Modules/clinic/mathmodule.c.h @@ -826,25 +826,59 @@ exit: } PyDoc_STRVAR(math_nextafter__doc__, -"nextafter($module, x, y, /)\n" +"nextafter($module, x, y, /, *, steps=None)\n" "--\n" "\n" -"Return the next floating-point value after x towards y."); +"Return the floating-point value the given number of steps after x towards y.\n" +"\n" +"If steps is not specified or is None, it defaults to 1.\n" +"\n" +"Raises a TypeError, if x or y is not a double, or if steps is not an integer.\n" +"Raises ValueError if steps is negative."); #define MATH_NEXTAFTER_METHODDEF \ - {"nextafter", _PyCFunction_CAST(math_nextafter), METH_FASTCALL, math_nextafter__doc__}, + {"nextafter", _PyCFunction_CAST(math_nextafter), METH_FASTCALL|METH_KEYWORDS, math_nextafter__doc__}, static PyObject * -math_nextafter_impl(PyObject *module, double x, double y); +math_nextafter_impl(PyObject *module, double x, double y, PyObject *steps); static PyObject * -math_nextafter(PyObject *module, PyObject *const *args, Py_ssize_t nargs) +math_nextafter(PyObject *module, PyObject *const *args, Py_ssize_t nargs, PyObject *kwnames) { PyObject *return_value = NULL; + #if defined(Py_BUILD_CORE) && !defined(Py_BUILD_CORE_MODULE) + + #define NUM_KEYWORDS 1 + static struct { + PyGC_Head _this_is_not_used; + PyObject_VAR_HEAD + PyObject *ob_item[NUM_KEYWORDS]; + } _kwtuple = { + .ob_base = PyVarObject_HEAD_INIT(&PyTuple_Type, NUM_KEYWORDS) + .ob_item = { &_Py_ID(steps), }, + }; + #undef NUM_KEYWORDS + #define KWTUPLE (&_kwtuple.ob_base.ob_base) + + #else // !Py_BUILD_CORE + # define KWTUPLE NULL + #endif // !Py_BUILD_CORE + + static const char * const _keywords[] = {"", "", "steps", NULL}; + static _PyArg_Parser _parser = { + .keywords = _keywords, + .fname = "nextafter", + .kwtuple = KWTUPLE, + }; + #undef KWTUPLE + PyObject *argsbuf[3]; + Py_ssize_t noptargs = nargs + (kwnames ? PyTuple_GET_SIZE(kwnames) : 0) - 2; double x; double y; + PyObject *steps = Py_None; - if (!_PyArg_CheckPositional("nextafter", nargs, 2, 2)) { + args = _PyArg_UnpackKeywords(args, nargs, NULL, kwnames, &_parser, 2, 2, 0, argsbuf); + if (!args) { goto exit; } if (PyFloat_CheckExact(args[0])) { @@ -867,7 +901,12 @@ math_nextafter(PyObject *module, PyObject *const *args, Py_ssize_t nargs) goto exit; } } - return_value = math_nextafter_impl(module, x, y); + if (!noptargs) { + goto skip_optional_kwonly; + } + steps = args[2]; +skip_optional_kwonly: + return_value = math_nextafter_impl(module, x, y, steps); exit: return return_value; @@ -911,4 +950,4 @@ math_ulp(PyObject *module, PyObject *arg) exit: return return_value; } -/*[clinic end generated code: output=a6437a3ba18c486a input=a9049054013a1b77]*/ +/*[clinic end generated code: output=91a0357265a2a553 input=a9049054013a1b77]*/ diff --git a/Modules/mathmodule.c b/Modules/mathmodule.c index f369b2c45ce..f26602d5871 100644 --- a/Modules/mathmodule.c +++ b/Modules/mathmodule.c @@ -3864,13 +3864,20 @@ math.nextafter x: double y: double / + * + steps: object = None -Return the next floating-point value after x towards y. +Return the floating-point value the given number of steps after x towards y. + +If steps is not specified or is None, it defaults to 1. + +Raises a TypeError, if x or y is not a double, or if steps is not an integer. +Raises ValueError if steps is negative. [clinic start generated code]*/ static PyObject * -math_nextafter_impl(PyObject *module, double x, double y) -/*[clinic end generated code: output=750c8266c1c540ce input=02b2d50cd1d9f9b6]*/ +math_nextafter_impl(PyObject *module, double x, double y, PyObject *steps) +/*[clinic end generated code: output=cc6511f02afc099e input=7f2a5842112af2b4]*/ { #if defined(_AIX) if (x == y) { @@ -3885,7 +3892,101 @@ math_nextafter_impl(PyObject *module, double x, double y) return PyFloat_FromDouble(y); } #endif - return PyFloat_FromDouble(nextafter(x, y)); + if (steps == Py_None) { + // fast path: we default to one step. + return PyFloat_FromDouble(nextafter(x, y)); + } + steps = PyNumber_Index(steps); + if (steps == NULL) { + return NULL; + } + assert(PyLong_CheckExact(steps)); + if (_PyLong_IsNegative((PyLongObject *)steps)) { + PyErr_SetString(PyExc_ValueError, + "steps must be a non-negative integer"); + Py_DECREF(steps); + return NULL; + } + + unsigned long long usteps_ull = PyLong_AsUnsignedLongLong(steps); + // Conveniently, uint64_t and double have the same number of bits + // on all the platforms we care about. + // So if an overflow occurs, we can just use UINT64_MAX. + Py_DECREF(steps); + if (usteps_ull >= UINT64_MAX) { + // This branch includes the case where an error occurred, since + // (unsigned long long)(-1) = ULLONG_MAX >= UINT64_MAX. Note that + // usteps_ull can be strictly larger than UINT64_MAX on a machine + // where unsigned long long has width > 64 bits. + if (PyErr_Occurred()) { + if (PyErr_ExceptionMatches(PyExc_OverflowError)) { + PyErr_Clear(); + } + else { + return NULL; + } + } + usteps_ull = UINT64_MAX; + } + assert(usteps_ull <= UINT64_MAX); + uint64_t usteps = (uint64_t)usteps_ull; + + if (usteps == 0) { + return PyFloat_FromDouble(x); + } + if (Py_IS_NAN(x)) { + return PyFloat_FromDouble(x); + } + if (Py_IS_NAN(y)) { + return PyFloat_FromDouble(y); + } + + // We assume that double and uint64_t have the same endianness. + // This is not guaranteed by the C-standard, but it is true for + // all platforms we care about. (The most likely form of violation + // would be a "mixed-endian" double.) + union pun {double f; uint64_t i;}; + union pun ux = {x}, uy = {y}; + if (ux.i == uy.i) { + return PyFloat_FromDouble(x); + } + + const uint64_t sign_bit = 1ULL<<63; + + uint64_t ax = ux.i & ~sign_bit; + uint64_t ay = uy.i & ~sign_bit; + + // opposite signs + if (((ux.i ^ uy.i) & sign_bit)) { + // NOTE: ax + ay can never overflow, because their most significant bit + // ain't set. + if (ax + ay <= usteps) { + return PyFloat_FromDouble(uy.f); + // This comparison has to use <, because <= would get +0.0 vs -0.0 + // wrong. + } else if (ax < usteps) { + union pun result = {.i = (uy.i & sign_bit) | (usteps - ax)}; + return PyFloat_FromDouble(result.f); + } else { + ux.i -= usteps; + return PyFloat_FromDouble(ux.f); + } + // same sign + } else if (ax > ay) { + if (ax - ay >= usteps) { + ux.i -= usteps; + return PyFloat_FromDouble(ux.f); + } else { + return PyFloat_FromDouble(uy.f); + } + } else { + if (ay - ax >= usteps) { + ux.i += usteps; + return PyFloat_FromDouble(ux.f); + } else { + return PyFloat_FromDouble(uy.f); + } + } }