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