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
14namespace 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  */
29template<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    * \deprecated Since Eigen 3.3, its usage is deprecated. Use Eigen::Index instead.
36    * \sa StorageIndex, \ref TopicPreprocessorDirectives.
37    */
38  typedef Eigen::Index Index;
39
40  // FIXME is it needed?
41  typedef typename internal::traits<Derived>::StorageKind StorageKind;
42
43  /** \returns a reference to the derived object */
44  EIGEN_DEVICE_FUNC
45  Derived& derived() { return *static_cast<Derived*>(this); }
46  /** \returns a const reference to the derived object */
47  EIGEN_DEVICE_FUNC
48  const Derived& derived() const { return *static_cast<const Derived*>(this); }
49
50  EIGEN_DEVICE_FUNC
51  inline Derived& const_cast_derived() const
52  { return *static_cast<Derived*>(const_cast<EigenBase*>(this)); }
53  EIGEN_DEVICE_FUNC
54  inline const Derived& const_derived() const
55  { return *static_cast<const Derived*>(this); }
56
57  /** \returns the number of rows. \sa cols(), RowsAtCompileTime */
58  EIGEN_DEVICE_FUNC
59  inline Index rows() const { return derived().rows(); }
60  /** \returns the number of columns. \sa rows(), ColsAtCompileTime*/
61  EIGEN_DEVICE_FUNC
62  inline Index cols() const { return derived().cols(); }
63  /** \returns the number of coefficients, which is rows()*cols().
64    * \sa rows(), cols(), SizeAtCompileTime. */
65  EIGEN_DEVICE_FUNC
66  inline Index size() const { return rows() * cols(); }
67
68  /** \internal Don't use it, but do the equivalent: \code dst = *this; \endcode */
69  template<typename Dest>
70  EIGEN_DEVICE_FUNC
71  inline void evalTo(Dest& dst) const
72  { derived().evalTo(dst); }
73
74  /** \internal Don't use it, but do the equivalent: \code dst += *this; \endcode */
75  template<typename Dest>
76  EIGEN_DEVICE_FUNC
77  inline void addTo(Dest& dst) const
78  {
79    // This is the default implementation,
80    // derived class can reimplement it in a more optimized way.
81    typename Dest::PlainObject res(rows(),cols());
82    evalTo(res);
83    dst += res;
84  }
85
86  /** \internal Don't use it, but do the equivalent: \code dst -= *this; \endcode */
87  template<typename Dest>
88  EIGEN_DEVICE_FUNC
89  inline void subTo(Dest& dst) const
90  {
91    // This is the default implementation,
92    // derived class can reimplement it in a more optimized way.
93    typename Dest::PlainObject res(rows(),cols());
94    evalTo(res);
95    dst -= res;
96  }
97
98  /** \internal Don't use it, but do the equivalent: \code dst.applyOnTheRight(*this); \endcode */
99  template<typename Dest>
100  EIGEN_DEVICE_FUNC inline void applyThisOnTheRight(Dest& dst) const
101  {
102    // This is the default implementation,
103    // derived class can reimplement it in a more optimized way.
104    dst = dst * this->derived();
105  }
106
107  /** \internal Don't use it, but do the equivalent: \code dst.applyOnTheLeft(*this); \endcode */
108  template<typename Dest>
109  EIGEN_DEVICE_FUNC inline void applyThisOnTheLeft(Dest& dst) const
110  {
111    // This is the default implementation,
112    // derived class can reimplement it in a more optimized way.
113    dst = this->derived() * dst;
114  }
115
116};
117
118/***************************************************************************
119* Implementation of matrix base methods
120***************************************************************************/
121
122/** \brief Copies the generic expression \a other into *this.
123  *
124  * \details The expression must provide a (templated) evalTo(Derived& dst) const
125  * function which does the actual job. In practice, this allows any user to write
126  * its own special matrix without having to modify MatrixBase
127  *
128  * \returns a reference to *this.
129  */
130template<typename Derived>
131template<typename OtherDerived>
132EIGEN_DEVICE_FUNC
133Derived& DenseBase<Derived>::operator=(const EigenBase<OtherDerived> &other)
134{
135  call_assignment(derived(), other.derived());
136  return derived();
137}
138
139template<typename Derived>
140template<typename OtherDerived>
141EIGEN_DEVICE_FUNC
142Derived& DenseBase<Derived>::operator+=(const EigenBase<OtherDerived> &other)
143{
144  call_assignment(derived(), other.derived(), internal::add_assign_op<Scalar,typename OtherDerived::Scalar>());
145  return derived();
146}
147
148template<typename Derived>
149template<typename OtherDerived>
150EIGEN_DEVICE_FUNC
151Derived& DenseBase<Derived>::operator-=(const EigenBase<OtherDerived> &other)
152{
153  call_assignment(derived(), other.derived(), internal::sub_assign_op<Scalar,typename OtherDerived::Scalar>());
154  return derived();
155}
156
157} // end namespace Eigen
158
159#endif // EIGEN_EIGENBASE_H
160