10#ifndef EIGEN_AUTODIFF_SCALAR_H
11#define EIGEN_AUTODIFF_SCALAR_H
17template<
typename A,
typename B>
18struct make_coherent_impl {
19 static void run(A&, B&) {}
23template<
typename A,
typename B>
24void make_coherent(
const A& a,
const B&b)
26 make_coherent_impl<A,B>::run(a.const_cast_derived(), b.const_cast_derived());
29template<
typename _DerType,
bool Enable>
struct auto_diff_special_op;
33template<
typename _DerType>
class AutoDiffScalar;
35template<
typename NewDerType>
36inline AutoDiffScalar<NewDerType> MakeAutoDiffScalar(
const typename NewDerType::Scalar& value,
const NewDerType &der) {
37 return AutoDiffScalar<NewDerType>(value,der);
66template<
typename _DerType>
68 :
public internal::auto_diff_special_op
69 <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
70 typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value>
73 typedef internal::auto_diff_special_op
74 <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
75 typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value> Base;
76 typedef typename internal::remove_all<_DerType>::type DerType;
77 typedef typename internal::traits<DerType>::Scalar Scalar;
78 typedef typename NumTraits<Scalar>::Real Real;
80 using Base::operator+;
81 using Base::operator*;
89 : m_value(value), m_derivatives(DerType::Zero(nbDer))
91 m_derivatives.coeffRef(derNumber) = Scalar(1);
99 if(m_derivatives.size()>0)
100 m_derivatives.setZero();
105 : m_value(value), m_derivatives(der)
108 template<
typename OtherDerType>
110#ifndef EIGEN_PARSED_BY_DOXYGEN
111 ,
typename internal::enable_if<internal::is_same<Scalar,
typename internal::traits<
typename internal::remove_all<OtherDerType>::type>::Scalar>::value,
void*>::type = 0
114 : m_value(other.value()), m_derivatives(other.derivatives())
117 friend std::ostream & operator << (std::ostream & s,
const AutoDiffScalar& a)
119 return s << a.value();
123 : m_value(other.value()), m_derivatives(other.derivatives())
126 template<
typename OtherDerType>
127 inline AutoDiffScalar& operator=(
const AutoDiffScalar<OtherDerType>& other)
129 m_value = other.value();
130 m_derivatives = other.derivatives();
136 m_value = other.value();
137 m_derivatives = other.derivatives();
144 if(m_derivatives.size()>0)
145 m_derivatives.setZero();
152 inline const Scalar& value()
const {
return m_value; }
153 inline Scalar& value() {
return m_value; }
155 inline const DerType& derivatives()
const {
return m_derivatives; }
156 inline DerType& derivatives() {
return m_derivatives; }
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 inline bool operator> (
const Scalar& other)
const {
return m_value > other; }
161 inline bool operator>=(
const Scalar& other)
const {
return m_value >= other; }
162 inline bool operator==(
const Scalar& other)
const {
return m_value == other; }
163 inline bool operator!=(
const Scalar& other)
const {
return m_value != other; }
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 friend inline bool operator> (
const Scalar& a,
const AutoDiffScalar& b) {
return a > b.value(); }
168 friend inline bool operator>=(
const Scalar& a,
const AutoDiffScalar& b) {
return a >= b.value(); }
169 friend inline bool operator==(
const Scalar& a,
const AutoDiffScalar& b) {
return a == b.value(); }
170 friend inline bool operator!=(
const Scalar& a,
const AutoDiffScalar& b) {
return a != 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 template<
typename OtherDerType>
inline bool operator> (
const AutoDiffScalar<OtherDerType>& b)
const {
return m_value > b.value(); }
175 template<
typename OtherDerType>
inline bool operator>=(
const AutoDiffScalar<OtherDerType>& b)
const {
return m_value >= b.value(); }
176 template<
typename OtherDerType>
inline bool operator==(
const AutoDiffScalar<OtherDerType>& b)
const {
return m_value == b.value(); }
177 template<
typename OtherDerType>
inline bool operator!=(
const AutoDiffScalar<OtherDerType>& b)
const {
return m_value != b.value(); }
179 inline const AutoDiffScalar<DerType&> operator+(
const Scalar& other)
const
181 return AutoDiffScalar<DerType&>(m_value + other, m_derivatives);
184 friend inline const AutoDiffScalar<DerType&> operator+(
const Scalar& a,
const AutoDiffScalar& b)
186 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
205 template<
typename OtherDerType>
206 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
const DerType,
const typename internal::remove_all<OtherDerType>::type> >
207 operator+(
const AutoDiffScalar<OtherDerType>& other)
const
209 internal::make_coherent(m_derivatives, other.derivatives());
210 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
const DerType,
const typename internal::remove_all<OtherDerType>::type> >(
211 m_value + other.value(),
212 m_derivatives + other.derivatives());
215 template<
typename OtherDerType>
217 operator+=(
const AutoDiffScalar<OtherDerType>& other)
219 (*this) = (*this) + other;
223 inline const AutoDiffScalar<DerType&> operator-(
const Scalar& b)
const
225 return AutoDiffScalar<DerType&>(m_value - b, m_derivatives);
228 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>,
const DerType> >
231 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>,
const DerType> >
232 (a - b.value(), -b.derivatives());
241 template<
typename OtherDerType>
242 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
const DerType,
const typename internal::remove_all<OtherDerType>::type> >
243 operator-(
const AutoDiffScalar<OtherDerType>& other)
const
245 internal::make_coherent(m_derivatives, other.derivatives());
246 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
const DerType,
const typename internal::remove_all<OtherDerType>::type> >(
247 m_value - other.value(),
248 m_derivatives - other.derivatives());
251 template<
typename OtherDerType>
253 operator-=(
const AutoDiffScalar<OtherDerType>& other)
255 *
this = *
this - other;
259 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>,
const DerType> >
262 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>,
const DerType> >(
267 inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) >
268 operator*(
const Scalar& other)
const
270 return MakeAutoDiffScalar(m_value * other, m_derivatives * other);
273 friend inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) >
276 return MakeAutoDiffScalar(a.value() * other, a.derivatives() * other);
295 inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) >
296 operator/(
const Scalar& other)
const
298 return MakeAutoDiffScalar(m_value / other, (m_derivatives * (Scalar(1)/other)));
301 friend inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) >
304 return MakeAutoDiffScalar(other / a.value(), a.derivatives() * (Scalar(-other) / (a.value()*a.value())));
323 template<
typename OtherDerType>
324 inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(
325 CwiseBinaryOp<internal::scalar_difference_op<Scalar> EIGEN_COMMA
326 const EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) EIGEN_COMMA
327 const EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(
typename internal::remove_all<OtherDerType>::type,Scalar,product) >,Scalar,product) >
328 operator/(
const AutoDiffScalar<OtherDerType>& other)
const
330 internal::make_coherent(m_derivatives, other.derivatives());
331 return MakeAutoDiffScalar(
332 m_value / other.value(),
333 ((m_derivatives * other.value()) - (other.derivatives() * m_value))
334 * (Scalar(1)/(other.value()*other.value())));
337 template<
typename OtherDerType>
338 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
339 const EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product),
340 const EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(
typename internal::remove_all<OtherDerType>::type,Scalar,product) > >
341 operator*(
const AutoDiffScalar<OtherDerType>& other)
const
343 internal::make_coherent(m_derivatives, other.derivatives());
344 return MakeAutoDiffScalar(
345 m_value * other.value(),
346 (m_derivatives * other.value()) + (other.derivatives() * m_value));
351 *
this = *
this * other;
355 template<
typename OtherDerType>
356 inline AutoDiffScalar& operator*=(
const AutoDiffScalar<OtherDerType>& other)
358 *
this = *
this * other;
364 *
this = *
this / other;
368 template<
typename OtherDerType>
369 inline AutoDiffScalar& operator/=(
const AutoDiffScalar<OtherDerType>& other)
371 *
this = *
this / other;
377 DerType m_derivatives;
383template<
typename _DerType>
384struct auto_diff_special_op<_DerType, true>
388 typedef typename remove_all<_DerType>::type DerType;
389 typedef typename traits<DerType>::Scalar Scalar;
390 typedef typename NumTraits<Scalar>::Real Real;
402 const AutoDiffScalar<_DerType>& derived()
const {
return *
static_cast<const AutoDiffScalar<_DerType>*
>(
this); }
403 AutoDiffScalar<_DerType>& derived() {
return *
static_cast<AutoDiffScalar<_DerType>*
>(
this); }
406 inline const AutoDiffScalar<DerType&> operator+(
const Real& other)
const
408 return AutoDiffScalar<DerType&>(derived().value() + other, derived().derivatives());
411 friend inline const AutoDiffScalar<DerType&> operator+(
const Real& a,
const AutoDiffScalar<_DerType>& b)
413 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
416 inline AutoDiffScalar<_DerType>& operator+=(
const Real& other)
418 derived().value() += other;
423 inline const AutoDiffScalar<typename CwiseUnaryOp<bind2nd_op<scalar_product_op<Scalar,Real> >, DerType>::Type >
424 operator*(
const Real& other)
const
426 return AutoDiffScalar<typename CwiseUnaryOp<bind2nd_op<scalar_product_op<Scalar,Real> >, DerType>::Type >(
427 derived().value() * other,
428 derived().derivatives() * other);
431 friend inline const AutoDiffScalar<typename CwiseUnaryOp<bind1st_op<scalar_product_op<Real,Scalar> >, DerType>::Type >
432 operator*(
const Real& other,
const AutoDiffScalar<_DerType>& a)
434 return AutoDiffScalar<typename CwiseUnaryOp<bind1st_op<scalar_product_op<Real,Scalar> >, DerType>::Type >(
436 a.derivatives() * other);
439 inline AutoDiffScalar<_DerType>& operator*=(
const Scalar& other)
441 *
this = *
this * other;
446template<
typename _DerType>
447struct auto_diff_special_op<_DerType, false>
449 void operator*()
const;
450 void operator-()
const;
451 void operator+()
const;
454template<
typename A_Scalar,
int A_Rows,
int A_Cols,
int A_Options,
int A_MaxRows,
int A_MaxCols,
typename B>
455struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>, B> {
456 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
457 static void run(A& a, B& b) {
458 if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
466template<
typename A,
typename B_Scalar,
int B_Rows,
int B_Cols,
int B_Options,
int B_MaxRows,
int B_MaxCols>
467struct make_coherent_impl<A, Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
468 typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
469 static void run(A& a, B& b) {
470 if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
478template<
typename A_Scalar,
int A_Rows,
int A_Cols,
int A_Options,
int A_MaxRows,
int A_MaxCols,
479 typename B_Scalar,
int B_Rows,
int B_Cols,
int B_Options,
int B_MaxRows,
int B_MaxCols>
480struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>,
481 Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
482 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
483 typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
484 static void run(A& a, B& b) {
485 if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
490 else if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
500template<
typename DerType,
typename BinOp>
501struct ScalarBinaryOpTraits<AutoDiffScalar<DerType>,typename DerType::Scalar,BinOp>
503 typedef AutoDiffScalar<DerType> ReturnType;
506template<
typename DerType,
typename BinOp>
507struct ScalarBinaryOpTraits<typename DerType::Scalar,AutoDiffScalar<DerType>, BinOp>
509 typedef AutoDiffScalar<DerType> ReturnType;
529#define EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(FUNC,CODE) \
530 template<typename DerType> \
531 inline const Eigen::AutoDiffScalar< \
532 EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(typename Eigen::internal::remove_all<DerType>::type, typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar, product) > \
533 FUNC(const Eigen::AutoDiffScalar<DerType>& x) { \
534 using namespace Eigen; \
535 EIGEN_UNUSED typedef typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar Scalar; \
539template<
typename DerType>
540inline const AutoDiffScalar<DerType>& conj(
const AutoDiffScalar<DerType>& x) {
return x; }
541template<
typename DerType>
542inline const AutoDiffScalar<DerType>& real(
const AutoDiffScalar<DerType>& x) {
return x; }
543template<
typename DerType>
544inline typename DerType::Scalar imag(
const AutoDiffScalar<DerType>&) {
return 0.; }
545template<
typename DerType,
typename T>
546inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (min)(
const AutoDiffScalar<DerType>& x,
const T& y) {
547 typedef AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> ADS;
548 return (x <= y ? ADS(x) : ADS(y));
550template<
typename DerType,
typename T>
551inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (max)(
const AutoDiffScalar<DerType>& x,
const T& y) {
552 typedef AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> ADS;
553 return (x >= y ? ADS(x) : ADS(y));
555template<
typename DerType,
typename T>
556inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (min)(
const T& x,
const AutoDiffScalar<DerType>& y) {
557 typedef AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> ADS;
558 return (x < y ? ADS(x) : ADS(y));
560template<
typename DerType,
typename T>
561inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (max)(
const T& x,
const AutoDiffScalar<DerType>& y) {
562 typedef AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> ADS;
563 return (x > y ? ADS(x) : ADS(y));
565template<
typename DerType>
566inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (min)(
const AutoDiffScalar<DerType>& x,
const AutoDiffScalar<DerType>& y) {
567 return (x.value() < y.value() ? x : y);
569template<
typename DerType>
570inline AutoDiffScalar<typename Eigen::internal::remove_all<DerType>::type::PlainObject> (max)(
const AutoDiffScalar<DerType>& x,
const AutoDiffScalar<DerType>& y) {
571 return (x.value() >= y.value() ? x : y);
575EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs,
577 return Eigen::MakeAutoDiffScalar(abs(x.value()), x.derivatives() * (x.value()<0 ? -1 : 1) );)
579EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs2,
581 return Eigen::MakeAutoDiffScalar(abs2(x.value()), x.derivatives() * (Scalar(2)*x.value()));)
583EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sqrt,
585 Scalar sqrtx = sqrt(x.value());
586 return Eigen::MakeAutoDiffScalar(sqrtx,x.derivatives() * (Scalar(0.5) / sqrtx));)
588EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cos,
591 return Eigen::MakeAutoDiffScalar(cos(x.value()), x.derivatives() * (-sin(x.value())));)
593EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sin,
596 return Eigen::MakeAutoDiffScalar(sin(x.value()),x.derivatives() * cos(x.value()));)
598EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(exp,
600 Scalar expx = exp(x.value());
601 return Eigen::MakeAutoDiffScalar(expx,x.derivatives() * expx);)
603EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(log,
605 return Eigen::MakeAutoDiffScalar(log(x.value()),x.derivatives() * (Scalar(1)/x.value()));)
607template<typename DerType>
609EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(
typename internal::remove_all<DerType>::type,
typename internal::traits<
typename internal::remove_all<DerType>::type>::Scalar,product) >
612 using namespace Eigen;
614 return Eigen::MakeAutoDiffScalar(pow(x.value(),y), x.derivatives() * (y * pow(x.value(),y-1)));
618template<
typename DerTypeA,
typename DerTypeB>
623 typedef typename internal::traits<typename internal::remove_all<DerTypeA>::type>::Scalar Scalar;
626 ret.value() = atan2(a.value(), b.value());
628 Scalar squared_hypot = a.value() * a.value() + b.value() * b.value();
631 ret.derivatives() = (a.derivatives() * b.value() - a.value() * b.derivatives()) / squared_hypot;
636EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tan,
639 return Eigen::MakeAutoDiffScalar(tan(x.value()),x.derivatives() * (Scalar(1)/numext::abs2(cos(x.value()))));)
641EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(asin,
644 return Eigen::MakeAutoDiffScalar(asin(x.value()),x.derivatives() * (Scalar(1)/sqrt(1-numext::abs2(x.value()))));)
646EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(acos,
649 return Eigen::MakeAutoDiffScalar(acos(x.value()),x.derivatives() * (Scalar(-1)/sqrt(1-numext::abs2(x.value()))));)
651EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tanh,
654 return Eigen::MakeAutoDiffScalar(tanh(x.value()),x.derivatives() * (Scalar(1)/numext::abs2(cosh(x.value()))));)
656EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sinh,
659 return Eigen::MakeAutoDiffScalar(sinh(x.value()),x.derivatives() * cosh(x.value()));)
661EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cosh,
664 return Eigen::MakeAutoDiffScalar(cosh(x.value()),x.derivatives() * sinh(x.value()));)
666#undef EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY
669 : NumTraits<
typename NumTraits<
typename internal::remove_all<DerType>::type::Scalar>::Real >
671 typedef typename internal::remove_all<DerType>::type DerTypeCleaned;
673 0, DerTypeCleaned::MaxRowsAtCompileTime, DerTypeCleaned::MaxColsAtCompileTime> > Real;
676 typedef typename NumTraits<typename DerTypeCleaned::Scalar>::Literal Literal;
678 RequireInitialization = 1
A scalar type replacement with automatic differentation capability.
Definition: AutoDiffScalar.h:71
AutoDiffScalar(const Scalar &value, const DerType &der)
Definition: AutoDiffScalar.h:104
AutoDiffScalar(const Real &value)
Definition: AutoDiffScalar.h:96
AutoDiffScalar(const Scalar &value, int nbDer, int derNumber)
Definition: AutoDiffScalar.h:88
AutoDiffScalar()
Definition: AutoDiffScalar.h:84
Namespace containing all symbols from the Eigen library.
Definition: AdolcForward:45