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 (at) 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 33 template<typename _DerType> class AutoDiffScalar; 34 35 template<typename NewDerType> 36 inline AutoDiffScalar<NewDerType> MakeAutoDiffScalar(const typename NewDerType::Scalar& value, const NewDerType &der) { 37 return AutoDiffScalar<NewDerType>(value,der); 38 } 39 40 /** \class AutoDiffScalar 41 * \brief A scalar type replacement with automatic differentation capability 42 * 43 * \param _DerType the vector type used to store/represent the derivatives. The base scalar type 44 * as well as the number of derivatives to compute are determined from this type. 45 * Typical choices include, e.g., \c Vector4f for 4 derivatives, or \c VectorXf 46 * if the number of derivatives is not known at compile time, and/or, the number 47 * of derivatives is large. 48 * Note that _DerType can also be a reference (e.g., \c VectorXf&) to wrap a 49 * existing vector into an AutoDiffScalar. 50 * Finally, _DerType can also be any Eigen compatible expression. 51 * 52 * This class represents a scalar value while tracking its respective derivatives using Eigen's expression 53 * template mechanism. 54 * 55 * It supports the following list of global math function: 56 * - std::abs, std::sqrt, std::pow, std::exp, std::log, std::sin, std::cos, 57 * - internal::abs, internal::sqrt, numext::pow, internal::exp, internal::log, internal::sin, internal::cos, 58 * - internal::conj, internal::real, internal::imag, numext::abs2. 59 * 60 * AutoDiffScalar can be used as the scalar type of an Eigen::Matrix object. However, 61 * in that case, the expression template mechanism only occurs at the top Matrix level, 62 * while derivatives are computed right away. 63 * 64 */ 65 66 template<typename _DerType> 67 class AutoDiffScalar 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> 71 { 72 public: 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; 79 80 using Base::operator+; 81 using Base::operator*; 82 83 /** Default constructor without any initialization. */ 84 AutoDiffScalar() {} 85 86 /** Constructs an active scalar from its \a value, 87 and initializes the \a nbDer derivatives such that it corresponds to the \a derNumber -th variable */ 88 AutoDiffScalar(const Scalar& value, int nbDer, int derNumber) 89 : m_value(value), m_derivatives(DerType::Zero(nbDer)) 90 { 91 m_derivatives.coeffRef(derNumber) = Scalar(1); 92 } 93 94 /** Conversion from a scalar constant to an active scalar. 95 * The derivatives are set to zero. */ 96 /*explicit*/ AutoDiffScalar(const Real& value) 97 : m_value(value) 98 { 99 if(m_derivatives.size()>0) 100 m_derivatives.setZero(); 101 } 102 103 /** Constructs an active scalar from its \a value and derivatives \a der */ 104 AutoDiffScalar(const Scalar& value, const DerType& der) 105 : m_value(value), m_derivatives(der) 106 {} 107 108 template<typename OtherDerType> 109 AutoDiffScalar(const AutoDiffScalar<OtherDerType>& other 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 112 #endif 113 ) 114 : m_value(other.value()), m_derivatives(other.derivatives()) 115 {} 116 117 friend std::ostream & operator << (std::ostream & s, const AutoDiffScalar& a) 118 { 119 return s << a.value(); 120 } 121 122 AutoDiffScalar(const AutoDiffScalar& other) 123 : m_value(other.value()), m_derivatives(other.derivatives()) 124 {} 125 126 template<typename OtherDerType> 127 inline AutoDiffScalar& operator=(const AutoDiffScalar<OtherDerType>& other) 128 { 129 m_value = other.value(); 130 m_derivatives = other.derivatives(); 131 return *this; 132 } 133 134 inline AutoDiffScalar& operator=(const AutoDiffScalar& other) 135 { 136 m_value = other.value(); 137 m_derivatives = other.derivatives(); 138 return *this; 139 } 140 141 inline AutoDiffScalar& operator=(const Scalar& other) 142 { 143 m_value = other; 144 if(m_derivatives.size()>0) 145 m_derivatives.setZero(); 146 return *this; 147 } 148 149 // inline operator const Scalar& () const { return m_value; } 150 // inline operator Scalar& () { return m_value; } 151 152 inline const Scalar& value() const { return m_value; } 153 inline Scalar& value() { return m_value; } 154 155 inline const DerType& derivatives() const { return m_derivatives; } 156 inline DerType& derivatives() { return m_derivatives; } 157 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; } 164 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(); } 171 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(); } 178 179 inline const AutoDiffScalar<DerType&> operator+(const Scalar& other) const 180 { 181 return AutoDiffScalar<DerType&>(m_value + other, m_derivatives); 182 } 183 184 friend inline const AutoDiffScalar<DerType&> operator+(const Scalar& a, const AutoDiffScalar& b) 185 { 186 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 187 } 188 189 // inline const AutoDiffScalar<DerType&> operator+(const Real& other) const 190 // { 191 // return AutoDiffScalar<DerType&>(m_value + other, m_derivatives); 192 // } 193 194 // friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar& b) 195 // { 196 // return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 197 // } 198 199 inline AutoDiffScalar& operator+=(const Scalar& other) 200 { 201 value() += other; 202 return *this; 203 } 204 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 208 { 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()); 213 } 214 215 template<typename OtherDerType> 216 inline AutoDiffScalar& 217 operator+=(const AutoDiffScalar<OtherDerType>& other) 218 { 219 (*this) = (*this) + other; 220 return *this; 221 } 222 223 inline const AutoDiffScalar<DerType&> operator-(const Scalar& b) const 224 { 225 return AutoDiffScalar<DerType&>(m_value - b, m_derivatives); 226 } 227 228 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 229 operator-(const Scalar& a, const AutoDiffScalar& b) 230 { 231 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 232 (a - b.value(), -b.derivatives()); 233 } 234 235 inline AutoDiffScalar& operator-=(const Scalar& other) 236 { 237 value() -= other; 238 return *this; 239 } 240 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 244 { 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()); 249 } 250 251 template<typename OtherDerType> 252 inline AutoDiffScalar& 253 operator-=(const AutoDiffScalar<OtherDerType>& other) 254 { 255 *this = *this - other; 256 return *this; 257 } 258 259 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 260 operator-() const 261 { 262 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >( 263 -m_value, 264 -m_derivatives); 265 } 266 267 inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) > 268 operator*(const Scalar& other) const 269 { 270 return MakeAutoDiffScalar(m_value * other, m_derivatives * other); 271 } 272 273 friend inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) > 274 operator*(const Scalar& other, const AutoDiffScalar& a) 275 { 276 return MakeAutoDiffScalar(a.value() * other, 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<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) > 296 operator/(const Scalar& other) const 297 { 298 return MakeAutoDiffScalar(m_value / other, (m_derivatives * (Scalar(1)/other))); 299 } 300 301 friend inline const AutoDiffScalar<EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(DerType,Scalar,product) > 302 operator/(const Scalar& other, const AutoDiffScalar& a) 303 { 304 return MakeAutoDiffScalar(other / a.value(), a.derivatives() * (Scalar(-other) / (a.value()*a.value()))); 305 } 306 307 // inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 308 // operator/(const Real& other) const 309 // { 310 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 311 // m_value / other, 312 // (m_derivatives * (Real(1)/other))); 313 // } 314 // 315 // friend inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 316 // operator/(const Real& other, const AutoDiffScalar& a) 317 // { 318 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 319 // other / a.value(), 320 // a.derivatives() * (-Real(1)/other)); 321 // } 322 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 329 { 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()))); 335 } 336 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 342 { 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)); 347 } 348 349 inline AutoDiffScalar& operator*=(const Scalar& other) 350 { 351 *this = *this * other; 352 return *this; 353 } 354 355 template<typename OtherDerType> 356 inline AutoDiffScalar& operator*=(const AutoDiffScalar<OtherDerType>& other) 357 { 358 *this = *this * other; 359 return *this; 360 } 361 362 inline AutoDiffScalar& operator/=(const Scalar& other) 363 { 364 *this = *this / other; 365 return *this; 366 } 367 368 template<typename OtherDerType> 369 inline AutoDiffScalar& operator/=(const AutoDiffScalar<OtherDerType>& other) 370 { 371 *this = *this / other; 372 return *this; 373 } 374 375 protected: 376 Scalar m_value; 377 DerType m_derivatives; 378 379 }; 380 381 namespace internal { 382 383 template<typename _DerType> 384 struct auto_diff_special_op<_DerType, true> 385 // : auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real, 386 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value> 387 { 388 typedef typename remove_all<_DerType>::type DerType; 389 typedef typename traits<DerType>::Scalar Scalar; 390 typedef typename NumTraits<Scalar>::Real Real; 391 392 // typedef auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real, 393 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value> Base; 394 395 // using Base::operator+; 396 // using Base::operator+=; 397 // using Base::operator-; 398 // using Base::operator-=; 399 // using Base::operator*; 400 // using Base::operator*=; 401 402 const AutoDiffScalar<_DerType>& derived() const { return *static_cast<const AutoDiffScalar<_DerType>*>(this); } 403 AutoDiffScalar<_DerType>& derived() { return *static_cast<AutoDiffScalar<_DerType>*>(this); } 404 405 406 inline const AutoDiffScalar<DerType&> operator+(const Real& other) const 407 { 408 return AutoDiffScalar<DerType&>(derived().value() + other, derived().derivatives()); 409 } 410 411 friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar<_DerType>& b) 412 { 413 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 414 } 415 416 inline AutoDiffScalar<_DerType>& operator+=(const Real& other) 417 { 418 derived().value() += other; 419 return derived(); 420 } 421 422 423 inline const AutoDiffScalar<typename CwiseUnaryOp<bind2nd_op<scalar_product_op<Scalar,Real> >, DerType>::Type > 424 operator*(const Real& other) const 425 { 426 return AutoDiffScalar<typename CwiseUnaryOp<bind2nd_op<scalar_product_op<Scalar,Real> >, DerType>::Type >( 427 derived().value() * other, 428 derived().derivatives() * other); 429 } 430 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) 433 { 434 return AutoDiffScalar<typename CwiseUnaryOp<bind1st_op<scalar_product_op<Real,Scalar> >, DerType>::Type >( 435 a.value() * other, 436 a.derivatives() * other); 437 } 438 439 inline AutoDiffScalar<_DerType>& operator*=(const Scalar& other) 440 { 441 *this = *this * other; 442 return derived(); 443 } 444 }; 445 446 template<typename _DerType> 447 struct auto_diff_special_op<_DerType, false> 448 { 449 void operator*() const; 450 void operator-() const; 451 void operator+() const; 452 }; 453 454 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols, typename B> 455 struct 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)) 459 { 460 a.resize(b.size()); 461 a.setZero(); 462 } 463 } 464 }; 465 466 template<typename A, typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols> 467 struct 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)) 471 { 472 b.resize(a.size()); 473 b.setZero(); 474 } 475 } 476 }; 477 478 template<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> 480 struct 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)) 486 { 487 a.resize(b.size()); 488 a.setZero(); 489 } 490 else if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0)) 491 { 492 b.resize(a.size()); 493 b.setZero(); 494 } 495 } 496 }; 497 498 } // end namespace internal 499 500 template<typename DerType, typename BinOp> 501 struct ScalarBinaryOpTraits<AutoDiffScalar<DerType>,typename DerType::Scalar,BinOp> 502 { 503 typedef AutoDiffScalar<DerType> ReturnType; 504 }; 505 506 template<typename DerType, typename BinOp> 507 struct ScalarBinaryOpTraits<typename DerType::Scalar,AutoDiffScalar<DerType>, BinOp> 508 { 509 typedef AutoDiffScalar<DerType> ReturnType; 510 }; 511 512 513 // The following is an attempt to let Eigen's known about expression template, but that's more tricky! 514 515 // template<typename DerType, typename BinOp> 516 // struct ScalarBinaryOpTraits<AutoDiffScalar<DerType>,AutoDiffScalar<DerType>, BinOp> 517 // { 518 // enum { Defined = 1 }; 519 // typedef AutoDiffScalar<typename DerType::PlainObject> ReturnType; 520 // }; 521 // 522 // template<typename DerType1,typename DerType2, typename BinOp> 523 // struct ScalarBinaryOpTraits<AutoDiffScalar<DerType1>,AutoDiffScalar<DerType2>, BinOp> 524 // { 525 // enum { Defined = 1 };//internal::is_same<typename DerType1::Scalar,typename DerType2::Scalar>::value }; 526 // typedef AutoDiffScalar<typename DerType1::PlainObject> ReturnType; 527 // }; 528 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; \ 536 CODE; \ 537 } 538 539 template<typename DerType> 540 inline const AutoDiffScalar<DerType>& conj(const AutoDiffScalar<DerType>& x) { return x; } 541 template<typename DerType> 542 inline const AutoDiffScalar<DerType>& real(const AutoDiffScalar<DerType>& x) { return x; } 543 template<typename DerType> 544 inline typename DerType::Scalar imag(const AutoDiffScalar<DerType>&) { return 0.; } 545 template<typename DerType, typename T> 546 inline 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)); 549 } 550 template<typename DerType, typename T> 551 inline 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)); 554 } 555 template<typename DerType, typename T> 556 inline 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)); 559 } 560 template<typename DerType, typename T> 561 inline 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)); 564 } 565 template<typename DerType> 566 inline 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); 568 } 569 template<typename DerType> 570 inline 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); 572 } 573 574 575 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs, 576 using std::abs; 577 return Eigen::MakeAutoDiffScalar(abs(x.value()), x.derivatives() * (x.value()<0 ? -1 : 1) );) 578 579 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs2, 580 using numext::abs2; 581 return Eigen::MakeAutoDiffScalar(abs2(x.value()), x.derivatives() * (Scalar(2)*x.value()));) 582 583 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sqrt, 584 using std::sqrt; 585 Scalar sqrtx = sqrt(x.value()); 586 return Eigen::MakeAutoDiffScalar(sqrtx,x.derivatives() * (Scalar(0.5) / sqrtx));) 587 588 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cos, 589 using std::cos; 590 using std::sin; 591 return Eigen::MakeAutoDiffScalar(cos(x.value()), x.derivatives() * (-sin(x.value())));) 592 593 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sin, 594 using std::sin; 595 using std::cos; 596 return Eigen::MakeAutoDiffScalar(sin(x.value()),x.derivatives() * cos(x.value()));) 597 598 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(exp, 599 using std::exp; 600 Scalar expx = exp(x.value()); 601 return Eigen::MakeAutoDiffScalar(expx,x.derivatives() * expx);) 602 603 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(log, 604 using std::log; 605 return Eigen::MakeAutoDiffScalar(log(x.value()),x.derivatives() * (Scalar(1)/x.value()));) 606 607 template<typename DerType> 608 inline const Eigen::AutoDiffScalar< 609 EIGEN_EXPR_BINARYOP_SCALAR_RETURN_TYPE(typename internal::remove_all<DerType>::type,typename internal::traits<typename internal::remove_all<DerType>::type>::Scalar,product) > 610 pow(const Eigen::AutoDiffScalar<DerType> &x, const typename internal::traits<typename internal::remove_all<DerType>::type>::Scalar &y) 611 { 612 using namespace Eigen; 613 using std::pow; 614 return Eigen::MakeAutoDiffScalar(pow(x.value(),y), x.derivatives() * (y * pow(x.value(),y-1))); 615 } 616 617 618 template<typename DerTypeA,typename DerTypeB> 619 inline const AutoDiffScalar<Matrix<typename internal::traits<typename internal::remove_all<DerTypeA>::type>::Scalar,Dynamic,1> > 620 atan2(const AutoDiffScalar<DerTypeA>& a, const AutoDiffScalar<DerTypeB>& b) 621 { 622 using std::atan2; 623 typedef typename internal::traits<typename internal::remove_all<DerTypeA>::type>::Scalar Scalar; 624 typedef AutoDiffScalar<Matrix<Scalar,Dynamic,1> > PlainADS; 625 PlainADS ret; 626 ret.value() = atan2(a.value(), b.value()); 627 628 Scalar squared_hypot = a.value() * a.value() + b.value() * b.value(); 629 630 // if (squared_hypot==0) the derivation is undefined and the following results in a NaN: 631 ret.derivatives() = (a.derivatives() * b.value() - a.value() * b.derivatives()) / squared_hypot; 632 633 return ret; 634 } 635 636 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tan, 637 using std::tan; 638 using std::cos; 639 return Eigen::MakeAutoDiffScalar(tan(x.value()),x.derivatives() * (Scalar(1)/numext::abs2(cos(x.value()))));) 640 641 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(asin, 642 using std::sqrt; 643 using std::asin; 644 return Eigen::MakeAutoDiffScalar(asin(x.value()),x.derivatives() * (Scalar(1)/sqrt(1-numext::abs2(x.value()))));) 645 646 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(acos, 647 using std::sqrt; 648 using std::acos; 649 return Eigen::MakeAutoDiffScalar(acos(x.value()),x.derivatives() * (Scalar(-1)/sqrt(1-numext::abs2(x.value()))));) 650 651 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tanh, 652 using std::cosh; 653 using std::tanh; 654 return Eigen::MakeAutoDiffScalar(tanh(x.value()),x.derivatives() * (Scalar(1)/numext::abs2(cosh(x.value()))));) 655 656 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sinh, 657 using std::sinh; 658 using std::cosh; 659 return Eigen::MakeAutoDiffScalar(sinh(x.value()),x.derivatives() * cosh(x.value()));) 660 661 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cosh, 662 using std::sinh; 663 using std::cosh; 664 return Eigen::MakeAutoDiffScalar(cosh(x.value()),x.derivatives() * sinh(x.value()));) 665 666 #undef EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY 667 668 template<typename DerType> struct NumTraits<AutoDiffScalar<DerType> > 669 : NumTraits< typename NumTraits<typename internal::remove_all<DerType>::type::Scalar>::Real > 670 { 671 typedef typename internal::remove_all<DerType>::type DerTypeCleaned; 672 typedef AutoDiffScalar<Matrix<typename NumTraits<typename DerTypeCleaned::Scalar>::Real,DerTypeCleaned::RowsAtCompileTime,DerTypeCleaned::ColsAtCompileTime, 673 0, DerTypeCleaned::MaxRowsAtCompileTime, DerTypeCleaned::MaxColsAtCompileTime> > Real; 674 typedef AutoDiffScalar<DerType> NonInteger; 675 typedef AutoDiffScalar<DerType> Nested; 676 typedef typename NumTraits<typename DerTypeCleaned::Scalar>::Literal Literal; 677 enum{ 678 RequireInitialization = 1 679 }; 680 }; 681 682 } 683 684 #endif // EIGEN_AUTODIFF_SCALAR_H 685