Change the underlying type of VectorXreal from double to float/long double

[zhanggiene]

I am using auto diff to write a numerical library that uses meta-programming so that not just double, but float/long double can also be used.

Currently, i can use eigen vector of type double to initialize ArrayXreal variable.

This example is working.

  Eigen::Matrix<double ,Eigen::Dynamic,1> xArrayF(3);
  autodiff::ArrayXreal a;
  xArrayF<<1,2,3,4;
  a=xArrayF;

This example is not working when the input vector is of float type.

  Eigen::Matrix<Float ,Eigen::Dynamic,1> xArrayF(3);
  autodiff::ArrayXreal a;
  xArrayF<<1,2,3,4;
  a=xArrayF;

This is the compilation error.

/usr/include/eigen3/Eigen/src/Core/AssignEvaluator.h: In instantiation of ‘void Eigen::internal::call_assignment_no_alias(Dst&, const Src&, const Func&) [with Dst = Eigen::Array<autodiff::detail::Real<1, double>, -1, 1>; Src = Eigen::Matrix<float, -1, 1>; Func = Eigen::internal::assign_op<autodiff::detail::Real<1, double>, float>]’:
/usr/include/eigen3/Eigen/src/Core/AssignEvaluator.h:804:27:   required from ‘void Eigen::internal::call_assignment(Dst&, const Src&, const Func&, typename Eigen::internal::enable_if<(! Eigen::internal::evaluator_assume_aliasing<Src>::value), void*>::type) [with Dst = Eigen::Array<autodiff::detail::Real<1, double>, -1, 1>; Src = Eigen::Matrix<float, -1, 1>; Func = Eigen::internal::assign_op<autodiff::detail::Real<1, double>, float>; typename Eigen::internal::enable_if<(! Eigen::internal::evaluator_assume_aliasing<Src>::value), void*>::type = void*]’
/usr/include/eigen3/Eigen/src/Core/AssignEvaluator.h:782:18:   required from ‘void Eigen::internal::call_assignment(Dst&, const Src&) [with Dst = Eigen::Array<autodiff::detail::Real<1, double>, -1, 1>; Src = Eigen::Matrix<float, -1, 1>]’
/usr/include/eigen3/Eigen/src/Core/PlainObjectBase.h:714:32:   required from ‘Derived& Eigen::PlainObjectBase<Derived>::_set(const Eigen::DenseBase<OtherDerived>&) [with OtherDerived = Eigen::Matrix<float, -1, 1>; Derived = Eigen::Array<autodiff::detail::Real<1, double>, -1, 1>]’
/usr/include/eigen3/Eigen/src/Core/Array.h:109:24:   required from ‘Eigen::Array<_Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols>& Eigen::Array<_Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols>::operator=(const Eigen::DenseBase<OtherDerived>&) [with OtherDerived = Eigen::Matrix<float, -1, 1>; _Scalar = autodiff::detail::Real<1, double>; int _Rows = -1; int _Cols = 1; int _Options = 0; int _MaxRows = -1; int _MaxCols = 1]’
/home/administrator/volesti/examples/Experiment/simple_hmc.cpp:283:5:   required from here
/usr/include/eigen3/Eigen/src/Core/util/XprHelper.h:819:96: error: static assertion failed: YOU_MIXED_DIFFERENT_NUMERIC_TYPES__YOU_NEED_TO_USE_THE_CAST_METHOD_OF_MATRIXBASE_TO_CAST_NUMERIC_TYPES_EXPLICITLY
  819 | gen::internal::has_ReturnType<ScalarBinaryOpTraits<LHS, RHS,BINOP> >::value), \
      | ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^~~~~~

/usr/include/eigen3/Eigen/src/Core/util/StaticAssert.h:33:54: note: in definition of macro ‘EIGEN_STATIC_ASSERT’
   33 |     #define EIGEN_STATIC_ASSERT(X,MSG) static_assert(X,#MSG);
      |                                                      ^
/usr/include/eigen3/Eigen/src/Core/AssignEvaluator.h:834:3: note: in expansion of macro ‘EIGEN_CHECK_BINARY_COMPATIBILIY’
  834 |   EIGEN_CHECK_BINARY_COMPATIBILIY(Func,typename ActualDstTypeCleaned::Scalar,typename Src::Scalar);
      |   ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
make[2]: *** [CMakeFiles/simple_hmc.dir/build.make:76: CMakeFiles/simple_hmc.dir/simple_hmc.o] Error 1
make[1]: *** [CMakeFiles/Makefile2:223: CMakeFiles/simple_hmc.dir/all] Error 2
make: *** [Makefile:91: all] Error 2

Ultimately i want to different wrt x of Eigen::Matrix<T ,Eigen::Dynamic,1> type where is T is a type supplied by user at compile time.

The current Autodiff::real only support double and there is no way to support other data types as well.

 using real1st = Real<1, double>;
using real = real1st;

[allanleal]

If you check the implementation of Real, you'll see that the underlying floating point type is a template type (T). So, you can achieve what you want with float (or any other type).

using realf = Real<1, float>;

You'll need to use Eigen types with float too. Maybe create new aliases for them.

[zhanggiene]

This requires modifying the source code. Is there any way for user to declare realf (as shown in the example above), without changing the source file?

[allanleal]

You don't need to modify the source code. You can define these new types in your own namespace:

namespace foo {

using realf = autodiff::detail::Real<1, float>;

}

[zhanggiene]

I have made examples using the "float" version of ArrayXreal. It might be helpful for other people. You can help to check it as well.

#include <iostream>

// autodiff include
#include <autodiff/forward/real.hpp>
#include <autodiff/forward/real/eigen.hpp>
#include <Eigen/Core>
#include <Eigen/Dense>



// The single-variable function for which derivatives are needed

namespace foo {

using realf = autodiff::detail::Real<1, float>;

}

autodiff::real f(const autodiff::ArrayXreal& x)
{
    return (x * x.exp()).sum(); // sum([xi * exp(xi) for i = 1:5])
}

foo::realf f2(const Eigen::Array<foo::realf,Eigen::Dynamic,1>& x)
{
    return (x * x.exp()).sum();
}

int main()
{

    using Eigen::VectorXd;

    autodiff::ArrayXreal x(5);                            // the input array x with 5 variables
    x << 1, 2, 3, 4, 5;                         // x = [1, 2, 3, 4, 5]

    autodiff::real u;                                     // the output scalar u = f(x) evaluated together with gradient below

    VectorXd g = autodiff::gradient(f, autodiff::wrt(x), autodiff::at(x), u); // evaluate the function value u and its gradient vector g = du/dx

    std::cout << "u = " << u << std::endl;      // print the evaluated output u
    std::cout << "g = \n" << g << std::endl;


    Eigen::Array<foo::realf,Eigen::Dynamic,1> x2(5);                            // the input array x with 5 variables
    x2 << 1, 2, 3, 4, 5;                         // x = [1, 2, 3, 4, 5]

    foo::realf u2;                                     // the output scalar u = f(x) evaluated together with gradient below

    auto g2 = autodiff::gradient(f2, autodiff::wrt(x2), autodiff::at(x2), u2); // evaluate the function value u and its gradient vector g = du/dx

    std::cout << "u2 = " << u2 << std::endl;      // print the evaluated output u
    std::cout << "g2 = \n" << g2 << std::endl;

   


}