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 /** \class AutoDiffScalar 34 * \brief A scalar type replacement with automatic differentation capability 35 * 36 * \param _DerType the vector type used to store/represent the derivatives. The base scalar type 37 * as well as the number of derivatives to compute are determined from this type. 38 * Typical choices include, e.g., \c Vector4f for 4 derivatives, or \c VectorXf 39 * if the number of derivatives is not known at compile time, and/or, the number 40 * of derivatives is large. 41 * Note that _DerType can also be a reference (e.g., \c VectorXf&) to wrap a 42 * existing vector into an AutoDiffScalar. 43 * Finally, _DerType can also be any Eigen compatible expression. 44 * 45 * This class represents a scalar value while tracking its respective derivatives using Eigen's expression 46 * template mechanism. 47 * 48 * It supports the following list of global math function: 49 * - std::abs, std::sqrt, std::pow, std::exp, std::log, std::sin, std::cos, 50 * - internal::abs, internal::sqrt, internal::pow, internal::exp, internal::log, internal::sin, internal::cos, 51 * - internal::conj, internal::real, internal::imag, internal::abs2. 52 * 53 * AutoDiffScalar can be used as the scalar type of an Eigen::Matrix object. However, 54 * in that case, the expression template mechanism only occurs at the top Matrix level, 55 * while derivatives are computed right away. 56 * 57 */ 58 59 template<typename _DerType> 60 class AutoDiffScalar 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 76 /** Default constructor without any initialization. */ 77 AutoDiffScalar() {} 78 79 /** Constructs an active scalar from its \a value, 80 and initializes the \a nbDer derivatives such that it corresponds to the \a derNumber -th variable */ 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 87 /** Conversion from a scalar constant to an active scalar. 88 * The derivatives are set to zero. */ 89 /*explicit*/ AutoDiffScalar(const Real& value) 90 : m_value(value) 91 { 92 if(m_derivatives.size()>0) 93 m_derivatives.setZero(); 94 } 95 96 /** Constructs an active scalar from its \a value and derivatives \a der */ 97 AutoDiffScalar(const Scalar& value, const DerType& der) 98 : m_value(value), m_derivatives(der) 99 {} 100 101 template<typename OtherDerType> 102 AutoDiffScalar(const AutoDiffScalar<OtherDerType>& other) 103 : m_value(other.value()), m_derivatives(other.derivatives()) 104 {} 105 106 friend std::ostream & operator << (std::ostream & s, const AutoDiffScalar& a) 107 { 108 return s << a.value(); 109 } 110 111 AutoDiffScalar(const AutoDiffScalar& other) 112 : m_value(other.value()), m_derivatives(other.derivatives()) 113 {} 114 115 template<typename OtherDerType> 116 inline AutoDiffScalar& operator=(const AutoDiffScalar<OtherDerType>& other) 117 { 118 m_value = other.value(); 119 m_derivatives = other.derivatives(); 120 return *this; 121 } 122 123 inline AutoDiffScalar& operator=(const AutoDiffScalar& other) 124 { 125 m_value = other.value(); 126 m_derivatives = other.derivatives(); 127 return *this; 128 } 129 130 // inline operator const Scalar& () const { return m_value; } 131 // inline operator Scalar& () { return m_value; } 132 133 inline const Scalar& value() const { return m_value; } 134 inline Scalar& value() { return m_value; } 135 136 inline const DerType& derivatives() const { return m_derivatives; } 137 inline DerType& derivatives() { return m_derivatives; } 138 139 inline bool operator< (const Scalar& other) const { return m_value < other; } 140 inline bool operator<=(const Scalar& other) const { return m_value <= other; } 141 inline bool operator> (const Scalar& other) const { return m_value > other; } 142 inline bool operator>=(const Scalar& other) const { return m_value >= other; } 143 inline bool operator==(const Scalar& other) const { return m_value == other; } 144 inline bool operator!=(const Scalar& other) const { return m_value != other; } 145 146 friend inline bool operator< (const Scalar& a, const AutoDiffScalar& b) { return a < b.value(); } 147 friend inline bool operator<=(const Scalar& a, const AutoDiffScalar& b) { return a <= b.value(); } 148 friend inline bool operator> (const Scalar& a, const AutoDiffScalar& b) { return a > b.value(); } 149 friend inline bool operator>=(const Scalar& a, const AutoDiffScalar& b) { return a >= b.value(); } 150 friend inline bool operator==(const Scalar& a, const AutoDiffScalar& b) { return a == b.value(); } 151 friend inline bool operator!=(const Scalar& a, const AutoDiffScalar& b) { return a != b.value(); } 152 153 template<typename OtherDerType> inline bool operator< (const AutoDiffScalar<OtherDerType>& b) const { return m_value < b.value(); } 154 template<typename OtherDerType> inline bool operator<=(const AutoDiffScalar<OtherDerType>& b) const { return m_value <= b.value(); } 155 template<typename OtherDerType> inline bool operator> (const AutoDiffScalar<OtherDerType>& b) const { return m_value > b.value(); } 156 template<typename OtherDerType> inline bool operator>=(const AutoDiffScalar<OtherDerType>& b) const { return m_value >= b.value(); } 157 template<typename OtherDerType> inline bool operator==(const AutoDiffScalar<OtherDerType>& b) const { return m_value == b.value(); } 158 template<typename OtherDerType> inline bool operator!=(const AutoDiffScalar<OtherDerType>& b) const { return m_value != b.value(); } 159 160 inline const AutoDiffScalar<DerType&> operator+(const Scalar& other) const 161 { 162 return AutoDiffScalar<DerType&>(m_value + other, m_derivatives); 163 } 164 165 friend inline const AutoDiffScalar<DerType&> operator+(const Scalar& a, const AutoDiffScalar& b) 166 { 167 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 168 } 169 170 // inline const AutoDiffScalar<DerType&> operator+(const Real& other) const 171 // { 172 // return AutoDiffScalar<DerType&>(m_value + other, m_derivatives); 173 // } 174 175 // friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar& b) 176 // { 177 // return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 178 // } 179 180 inline AutoDiffScalar& operator+=(const Scalar& other) 181 { 182 value() += other; 183 return *this; 184 } 185 186 template<typename OtherDerType> 187 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> > 188 operator+(const AutoDiffScalar<OtherDerType>& other) const 189 { 190 internal::make_coherent(m_derivatives, other.derivatives()); 191 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> >( 192 m_value + other.value(), 193 m_derivatives + other.derivatives()); 194 } 195 196 template<typename OtherDerType> 197 inline AutoDiffScalar& 198 operator+=(const AutoDiffScalar<OtherDerType>& other) 199 { 200 (*this) = (*this) + other; 201 return *this; 202 } 203 204 inline const AutoDiffScalar<DerType&> operator-(const Scalar& b) const 205 { 206 return AutoDiffScalar<DerType&>(m_value - b, m_derivatives); 207 } 208 209 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 210 operator-(const Scalar& a, const AutoDiffScalar& b) 211 { 212 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 213 (a - b.value(), -b.derivatives()); 214 } 215 216 inline AutoDiffScalar& operator-=(const Scalar& other) 217 { 218 value() -= other; 219 return *this; 220 } 221 222 template<typename OtherDerType> 223 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> > 224 operator-(const AutoDiffScalar<OtherDerType>& other) const 225 { 226 internal::make_coherent(m_derivatives, other.derivatives()); 227 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> >( 228 m_value - other.value(), 229 m_derivatives - other.derivatives()); 230 } 231 232 template<typename OtherDerType> 233 inline AutoDiffScalar& 234 operator-=(const AutoDiffScalar<OtherDerType>& other) 235 { 236 *this = *this - other; 237 return *this; 238 } 239 240 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> > 241 operator-() const 242 { 243 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >( 244 -m_value, 245 -m_derivatives); 246 } 247 248 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> > 249 operator*(const Scalar& other) const 250 { 251 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >( 252 m_value * other, 253 (m_derivatives * other)); 254 } 255 256 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> > 257 operator*(const Scalar& other, const AutoDiffScalar& a) 258 { 259 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >( 260 a.value() * other, 261 a.derivatives() * other); 262 } 263 264 // inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 265 // operator*(const Real& other) const 266 // { 267 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 268 // m_value * other, 269 // (m_derivatives * other)); 270 // } 271 // 272 // friend inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 273 // operator*(const Real& other, const AutoDiffScalar& a) 274 // { 275 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 276 // a.value() * other, 277 // a.derivatives() * other); 278 // } 279 280 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> > 281 operator/(const Scalar& other) const 282 { 283 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >( 284 m_value / other, 285 (m_derivatives * (Scalar(1)/other))); 286 } 287 288 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> > 289 operator/(const Scalar& other, const AutoDiffScalar& a) 290 { 291 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >( 292 other / a.value(), 293 a.derivatives() * (Scalar(-other) / (a.value()*a.value()))); 294 } 295 296 // inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 297 // operator/(const Real& other) const 298 // { 299 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 300 // m_value / other, 301 // (m_derivatives * (Real(1)/other))); 302 // } 303 // 304 // friend inline const AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type > 305 // operator/(const Real& other, const AutoDiffScalar& a) 306 // { 307 // return AutoDiffScalar<typename CwiseUnaryOp<internal::scalar_multiple_op<Real>, DerType>::Type >( 308 // other / a.value(), 309 // a.derivatives() * (-Real(1)/other)); 310 // } 311 312 template<typename OtherDerType> 313 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, 314 const CwiseBinaryOp<internal::scalar_difference_op<Scalar>, 315 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>, 316 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > > 317 operator/(const AutoDiffScalar<OtherDerType>& other) const 318 { 319 internal::make_coherent(m_derivatives, other.derivatives()); 320 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, 321 const CwiseBinaryOp<internal::scalar_difference_op<Scalar>, 322 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>, 323 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > >( 324 m_value / other.value(), 325 ((m_derivatives * other.value()) - (m_value * other.derivatives())) 326 * (Scalar(1)/(other.value()*other.value()))); 327 } 328 329 template<typename OtherDerType> 330 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>, 331 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>, 332 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type> > > 333 operator*(const AutoDiffScalar<OtherDerType>& other) const 334 { 335 internal::make_coherent(m_derivatives, other.derivatives()); 336 return AutoDiffScalar<const CwiseBinaryOp<internal::scalar_sum_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 } 342 343 inline AutoDiffScalar& operator*=(const Scalar& other) 344 { 345 *this = *this * other; 346 return *this; 347 } 348 349 template<typename OtherDerType> 350 inline AutoDiffScalar& operator*=(const AutoDiffScalar<OtherDerType>& other) 351 { 352 *this = *this * other; 353 return *this; 354 } 355 356 inline AutoDiffScalar& operator/=(const Scalar& other) 357 { 358 *this = *this / other; 359 return *this; 360 } 361 362 template<typename OtherDerType> 363 inline AutoDiffScalar& operator/=(const AutoDiffScalar<OtherDerType>& other) 364 { 365 *this = *this / other; 366 return *this; 367 } 368 369 protected: 370 Scalar m_value; 371 DerType m_derivatives; 372 373 }; 374 375 namespace internal { 376 377 template<typename _DerType> 378 struct auto_diff_special_op<_DerType, true> 379 // : auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real, 380 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value> 381 { 382 typedef typename remove_all<_DerType>::type DerType; 383 typedef typename traits<DerType>::Scalar Scalar; 384 typedef typename NumTraits<Scalar>::Real Real; 385 386 // typedef auto_diff_scalar_op<_DerType, typename NumTraits<Scalar>::Real, 387 // is_same<Scalar,typename NumTraits<Scalar>::Real>::value> Base; 388 389 // using Base::operator+; 390 // using Base::operator+=; 391 // using Base::operator-; 392 // using Base::operator-=; 393 // using Base::operator*; 394 // using Base::operator*=; 395 396 const AutoDiffScalar<_DerType>& derived() const { return *static_cast<const AutoDiffScalar<_DerType>*>(this); } 397 AutoDiffScalar<_DerType>& derived() { return *static_cast<AutoDiffScalar<_DerType>*>(this); } 398 399 400 inline const AutoDiffScalar<DerType&> operator+(const Real& other) const 401 { 402 return AutoDiffScalar<DerType&>(derived().value() + other, derived().derivatives()); 403 } 404 405 friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar<_DerType>& b) 406 { 407 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives()); 408 } 409 410 inline AutoDiffScalar<_DerType>& operator+=(const Real& other) 411 { 412 derived().value() += other; 413 return derived(); 414 } 415 416 417 inline const AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type > 418 operator*(const Real& other) const 419 { 420 return AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >( 421 derived().value() * other, 422 derived().derivatives() * other); 423 } 424 425 friend inline const AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type > 426 operator*(const Real& other, const AutoDiffScalar<_DerType>& a) 427 { 428 return AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >( 429 a.value() * other, 430 a.derivatives() * other); 431 } 432 433 inline AutoDiffScalar<_DerType>& operator*=(const Scalar& other) 434 { 435 *this = *this * other; 436 return derived(); 437 } 438 }; 439 440 template<typename _DerType> 441 struct auto_diff_special_op<_DerType, false> 442 { 443 void operator*() const; 444 void operator-() const; 445 void operator+() const; 446 }; 447 448 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols, typename B> 449 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>, B> { 450 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A; 451 static void run(A& a, B& b) { 452 if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0)) 453 { 454 a.resize(b.size()); 455 a.setZero(); 456 } 457 } 458 }; 459 460 template<typename A, typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols> 461 struct make_coherent_impl<A, Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > { 462 typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B; 463 static void run(A& a, B& b) { 464 if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0)) 465 { 466 b.resize(a.size()); 467 b.setZero(); 468 } 469 } 470 }; 471 472 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols, 473 typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols> 474 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>, 475 Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > { 476 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A; 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((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0)) 480 { 481 a.resize(b.size()); 482 a.setZero(); 483 } 484 else if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0)) 485 { 486 b.resize(a.size()); 487 b.setZero(); 488 } 489 } 490 }; 491 492 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> 493 { 494 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType; 495 }; 496 497 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> > 498 { 499 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType; 500 }; 501 502 template<typename DerType> 503 struct scalar_product_traits<AutoDiffScalar<DerType>,typename DerType::Scalar> 504 { 505 typedef AutoDiffScalar<DerType> ReturnType; 506 }; 507 508 } // end namespace internal 509 510 #define EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(FUNC,CODE) \ 511 template<typename DerType> \ 512 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> > \ 513 FUNC(const Eigen::AutoDiffScalar<DerType>& x) { \ 514 using namespace Eigen; \ 515 typedef typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar Scalar; \ 516 typedef AutoDiffScalar<CwiseUnaryOp<Eigen::internal::scalar_multiple_op<Scalar>, const typename Eigen::internal::remove_all<DerType>::type> > ReturnType; \ 517 CODE; \ 518 } 519 520 template<typename DerType> 521 inline const AutoDiffScalar<DerType>& conj(const AutoDiffScalar<DerType>& x) { return x; } 522 template<typename DerType> 523 inline const AutoDiffScalar<DerType>& real(const AutoDiffScalar<DerType>& x) { return x; } 524 template<typename DerType> 525 inline typename DerType::Scalar imag(const AutoDiffScalar<DerType>&) { return 0.; } 526 template<typename DerType, typename T> 527 inline AutoDiffScalar<DerType> (min)(const AutoDiffScalar<DerType>& x, const T& y) { return (x <= y ? x : y); } 528 template<typename DerType, typename T> 529 inline AutoDiffScalar<DerType> (max)(const AutoDiffScalar<DerType>& x, const T& y) { return (x >= y ? x : y); } 530 template<typename DerType, typename T> 531 inline AutoDiffScalar<DerType> (min)(const T& x, const AutoDiffScalar<DerType>& y) { return (x < y ? x : y); } 532 template<typename DerType, typename T> 533 inline AutoDiffScalar<DerType> (max)(const T& x, const AutoDiffScalar<DerType>& y) { return (x > y ? x : y); } 534 535 #define sign(x) x >= 0 ? 1 : -1 // required for abs function below 536 537 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs, 538 using std::abs; 539 return ReturnType(abs(x.value()), x.derivatives() * (sign(x.value())));) 540 541 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs2, 542 using internal::abs2; 543 return ReturnType(abs2(x.value()), x.derivatives() * (Scalar(2)*x.value()));) 544 545 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sqrt, 546 using std::sqrt; 547 Scalar sqrtx = sqrt(x.value()); 548 return ReturnType(sqrtx,x.derivatives() * (Scalar(0.5) / sqrtx));) 549 550 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cos, 551 using std::cos; 552 using std::sin; 553 return ReturnType(cos(x.value()), x.derivatives() * (-sin(x.value())));) 554 555 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sin, 556 using std::sin; 557 using std::cos; 558 return ReturnType(sin(x.value()),x.derivatives() * cos(x.value()));) 559 560 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(exp, 561 using std::exp; 562 Scalar expx = exp(x.value()); 563 return ReturnType(expx,x.derivatives() * expx);) 564 565 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(log, 566 using std::log; 567 return ReturnType(log(x.value()),x.derivatives() * (Scalar(1)/x.value()));) 568 569 template<typename DerType> 570 inline const Eigen::AutoDiffScalar<Eigen::CwiseUnaryOp<Eigen::internal::scalar_multiple_op<typename Eigen::internal::traits<DerType>::Scalar>, const DerType> > 571 pow(const Eigen::AutoDiffScalar<DerType>& x, typename Eigen::internal::traits<DerType>::Scalar y) 572 { 573 using namespace Eigen; 574 typedef typename Eigen::internal::traits<DerType>::Scalar Scalar; 575 return AutoDiffScalar<CwiseUnaryOp<Eigen::internal::scalar_multiple_op<Scalar>, const DerType> >( 576 std::pow(x.value(),y), 577 x.derivatives() * (y * std::pow(x.value(),y-1))); 578 } 579 580 581 template<typename DerTypeA,typename DerTypeB> 582 inline const AutoDiffScalar<Matrix<typename internal::traits<DerTypeA>::Scalar,Dynamic,1> > 583 atan2(const AutoDiffScalar<DerTypeA>& a, const AutoDiffScalar<DerTypeB>& b) 584 { 585 using std::atan2; 586 using std::max; 587 typedef typename internal::traits<DerTypeA>::Scalar Scalar; 588 typedef AutoDiffScalar<Matrix<Scalar,Dynamic,1> > PlainADS; 589 PlainADS ret; 590 ret.value() = atan2(a.value(), b.value()); 591 592 Scalar tmp2 = a.value() * a.value(); 593 Scalar tmp3 = b.value() * b.value(); 594 Scalar tmp4 = tmp3/(tmp2+tmp3); 595 596 if (tmp4!=0) 597 ret.derivatives() = (a.derivatives() * b.value() - a.value() * b.derivatives()) * (tmp2+tmp3); 598 599 return ret; 600 } 601 602 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tan, 603 using std::tan; 604 using std::cos; 605 return ReturnType(tan(x.value()),x.derivatives() * (Scalar(1)/internal::abs2(cos(x.value()))));) 606 607 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(asin, 608 using std::sqrt; 609 using std::asin; 610 return ReturnType(asin(x.value()),x.derivatives() * (Scalar(1)/sqrt(1-internal::abs2(x.value()))));) 611 612 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(acos, 613 using std::sqrt; 614 using std::acos; 615 return ReturnType(acos(x.value()),x.derivatives() * (Scalar(-1)/sqrt(1-internal::abs2(x.value()))));) 616 617 #undef EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY 618 619 template<typename DerType> struct NumTraits<AutoDiffScalar<DerType> > 620 : NumTraits< typename NumTraits<typename DerType::Scalar>::Real > 621 { 622 typedef AutoDiffScalar<Matrix<typename NumTraits<typename DerType::Scalar>::Real,DerType::RowsAtCompileTime,DerType::ColsAtCompileTime> > Real; 623 typedef AutoDiffScalar<DerType> NonInteger; 624 typedef AutoDiffScalar<DerType>& Nested; 625 enum{ 626 RequireInitialization = 1 627 }; 628 }; 629 630 } 631 632 #endif // EIGEN_AUTODIFF_SCALAR_H 633