Re: [eigen] Using Eigen with CppAD scalar types

[Brad Bell <bradbell@xxxxxxxxxx>]

Sorry for taking so long to get back to this. Please see my questions below.

On 06/18/2012 01:00 PM, Gael Guennebaud wrote:
> Hi,
>
> thanks for sharing this information. This reminded me that the doc on
> custom scalar types, as well as the AdolcForward unsupported module
> were both deprecated.So in theory, your support code should be rather
> written like:
>
> namespace Eigen {
> 	template<class Base>  struct NumTraits<  CppAD::AD<Base>  >
>           : NumTraits<Base>
> 	{	typedef CppAD::AD<Base>      Real;
> 		typedef CppAD::AD<Base>      NonInteger;
> 		typedef CppAD::AD<Base>&    Nested;
>
> 		enum {
> 			IsComplex             = 0 ,
> 			IsInteger             = 0 ,
> 			IsSigned              = 1 ,
> 			RequireInitialization = 1 ,
> 			ReadCost              = 1 ,
> 			AddCost               = 2 ,
> 			MulCost               = 2
> 		};
> 	};
> }
Would the following be OK ?

namespace Eigen {
     template<class Base>  struct NumTraits<  CppAD::AD<Base>  >
     {    typedef CppAD::AD<Base>    Real;
          typedef CppAD::AD<Base>    NonInteger;
          typedef CppAD::AD<Base>    Nested;

          enum {
               IsComplex             = 0 ,
               IsInteger             = 0 ,
               IsSigned              = 1 ,
               RequireInitialization = 1 ,
               ReadCost              = 1 ,
               AddCost               = 2 ,
               MulCost               = 2
          };
     };
     // machine epsilon with type of real part of x
     static inline Base epsilon(const AD<Base>&  x)
     {    return CppAD::numeric_limits<Base>::epsilon(); }
     // relaxed version of machine epsilon for comparison of different
     // operations that should result in the same value
     static inline Base dummy_precision(const AD<Base>&  x)
     {    return 100. * CppAD::numeric_limits<Base>::epsilon(); }
     // minimum normalized positive value
     static inline Base lowest(void)
     {    return CppAD::numeric_limits<  AD<Base>  >::min(); }
     // maximum finite value
     static inline Base highest(void)
     {    return CppAD::numeric_limits<  AD<Base>  >::max(); }
}


Note that the return type is Base, instead of AD<Base>, so it has a 
different signature than the template class
    template<typename T> struct GenericNumTraits
in the file NumTraits.h. All the standard numerical operations are 
defined between Base and AD<Base> objects (no matter what the type of 
Base; e.g., Base could be AD<double> which is used when one wants to 
differentiate just to define the function that needs to be differentiated).

> Making it inherits from NumTraits<Base>  permits to get the epsilon,
> dummy_precision, lowest and highest methods defined automatically from
> the one of the Base type.
>
> The Nested type should be reference to limit copies.
>
> Then you should only have to implement the few Eigen specific
> functions in the CppAD namespace:
>
> 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
> 	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; }
> }
>
> that's all you should need.

Does this mean that the CppAD::internal namespace is not used any more ?

> gael
> On Mon, Jun 18, 2012 at 4:42 PM, Brad Bell<bradbell@xxxxxxxxxx>  wrote:
> > On 06/18/2012 07:36 AM, Brad Bell wrote:
> > > I think that I have correctly connected CppAD::AD<Base>  types to Eigen
> > > (for the case where Base is a real type).
> > >
> > > Here is the Eigen custom scalar type interface
> > >      http://www.coin-or.org/CppAD/Doc/eigen_det.cpp.xml
> > Sorry for the second posting, the correct link for the scalar type interface
> > is
> >     http://www.coin-or.org/CppAD/Doc/cppad_eigen.hpp.xml
> >
> > > Here is an example usage
> > >     http://www.coin-or.org/CppAD/Doc/eigen_det.cpp.xml
> > >
> > > Here is the corresponding C++ source code
> > >
> > >   https://projects.coin-or.org/CppAD/browser/trunk/cppad/example/cppad_eigen.hpp
> > >     https://projects.coin-or.org/CppAD/browser/trunk/example/eigen_det.cpp