From 0106520d80be2f8b95f6ed7a9ada348eadafcca5 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Vladim=C3=ADr=20Vondru=C5=A1?= <mosra@centrum.cz>
Date: Tue, 5 Mar 2024 13:11:04 +0100
Subject: [PATCH] Math: tighten up lerp() behavior with boundary preservation
 tests.

Got a suggestion that lerp() could be optimized to be one arithmetic
operation less. While valid in certain cases, it would break in case the
endpoints have wildly different magnitudes. Unfortunately that was only
my personal knowledge, not backed by regression tests. Now it is.
---
 src/Magnum/Math/Test/FunctionsTest.cpp | 42 ++++++++++++++++++++++++++
 src/Magnum/Math/Vector.h               |  3 ++
 2 files changed, 45 insertions(+)

diff --git a/src/Magnum/Math/Test/FunctionsTest.cpp b/src/Magnum/Math/Test/FunctionsTest.cpp
index e7834c7bb..114068016 100644
--- a/src/Magnum/Math/Test/FunctionsTest.cpp
+++ b/src/Magnum/Math/Test/FunctionsTest.cpp
@@ -63,6 +63,8 @@ struct FunctionsTest: TestSuite::Tester {
     void sqrt();
     void sqrtInverted();
     void lerp();
+    void lerpLimits();
+    void lerpInfinity();
     void lerpBool();
     void lerpInverted();
     void select();
@@ -128,6 +130,8 @@ FunctionsTest::FunctionsTest() {
               &FunctionsTest::sqrt,
               &FunctionsTest::sqrtInverted,
               &FunctionsTest::lerp,
+              &FunctionsTest::lerpLimits,
+              &FunctionsTest::lerpInfinity,
               &FunctionsTest::lerpBool,
               &FunctionsTest::lerpInverted,
               &FunctionsTest::select,
@@ -403,6 +407,44 @@ void FunctionsTest::lerp() {
     CORRADE_COMPARE(Math::lerp(2.0_usec, 5.0_usec, 0.5f), 3.5_usec);
 }
 
+template<class T, class U> T lerpOptimized(const T& a, const T& b, U t) {
+    /* One multiplication and two additions, while `T((U(1) - t)*a + t*b)` is
+       two multiplications, addition and subtraction. Doesn't correctly
+       preserve boundary values. */
+    return t*(b - a) + a;
+}
+
+void FunctionsTest::lerpLimits() {
+    CORRADE_COMPARE(Math::lerp(1.0e10f, 1.0e-5f, 0.0f), 1.0e10f);
+    CORRADE_COMPARE(Math::lerp(1.0e10f, 1.0e-5f, 1.0f), 1.0e-5f);
+    CORRADE_COMPARE(Math::lerp(1.0e-5f, 1.0e10f, 0.0f), 1.0e-5f);
+    CORRADE_COMPARE(Math::lerp(1.0e-5f, 1.0e10f, 1.0f), 1.0e10f);
+
+    CORRADE_COMPARE(lerpOptimized(1.0e10f, 1.0e-5f, 0.0f), 1.0e10f);
+    {
+        CORRADE_EXPECT_FAIL("\"Optimized\" version of a lerp doesn't correctly preserve boundary values with wildly different magnitudes.");
+        CORRADE_COMPARE(lerpOptimized(1.0e10f, 1.0e-5f, 1.0f), 1.0e-5f);
+    }
+}
+
+void FunctionsTest::lerpInfinity() {
+    CORRADE_COMPARE(Math::lerp(Constants::inf(), 0.0f, 0.0f), Constants::inf());
+    CORRADE_COMPARE(Math::lerp(0.0f, Constants::inf(), 1.0f), Constants::inf());
+    {
+        CORRADE_EXPECT_FAIL("Lerp with infinity doesn't correctly preserve the other boundary value.");
+        CORRADE_COMPARE(Math::lerp(Constants::inf(), 0.0f, 1.0f), 0.0f);
+        CORRADE_COMPARE(Math::lerp(0.0f, Constants::inf(), 0.0f), 0.0f);
+    }
+
+    CORRADE_COMPARE(lerpOptimized(0.0f, Constants::inf(), 1.0f), Constants::inf());
+    {
+        CORRADE_EXPECT_FAIL("\"Optimized\" version of a lerp doesn't correctly preserve boundary values if an infinity is present.");
+        CORRADE_COMPARE(lerpOptimized(Constants::inf(), 0.0f, 0.0f), Constants::inf());
+        CORRADE_COMPARE(lerpOptimized(Constants::inf(), 0.0f, 1.0f), 0.0f);
+        CORRADE_COMPARE(lerpOptimized(0.0f, Constants::inf(), 0.0f), 0.0f);
+    }
+}
+
 void FunctionsTest::lerpBool() {
     /* Scalar interpolation phase */
     CORRADE_COMPARE(Math::lerp(Vector3i{1, 2, 3}, Vector3i{5, 6, 7}, true), (Vector3i{5, 6, 7}));
diff --git a/src/Magnum/Math/Vector.h b/src/Magnum/Math/Vector.h
index 8aee3ba35..21aaa3697 100644
--- a/src/Magnum/Math/Vector.h
+++ b/src/Magnum/Math/Vector.h
@@ -74,6 +74,9 @@ namespace Implementation {
     template<std::size_t, class, class> struct VectorConverter;
     /* Needed by DualQuaternion and Functions.h (to avoid dependency between them) */
     template<class T, class U> T lerp(const T& a, const T& b, U t) {
+        /* While `t*(b - a) + a` is one ALU op less, the following is
+           guaranteed to correctly preserves exact boundary values with t being
+           0 or 1. See FunctionsTest::lerpLimits() for details. */
         return T((U(1) - t)*a + t*b);
     }