Enable Use of Eigen Linear Algebra Package with CppAD

@(@\newcommand{\W}[1]{ \; #1 \; } \newcommand{\R}[1]{ {\rm #1} } \newcommand{\B}[1]{ {\bf #1} } \newcommand{\D}[2]{ \frac{\partial #1}{\partial #2} } \newcommand{\DD}[3]{ \frac{\partial^2 #1}{\partial #2 \partial #3} } \newcommand{\Dpow}[2]{ \frac{\partial^{#1}}{\partial {#2}^{#1}} } \newcommand{\dpow}[2]{ \frac{ {\rm d}^{#1}}{{\rm d}\, {#2}^{#1}} }@)@Enable Use of Eigen Linear Algebra Package with CppAD
Syntax
# include <cppad/example/cppad_eigen.hpp>

Purpose
Enables the use of the eigen linear algebra package with the type AD<Base> ; see custom scalar types .

Example
The files eigen_array.cpp and eigen_det.cpp contain an example and test of this include file. They return true if they succeed and false otherwise.

Include Files
The file cppad_eigen.hpp includes both <cppad/cppad.hpp> and <Eigen/Core>. The file eigen_plugin.hpp defines value_type in the Eigen matrix class so its vectors are simple vectors (not necessary for eigen-3.3.3 and later).


# define EIGEN_MATRIXBASE_PLUGIN <cppad/example/eigen_plugin.hpp>
# include <Eigen/Core>
# include <cppad/cppad.hpp>

Eigen NumTraits
Eigen needs the following definitions to work properly with AD<Base> scalars:

namespace Eigen {
     template <class Base> struct NumTraits< CppAD::AD<Base> >
     {     // type that corresponds to the real part of an AD<Base> value
          typedef CppAD::AD<Base>   Real;
          // type for AD<Base> operations that result in non-integer values
          typedef CppAD::AD<Base>   NonInteger;
          //  type to use for numeric literals such as "2" or "0.5".
          typedef CppAD::AD<Base>   Literal;
          // type for nested value inside an AD<Base> expression tree
          typedef CppAD::AD<Base>   Nested;

          enum {
               // does not support complex Base types
               IsComplex             = 0 ,
               // does not support integer Base types
               IsInteger             = 0 ,
               // only support signed Base types
               IsSigned              = 1 ,
               // must initialize an AD<Base> object
               RequireInitialization = 1 ,
               // computational cost of the corresponding operations
               ReadCost              = 1 ,
               AddCost               = 2 ,
               MulCost               = 2
          };

          // machine epsilon with type of real part of x
          // (use assumption that Base is not complex)
          static CppAD::AD<Base> epsilon(void)
          {     return CppAD::numeric_limits< CppAD::AD<Base> >::epsilon(); }

          // relaxed version of machine epsilon for comparison of different
          // operations that should result in the same value
          static CppAD::AD<Base> dummy_precision(void)
          {     return 100. *
                    CppAD::numeric_limits< CppAD::AD<Base> >::epsilon();
          }

          // minimum normalized positive value
          static CppAD::AD<Base> lowest(void)
          {     return CppAD::numeric_limits< CppAD::AD<Base> >::min(); }

          // maximum finite value
          static CppAD::AD<Base> highest(void)
          {     return CppAD::numeric_limits< CppAD::AD<Base> >::max(); }

          // number of decimal digits that can be represented without change.
          static int digits10(void)
          {     return CppAD::numeric_limits< CppAD::AD<Base> >::digits10; }
     };
}

CppAD Namespace
Eigen also needs the following definitions to work properly with AD<Base> scalars:

namespace CppAD {
          // functions that return references
          template <class Base> const AD<Base>& conj(const AD<Base>& x)
          {     return x; }
          template <class Base> const AD<Base>& real(const AD<Base>& x)
          {     return x; }

          // functions that return values (note abs is defined by cppad.hpp)
          template <class Base> AD<Base> imag(const AD<Base>& x)
          {     return CppAD::AD<Base>(0.); }
          template <class Base> AD<Base> abs2(const AD<Base>& x)
          {     return x * x; }
}

Input File: cppad/example/cppad_eigen.hpp