AutoDiffVector.h
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // Eigen is free software; you can redistribute it and/or
7 // modify it under the terms of the GNU Lesser General Public
8 // License as published by the Free Software Foundation; either
9 // version 3 of the License, or (at your option) any later version.
10 //
11 // Alternatively, you can redistribute it and/or
12 // modify it under the terms of the GNU General Public License as
13 // published by the Free Software Foundation; either version 2 of
14 // the License, or (at your option) any later version.
15 //
16 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19 // GNU General Public License for more details.
20 //
21 // You should have received a copy of the GNU Lesser General Public
22 // License and a copy of the GNU General Public License along with
23 // Eigen. If not, see <http://www.gnu.org/licenses/>.
24 
25 #ifndef EIGEN_AUTODIFF_VECTOR_H
26 #define EIGEN_AUTODIFF_VECTOR_H
27 
28 namespace Eigen {
29 
30 /* \class AutoDiffScalar
31  * \brief A scalar type replacement with automatic differentation capability
32  *
33  * \param DerType the vector type used to store/represent the derivatives (e.g. Vector3f)
34  *
35  * This class represents a scalar value while tracking its respective derivatives.
36  *
37  * It supports the following list of global math function:
38  * - std::abs, std::sqrt, std::pow, std::exp, std::log, std::sin, std::cos,
39  * - internal::abs, internal::sqrt, internal::pow, internal::exp, internal::log, internal::sin, internal::cos,
40  * - internal::conj, internal::real, internal::imag, internal::abs2.
41  *
42  * AutoDiffScalar can be used as the scalar type of an Eigen::Matrix object. However,
43  * in that case, the expression template mechanism only occurs at the top Matrix level,
44  * while derivatives are computed right away.
45  *
46  */
47 template<typename ValueType, typename JacobianType>
48 class AutoDiffVector
49 {
50  public:
51  //typedef typename internal::traits<ValueType>::Scalar Scalar;
52  typedef typename internal::traits<ValueType>::Scalar BaseScalar;
53  typedef AutoDiffScalar<Matrix<BaseScalar,JacobianType::RowsAtCompileTime,1> > ActiveScalar;
54  typedef ActiveScalar Scalar;
55  typedef AutoDiffScalar<typename JacobianType::ColXpr> CoeffType;
56  typedef typename JacobianType::Index Index;
57 
58  inline AutoDiffVector() {}
59 
60  inline AutoDiffVector(const ValueType& values)
61  : m_values(values)
62  {
63  m_jacobian.setZero();
64  }
65 
66 
67  CoeffType operator[] (Index i) { return CoeffType(m_values[i], m_jacobian.col(i)); }
68  const CoeffType operator[] (Index i) const { return CoeffType(m_values[i], m_jacobian.col(i)); }
69 
70  CoeffType operator() (Index i) { return CoeffType(m_values[i], m_jacobian.col(i)); }
71  const CoeffType operator() (Index i) const { return CoeffType(m_values[i], m_jacobian.col(i)); }
72 
73  CoeffType coeffRef(Index i) { return CoeffType(m_values[i], m_jacobian.col(i)); }
74  const CoeffType coeffRef(Index i) const { return CoeffType(m_values[i], m_jacobian.col(i)); }
75 
76  Index size() const { return m_values.size(); }
77 
78  // FIXME here we could return an expression of the sum
79  Scalar sum() const { /*std::cerr << "sum \n\n";*/ /*std::cerr << m_jacobian.rowwise().sum() << "\n\n";*/ return Scalar(m_values.sum(), m_jacobian.rowwise().sum()); }
80 
81 
82  inline AutoDiffVector(const ValueType& values, const JacobianType& jac)
83  : m_values(values), m_jacobian(jac)
84  {}
85 
86  template<typename OtherValueType, typename OtherJacobianType>
87  inline AutoDiffVector(const AutoDiffVector<OtherValueType, OtherJacobianType>& other)
88  : m_values(other.values()), m_jacobian(other.jacobian())
89  {}
90 
91  inline AutoDiffVector(const AutoDiffVector& other)
92  : m_values(other.values()), m_jacobian(other.jacobian())
93  {}
94 
95  template<typename OtherValueType, typename OtherJacobianType>
96  inline AutoDiffVector& operator=(const AutoDiffVector<OtherValueType, OtherJacobianType>& other)
97  {
98  m_values = other.values();
99  m_jacobian = other.jacobian();
100  return *this;
101  }
102 
103  inline AutoDiffVector& operator=(const AutoDiffVector& other)
104  {
105  m_values = other.values();
106  m_jacobian = other.jacobian();
107  return *this;
108  }
109 
110  inline const ValueType& values() const { return m_values; }
111  inline ValueType& values() { return m_values; }
112 
113  inline const JacobianType& jacobian() const { return m_jacobian; }
114  inline JacobianType& jacobian() { return m_jacobian; }
115 
116  template<typename OtherValueType,typename OtherJacobianType>
117  inline const AutoDiffVector<
118  typename MakeCwiseBinaryOp<internal::scalar_sum_op<BaseScalar>,ValueType,OtherValueType>::Type,
119  typename MakeCwiseBinaryOp<internal::scalar_sum_op<BaseScalar>,JacobianType,OtherJacobianType>::Type >
120  operator+(const AutoDiffVector<OtherValueType,OtherJacobianType>& other) const
121  {
122  return AutoDiffVector<
123  typename MakeCwiseBinaryOp<internal::scalar_sum_op<BaseScalar>,ValueType,OtherValueType>::Type,
124  typename MakeCwiseBinaryOp<internal::scalar_sum_op<BaseScalar>,JacobianType,OtherJacobianType>::Type >(
125  m_values + other.values(),
126  m_jacobian + other.jacobian());
127  }
128 
129  template<typename OtherValueType, typename OtherJacobianType>
130  inline AutoDiffVector&
131  operator+=(const AutoDiffVector<OtherValueType,OtherJacobianType>& other)
132  {
133  m_values += other.values();
134  m_jacobian += other.jacobian();
135  return *this;
136  }
137 
138  template<typename OtherValueType,typename OtherJacobianType>
139  inline const AutoDiffVector<
140  typename MakeCwiseBinaryOp<internal::scalar_difference_op<Scalar>,ValueType,OtherValueType>::Type,
141  typename MakeCwiseBinaryOp<internal::scalar_difference_op<Scalar>,JacobianType,OtherJacobianType>::Type >
142  operator-(const AutoDiffVector<OtherValueType,OtherJacobianType>& other) const
143  {
144  return AutoDiffVector<
145  typename MakeCwiseBinaryOp<internal::scalar_difference_op<Scalar>,ValueType,OtherValueType>::Type,
146  typename MakeCwiseBinaryOp<internal::scalar_difference_op<Scalar>,JacobianType,OtherJacobianType>::Type >(
147  m_values - other.values(),
148  m_jacobian - other.jacobian());
149  }
150 
151  template<typename OtherValueType, typename OtherJacobianType>
152  inline AutoDiffVector&
153  operator-=(const AutoDiffVector<OtherValueType,OtherJacobianType>& other)
154  {
155  m_values -= other.values();
156  m_jacobian -= other.jacobian();
157  return *this;
158  }
159 
160  inline const AutoDiffVector<
161  typename MakeCwiseUnaryOp<internal::scalar_opposite_op<Scalar>, ValueType>::Type,
162  typename MakeCwiseUnaryOp<internal::scalar_opposite_op<Scalar>, JacobianType>::Type >
163  operator-() const
164  {
165  return AutoDiffVector<
166  typename MakeCwiseUnaryOp<internal::scalar_opposite_op<Scalar>, ValueType>::Type,
167  typename MakeCwiseUnaryOp<internal::scalar_opposite_op<Scalar>, JacobianType>::Type >(
168  -m_values,
169  -m_jacobian);
170  }
171 
172  inline const AutoDiffVector<
173  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, ValueType>::Type,
174  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>::Type>
175  operator*(const BaseScalar& other) const
176  {
177  return AutoDiffVector<
178  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, ValueType>::Type,
179  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>::Type >(
180  m_values * other,
181  m_jacobian * other);
182  }
183 
184  friend inline const AutoDiffVector<
185  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, ValueType>::Type,
186  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>::Type >
187  operator*(const Scalar& other, const AutoDiffVector& v)
188  {
189  return AutoDiffVector<
190  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, ValueType>::Type,
191  typename MakeCwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>::Type >(
192  v.values() * other,
193  v.jacobian() * other);
194  }
195 
196 // template<typename OtherValueType,typename OtherJacobianType>
197 // inline const AutoDiffVector<
198 // CwiseBinaryOp<internal::scalar_multiple_op<Scalar>, ValueType, OtherValueType>
199 // CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
200 // CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>,
201 // CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, OtherJacobianType> > >
202 // operator*(const AutoDiffVector<OtherValueType,OtherJacobianType>& other) const
203 // {
204 // return AutoDiffVector<
205 // CwiseBinaryOp<internal::scalar_multiple_op<Scalar>, ValueType, OtherValueType>
206 // CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
207 // CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, JacobianType>,
208 // CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, OtherJacobianType> > >(
209 // m_values.cwise() * other.values(),
210 // (m_jacobian * other.values()) + (m_values * other.jacobian()));
211 // }
212 
213  inline AutoDiffVector& operator*=(const Scalar& other)
214  {
215  m_values *= other;
216  m_jacobian *= other;
217  return *this;
218  }
219 
220  template<typename OtherValueType,typename OtherJacobianType>
221  inline AutoDiffVector& operator*=(const AutoDiffVector<OtherValueType,OtherJacobianType>& other)
222  {
223  *this = *this * other;
224  return *this;
225  }
226 
227  protected:
228  ValueType m_values;
229  JacobianType m_jacobian;
230 
231 };
232 
233 }
234 
235 #endif // EIGEN_AUTODIFF_VECTOR_H