Math: added "shortest path" alternatives to lerp(), slerp() and sclerp().

Before neither of the lerp(), slerp() had the shortest path check, while
sclerp() had it. Now, to be consistent, none of them has it and there
are lerpShortestPath(), slerpShortestPath() and sclerpShortestPath()
functions that have the shortest path check.

This is different from other engines, where there's usually only the
shortest path interpolation by default and either an optional
"non-shortest-path" interpolation or no alternative at all. I like to
give the users a choice, so there's both versions and the
non-shortest-path version is the default, because -- at least in case of
lerp() -- this results in a quite significant perf difference (15%
faster), so why not have it. Preprocess your data instead ;)

doc/changelog.dox

@@ -72,6 +72,11 @@ See also:
 -   Added @ref Math::lerp(const Complex<T>&, const Complex<T>&, T) and
     @ref Math::slerp(const Complex<T>&, const Complex<T>&, T) for feature
     parity with @ref Math::Quaternion
+-   Added @ref Math::lerpShortestPath(const Quaternion<T>&, const Quaternion<T>&, T),
+    @ref Math::slerpShortestPath(const Quaternion<T>&, const Quaternion<T>&, T)
+    and @ref Math::sclerpShortestPath(const DualQuaternion<T>&, const DualQuaternion<T>&, T)
+    variants; clarified that the original versions are explicitly *not*
+    shortest-path
 -   Added @ref Math::Range2D::x(), @ref Math::Range3D::x(),
     @ref Math::Range2D::y(), @ref Math::Range3D::y(), @ref Math::Range3D::z()
     and @ref Math::Range3D::y() for slicing ranges into lower dimensions

src/Magnum/Math/DualQuaternion.h

@@ -48,26 +48,107 @@ namespace Implementation {
 @param normalizedB  Second dual quaternion
 @param t            Interpolation phase (from range @f$ [0; 1] @f$)
 
-Expects that both dual quaternions are normalized. @f[
+Expects that both dual quaternions are normalized. If the real parts are the
-\begin{array}{rcl}
+same or one is a negation of the other, returns the @ref DualQuaternion::rotation()
-    l + \epsilon m & = & \hat q_A^* \hat q_B \\
+(real) part combined with interpolated @ref DualQuaternion::translation(): @f[
-    \frac{\hat a} 2 & = & acos \left( l_S \right) - \epsilon m_S \frac 1 {|l_V|} \\
+    \begin{array}{rcl}
-    \hat {\boldsymbol n} & = & \boldsymbol n_0 + \epsilon \boldsymbol n_\epsilon
+        d & = & q_{A_0} \cdot q_{B_0} \\[5pt]
-    ~~~~~~~~ \boldsymbol n_0 = l_V \frac 1 {|l_V|}
+        {\hat q}_{ScLERP} & = & 2 \left[(1 - t)(q_{A_\epsilon} q_{A_0}^*)_V + t (q_{B_\epsilon} q_{B_0}^*)_V \right] q_A, ~ {\color{m-primary} \text{if} ~ d \ge 1}
-    ~~~~~~~~ \boldsymbol n_\epsilon = \left( m_V - {\boldsymbol n}_0 \frac {a_\epsilon} 2 l_S \right)\frac 1 {|l_V|} \\
+    \end{array}
-    {\hat q}_{ScLERP} & = & \hat q_A (\hat q_A^* \hat q_B)^t =
-        \hat q_A \left[ \hat {\boldsymbol n} sin \left( t \frac {\hat a} 2 \right), cos \left( t \frac {\hat a} 2 \right) \right] \\
-\end{array}
 @f]
-@see @ref DualQuaternion::isNormalized(),
+
-    @ref slerp(const Quaternion<T>&, const Quaternion<T>&, T),
+@m_class{m-noindent}
-    @ref lerp(const T&, const T&, U)
+
+otherwise, the interpolation is performed as: @f[
+    \begin{array}{rcl}
+        l + \epsilon m & = & \hat q_A^* \hat q_B \\[5pt]
+        \frac{\hat a} 2 & = & \arccos \left( l_S \right) - \epsilon m_S \frac 1 {|\boldsymbol{l}_V|} \\[5pt]
+        \hat {\boldsymbol n} & = & \boldsymbol n_0 + \epsilon \boldsymbol n_\epsilon,
+        ~~~~~~~~ \boldsymbol n_0 = \boldsymbol{l}_V \frac 1 {|\boldsymbol{l}_V|},
+        ~~~~~~~~ \boldsymbol n_\epsilon = \left(\boldsymbol{m}_V - {\boldsymbol n}_0 \frac {a_\epsilon} 2 l_S \right)\frac 1 {|\boldsymbol{l}_V|} \\[5pt]
+        {\hat q}_{ScLERP} & = & \hat q_A (\hat q_A^* \hat q_B)^t =
+            \hat q_A \left[ \hat {\boldsymbol n} \sin \left( t \frac {\hat a} 2 \right), \cos \left( t \frac {\hat a} 2 \right) \right]
+    \end{array}
+@f]
+
+Note that this function does not check for shortest path interpolation, see
+@ref sclerpShortestPath() for an alternative.
+@see @ref DualQuaternion::isNormalized(), @ref DualQuaternion::quaternionConjugated(),
+    @ref lerp(const Quaternion<T>&, const Quaternion<T>&, T),
+    @ref slerp(const Quaternion<T>&, const Quaternion<T>&, T)
 */
 template<class T> inline DualQuaternion<T> sclerp(const DualQuaternion<T>& normalizedA, const DualQuaternion<T>& normalizedB, const T t) {
     CORRADE_ASSERT(normalizedA.isNormalized() && normalizedB.isNormalized(),
         "Math::sclerp(): dual quaternions must be normalized", {});
     const T cosHalfAngle = dot(normalizedA.real(), normalizedB.real());
 
+    /* Avoid division by zero: interpolate just the translation part */
+    /** @todo could this be optimized somehow? */
+    if(std::abs(cosHalfAngle) >= T(1))
+        return DualQuaternion<T>::translation(Implementation::lerp(normalizedA.translation(), normalizedB.translation(), t))*DualQuaternion<T>{normalizedA.real()};
+
+    /* l + εm = q_A^**q_B */
+    const DualQuaternion<T> diff = normalizedA.quaternionConjugated()*normalizedB;
+    const Quaternion<T>& l = diff.real();
+    const Quaternion<T>& m = diff.dual();
+
+    /* a/2 = acos(l_S) - εm_S/|l_V| */
+    const T invr = l.vector().lengthInverted();
+    const Dual<T> aHalf{std::acos(l.scalar()), -m.scalar()*invr};
+
+    /* direction = n_0 = l_V/|l_V|
+       moment = n_ε = (m_V - n_0*(a_ε/2)*l_S)/|l_V| */
+    const Vector3<T> direction = l.vector()*invr;
+    const Vector3<T> moment = (m.vector() - direction*(aHalf.dual()*l.scalar()))*invr;
+    const Dual<Vector3<T>> n{direction, moment};
+
+    /* q_ScLERP = q_A*(cos(t*a/2) + n*sin(t*a/2)) */
+    Dual<T> sin, cos;
+    std::tie(sin, cos) = Math::sincos(t*Dual<Rad<T>>(aHalf));
+    return normalizedA*DualQuaternion<T>{n*sin, cos};
+}
+
+/** @relatesalso DualQuaternion
+@brief Screw linear shortest-path interpolation of two dual quaternions
+@param normalizedA  First dual quaternion
+@param normalizedB  Second dual quaternion
+@param t            Interpolation phase (from range @f$ [0; 1] @f$)
+
+Unlike @ref sclerp(const DualQuaternion<T>&, const DualQuaternion<T>&, T) this
+function interpolates on the shortest path. Expects that both dual quaternions
+are normalized. If the real parts are the same or one is a negation of the
+other, returns the @ref DualQuaternion::rotation() (real) part combined with
+interpolated @ref DualQuaternion::translation(): @f[
+    \begin{array}{rcl}
+        d & = & q_{A_0} \cdot q_{B_0} \\[5pt]
+        {\hat q}_{ScLERP} & = & 2 \left((1 - t)(q_{A_\epsilon} q_{A_0}^*)_V + t (q_{B_\epsilon} q_{B_0}^*)_V \right) (q_{A_0} + \epsilon [\boldsymbol{0}, 0]), ~ {\color{m-primary} \text{if} ~ d \ge 1}
+    \end{array}
+@f]
+
+@m_class{m-noindent}
+
+otherwise, the interpolation is performed as: @f[
+    \begin{array}{rcl}
+        l + \epsilon m & = & \begin{cases}
+            \phantom{-}\hat q_A^* \hat q_B, & d \ge 0 \\
+            -\hat q_A^* \hat q_B, & d < 0 \\
+        \end{cases} \\[15pt]
+        \frac{\hat a} 2 & = & \arccos \left( l_S \right) - \epsilon m_S \frac 1 {|\boldsymbol{l}_V|} \\[5pt]
+        \hat {\boldsymbol n} & = & \boldsymbol n_0 + \epsilon \boldsymbol n_\epsilon,
+        ~~~~~~~~ \boldsymbol n_0 = \boldsymbol{l}_V \frac 1 {|\boldsymbol{l}_V|},
+        ~~~~~~~~ \boldsymbol n_\epsilon = \left(\boldsymbol{m}_V - {\boldsymbol n}_0 \frac {a_\epsilon} 2 l_S \right)\frac 1 {|\boldsymbol{l}_V|} \\[5pt]
+        {\hat q}_{ScLERP} & = & \hat q_A (\hat q_A^* \hat q_B)^t =
+            \hat q_A \left[ \hat {\boldsymbol n} \sin \left( t \frac {\hat a} 2 \right), \cos \left( t \frac {\hat a} 2 \right) \right]
+    \end{array}
+@f]
+@see @ref DualQuaternion::isNormalized(), @ref lerpShortestPath(),
+    @ref slerpShortestPath()
+*/
+template<class T> inline DualQuaternion<T> sclerpShortestPath(const DualQuaternion<T>& normalizedA, const DualQuaternion<T>& normalizedB, const T t) {
+    CORRADE_ASSERT(normalizedA.isNormalized() && normalizedB.isNormalized(),
+        "Math::sclerp(): dual quaternions must be normalized", {});
+    const T cosHalfAngle = dot(normalizedA.real(), normalizedB.real());
+
     /* Avoid division by zero: interpolate just the translation part */
     /** @todo could this be optimized somehow? */
     if(std::abs(cosHalfAngle) >= T(1))

src/Magnum/Math/Quaternion.h

@@ -87,6 +87,9 @@ template<class T> inline Rad<T> angle(const Quaternion<T>& normalizedA, const Qu
 Expects that both quaternions are normalized. @f[
     q_{LERP} = \frac{(1 - t) q_A + t q_B}{|(1 - t) q_A + t q_B|}
 @f]
+
+Note that this function does not check for shortest path interpolation, see
+@ref lerpShortestPath() for an alternative.
 @see @ref Quaternion::isNormalized(),
     @ref slerp(const Quaternion<T>&, const Quaternion<T>&, T), @ref sclerp(),
     @ref lerp(const T&, const T&, U),
@@ -101,6 +104,31 @@ template<class T> inline Quaternion<T> lerp(const Quaternion<T>& normalizedA, co
     return ((T(1) - t)*normalizedA + t*normalizedB).normalized();
 }
 
+/** @relatesalso Quaternion
+@brief Linear shortest-path interpolation of two quaternions
+@param normalizedA  First quaternion
+@param normalizedB  Second quaternion
+@param t            Interpolation phase (from range @f$ [0; 1] @f$)
+
+Unlike @ref lerp(const Quaternion<T>&, const Quaternion<T>&, T), this
+interpolates on the shortest path at some performance expense. Expects that
+both quaternions are normalized. @f[
+    \begin{array}{rcl}
+        d & = & q_A \cdot q_B \\[5pt]
+        q'_A & = & \begin{cases}
+                \phantom{-}q_A, & d \ge 0 \\
+                -q_A, & d < 0
+            \end{cases} \\[15pt]
+        q_{LERP} & = & \cfrac{(1 - t) q'_A + t q_B}{|(1 - t) q'_A + t q_B|}
+    \end{array}
+@f]
+@see @ref Quaternion::isNormalized(), @ref slerpShortestPath(),
+    @ref sclerpShortestPath()
+*/
+template<class T> inline Quaternion<T> lerpShortestPath(const Quaternion<T>& normalizedA, const Quaternion<T>& normalizedB, T t) {
+    return lerp(dot(normalizedA, normalizedB) < T(0) ? -normalizedA : normalizedA, normalizedB, t);
+}
+
 /** @relatesalso Quaternion
 @brief Spherical linear interpolation of two quaternions
 @param normalizedA  First quaternion
@@ -108,11 +136,24 @@ template<class T> inline Quaternion<T> lerp(const Quaternion<T>& normalizedA, co
 @param t            Interpolation phase (from range @f$ [0; 1] @f$)
 
 Expects that both quaternions are normalized. If the quaternions are the same
-or one is a negation of the other, returns the first argument. @f[
+or one is a negation of the other, it just returns the first argument: @f[
-    q_{SLERP} = \frac{sin((1 - t) \theta) q_A + sin(t \theta) q_B}{sin \theta}
+    \begin{array}{rcl}
-    ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~
+        d & = & q_A \cdot q_B \\[5pt]
-    \theta = acos \left( \frac{q_A \cdot q_B}{|q_A| \cdot |q_B|} \right) = acos(q_A \cdot q_B)
+        q_{SLERP} & = & q_A, ~ {\color{m-primary} \text{if} ~ d \ge 1}
+    \end{array}
 @f]
+
+@m_class{m-noindent}
+
+otherwise, the interpolation is performed as: @f[
+    \begin{array}{rcl}
+        \theta & = & \arccos \left( \frac{q_A \cdot q_B}{|q_A| |q_B|} \right) = \arccos(q_A \cdot q_B) = \arccos(d) \\[5pt]
+        q_{SLERP} & = & \cfrac{\sin((1 - t) \theta) q_A + \sin(t \theta) q_B}{\sin(\theta)}
+    \end{array}
+@f]
+
+Note that this function does not check for shortest path interpolation, see
+@ref slerpShortestPath() for an alternative.
 @see @ref Quaternion::isNormalized(), @ref lerp(const Quaternion<T>&, const Quaternion<T>&, T),
     @ref slerp(const Complex<T>&, const Complex<T>&, T), @ref sclerp()
 */
@@ -128,11 +169,62 @@ template<class T> inline Quaternion<T> slerp(const Quaternion<T>& normalizedA, c
     return (std::sin((T(1) - t)*a)*normalizedA + std::sin(t*a)*normalizedB)/std::sin(a);
 }
 
+/** @relatesalso Quaternion
+@brief Spherical linear shortest-path interpolation of two quaternions
+@param normalizedA  First quaternion
+@param normalizedB  Second quaternion
+@param t            Interpolation phase (from range @f$ [0; 1] @f$)
+
+Unlike @ref slerp(const Quaternion<T>&, const Quaternion<T>&, T) this function
+interpolates on the shortest path. Expects that both quaternions are
+normalized. If the quaternions are the same or one is a negation of the other,
+it just returns the first argument: @f[
+    \begin{array}{rcl}
+        d & = & q_A \cdot q_B \\
+        q_{SLERP} & = & q_A, ~ {\color{m-primary} \text{if} ~ d \ge 1}
+    \end{array}
+@f]
+
+@m_class{m-noindent}
+
+otherwise, the interpolation is performed as: @f[
+    \begin{array}{rcl}
+        q'_A & = & \begin{cases}
+                \phantom{-}q_A, & d \ge 0 \\
+                -q_A, & d < 0
+            \end{cases} \\[15pt]
+        \theta & = & \arccos \left( \frac{|q'_A \cdot q_B|}{|q'_A| |q_B|} \right) = \arccos(|q'_A \cdot q_B|) = \arccos(|d|) \\[5pt]
+        q_{SLERP} & = & \cfrac{\sin((1 - t) \theta) q'_A + \sin(t \theta) q_B}{\sin(\theta)}
+    \end{array}
+@f]
+@see @ref Quaternion::isNormalized(), @ref lerpShortestPath(),
+        @ref sclerpShortestPath()
+*/
+template<class T> inline Quaternion<T> slerpShortestPath(const Quaternion<T>& normalizedA, const Quaternion<T>& normalizedB, T t) {
+    CORRADE_ASSERT(normalizedA.isNormalized() && normalizedB.isNormalized(),
+        "Math::slerp(): quaternions must be normalized", {});
+    const T cosHalfAngle = dot(normalizedA, normalizedB);
+
+    /* Avoid division by zero */
+    if(std::abs(cosHalfAngle) >= T(1)) return Quaternion<T>{normalizedA};
+
+    const Quaternion<T> shortestNormalizedA = cosHalfAngle < 0 ? -normalizedA : normalizedA;
+
+    const T a = std::acos(std::abs(cosHalfAngle));
+    return (std::sin((T(1) - t)*a)*shortestNormalizedA + std::sin(t*a)*normalizedB)/std::sin(a);
+}
+
 /**
 @brief Quaternion
 @tparam T   Underlying data type
 
-Represents 3D rotation. See @ref transformations for brief introduction.
+Represents 3D rotation. Usually denoted as the following in equations, with
+@f$ \boldsymbol{q}_V @f$ being the @ref vector() part and @f$ q_S @f$ being the
+@ref scalar() part: @f[
+    q = [\boldsymbol{q}_V, q_S]
+@f]
+
+See @ref transformations for a brief introduction.
 @see @ref Magnum::Quaternion, @ref Magnum::Quaterniond, @ref DualQuaternion,
     @ref Matrix4
 */

src/Magnum/Math/Test/DualQuaternionTest.cpp

@@ -90,6 +90,7 @@ struct DualQuaternionTest: Corrade::TestSuite::Tester {
     void transformPointNormalized();
 
     void sclerp();
+    void sclerpShortestPath();
 
     void debug();
     void configuration();
@@ -144,6 +145,7 @@ DualQuaternionTest::DualQuaternionTest() {
               &DualQuaternionTest::transformPointNormalized,
 
               &DualQuaternionTest::sclerp,
+              &DualQuaternionTest::sclerpShortestPath,
 
               &DualQuaternionTest::debug,
               &DualQuaternionTest::configuration});
@@ -480,40 +482,98 @@ void DualQuaternionTest::transformPointNormalized() {
 
 void DualQuaternionTest::sclerp() {
     auto from = DualQuaternion::translation({20.0f, 0.0f, 0.0f})*
-        DualQuaternion::rotation(180.0_degf, Vector3::yAxis());
+        DualQuaternion::rotation(65.0_degf, Vector3::yAxis());
     auto to = DualQuaternion::translation({42.0f, 42.0f, 42.0f})*
         DualQuaternion::rotation(75.0_degf, Vector3::xAxis());
 
+    const DualQuaternion begin = Math::sclerp(from, to, 0.0f);
+    const DualQuaternion beginShortestPath = Math::sclerpShortestPath(from, to, 0.0f);
+    const DualQuaternion end = Math::sclerp(from, to, 1.0f);
+    const DualQuaternion endShortestPath = Math::sclerpShortestPath(from, to, 1.0f);
+    CORRADE_COMPARE(begin, from);
+    CORRADE_COMPARE(beginShortestPath, from);
+    CORRADE_COMPARE(end, to);
+    CORRADE_COMPARE(endShortestPath, to);
+
     DualQuaternion expected1{
-        {{0.23296291f, 0.92387953f, 0.0f}, 0.30360317f},
+        {{0.170316f, 0.424975f, 0.0f}, 0.889038f},
-        {{2.23561910f, 2.81697198f, 10.72224091f}, -10.28763633f}};
+        {{10.689f, 7.47059f, 5.33428f}, -5.61881f}};
     DualQuaternion expected2{
-        {{0.44376798f, 0.68454710f, 0.0f}, 0.57832969f},
+        {{0.34568f, 0.282968f, 0.0f}, 0.89467f},
-        {{5.76439487f, 11.16130665f, 9.67126701f}, -17.63439459f}};
+        {{12.8764f, 15.8357f, 5.03088f}, -9.98371f}};
     DualQuaternion expected3{
-        {{0.59797859f, 0.18738131f, 0.0f}, 0.77930087f},
+        {{0.550678f, 0.072563f, 0.0f}, 0.831558f},
-        {{13.40962790f, 25.45212445f, 5.68158104f}, -16.40948111f}};
+        {{15.6916f, 26.3477f, 4.23219f}, -12.6905f}};
 
     const DualQuaternion interp1 = Math::sclerp(from, to, 0.25f);
+    const DualQuaternion interp1ShortestPath = Math::sclerpShortestPath(from, to, 0.25f);
     const DualQuaternion interp2 = Math::sclerp(from, to, 0.52f);
+    const DualQuaternion interp2ShortestPath = Math::sclerpShortestPath(from, to, 0.52f);
     const DualQuaternion interp3 = Math::sclerp(from, to, 0.88f);
+    const DualQuaternion interp3ShortestPath = Math::sclerpShortestPath(from, to, 0.88f);
 
     CORRADE_COMPARE(interp1, expected1);
+    CORRADE_COMPARE(interp1ShortestPath, expected1);
     CORRADE_COMPARE(interp2, expected2);
+    CORRADE_COMPARE(interp2ShortestPath, expected2);
     CORRADE_COMPARE(interp3, expected3);
+    CORRADE_COMPARE(interp3ShortestPath, expected3);
 
     /* Edge cases: */
 
     /* Dual quaternions with identical rotation */
     CORRADE_COMPARE(Math::sclerp(from, from, 0.42f), from);
+    CORRADE_COMPARE(Math::sclerpShortestPath(from, from, 0.42f), from);
     CORRADE_COMPARE(Math::sclerp(from, -from, 0.42f), from);
+    CORRADE_COMPARE(Math::sclerpShortestPath(from, -from, 0.42f), from);
 
     /* No difference in rotation, but in translation */
-    auto rotation = DualQuaternion::rotation(35.0_degf, Vector3{0.3f, 0.2f, 0.1f}.normalized());
+    {
-    auto interpolateTranslation = Math::sclerp(
+        auto rotation = DualQuaternion::rotation(35.0_degf, Vector3{0.3f, 0.2f, 0.1f}.normalized());
-        DualQuaternion::translation({1.0f, 2.0f, 4.0f})*rotation, DualQuaternion::translation({5.0f, -6.0f, 2.0f})*rotation, 0.25f);
+        auto a = DualQuaternion::translation({1.0f, 2.0f, 4.0f})*rotation;
-    CORRADE_VERIFY(interpolateTranslation.isNormalized());
+        auto b = DualQuaternion::translation({5.0f, -6.0f, 2.0f})*rotation;
-    CORRADE_COMPARE(interpolateTranslation, DualQuaternion::translation({2.0f, 0.0f, 3.5f})*rotation);
+        auto expected = DualQuaternion::translation({2.0f, 0.0f, 3.5f})*rotation;
+
+        auto interpolateTranslation = Math::sclerp(a, b, 0.25f);
+        auto interpolateTranslationShortestPath = Math::sclerpShortestPath(a, b, 0.25f);
+        CORRADE_VERIFY(interpolateTranslation.isNormalized());
+        CORRADE_VERIFY(interpolateTranslationShortestPath.isNormalized());
+        CORRADE_COMPARE(interpolateTranslation, expected);
+        CORRADE_COMPARE(interpolateTranslationShortestPath, expected);
+    }
+}
+
+void DualQuaternionTest::sclerpShortestPath() {
+    DualQuaternion a = DualQuaternion::translation({1.5f, 0.3f, 0.0f})*
+        DualQuaternion::rotation(0.0_degf, Vector3::zAxis());
+    DualQuaternion b = DualQuaternion::translation({3.5f, 0.3f, 1.0f})*
+        DualQuaternion::rotation(225.0_degf, Vector3::zAxis());
+
+    DualQuaternion sclerp = Math::sclerp(a, b, 0.25f);
+    DualQuaternion sclerpShortestPath = Math::sclerpShortestPath(a, b, 0.25f);
+
+    CORRADE_VERIFY(sclerp.isNormalized());
+    CORRADE_VERIFY(sclerpShortestPath.isNormalized());
+    CORRADE_COMPARE(sclerp.rotation().axis(), Vector3::zAxis());
+    /** @todo why is this inverted compared to QuaternionTest::slerpShortestPath()? */
+    CORRADE_COMPARE(sclerpShortestPath.rotation().axis(), -Vector3::zAxis());
+    CORRADE_COMPARE(sclerp.rotation().angle(), 56.25_degf);
+    /* Because the axis is inverted, this is also inverted compared to
+       QuaternionTest::slerpShortestPath() */
+    CORRADE_COMPARE(sclerpShortestPath.rotation().angle(), 360.0_degf - 326.25_degf);
+
+    CORRADE_COMPARE(sclerp, (DualQuaternion{
+        {{0.0f, 0.0f, 0.471397f}, 0.881921f},
+        {{0.536892f, -0.692656f, 0.11024f}, -0.0589246f}}));
+    /* Also inverted compared to QuaternionTest::slerpShortestPath() */
+    CORRADE_COMPARE(sclerpShortestPath, (DualQuaternion{
+        {{0.0f, 0.0f, -0.290285f}, 0.95694f},
+        {{0.794402f, 0.651539f, 0.119618f}, 0.0362856f}}));
+
+    /* Translation along Z should be the same in both, in 25% of the way.
+       Translation in the XY plane is along a screw, so that's different. */
+    CORRADE_COMPARE(sclerpShortestPath.translation().z(), 0.25f);
+    CORRADE_COMPARE(sclerpShortestPath.translation().z(), 0.25f);
 }
 
 void DualQuaternionTest::debug() {

src/Magnum/Math/Test/QuaternionTest.cpp

@@ -89,12 +89,16 @@ struct QuaternionTest: Corrade::TestSuite::Tester {
     void rotation();
     void angle();
     void matrix();
+
     void lerp();
+    void lerpShortestPath();
     void lerp2D();
     void lerpNotNormalized();
     void slerp();
+    void slerpShortestPath();
     void slerp2D();
     void slerpNotNormalized();
+
     void transformVector();
     void transformVectorNormalized();
 
@@ -150,12 +154,16 @@ QuaternionTest::QuaternionTest() {
               &QuaternionTest::rotation,
               &QuaternionTest::angle,
               &QuaternionTest::matrix,
+
               &QuaternionTest::lerp,
+              &QuaternionTest::lerpShortestPath,
               &QuaternionTest::lerp2D,
               &QuaternionTest::lerpNotNormalized,
               &QuaternionTest::slerp,
+              &QuaternionTest::slerpShortestPath,
               &QuaternionTest::slerp2D,
               &QuaternionTest::slerpNotNormalized,
+
               &QuaternionTest::transformVector,
               &QuaternionTest::transformVectorNormalized,
 
@@ -481,10 +489,34 @@ void QuaternionTest::matrix() {
 void QuaternionTest::lerp() {
     Quaternion a = Quaternion::rotation(15.0_degf, Vector3(1.0f/Constants<Float>::sqrt3()));
     Quaternion b = Quaternion::rotation(23.0_degf, Vector3::xAxis());
+
     Quaternion lerp = Math::lerp(a, b, 0.35f);
+    Quaternion lerpShortestPath = Math::lerpShortestPath(a, b, 0.35f);
+    Quaternion expected{{0.119127f, 0.049134f, 0.049134f}, 0.990445f};
 
+    /* Both should give the same result */
     CORRADE_VERIFY(lerp.isNormalized());
-    CORRADE_COMPARE(lerp, Quaternion({0.119127f, 0.049134f, 0.049134f}, 0.990445f));
+    CORRADE_VERIFY(lerpShortestPath.isNormalized());
+    CORRADE_COMPARE(lerp, expected);
+    CORRADE_COMPARE(lerpShortestPath, expected);
+}
+
+void QuaternionTest::lerpShortestPath() {
+    Quaternion a = Quaternion::rotation(0.0_degf, Vector3::zAxis());
+    Quaternion b = Quaternion::rotation(225.0_degf, Vector3::zAxis());
+
+    Quaternion slerp = Math::lerp(a, b, 0.25f);
+    Quaternion slerpShortestPath = Math::lerpShortestPath(a, b, 0.25f);
+
+    CORRADE_VERIFY(slerp.isNormalized());
+    CORRADE_VERIFY(slerpShortestPath.isNormalized());
+    CORRADE_COMPARE(slerp.axis(), Vector3::zAxis());
+    CORRADE_COMPARE(slerpShortestPath.axis(), Vector3::zAxis());
+    CORRADE_COMPARE(slerp.angle(), 38.8848_degf);
+    CORRADE_COMPARE(slerpShortestPath.angle(), 329.448_degf);
+
+    CORRADE_COMPARE(slerp, (Quaternion{{0.0f, 0.0f, 0.332859f}, 0.942977f}));
+    CORRADE_COMPARE(slerpShortestPath, (Quaternion{{0.0f, 0.0f, 0.26347f}, -0.964667f}));
 }
 
 void QuaternionTest::lerp2D() {
@@ -513,14 +545,40 @@ void QuaternionTest::lerpNotNormalized() {
 void QuaternionTest::slerp() {
     Quaternion a = Quaternion::rotation(15.0_degf, Vector3(1.0f/Constants<Float>::sqrt3()));
     Quaternion b = Quaternion::rotation(23.0_degf, Vector3::xAxis());
+
     Quaternion slerp = Math::slerp(a, b, 0.35f);
+    Quaternion slerpShortestPath = Math::slerpShortestPath(a, b, 0.35f);
+    Quaternion expected{{0.1191653f, 0.0491109f, 0.0491109f}, 0.9904423f};
 
+    /* Both should give the same result */
     CORRADE_VERIFY(slerp.isNormalized());
-    CORRADE_COMPARE(slerp, Quaternion({0.1191653f, 0.0491109f, 0.0491109f}, 0.9904423f));
+    CORRADE_COMPARE(slerp, expected);
+    CORRADE_VERIFY(slerpShortestPath.isNormalized());
+    CORRADE_COMPARE(slerpShortestPath, expected);
 
     /* Avoid division by zero */
     CORRADE_COMPARE(Math::slerp(a, a, 0.25f), a);
     CORRADE_COMPARE(Math::slerp(a, -a, 0.42f), a);
+    CORRADE_COMPARE(Math::slerpShortestPath(a, a, 0.25f), a);
+    CORRADE_COMPARE(Math::slerpShortestPath(a, -a, 0.25f), a);
+}
+
+void QuaternionTest::slerpShortestPath() {
+    Quaternion a = Quaternion::rotation(0.0_degf, Vector3::zAxis());
+    Quaternion b = Quaternion::rotation(225.0_degf, Vector3::zAxis());
+
+    Quaternion slerp = Math::slerp(a, b, 0.25f);
+    Quaternion slerpShortestPath = Math::slerpShortestPath(a, b, 0.25f);
+
+    CORRADE_VERIFY(slerp.isNormalized());
+    CORRADE_VERIFY(slerpShortestPath.isNormalized());
+    CORRADE_COMPARE(slerp.axis(), Vector3::zAxis());
+    CORRADE_COMPARE(slerpShortestPath.axis(), Vector3::zAxis());
+    CORRADE_COMPARE(slerp.angle(), 56.25_degf);
+    CORRADE_COMPARE(slerpShortestPath.angle(), 326.25_degf);
+
+    CORRADE_COMPARE(slerp, (Quaternion{{0.0f, 0.0f, 0.471397f}, 0.881921f}));
+    CORRADE_COMPARE(slerpShortestPath, (Quaternion{{0.0f, 0.0f, 0.290285f}, -0.95694f}));
 }
 
 void QuaternionTest::slerp2D() {