1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com> 5 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #ifndef EIGEN_EIGENBASE_H 12 #define EIGEN_EIGENBASE_H 13 14 namespace Eigen { 15 16 /** \class EigenBase 17 * \ingroup Core_Module 18 * 19 * Common base class for all classes T such that MatrixBase has an operator=(T) and a constructor MatrixBase(T). 20 * 21 * In other words, an EigenBase object is an object that can be copied into a MatrixBase. 22 * 23 * Besides MatrixBase-derived classes, this also includes special matrix classes such as diagonal matrices, etc. 24 * 25 * Notice that this class is trivial, it is only used to disambiguate overloaded functions. 26 * 27 * \sa \blank \ref TopicClassHierarchy 28 */ 29 template<typename Derived> struct EigenBase 30 { 31 // typedef typename internal::plain_matrix_type<Derived>::type PlainObject; 32 33 /** \brief The interface type of indices 34 * \details To change this, \c \#define the preprocessor symbol \c EIGEN_DEFAULT_DENSE_INDEX_TYPE. 35 * \sa StorageIndex, \ref TopicPreprocessorDirectives. 36 * DEPRECATED: Since Eigen 3.3, its usage is deprecated. Use Eigen::Index instead. 37 * Deprecation is not marked with a doxygen comment because there are too many existing usages to add the deprecation attribute. 38 */ 39 typedef Eigen::Index Index; 40 41 // FIXME is it needed? 42 typedef typename internal::traits<Derived>::StorageKind StorageKind; 43 44 /** \returns a reference to the derived object */ 45 EIGEN_DEVICE_FUNC derivedEigenBase46 Derived& derived() { return *static_cast<Derived*>(this); } 47 /** \returns a const reference to the derived object */ 48 EIGEN_DEVICE_FUNC derivedEigenBase49 const Derived& derived() const { return *static_cast<const Derived*>(this); } 50 51 EIGEN_DEVICE_FUNC const_cast_derivedEigenBase52 inline Derived& const_cast_derived() const 53 { return *static_cast<Derived*>(const_cast<EigenBase*>(this)); } 54 EIGEN_DEVICE_FUNC const_derivedEigenBase55 inline const Derived& const_derived() const 56 { return *static_cast<const Derived*>(this); } 57 58 /** \returns the number of rows. \sa cols(), RowsAtCompileTime */ 59 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR rowsEigenBase60 inline Index rows() const EIGEN_NOEXCEPT { return derived().rows(); } 61 /** \returns the number of columns. \sa rows(), ColsAtCompileTime*/ 62 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR colsEigenBase63 inline Index cols() const EIGEN_NOEXCEPT { return derived().cols(); } 64 /** \returns the number of coefficients, which is rows()*cols(). 65 * \sa rows(), cols(), SizeAtCompileTime. */ 66 EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR sizeEigenBase67 inline Index size() const EIGEN_NOEXCEPT { return rows() * cols(); } 68 69 /** \internal Don't use it, but do the equivalent: \code dst = *this; \endcode */ 70 template<typename Dest> 71 EIGEN_DEVICE_FUNC evalToEigenBase72 inline void evalTo(Dest& dst) const 73 { derived().evalTo(dst); } 74 75 /** \internal Don't use it, but do the equivalent: \code dst += *this; \endcode */ 76 template<typename Dest> 77 EIGEN_DEVICE_FUNC addToEigenBase78 inline void addTo(Dest& dst) const 79 { 80 // This is the default implementation, 81 // derived class can reimplement it in a more optimized way. 82 typename Dest::PlainObject res(rows(),cols()); 83 evalTo(res); 84 dst += res; 85 } 86 87 /** \internal Don't use it, but do the equivalent: \code dst -= *this; \endcode */ 88 template<typename Dest> 89 EIGEN_DEVICE_FUNC subToEigenBase90 inline void subTo(Dest& dst) const 91 { 92 // This is the default implementation, 93 // derived class can reimplement it in a more optimized way. 94 typename Dest::PlainObject res(rows(),cols()); 95 evalTo(res); 96 dst -= res; 97 } 98 99 /** \internal Don't use it, but do the equivalent: \code dst.applyOnTheRight(*this); \endcode */ 100 template<typename Dest> applyThisOnTheRightEigenBase101 EIGEN_DEVICE_FUNC inline void applyThisOnTheRight(Dest& dst) const 102 { 103 // This is the default implementation, 104 // derived class can reimplement it in a more optimized way. 105 dst = dst * this->derived(); 106 } 107 108 /** \internal Don't use it, but do the equivalent: \code dst.applyOnTheLeft(*this); \endcode */ 109 template<typename Dest> applyThisOnTheLeftEigenBase110 EIGEN_DEVICE_FUNC inline void applyThisOnTheLeft(Dest& dst) const 111 { 112 // This is the default implementation, 113 // derived class can reimplement it in a more optimized way. 114 dst = this->derived() * dst; 115 } 116 117 }; 118 119 /*************************************************************************** 120 * Implementation of matrix base methods 121 ***************************************************************************/ 122 123 /** \brief Copies the generic expression \a other into *this. 124 * 125 * \details The expression must provide a (templated) evalTo(Derived& dst) const 126 * function which does the actual job. In practice, this allows any user to write 127 * its own special matrix without having to modify MatrixBase 128 * 129 * \returns a reference to *this. 130 */ 131 template<typename Derived> 132 template<typename OtherDerived> 133 EIGEN_DEVICE_FUNC 134 Derived& DenseBase<Derived>::operator=(const EigenBase<OtherDerived> &other) 135 { 136 call_assignment(derived(), other.derived()); 137 return derived(); 138 } 139 140 template<typename Derived> 141 template<typename OtherDerived> 142 EIGEN_DEVICE_FUNC 143 Derived& DenseBase<Derived>::operator+=(const EigenBase<OtherDerived> &other) 144 { 145 call_assignment(derived(), other.derived(), internal::add_assign_op<Scalar,typename OtherDerived::Scalar>()); 146 return derived(); 147 } 148 149 template<typename Derived> 150 template<typename OtherDerived> 151 EIGEN_DEVICE_FUNC 152 Derived& DenseBase<Derived>::operator-=(const EigenBase<OtherDerived> &other) 153 { 154 call_assignment(derived(), other.derived(), internal::sub_assign_op<Scalar,typename OtherDerived::Scalar>()); 155 return derived(); 156 } 157 158 } // end namespace Eigen 159 160 #endif // EIGEN_EIGENBASE_H 161