Browse Source

Math: ability to multiply and divide dual and real number.

pull/118/head
Vladimír Vondruš 11 years ago
parent
commit
ce9be12708
  1. 56
      src/Magnum/Math/Dual.h
  2. 20
      src/Magnum/Math/Test/DualTest.cpp

56
src/Magnum/Math/Dual.h

@ -37,6 +37,10 @@
namespace Magnum { namespace Math {
namespace Implementation {
CORRADE_HAS_TYPE(IsDual, decltype(std::declval<const T>().dual()));
}
/**
@brief Dual number
@tparam T Underlying data type
@ -165,22 +169,52 @@ template<class T> class Dual {
* @f[
* \hat a \hat b = a_0 b_0 + \epsilon (a_0 b_\epsilon + a_\epsilon b_0)
* @f]
* @see @ref operator*(const U&) const,
* @ref operator*(const T&, const Dual<U>&)
*/
template<class U> auto operator*(const Dual<U>& other) const -> Dual<decltype(std::declval<T>()*std::declval<U>())> {
return {_real*other._real, _real*other._dual + _dual*other._real};
}
/**
* @brief Multiply by real number
*
* Equivalent to the above assuming that @f$ b_\epsilon = 0 @f$.
* @f[
* \hat a \hat b = a_0 b_0 + \epsilon (a_0 b_\epsilon + a_\epsilon b_0) = a_0 b_0 + \epsilon a_\epsilon b_0
* @f]
* @see @ref operator*(const Dual<U>&) const,
* @ref operator*(const T&, const Dual<U>&)
*/
template<class U, class V = typename std::enable_if<!Implementation::IsDual<U>::value, void>::type> Dual<decltype(std::declval<T>()*std::declval<U>())> operator*(const U& other) const {
return {_real*other, _dual*other};
}
/**
* @brief Divide by dual number
*
* @f[
* \frac{\hat a}{\hat b} = \frac{a_0}{b_0} + \epsilon \frac{a_\epsilon b_0 - a_0 b_\epsilon}{b_0^2}
* @f]
* @see @ref operator/(const U&) const
*/
template<class U> auto operator/(const Dual<U>& other) const -> Dual<decltype(std::declval<T>()/std::declval<U>())> {
return {_real/other._real, (_dual*other._real - _real*other._dual)/(other._real*other._real)};
}
/**
* @brief Divide by real number
*
* Equivalent to the above assuming that @f$ b_\epsilon = 0 @f$.
* @f[
* \frac{\hat a}{\hat b} = \frac{a_0}{b_0} + \epsilon \frac{a_\epsilon b_0 - a_0 b_\epsilon}{b_0^2} = \frac{a_0}{b_0} + \epsilon \frac{a_\epsilon}{b_0}
* @f]
* @see @ref operator/(const Dual<U>&) const
*/
template<class U, class V = typename std::enable_if<!Implementation::IsDual<U>::value, Dual<decltype(std::declval<T>()/std::declval<U>())>>::type> V operator/(const U& other) const {
return {_real/other, _dual/other};
}
/**
* @brief Conjugated dual number
*
@ -196,6 +230,19 @@ template<class T> class Dual {
T _real, _dual;
};
/** @relates Dual
@brief Multiply real number by dual number
Equivalent to @ref Dual::operator*(const Dual<U>&) const assuming that
@f$ a_\epsilon = 0 @f$.
@f[
\hat a \hat b = a_0 b_0 + \epsilon (a_0 b_\epsilon + a_\epsilon b_0) = a_0 b_0 + \epsilon a_0 b_\epsilon
@f]
*/
template<class T, class U, class V = typename std::enable_if<!Implementation::IsDual<T>::value, Dual<decltype(std::declval<T>()*std::declval<U>())>>::type> inline V operator*(const T& a, const Dual<U>& b) {
return {a*b.real(), a*b.dual()};
}
#ifndef DOXYGEN_GENERATING_OUTPUT
#define MAGNUM_DUAL_SUBCLASS_IMPLEMENTATION(Type, Underlying, Multiplicable) \
Type<T> operator-() const { \
@ -218,8 +265,14 @@ template<class T> class Dual {
Type<T> operator*(const Math::Dual<Multiplicable>& other) const { \
return Math::Dual<Underlying<T>>::operator*(other); \
} \
Type<T> operator*(const Multiplicable& other) const { \
return Math::Dual<Underlying<T>>::operator*(other); \
} \
Type<T> operator/(const Math::Dual<Multiplicable>& other) const { \
return Math::Dual<Underlying<T>>::operator/(other); \
} \
Type<T> operator/(const Multiplicable& other) const { \
return Math::Dual<Underlying<T>>::operator/(other); \
}
/* DualComplex needs its own special implementation of multiplication/division */
@ -241,6 +294,9 @@ template<class T> class Dual {
template<class T> inline Type<T> operator*(const Math::Dual<Multiplicable>& a, const Type<T>& b) { \
return a*static_cast<const Math::Dual<Underlying<T>>&>(b); \
} \
template<class T> inline Type<T> operator*(const Multiplicable& a, const Type<T>& b) { \
return a*static_cast<const Math::Dual<Underlying<T>>&>(b); \
} \
template<class T> inline Type<T> operator/(const Math::Dual<Multiplicable>& a, const Type<T>& b) { \
return a/static_cast<const Math::Dual<Underlying<T>>&>(b); \
}

20
src/Magnum/Math/Test/DualTest.cpp

@ -45,6 +45,7 @@ struct DualTest: Corrade::TestSuite::Tester {
void addSubtract();
void negated();
void multiplyDivide();
void multiplyDivideScalar();
void multiplyDivideDifferentType();
void conjugated();
@ -72,6 +73,7 @@ DualTest::DualTest() {
&DualTest::addSubtract,
&DualTest::negated,
&DualTest::multiplyDivide,
&DualTest::multiplyDivideScalar,
&DualTest::multiplyDivideDifferentType,
&DualTest::conjugated,
@ -157,6 +159,15 @@ void DualTest::multiplyDivide() {
CORRADE_COMPARE(c/b, a);
}
void DualTest::multiplyDivideScalar() {
Dual a{1.5f, -4.0f};
Dual b{-3.0f, 8.0f};
CORRADE_COMPARE(a*-2.0f, b);
CORRADE_COMPARE(-2.0f*a, b);
CORRADE_COMPARE(b/-2.0f, a);
}
void DualTest::multiplyDivideDifferentType() {
DualVector2 a{{1.5f, 2.0f}, {-4.0f, 1.3f}};
Dual b{-2.0f, 0.5f};
@ -211,6 +222,10 @@ void DualTest::subclassTypes() {
CORRADE_VERIFY((std::is_same<decltype(c*a), DualVec2>::value));
CORRADE_VERIFY((std::is_same<decltype(a/c), DualVec2>::value));
CORRADE_VERIFY((std::is_same<decltype(c/a), DualVec2>::value));
CORRADE_VERIFY((std::is_same<decltype(a*5.0f), DualVec2>::value));
CORRADE_VERIFY((std::is_same<decltype(5.0f*a), DualVec2>::value));
CORRADE_VERIFY((std::is_same<decltype(a/5.0f), DualVec2>::value));
}
void DualTest::subclass() {
@ -235,6 +250,11 @@ void DualTest::subclass() {
CORRADE_COMPARE(e*a, f);
CORRADE_COMPARE(f/e, a);
CORRADE_COMPARE(e/a, g);
const DualVec2 h{Vector2{-3.0f, -4.0f}, Vector2{8.0f, -2.6f}};
CORRADE_COMPARE(a*-2.0f, h);
CORRADE_COMPARE(-2.0f*a, h);
CORRADE_COMPARE(h/-2.0f, a);
}
void DualTest::debug() {

Loading…
Cancel
Save