diff --git a/Doc/library/math.rst b/Doc/library/math.rst index efe411e5a43f27d..dfc13742e47aeee 100644 --- a/Doc/library/math.rst +++ b/Doc/library/math.rst @@ -36,6 +36,7 @@ noted otherwise, all return values are floats. :func:`fmax(x, y) ` Maximum of two floating-point values :func:`fmin(x, y) ` Minimum of two floating-point values :func:`fmod(x, y) ` Remainder of division ``x / y`` +:func:`lerp(a, b, t) ` Linear interpolation between *a* and *b* using parameter *t* :func:`modf(x) ` Fractional and integer parts of *x* :func:`remainder(x, y) ` Remainder of *x* with respect to *y* :func:`trunc(x) ` Integer part of *x* @@ -351,6 +352,31 @@ Floating point manipulation functions :func:`frexp`. +.. function:: lerp(a, b, t) + + Linear interpolation between *a* and *b*, using *t* as the interpolation + parameter. Returns ``a + t * (b - a)``, computed in a way that is + guaranteed to be numerically well-behaved: + + * ``lerp(a, b, 0.0)`` returns *a* exactly. + * ``lerp(a, b, 1.0)`` returns *b* exactly. + * ``lerp(a, b, t)`` is monotonic in *t*, for finite *a* and *b* with + ``a <= b`` and finite *t*. + * *t* is not restricted to ``[0.0, 1.0]``; values outside that range + perform extrapolation. + + The naive formula ``a + t * (b - a)`` does not guarantee any of the + above; it can, for example, fail to return *b* exactly when *t* is + ``1.0`` due to floating-point rounding. :func:`lerp` uses the + algorithm standardized for C++20's ``std::lerp``, which fixes these + issues. As a consequence, results involving infinities or NaN follow + from ordinary floating-point arithmetic rules on each branch of that + algorithm, rather than a general guarantee that NaN inputs always + produce a NaN result. + + .. versionadded:: 3.16 + + .. function:: nextafter(x, y, steps=1) Return the floating-point value *steps* steps after *x* towards *y*. diff --git a/Doc/whatsnew/3.16.rst b/Doc/whatsnew/3.16.rst index 662defa709a246c..3d289e59cb0e25d 100644 --- a/Doc/whatsnew/3.16.rst +++ b/Doc/whatsnew/3.16.rst @@ -371,6 +371,11 @@ math 754-2019 and standardized in C23. (Contributed by Jeff Epler in :gh:`150534`.) +* Added :func:`math.lerp` for linear interpolation, using the algorithm + standardized for C++20's ``std::lerp``, which guarantees exact + endpoints and monotonicity. + (Contributed by PhysicistJohn in :gh:`154573`.) + os -- diff --git a/Lib/test/test_math.py b/Lib/test/test_math.py index cc37697189af3e5..89cea77d86e1e61 100644 --- a/Lib/test/test_math.py +++ b/Lib/test/test_math.py @@ -2847,6 +2847,112 @@ def assertIsNegativeZero(self, value): ) +class LerpTests(unittest.TestCase): + """ Tests for math.lerp. """ + + def test_basic(self): + self.assertEqual(math.lerp(0.0, 10.0, 0.5), 5.0) + self.assertEqual(math.lerp(1.0, 3.0, 2.0), 5.0) # extrapolation + self.assertEqual(math.lerp(1.0, 3.0, -1.0), -1.0) # extrapolation + self.assertEqual(math.lerp(-5.0, 5.0, 0.5), 0.0) + + def test_int_arguments_are_accepted(self): + # Clinic converts ints to float, like other math functions. + self.assertEqual(math.lerp(0, 10, 0), 0.0) + self.assertEqual(math.lerp(0, 10, 1), 10.0) + + def test_endpoints_are_exact(self): + # t == 0.0 and t == 1.0 must return the endpoints exactly, even + # where the naive "a + t * (b - a)" formula would round off. + cases = [ + (0.0, 10.0), (10.0, 0.0), (-5.0, 5.0), (1e300, -1e300), + (524560164915884.0, -995787893297778.6), + ] + for a, b in cases: + with self.subTest(a=a, b=b): + self.assertEqual(math.lerp(a, b, 0.0), a) + self.assertEqual(math.lerp(a, b, 1.0), b) + + def test_naive_formula_would_be_inexact_at_t1(self): + # Concrete demonstration of why lerp() special-cases t == 1.0 + # rather than using "a + t * (b - a)" directly. + a = 524560164915884.0 + b = -995787893297778.6 + naive = a + 1.0 * (b - a) + self.assertNotEqual(naive, b) + self.assertEqual(math.lerp(a, b, 1.0), b) + + def test_a_equals_b(self): + # lerp(a, a, t) == a for any finite t, including outside [0, 1]. + for t in [-100.0, -1.0, 0.0, 0.5, 1.0, 2.0, 100.0]: + with self.subTest(t=t): + self.assertEqual(math.lerp(3.5, 3.5, t), 3.5) + + def test_monotonic(self): + random.seed(163875206) + for _ in range(1000): + a = random.uniform(-1e10, 1e10) + b = random.uniform(-1e10, 1e10) + lo, hi = min(a, b), max(a, b) + ts = sorted(random.uniform(-2.0, 2.0) for _ in range(8)) + values = [math.lerp(lo, hi, t) for t in ts] + self.assertEqual(values, sorted(values)) + + def test_matches_cpp20_stdlerp_reference_values(self): + # Bit-for-bit reference values captured from libc++'s std::lerp + # (LLVM main, clang -std=c++20), including edge cases involving + # zero, signed infinities, and NaN. math.lerp intentionally + # implements the same standardized algorithm, so it should match + # exactly, including cases where NaN in 'a' does not propagate + # (a documented characteristic of the algorithm itself, not a + # special case added here). + nan = math.nan + inf = math.inf + reference = [ + ((0.0, 10.0, 0.5), 5.0), + ((10.0, 0.0, 0.0), 10.0), + ((10.0, 0.0, 1.0), 0.0), + ((1.0, 1.0, 100.0), 1.0), + ((1.0, 1.0, -100.0), 1.0), + ((-0.0, 10.0, 0.5), 5.0), + ((1e300, 1e300, 5.0), 1e300), + ((1e308, -1e308, 3.0), -inf), + ((-1e308, 1e308, 3.0), inf), + ((5.0, -5.0, 0.5), 0.0), + ] + for (a, b, t), expected in reference: + with self.subTest(a=a, b=b, t=t): + self.assertEqual(math.lerp(a, b, t), expected) + + nan_cases = [ + (5.0, nan, 0.5), + (0.0, 10.0, nan), + (-inf, inf, 0.5), + ] + for a, b, t in nan_cases: + with self.subTest(a=a, b=b, t=t): + self.assertTrue(math.isnan(math.lerp(a, b, t))) + + # These do NOT produce NaN, even though 'a' or 't' is infinite/NaN + # in one operand: this matches std::lerp's documented behavior, + # where the t == 1 and clamping branches can return a or b + # without evaluating an expression that would propagate NaN. + self.assertEqual(math.lerp(nan, 5.0, 0.5), 5.0) + self.assertEqual(math.lerp(inf, 5.0, 0.5), 5.0) + self.assertEqual(math.lerp(5.0, inf, 0.5), inf) + + def test_wrong_number_of_arguments(self): + self.assertRaises(TypeError, math.lerp) + self.assertRaises(TypeError, math.lerp, 1.0) + self.assertRaises(TypeError, math.lerp, 1.0, 2.0) + self.assertRaises(TypeError, math.lerp, 1.0, 2.0, 3.0, 4.0) + + def test_wrong_argument_type(self): + self.assertRaises(TypeError, math.lerp, 'a', 2.0, 0.5) + self.assertRaises(TypeError, math.lerp, 1.0, 'b', 0.5) + self.assertRaises(TypeError, math.lerp, 1.0, 2.0, 't') + + def load_tests(loader, tests, pattern): from doctest import DocFileSuite tests.addTest(DocFileSuite(os.path.join("mathdata", "ieee754.txt"))) diff --git a/Misc/NEWS.d/next/Library/2026-07-23-23-00-00.gh-issue-154573.math-lerp.rst b/Misc/NEWS.d/next/Library/2026-07-23-23-00-00.gh-issue-154573.math-lerp.rst new file mode 100644 index 000000000000000..f2176bfa5300289 --- /dev/null +++ b/Misc/NEWS.d/next/Library/2026-07-23-23-00-00.gh-issue-154573.math-lerp.rst @@ -0,0 +1,5 @@ +Add :func:`math.lerp`, computing linear interpolation between two +values using the algorithm standardized for C++20's ``std::lerp``, +which guarantees exact endpoints (``lerp(a, b, 0.0) == a`` and +``lerp(a, b, 1.0) == b``) and monotonicity, properties the naive +``a + t * (b - a)`` formula does not have. Patch by PhysicistJohn. diff --git a/Modules/clinic/mathmodule.c.h b/Modules/clinic/mathmodule.c.h index a5e583c180defee..6ebef5767bd1913 100644 --- a/Modules/clinic/mathmodule.c.h +++ b/Modules/clinic/mathmodule.c.h @@ -271,6 +271,73 @@ math_ldexp(PyObject *module, PyObject *const *args, Py_ssize_t nargs) return return_value; } +PyDoc_STRVAR(math_lerp__doc__, +"lerp($module, a, b, t, /)\n" +"--\n" +"\n" +"Linear interpolation between a and b using parameter t.\n" +"\n" +"Returns a + t * (b - a), computed in a way that guarantees intended\n" +"mathematical behavior:\n" +"\n" +"* lerp(a, b, 0.0) == a\n" +"* lerp(a, b, 1.0) == b\n" +"* lerp(a, b, t) is monotonic in t, when a <= b\n" +"* t outside [0.0, 1.0] is allowed, for extrapolation"); + +#define MATH_LERP_METHODDEF \ + {"lerp", _PyCFunction_CAST(math_lerp), METH_FASTCALL, math_lerp__doc__}, + +static PyObject * +math_lerp_impl(PyObject *module, double a, double b, double t); + +static PyObject * +math_lerp(PyObject *module, PyObject *const *args, Py_ssize_t nargs) +{ + PyObject *return_value = NULL; + double a; + double b; + double t; + + if (!_PyArg_CheckPositional("lerp", nargs, 3, 3)) { + goto exit; + } + if (PyFloat_CheckExact(args[0])) { + a = PyFloat_AS_DOUBLE(args[0]); + } + else + { + a = PyFloat_AsDouble(args[0]); + if (a == -1.0 && PyErr_Occurred()) { + goto exit; + } + } + if (PyFloat_CheckExact(args[1])) { + b = PyFloat_AS_DOUBLE(args[1]); + } + else + { + b = PyFloat_AsDouble(args[1]); + if (b == -1.0 && PyErr_Occurred()) { + goto exit; + } + } + if (PyFloat_CheckExact(args[2])) { + t = PyFloat_AS_DOUBLE(args[2]); + } + else + { + t = PyFloat_AsDouble(args[2]); + if (t == -1.0 && PyErr_Occurred()) { + goto exit; + } + } + return_value = math_lerp_impl(module, a, b, t); + +exit: + return return_value; +} + PyDoc_STRVAR(math_modf__doc__, "modf($module, x, /)\n" "--\n" @@ -1164,4 +1231,4 @@ math_ulp(PyObject *module, PyObject *arg) exit: return return_value; } -/*[clinic end generated code: output=3452ce8caa2d1bd7 input=a9049054013a1b77]*/ +/*[clinic end generated code: output=3bd652589e637390 input=a9049054013a1b77]*/ diff --git a/Modules/mathmodule.c b/Modules/mathmodule.c index eaa1850b8aee1ac..60b37b5b59d899a 100644 --- a/Modules/mathmodule.c +++ b/Modules/mathmodule.c @@ -1719,6 +1719,53 @@ math_ldexp_impl(PyObject *module, double x, PyObject *i) } +/*[clinic input] +math.lerp + + a: double + b: double + t: double + / + +Linear interpolation between a and b using parameter t. + +Returns a + t * (b - a), computed in a way that guarantees intended +mathematical behavior: + +* lerp(a, b, 0.0) == a +* lerp(a, b, 1.0) == b +* lerp(a, b, t) is monotonic in t, when a <= b +* t outside [0.0, 1.0] is allowed, for extrapolation +[clinic start generated code]*/ + +static PyObject * +math_lerp_impl(PyObject *module, double a, double b, double t) +/*[clinic end generated code: output=7ca02ba2e81fd287 input=978d065447b6b7e0]*/ +{ + /* This is the algorithm standardized for std::lerp in C++20 + * ([c.math.lerp], from WG21 P0811R3 "Well-behaved interpolation + * for numbers and pointers" by S. Davis Herring). It guarantees + * properties the naive "a + t * (b - a)" formula does not: exact + * endpoints, monotonicity, and lerp(a, a, t) == a for finite t. + */ + if ((a <= 0.0 && b >= 0.0) || (a >= 0.0 && b <= 0.0)) { + return PyFloat_FromDouble(t * b + (1.0 - t) * a); + } + + if (t == 1.0) { + return PyFloat_FromDouble(b); + } + + double x = a + t * (b - a); + if ((t > 1.0) == (b > a)) { + return PyFloat_FromDouble(b < x ? x : b); + } + else { + return PyFloat_FromDouble(x < b ? x : b); + } +} + + /*[clinic input] math.modf @@ -3264,6 +3311,7 @@ static PyMethodDef math_methods[] = { MATH_ISINF_METHODDEF MATH_ISNAN_METHODDEF MATH_LDEXP_METHODDEF + MATH_LERP_METHODDEF {"lgamma", math_lgamma, METH_O, math_lgamma_doc}, {"log", _PyCFunction_CAST(math_log), METH_FASTCALL, math_log_doc}, {"log1p", math_log1p, METH_O, math_log1p_doc},