![]() |
Eigen-unsupported
3.3.3
|
00001 // This file is part of Eigen, a lightweight C++ template library 00002 // for linear algebra. 00003 // 00004 // Copyright (C) 2009 Gael Guennebaud <gael.guennebaud@inria.fr> 00005 // 00006 // This Source Code Form is subject to the terms of the Mozilla 00007 // Public License v. 2.0. If a copy of the MPL was not distributed 00008 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 00009 00010 #ifndef EIGEN_AUTODIFF_JACOBIAN_H 00011 #define EIGEN_AUTODIFF_JACOBIAN_H 00012 00013 namespace Eigen 00014 { 00015 00016 template<typename Functor> class AutoDiffJacobian : public Functor 00017 { 00018 public: 00019 AutoDiffJacobian() : Functor() {} 00020 AutoDiffJacobian(const Functor& f) : Functor(f) {} 00021 00022 // forward constructors 00023 #if EIGEN_HAS_VARIADIC_TEMPLATES 00024 template<typename... T> 00025 AutoDiffJacobian(const T& ...Values) : Functor(Values...) {} 00026 #else 00027 template<typename T0> 00028 AutoDiffJacobian(const T0& a0) : Functor(a0) {} 00029 template<typename T0, typename T1> 00030 AutoDiffJacobian(const T0& a0, const T1& a1) : Functor(a0, a1) {} 00031 template<typename T0, typename T1, typename T2> 00032 AutoDiffJacobian(const T0& a0, const T1& a1, const T2& a2) : Functor(a0, a1, a2) {} 00033 #endif 00034 00035 typedef typename Functor::InputType InputType; 00036 typedef typename Functor::ValueType ValueType; 00037 typedef typename ValueType::Scalar Scalar; 00038 00039 enum { 00040 InputsAtCompileTime = InputType::RowsAtCompileTime, 00041 ValuesAtCompileTime = ValueType::RowsAtCompileTime 00042 }; 00043 00044 typedef Matrix<Scalar, ValuesAtCompileTime, InputsAtCompileTime> JacobianType; 00045 typedef typename JacobianType::Index Index; 00046 00047 typedef Matrix<Scalar, InputsAtCompileTime, 1> DerivativeType; 00048 typedef AutoDiffScalar<DerivativeType> ActiveScalar; 00049 00050 typedef Matrix<ActiveScalar, InputsAtCompileTime, 1> ActiveInput; 00051 typedef Matrix<ActiveScalar, ValuesAtCompileTime, 1> ActiveValue; 00052 00053 #if EIGEN_HAS_VARIADIC_TEMPLATES 00054 // Some compilers don't accept variadic parameters after a default parameter, 00055 // i.e., we can't just write _jac=0 but we need to overload operator(): 00056 EIGEN_STRONG_INLINE 00057 void operator() (const InputType& x, ValueType* v) const 00058 { 00059 this->operator()(x, v, 0); 00060 } 00061 template<typename... ParamsType> 00062 void operator() (const InputType& x, ValueType* v, JacobianType* _jac, 00063 const ParamsType&... Params) const 00064 #else 00065 void operator() (const InputType& x, ValueType* v, JacobianType* _jac=0) const 00066 #endif 00067 { 00068 eigen_assert(v!=0); 00069 00070 if (!_jac) 00071 { 00072 #if EIGEN_HAS_VARIADIC_TEMPLATES 00073 Functor::operator()(x, v, Params...); 00074 #else 00075 Functor::operator()(x, v); 00076 #endif 00077 return; 00078 } 00079 00080 JacobianType& jac = *_jac; 00081 00082 ActiveInput ax = x.template cast<ActiveScalar>(); 00083 ActiveValue av(jac.rows()); 00084 00085 if(InputsAtCompileTime==Dynamic) 00086 for (Index j=0; j<jac.rows(); j++) 00087 av[j].derivatives().resize(x.rows()); 00088 00089 for (Index i=0; i<jac.cols(); i++) 00090 ax[i].derivatives() = DerivativeType::Unit(x.rows(),i); 00091 00092 #if EIGEN_HAS_VARIADIC_TEMPLATES 00093 Functor::operator()(ax, &av, Params...); 00094 #else 00095 Functor::operator()(ax, &av); 00096 #endif 00097 00098 for (Index i=0; i<jac.rows(); i++) 00099 { 00100 (*v)[i] = av[i].value(); 00101 jac.row(i) = av[i].derivatives(); 00102 } 00103 } 00104 }; 00105 00106 } 00107 00108 #endif // EIGEN_AUTODIFF_JACOBIAN_H