00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025 #ifndef EIGEN_AUTODIFF_SCALAR_H
00026 #define EIGEN_AUTODIFF_SCALAR_H
00027
00028 namespace Eigen {
00029
00030 namespace internal {
00031
00032 template<typename A, typename B>
00033 struct make_coherent_impl {
00034 static void run(A&, B&) {}
00035 };
00036
00037
00038 template<typename A, typename B>
00039 void make_coherent(const A& a, const B&b)
00040 {
00041 make_coherent_impl<A,B>::run(a.const_cast_derived(), b.const_cast_derived());
00042 }
00043
00044 template<typename _DerType, bool Enable> struct auto_diff_special_op;
00045
00046 }
00047
00074 template<typename _DerType>
00075 class AutoDiffScalar
00076 : public internal::auto_diff_special_op
00077 <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
00078 typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value>
00079 {
00080 public:
00081 typedef internal::auto_diff_special_op
00082 <_DerType, !internal::is_same<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar,
00083 typename NumTraits<typename internal::traits<typename internal::remove_all<_DerType>::type>::Scalar>::Real>::value> Base;
00084 typedef typename internal::remove_all<_DerType>::type DerType;
00085 typedef typename internal::traits<DerType>::Scalar Scalar;
00086 typedef typename NumTraits<Scalar>::Real Real;
00087
00088 using Base::operator+;
00089 using Base::operator*;
00090
00092 AutoDiffScalar() {}
00093
00096 AutoDiffScalar(const Scalar& value, int nbDer, int derNumber)
00097 : m_value(value), m_derivatives(DerType::Zero(nbDer))
00098 {
00099 m_derivatives.coeffRef(derNumber) = Scalar(1);
00100 }
00101
00104 AutoDiffScalar(const Real& value)
00105 : m_value(value)
00106 {
00107 if(m_derivatives.size()>0)
00108 m_derivatives.setZero();
00109 }
00110
00112 AutoDiffScalar(const Scalar& value, const DerType& der)
00113 : m_value(value), m_derivatives(der)
00114 {}
00115
00116 template<typename OtherDerType>
00117 AutoDiffScalar(const AutoDiffScalar<OtherDerType>& other)
00118 : m_value(other.value()), m_derivatives(other.derivatives())
00119 {}
00120
00121 friend std::ostream & operator << (std::ostream & s, const AutoDiffScalar& a)
00122 {
00123 return s << a.value();
00124 }
00125
00126 AutoDiffScalar(const AutoDiffScalar& other)
00127 : m_value(other.value()), m_derivatives(other.derivatives())
00128 {}
00129
00130 template<typename OtherDerType>
00131 inline AutoDiffScalar& operator=(const AutoDiffScalar<OtherDerType>& other)
00132 {
00133 m_value = other.value();
00134 m_derivatives = other.derivatives();
00135 return *this;
00136 }
00137
00138 inline AutoDiffScalar& operator=(const AutoDiffScalar& other)
00139 {
00140 m_value = other.value();
00141 m_derivatives = other.derivatives();
00142 return *this;
00143 }
00144
00145
00146
00147
00148 inline const Scalar& value() const { return m_value; }
00149 inline Scalar& value() { return m_value; }
00150
00151 inline const DerType& derivatives() const { return m_derivatives; }
00152 inline DerType& derivatives() { return m_derivatives; }
00153
00154 inline bool operator< (const Scalar& other) const { return m_value < other; }
00155 inline bool operator<=(const Scalar& other) const { return m_value <= other; }
00156 inline bool operator> (const Scalar& other) const { return m_value > other; }
00157 inline bool operator>=(const Scalar& other) const { return m_value >= other; }
00158 inline bool operator==(const Scalar& other) const { return m_value == other; }
00159 inline bool operator!=(const Scalar& other) const { return m_value != other; }
00160
00161 friend inline bool operator< (const Scalar& a, const AutoDiffScalar& b) { return a < b.value(); }
00162 friend inline bool operator<=(const Scalar& a, const AutoDiffScalar& b) { return a <= b.value(); }
00163 friend inline bool operator> (const Scalar& a, const AutoDiffScalar& b) { return a > b.value(); }
00164 friend inline bool operator>=(const Scalar& a, const AutoDiffScalar& b) { return a >= b.value(); }
00165 friend inline bool operator==(const Scalar& a, const AutoDiffScalar& b) { return a == b.value(); }
00166 friend inline bool operator!=(const Scalar& a, const AutoDiffScalar& b) { return a != b.value(); }
00167
00168 template<typename OtherDerType> inline bool operator< (const AutoDiffScalar<OtherDerType>& b) const { return m_value < b.value(); }
00169 template<typename OtherDerType> inline bool operator<=(const AutoDiffScalar<OtherDerType>& b) const { return m_value <= b.value(); }
00170 template<typename OtherDerType> inline bool operator> (const AutoDiffScalar<OtherDerType>& b) const { return m_value > b.value(); }
00171 template<typename OtherDerType> inline bool operator>=(const AutoDiffScalar<OtherDerType>& b) const { return m_value >= b.value(); }
00172 template<typename OtherDerType> inline bool operator==(const AutoDiffScalar<OtherDerType>& b) const { return m_value == b.value(); }
00173 template<typename OtherDerType> inline bool operator!=(const AutoDiffScalar<OtherDerType>& b) const { return m_value != b.value(); }
00174
00175 inline const AutoDiffScalar<DerType&> operator+(const Scalar& other) const
00176 {
00177 return AutoDiffScalar<DerType&>(m_value + other, m_derivatives);
00178 }
00179
00180 friend inline const AutoDiffScalar<DerType&> operator+(const Scalar& a, const AutoDiffScalar& b)
00181 {
00182 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
00183 }
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194
00195 inline AutoDiffScalar& operator+=(const Scalar& other)
00196 {
00197 value() += other;
00198 return *this;
00199 }
00200
00201 template<typename OtherDerType>
00202 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> >
00203 operator+(const AutoDiffScalar<OtherDerType>& other) const
00204 {
00205 internal::make_coherent(m_derivatives, other.derivatives());
00206 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,const DerType,const typename internal::remove_all<OtherDerType>::type> >(
00207 m_value + other.value(),
00208 m_derivatives + other.derivatives());
00209 }
00210
00211 template<typename OtherDerType>
00212 inline AutoDiffScalar&
00213 operator+=(const AutoDiffScalar<OtherDerType>& other)
00214 {
00215 (*this) = (*this) + other;
00216 return *this;
00217 }
00218
00219 inline const AutoDiffScalar<DerType&> operator-(const Scalar& b) const
00220 {
00221 return AutoDiffScalar<DerType&>(m_value - b, m_derivatives);
00222 }
00223
00224 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
00225 operator-(const Scalar& a, const AutoDiffScalar& b)
00226 {
00227 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
00228 (a - b.value(), -b.derivatives());
00229 }
00230
00231 inline AutoDiffScalar& operator-=(const Scalar& other)
00232 {
00233 value() -= other;
00234 return *this;
00235 }
00236
00237 template<typename OtherDerType>
00238 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> >
00239 operator-(const AutoDiffScalar<OtherDerType>& other) const
00240 {
00241 internal::make_coherent(m_derivatives, other.derivatives());
00242 return AutoDiffScalar<CwiseBinaryOp<internal::scalar_difference_op<Scalar>, const DerType,const typename internal::remove_all<OtherDerType>::type> >(
00243 m_value - other.value(),
00244 m_derivatives - other.derivatives());
00245 }
00246
00247 template<typename OtherDerType>
00248 inline AutoDiffScalar&
00249 operator-=(const AutoDiffScalar<OtherDerType>& other)
00250 {
00251 *this = *this - other;
00252 return *this;
00253 }
00254
00255 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >
00256 operator-() const
00257 {
00258 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_opposite_op<Scalar>, const DerType> >(
00259 -m_value,
00260 -m_derivatives);
00261 }
00262
00263 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
00264 operator*(const Scalar& other) const
00265 {
00266 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
00267 m_value * other,
00268 (m_derivatives * other));
00269 }
00270
00271 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
00272 operator*(const Scalar& other, const AutoDiffScalar& a)
00273 {
00274 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
00275 a.value() * other,
00276 a.derivatives() * other);
00277 }
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
00296 operator/(const Scalar& other) const
00297 {
00298 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
00299 m_value / other,
00300 (m_derivatives * (Scalar(1)/other)));
00301 }
00302
00303 friend inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >
00304 operator/(const Scalar& other, const AutoDiffScalar& a)
00305 {
00306 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType> >(
00307 other / a.value(),
00308 a.derivatives() * (Scalar(-other) / (a.value()*a.value())));
00309 }
00310
00311
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326
00327 template<typename OtherDerType>
00328 inline const AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>,
00329 const CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
00330 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
00331 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > >
00332 operator/(const AutoDiffScalar<OtherDerType>& other) const
00333 {
00334 internal::make_coherent(m_derivatives, other.derivatives());
00335 return AutoDiffScalar<CwiseUnaryOp<internal::scalar_multiple_op<Scalar>,
00336 const CwiseBinaryOp<internal::scalar_difference_op<Scalar>,
00337 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
00338 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > > >(
00339 m_value / other.value(),
00340 ((m_derivatives * other.value()) - (m_value * other.derivatives()))
00341 * (Scalar(1)/(other.value()*other.value())));
00342 }
00343
00344 template<typename OtherDerType>
00345 inline const AutoDiffScalar<CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
00346 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
00347 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type> > >
00348 operator*(const AutoDiffScalar<OtherDerType>& other) const
00349 {
00350 internal::make_coherent(m_derivatives, other.derivatives());
00351 return AutoDiffScalar<const CwiseBinaryOp<internal::scalar_sum_op<Scalar>,
00352 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const DerType>,
00353 const CwiseUnaryOp<internal::scalar_multiple_op<Scalar>, const typename internal::remove_all<OtherDerType>::type > > >(
00354 m_value * other.value(),
00355 (m_derivatives * other.value()) + (m_value * other.derivatives()));
00356 }
00357
00358 inline AutoDiffScalar& operator*=(const Scalar& other)
00359 {
00360 *this = *this * other;
00361 return *this;
00362 }
00363
00364 template<typename OtherDerType>
00365 inline AutoDiffScalar& operator*=(const AutoDiffScalar<OtherDerType>& other)
00366 {
00367 *this = *this * other;
00368 return *this;
00369 }
00370
00371 inline AutoDiffScalar& operator/=(const Scalar& other)
00372 {
00373 *this = *this / other;
00374 return *this;
00375 }
00376
00377 template<typename OtherDerType>
00378 inline AutoDiffScalar& operator/=(const AutoDiffScalar<OtherDerType>& other)
00379 {
00380 *this = *this / other;
00381 return *this;
00382 }
00383
00384 protected:
00385 Scalar m_value;
00386 DerType m_derivatives;
00387
00388 };
00389
00390 namespace internal {
00391
00392 template<typename _DerType>
00393 struct auto_diff_special_op<_DerType, true>
00394
00395
00396 {
00397 typedef typename remove_all<_DerType>::type DerType;
00398 typedef typename traits<DerType>::Scalar Scalar;
00399 typedef typename NumTraits<Scalar>::Real Real;
00400
00401
00402
00403
00404
00405
00406
00407
00408
00409
00410
00411 const AutoDiffScalar<_DerType>& derived() const { return *static_cast<const AutoDiffScalar<_DerType>*>(this); }
00412 AutoDiffScalar<_DerType>& derived() { return *static_cast<AutoDiffScalar<_DerType>*>(this); }
00413
00414
00415 inline const AutoDiffScalar<DerType&> operator+(const Real& other) const
00416 {
00417 return AutoDiffScalar<DerType&>(derived().value() + other, derived().derivatives());
00418 }
00419
00420 friend inline const AutoDiffScalar<DerType&> operator+(const Real& a, const AutoDiffScalar<_DerType>& b)
00421 {
00422 return AutoDiffScalar<DerType&>(a + b.value(), b.derivatives());
00423 }
00424
00425 inline AutoDiffScalar<_DerType>& operator+=(const Real& other)
00426 {
00427 derived().value() += other;
00428 return derived();
00429 }
00430
00431
00432 inline const AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >
00433 operator*(const Real& other) const
00434 {
00435 return AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >(
00436 derived().value() * other,
00437 derived().derivatives() * other);
00438 }
00439
00440 friend inline const AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >
00441 operator*(const Real& other, const AutoDiffScalar<_DerType>& a)
00442 {
00443 return AutoDiffScalar<typename CwiseUnaryOp<scalar_multiple2_op<Scalar,Real>, DerType>::Type >(
00444 a.value() * other,
00445 a.derivatives() * other);
00446 }
00447
00448 inline AutoDiffScalar<_DerType>& operator*=(const Scalar& other)
00449 {
00450 *this = *this * other;
00451 return derived();
00452 }
00453 };
00454
00455 template<typename _DerType>
00456 struct auto_diff_special_op<_DerType, false>
00457 {
00458 void operator*() const;
00459 void operator-() const;
00460 void operator+() const;
00461 };
00462
00463 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols, typename B>
00464 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>, B> {
00465 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
00466 static void run(A& a, B& b) {
00467 if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
00468 {
00469 a.resize(b.size());
00470 a.setZero();
00471 }
00472 }
00473 };
00474
00475 template<typename A, typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols>
00476 struct make_coherent_impl<A, Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
00477 typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
00478 static void run(A& a, B& b) {
00479 if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
00480 {
00481 b.resize(a.size());
00482 b.setZero();
00483 }
00484 }
00485 };
00486
00487 template<typename A_Scalar, int A_Rows, int A_Cols, int A_Options, int A_MaxRows, int A_MaxCols,
00488 typename B_Scalar, int B_Rows, int B_Cols, int B_Options, int B_MaxRows, int B_MaxCols>
00489 struct make_coherent_impl<Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols>,
00490 Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> > {
00491 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> A;
00492 typedef Matrix<B_Scalar, B_Rows, B_Cols, B_Options, B_MaxRows, B_MaxCols> B;
00493 static void run(A& a, B& b) {
00494 if((A_Rows==Dynamic || A_Cols==Dynamic) && (a.size()==0))
00495 {
00496 a.resize(b.size());
00497 a.setZero();
00498 }
00499 else if((B_Rows==Dynamic || B_Cols==Dynamic) && (b.size()==0))
00500 {
00501 b.resize(a.size());
00502 b.setZero();
00503 }
00504 }
00505 };
00506
00507 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>
00508 {
00509 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType;
00510 };
00511
00512 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> >
00513 {
00514 typedef Matrix<A_Scalar, A_Rows, A_Cols, A_Options, A_MaxRows, A_MaxCols> ReturnType;
00515 };
00516
00517 template<typename DerType>
00518 struct scalar_product_traits<AutoDiffScalar<DerType>,typename DerType::Scalar>
00519 {
00520 typedef AutoDiffScalar<DerType> ReturnType;
00521 };
00522
00523 }
00524
00525 #define EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(FUNC,CODE) \
00526 template<typename DerType> \
00527 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> > \
00528 FUNC(const Eigen::AutoDiffScalar<DerType>& x) { \
00529 using namespace Eigen; \
00530 typedef typename Eigen::internal::traits<typename Eigen::internal::remove_all<DerType>::type>::Scalar Scalar; \
00531 typedef AutoDiffScalar<CwiseUnaryOp<Eigen::internal::scalar_multiple_op<Scalar>, const typename Eigen::internal::remove_all<DerType>::type> > ReturnType; \
00532 CODE; \
00533 }
00534
00535 template<typename DerType>
00536 inline const AutoDiffScalar<DerType>& conj(const AutoDiffScalar<DerType>& x) { return x; }
00537 template<typename DerType>
00538 inline const AutoDiffScalar<DerType>& real(const AutoDiffScalar<DerType>& x) { return x; }
00539 template<typename DerType>
00540 inline typename DerType::Scalar imag(const AutoDiffScalar<DerType>&) { return 0.; }
00541 template<typename DerType, typename T>
00542 inline AutoDiffScalar<DerType> (min)(const AutoDiffScalar<DerType>& x, const T& y) { return (x <= y ? x : y); }
00543 template<typename DerType, typename T>
00544 inline AutoDiffScalar<DerType> (max)(const AutoDiffScalar<DerType>& x, const T& y) { return (x >= y ? x : y); }
00545 template<typename DerType, typename T>
00546 inline AutoDiffScalar<DerType> (min)(const T& x, const AutoDiffScalar<DerType>& y) { return (x < y ? x : y); }
00547 template<typename DerType, typename T>
00548 inline AutoDiffScalar<DerType> (max)(const T& x, const AutoDiffScalar<DerType>& y) { return (x > y ? x : y); }
00549
00550 #define sign(x) x >= 0 ? 1 : -1 // required for abs function below
00551
00552 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs,
00553 using std::abs;
00554 return ReturnType(abs(x.value()), x.derivatives() * (sign(x.value())));)
00555
00556 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(abs2,
00557 using internal::abs2;
00558 return ReturnType(abs2(x.value()), x.derivatives() * (Scalar(2)*x.value()));)
00559
00560 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sqrt,
00561 using std::sqrt;
00562 Scalar sqrtx = sqrt(x.value());
00563 return ReturnType(sqrtx,x.derivatives() * (Scalar(0.5) / sqrtx));)
00564
00565 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(cos,
00566 using std::cos;
00567 using std::sin;
00568 return ReturnType(cos(x.value()), x.derivatives() * (-sin(x.value())));)
00569
00570 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(sin,
00571 using std::sin;
00572 using std::cos;
00573 return ReturnType(sin(x.value()),x.derivatives() * cos(x.value()));)
00574
00575 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(exp,
00576 using std::exp;
00577 Scalar expx = exp(x.value());
00578 return ReturnType(expx,x.derivatives() * expx);)
00579
00580 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(log,
00581 using std::log;
00582 return ReturnType(log(x.value()),x.derivatives() * (Scalar(1)/x.value()));)
00583
00584 template<typename DerType>
00585 inline const Eigen::AutoDiffScalar<Eigen::CwiseUnaryOp<Eigen::internal::scalar_multiple_op<typename Eigen::internal::traits<DerType>::Scalar>, const DerType> >
00586 pow(const Eigen::AutoDiffScalar<DerType>& x, typename Eigen::internal::traits<DerType>::Scalar y)
00587 {
00588 using namespace Eigen;
00589 typedef typename Eigen::internal::traits<DerType>::Scalar Scalar;
00590 return AutoDiffScalar<CwiseUnaryOp<Eigen::internal::scalar_multiple_op<Scalar>, const DerType> >(
00591 std::pow(x.value(),y),
00592 x.derivatives() * (y * std::pow(x.value(),y-1)));
00593 }
00594
00595
00596 template<typename DerTypeA,typename DerTypeB>
00597 inline const AutoDiffScalar<Matrix<typename internal::traits<DerTypeA>::Scalar,Dynamic,1> >
00598 atan2(const AutoDiffScalar<DerTypeA>& a, const AutoDiffScalar<DerTypeB>& b)
00599 {
00600 using std::atan2;
00601 using std::max;
00602 typedef typename internal::traits<DerTypeA>::Scalar Scalar;
00603 typedef AutoDiffScalar<Matrix<Scalar,Dynamic,1> > PlainADS;
00604 PlainADS ret;
00605 ret.value() = atan2(a.value(), b.value());
00606
00607 Scalar tmp2 = a.value() * a.value();
00608 Scalar tmp3 = b.value() * b.value();
00609 Scalar tmp4 = tmp3/(tmp2+tmp3);
00610
00611 if (tmp4!=0)
00612 ret.derivatives() = (a.derivatives() * b.value() - a.value() * b.derivatives()) * (tmp2+tmp3);
00613
00614 return ret;
00615 }
00616
00617 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(tan,
00618 using std::tan;
00619 using std::cos;
00620 return ReturnType(tan(x.value()),x.derivatives() * (Scalar(1)/internal::abs2(cos(x.value()))));)
00621
00622 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(asin,
00623 using std::sqrt;
00624 using std::asin;
00625 return ReturnType(asin(x.value()),x.derivatives() * (Scalar(1)/sqrt(1-internal::abs2(x.value()))));)
00626
00627 EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY(acos,
00628 using std::sqrt;
00629 using std::acos;
00630 return ReturnType(acos(x.value()),x.derivatives() * (Scalar(-1)/sqrt(1-internal::abs2(x.value()))));)
00631
00632 #undef EIGEN_AUTODIFF_DECLARE_GLOBAL_UNARY
00633
00634 template<typename DerType> struct NumTraits<AutoDiffScalar<DerType> >
00635 : NumTraits< typename NumTraits<typename DerType::Scalar>::Real >
00636 {
00637 typedef AutoDiffScalar<Matrix<typename NumTraits<typename DerType::Scalar>::Real,DerType::RowsAtCompileTime,DerType::ColsAtCompileTime> > Real;
00638 typedef AutoDiffScalar<DerType> NonInteger;
00639 typedef AutoDiffScalar<DerType>& Nested;
00640 enum{
00641 RequireInitialization = 1
00642 };
00643 };
00644
00645 }
00646
00647 #endif // EIGEN_AUTODIFF_SCALAR_H