AutoDiffScalar.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 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #ifndef EIGEN_AUTODIFF_SCALAR_H
11 #define EIGEN_AUTODIFF_SCALAR_H
12 
13 namespace Eigen {
14 
15 namespace internal {
16 
17 template<typename A, typename B>
18 struct make_coherent_impl {
19  static void run(A&, B&) {}
20 };
21 
22 // resize a to match b is a.size()==0, and conversely.
23 template<typename A, typename B>
24 void make_coherent(const A& a, const B&b)
25 {
26  make_coherent_impl<A,B>::run(a.const_cast_derived(), b.const_cast_derived());
27 }
28 
29 template<typename _DerType, bool Enable> struct auto_diff_special_op;
30 
31 } // end namespace internal
32 
59 template<typename _DerType>
61  : public internal::auto_diff_special_op
62  <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
63  typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value>
64 {
65  public:
66  typedef internal::auto_diff_special_op
67  <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
68  typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value> Base;
69  typedef typename internal::remove_all<_DerType>::type DerType;
70  typedef typename internal::traits<DerType>::Scalar Scalar;
71  typedef typename NumTraits<Scalar>::Real Real;
72 
73  using Base::operator+;
74  using Base::operator*;
75 
78 
81  AutoDiffScalar(const Scalar& value, int nbDer, int derNumber)
82  : m_value(value), m_derivatives(DerType::Zero(nbDer))
83  {
84  m_derivatives.coeffRef(derNumber) = Scalar(1);
85  }
86 
89  /*explicit*/ AutoDiffScalar(const Real& value)
90  : m_value(value)
91  {
92  if(m_derivatives.size()>0)
93  m_derivatives.setZero();
94  }
95 
97  AutoDiffScalar(const Scalar& value, const DerType& der)
98  : m_value(value), m_derivatives(der)
99  {}
100 
101  template<typename OtherDerType>
103 #ifndef EIGEN_PARSED_BY_DOXYGEN
104  , typename internal::enable_if<internal::is_same<Scalar,typename OtherDerType::Scalar>::value,void*>::type = 0
105 #endif
106  )
107  : m_value(other.value()), m_derivatives(other.derivatives())
108  {}
109 
110  friend std::ostream & operator << (std::ostream & s, const AutoDiffScalar& a)
111  {
112  return s << a.value();
113  }
114 
115  AutoDiffScalar(const AutoDiffScalar& other)
116  : m_value(other.value()), m_derivatives(other.derivatives())
117  {}
118 
119  template<typename OtherDerType>
120  inline AutoDiffScalar& operator=(const AutoDiffScalar<OtherDerType>& other)
121  {
122  m_value = other.value();
123  m_derivatives = other.derivatives();
124  return *this;
125  }
126 
127  inline AutoDiffScalar& operator=(const AutoDiffScalar& other)
128  {
129  m_value = other.value();
130  m_derivatives = other.derivatives();
131  return *this;
132  }
133 
134  inline AutoDiffScalar& operator=(const Scalar& other)
135  {
136  m_value = other;
137  if(m_derivatives.size()>0)
138  m_derivatives.setZero();
139  return *this;
140  }
141 
142 // inline operator const Scalar& () const { return m_value; }
143 // inline operator Scalar& () { return m_value; }
144 
145  inline const Scalar& value() const { return m_value; }
146  inline Scalar& value() { return m_value; }
147 
148  inline const DerType& derivatives() const { return m_derivatives; }
149  inline DerType& derivatives() { return m_derivatives; }
150 
151  inline bool operator< (const Scalar& other) const { return m_value < other; }
152  inline bool operator<=(const Scalar& other) const { return m_value <= other; }
153  inline bool operator> (const Scalar& other) const { return m_value > other; }
154  inline bool operator>=(const Scalar& other) const { return m_value >= other; }
155  inline bool operator==(const Scalar& other) const { return m_value == other; }
156  inline bool operator!=(const Scalar& other) const { return m_value != other; }
157 
158  friend inline bool operator< (const Scalar& a, const AutoDiffScalar& b) { return a < b.value(); }
159  friend inline bool operator<=(const Scalar& a, const AutoDiffScalar& b) { return a <= b.value(); }
160  friend inline bool operator> (const Scalar& a, const AutoDiffScalar& b) { return a > b.value(); }
161  friend inline bool operator>=(const Scalar& a, const AutoDiffScalar& b) { return a >= b.value(); }
162  friend inline bool operator==(const Scalar& a, const AutoDiffScalar& b) { return a == b.value(); }
163  friend inline bool operator!=(const Scalar& a, const AutoDiffScalar& b) { return a != b.value(); }
164 
165  template<typename OtherDerType> inline bool operator< (const AutoDiffScalar<OtherDerType>& b) const { return m_value < b.value(); }
166  template<typename OtherDerType> inline bool operator<=(const AutoDiffScalar<OtherDerType>& b) const { return m_value <= b.value(); }
167  template<typename OtherDerType> inline bool operator> (const AutoDiffScalar<OtherDerType>& b) const { return m_value > b.value(); }
168  template<typename OtherDerType> inline bool operator>=(const AutoDiffScalar<OtherDerType>& b) const { return m_value >= b.value(); }
169  template<typename OtherDerType> inline bool operator==(const AutoDiffScalar<OtherDerType>& b) const { return m_value == b.value(); }
170  template<typename OtherDerType> inline bool operator!=(const AutoDiffScalar<OtherDerType>& b) const { return m_value != b.value(); }
171 
172  inline const AutoDiffScalar<DerType&> operator+(const Scalar& other) const
173  {
174  return AutoDiffScalar<DerType&>(m_value + other, m_derivatives);
175  }
176 
177  friend inline const AutoDiffScalar<DerType&> operator+(const Scalar& a, const AutoDiffScalar& b)
178  {
179  return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
180  }
181 
182 // inline const AutoDiffScalar<DerType&> operator+(const Real& other) const
183 // {
184 // return AutoDiffScalar<DerType&>(m_value + other, m_derivatives);
185 // }
186 
187 // friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar& b)
188 // {
189 // return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
190 // }
191 
192  inline AutoDiffScalar& operator+=(const Scalar& other)
193  {
194  value() += other;
195  return *this;
196  }
197 
198  template<typename OtherDerType>
199  inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> >
200  operator+(const AutoDiffScalar<OtherDerType>& other) const
201  {
202  internal::make_coherent(m_derivatives, other.derivatives());
203  return AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> >(
204  m_value + other.value(),
205  m_derivatives + other.derivatives());
206  }
207 
208  template<typename OtherDerType>
209  inline AutoDiffScalar&
210  operator+=(const AutoDiffScalar<OtherDerType>& other)
211  {
212  (*this) = (*this) + other;
213  return *this;
214  }
215 
216  inline const AutoDiffScalar<DerType&> operator-(const Scalar& b) const
217  {
218  return AutoDiffScalar<DerType&>(m_value - b, m_derivatives);
219  }
220 
221  friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
222  operator-(const Scalar& a, const AutoDiffScalar& b)
223  {
224  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
225  (a - b.value(), -b.derivatives());
226  }
227 
228  inline AutoDiffScalar& operator-=(const Scalar& other)
229  {
230  value() -= other;
231  return *this;
232  }
233 
234  template<typename OtherDerType>
235  inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> >
236  operator-(const AutoDiffScalar<OtherDerType>& other) const
237  {
238  internal::make_coherent(m_derivatives, other.derivatives());
239  return AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> >(
240  m_value - other.value(),
241  m_derivatives - other.derivatives());
242  }
243 
244  template<typename OtherDerType>
245  inline AutoDiffScalar&
246  operator-=(const AutoDiffScalar<OtherDerType>& other)
247  {
248  *this = *this - other;
249  return *this;
250  }
251 
252  inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
253  operator-() const
254  {
255  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >(
256  -m_value,
257  -m_derivatives);
258  }
259 
261  operator*(const Scalar& other) const
262  {
263  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
264  m_value * other,
265  (m_derivatives * other));
266  }
267 
268  friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
269  operator*(const Scalar& other, const AutoDiffScalar& a)
270  {
271  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
272  a.value() * other,
273  a.derivatives() * other);
274  }
275 
276 // inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >
277 // operator*(const Real& other) const
278 // {
279 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >(
280 // m_value * other,
281 // (m_derivatives * other));
282 // }
283 //
284 // friend inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >
285 // operator*(const Real& other, const AutoDiffScalar& a)
286 // {
287 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >(
288 // a.value() * other,
289 // a.derivatives() * other);
290 // }
291 
292  inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
293  operator/(const Scalar& other) const
294  {
295  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
296  m_value / other,
297  (m_derivatives * (Scalar(1)/other)));
298  }
299 
300  friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
301  operator/(const Scalar& other, const AutoDiffScalar& a)
302  {
303  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
304  other / a.value(),
305  a.derivatives() * (Scalar(-other) / (a.value()*a.value())));
306  }
307 
308 // inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >
309 // operator/(const Real& other) const
310 // {
311 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >(
312 // m_value / other,
313 // (m_derivatives * (Real(1)/other)));
314 // }
315 //
316 // friend inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >
317 // operator/(const Real& other, const AutoDiffScalar& a)
318 // {
319 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >(
320 // other / a.value(),
321 // a.derivatives() * (-Real(1)/other));
322 // }
323 
324  template<typename OtherDerType>
325  inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>,
326  const CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
327  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
328  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > >
329  operator/(const AutoDiffScalar<OtherDerType>& other) const
330  {
331  internal::make_coherent(m_derivatives, other.derivatives());
332  return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>,
333  const CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
334  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
335  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > >(
336  m_value / other.value(),
337  ((m_derivatives * other.value()) - (m_value * other.derivatives()))
338  * (Scalar(1)/(other.value()*other.value())));
339  }
340 
341  template<typename OtherDerType>
343  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
344  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type> > >
345  operator*(const AutoDiffScalar<OtherDerType>& other) const
346  {
347  internal::make_coherent(m_derivatives, other.derivatives());
349  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
350  const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > >(
351  m_value * other.value(),
352  (m_derivatives * other.value()) + (m_value * other.derivatives()));
353  }
354 
355  inline AutoDiffScalar& operator*=(const Scalar& other)
356  {
357  *this = *this * other;
358  return *this;
359  }
360 
361  template<typename OtherDerType>
362  inline AutoDiffScalar& operator*=(const AutoDiffScalar<OtherDerType>& other)
363  {
364  *this = *this * other;
365  return *this;
366  }
367 
368  inline AutoDiffScalar& operator/=(const Scalar& other)
369  {
370  *this = *this / other;
371  return *this;
372  }
373 
374  template<typename OtherDerType>
375  inline AutoDiffScalar& operator/=(const AutoDiffScalar<OtherDerType>& other)
376  {
377  *this = *this / other;
378  return *this;
379  }
380 
381  protected:
382  Scalar m_value;
383  DerType m_derivatives;
384 
385 };
386 
387 namespace internal {
388 
389 template<typename _DerType>
390 struct auto_diff_special_op<_DerType, true>
391 // : auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real,
392 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value>
393 {
394  typedef typename remove_all<_DerType>::type DerType;
395  typedef typename traits<DerType>::Scalar Scalar;
396  typedef typename NumTraits<Scalar>::Real Real;
397 
398 // typedef auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real,
399 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value> Base;
400 
401 // using Base::operator+;
402 // using Base::operator+=;
403 // using Base::operator-;
404 // using Base::operator-=;
405 // using Base::operator*;
406 // using Base::operator*=;
407 
408  const AutoDiffScalar<_DerType>& derived() const { return *static_cast<const AutoDiffScalar<_DerType>*>(this); }
409  AutoDiffScalar<_DerType>& derived() { return *static_cast<AutoDiffScalar<_DerType>*>(this); }
410 
411 
412  inline const AutoDiffScalar<DerType&> operator+(const Real& other) const
413  {
414  return AutoDiffScalar<DerType&>(derived().value() + other, derived().derivatives());
415  }
416 
417  friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar<_DerType>& b)
418  {
419  return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
420  }
421 
422  inline AutoDiffScalar<_DerType>& operator+=(const Real& other)
423  {
424  derived().value() += other;
425  return derived();
426  }
427 
428 
430  operator*(const Real& other) const
431  {
433  derived().value() * other,
434  derived().derivatives() * other);
435  }
436 
438  operator*(const Real& other, const AutoDiffScalar<_DerType>& a)
439  {
440  return AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >(
441  a.value() * other,
442  a.derivatives() * other);
443  }
444 
445  inline AutoDiffScalar<_DerType>& operator*=(const Scalar& other)
446  {
447  *this = *this * other;
448  return derived();
449  }
450 };
451 
452 template<typename _DerType>
453 struct auto_diff_special_op<_DerType, false>
454 {
455  void operator*() const;
456  void operator-() const;
457  void operator+() const;
458 };
459 
460 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols, typename B>
461 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>, B> {
462  typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
463  static void run(A& a, B& b) {
464  if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
465  {
466  a.resize(b.size());
467  a.setZero();
468  }
469  }
470 };
471 
472 template<typename A, typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols>
473 struct make_coherent_impl<A, Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
474  typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
475  static void run(A& a, B& b) {
476  if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
477  {
478  b.resize(a.size());
479  b.setZero();
480  }
481  }
482 };
483 
484 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols,
485  typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols>
486 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>,
487  Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
488  typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
489  typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
490  static void run(A& a, B& b) {
491  if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
492  {
493  a.resize(b.size());
494  a.setZero();
495  }
496  else if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
497  {
498  b.resize(a.size());
499  b.setZero();
500  }
501  }
502 };
503 
504 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols>
505 struct scalar_product_traits<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>,A_Scalar>
506 {
507  enum { Defined = 1 };
508  typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType;
509 };
510 
511 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols>
512 struct scalar_product_traits<A_Scalar, Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> >
513 {
514  enum { Defined = 1 };
515  typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType;
516 };
517 
518 template<typename DerType>
519 struct scalar_product_traits<AutoDiffScalar<DerType>,typename DerType::Scalar>
520 {
521  enum { Defined = 1 };
522  typedef AutoDiffScalar<DerType> ReturnType;
523 };
524 
525 template<typename DerType>
526 struct scalar_product_traits<typename DerType::Scalar,AutoDiffScalar<DerType> >
527 {
528  enum { Defined = 1 };
529  typedef AutoDiffScalar<DerType> ReturnType;
530 };
531 
532 } // end namespace internal
533 
534 #define EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(FUNC,CODE) \
535  template<typename DerType> \
536  inline const Eigen::AutoDiffScalar<Eigen::CwiseUnaryOp<Eigen::internal::scalar_multiple_op<typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar>, const typename Eigen::internal::remove_all<DerType>::type> > \
537  FUNC(const Eigen::AutoDiffScalar<DerType>& x) { \
538  using namespace Eigen; \
539  typedef typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar Scalar; \
540  typedef AutoDiffScalar<CwiseUnaryOp<Eigen::internal::scalar_multiple_op<Scalar>, const typename Eigen::internal::remove_all<DerType>::type> > ReturnType; \
541  CODE; \
542  }
543 
544 template<typename DerType>
545 inline const AutoDiffScalar<DerType>& conj(const AutoDiffScalar<DerType>& x) { return x; }
546 template<typename DerType>
547 inline const AutoDiffScalar<DerType>& real(const AutoDiffScalar<DerType>& x) { return x; }
548 template<typename DerType>
549 inline typename DerType::Scalar imag(const AutoDiffScalar<DerType>&) { return 0.; }
550 template<typename DerType, typename T>
551 inline AutoDiffScalar<DerType> (min)(const AutoDiffScalar<DerType>& x, const T& y) { return (x <= y ? x : y); }
552 template<typename DerType, typename T>
553 inline AutoDiffScalar<DerType> (max)(const AutoDiffScalar<DerType>& x, const T& y) { return (x >= y ? x : y); }
554 template<typename DerType, typename T>
555 inline AutoDiffScalar<DerType> (min)(const T& x, const AutoDiffScalar<DerType>& y) { return (x < y ? x : y); }
556 template<typename DerType, typename T>
557 inline AutoDiffScalar<DerType> (max)(const T& x, const AutoDiffScalar<DerType>& y) { return (x > y ? x : y); }
558 
559 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs,
560  using std::abs;
561  return ReturnType(abs(x.value()), x.derivatives() * (x.value()<0 ? -1 : 1) );)
562 
563 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs2,
564  using numext::abs2;
565  return ReturnType(abs2(x.value()), x.derivatives() * (Scalar(2)*x.value()));)
566 
567 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sqrt,
568  using std::sqrt;
569  Scalar sqrtx = sqrt(x.value());
570  return ReturnType(sqrtx,x.derivatives() * (Scalar(0.5) / sqrtx));)
571 
572 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cos,
573  using std::cos;
574  using std::sin;
575  return ReturnType(cos(x.value()), x.derivatives() * (-sin(x.value())));)
576 
577 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sin,
578  using std::sin;
579  using std::cos;
580  return ReturnType(sin(x.value()),x.derivatives() * cos(x.value()));)
581 
582 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(exp,
583  using std::exp;
584  Scalar expx = exp(x.value());
585  return ReturnType(expx,x.derivatives() * expx);)
586 
587 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(log,
588  using std::log;
589  return ReturnType(log(x.value()),x.derivatives() * (Scalar(1)/x.value()));)
590 
591 template<typename DerType>
592 inline const Eigen::AutoDiffScalar<Eigen::CwiseUnaryOp<Eigen::internal::scalar_multiple_op<typename Eigen::internal::traits<DerType>::Scalar>, const DerType> >
593 pow(const Eigen::AutoDiffScalar<DerType>& x, typename Eigen::internal::traits<DerType>::Scalar y)
594 {
595  using namespace Eigen;
596  typedef typename Eigen::internal::traits<DerType>::Scalar Scalar;
598  std::pow(x.value(),y),
599  x.derivatives() * (y * std::pow(x.value(),y-1)));
600 }
601 
602 
603 template<typename DerTypeA,typename DerTypeB>
605 atan2(const AutoDiffScalar<DerTypeA>& a, const AutoDiffScalar<DerTypeB>& b)
606 {
607  using std::atan2;
608  typedef typename internal::traits<DerTypeA>::Scalar Scalar;
609  typedef AutoDiffScalar<Matrix<Scalar,Dynamic,1> > PlainADS;
610  PlainADS ret;
611  ret.value() = atan2(a.value(), b.value());
612 
613  Scalar squared_hypot = a.value() * a.value() + b.value() * b.value();
614 
615  // if (squared_hypot==0) the derivation is undefined and the following results in a NaN:
616  ret.derivatives() = (a.derivatives() * b.value() - a.value() * b.derivatives()) / squared_hypot;
617 
618  return ret;
619 }
620 
621 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tan,
622  using std::tan;
623  using std::cos;
624  return ReturnType(tan(x.value()),x.derivatives() * (Scalar(1)/numext::abs2(cos(x.value()))));)
625 
626 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(asin,
627  using std::sqrt;
628  using std::asin;
629  return ReturnType(asin(x.value()),x.derivatives() * (Scalar(1)/sqrt(1-numext::abs2(x.value()))));)
630 
631 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(acos,
632  using std::sqrt;
633  using std::acos;
634  return ReturnType(acos(x.value()),x.derivatives() * (Scalar(-1)/sqrt(1-numext::abs2(x.value()))));)
635 
636 #undef EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY
637 
638 template<typename DerType> struct NumTraits<AutoDiffScalar<DerType> >
639  : NumTraits< typename NumTraits<typename DerType::Scalar>::Real >
640 {
641  typedef AutoDiffScalar<Matrix<typename NumTraits<typename DerType::Scalar>::Real,DerType::RowsAtCompileTime,DerType::ColsAtCompileTime,
642  DerType::Options, DerType::MaxRowsAtCompileTime, DerType::MaxColsAtCompileTime> > Real;
643  typedef AutoDiffScalar<DerType> NonInteger;
644  typedef AutoDiffScalar<DerType> Nested;
645  enum{
646  RequireInitialization = 1
647  };
648 };
649 
650 }
651 
652 #endif // EIGEN_AUTODIFF_SCALAR_H
AutoDiffScalar(const Scalar &value, int nbDer, int derNumber)
Definition: AutoDiffScalar.h:81
A scalar type replacement with automatic differentation capability.
Definition: AutoDiffScalar.h:60
AutoDiffScalar()
Definition: AutoDiffScalar.h:77
Namespace containing all symbols from the Eigen library.
Definition: CXX11Meta.h:13
AutoDiffScalar(const Real &value)
Definition: AutoDiffScalar.h:89
AutoDiffScalar(const Scalar &value, const DerType &der)
Definition: AutoDiffScalar.h:97