autodiff / autodiff

automatic differentiation made easier for C++
https://autodiff.github.io
MIT License
1.63k stars 172 forks source link

Can't add to VectorXdual #60

Open ibell opened 5 years ago

ibell commented 5 years ago

On MSVC 2019, this code

// C++ includes
#include <iostream>

// autodiff include
#include <autodiff/forward.hpp>
#include "Eigen/Dense"
#include <autodiff/forward/eigen.hpp>
using namespace autodiff;

int main()
{
    Eigen::VectorXdual x(2);  // the input vector x
    x << -0.5, 0.5;
    Eigen::VectorXdual x2 = x + 0.5 * 0.001;
    std::cout << x2 << "\n";
}

yields the compile error

1>D:\Code\helmderiv\snippets\bugdualvec.cpp(14,31): error C2676:  binary '+': 'Eigen::VectorXdual' does not define this operator or a conversion to a type acceptable to the predefined operator
ibell commented 5 years ago

Sorry, my fault (sort of). Eigen::Vector don't allow for coefficient-wise addition, but Eigen::Array do. Could we add Eigen::Array typedefs too for generality? Of course also possible to do two conversions, one to array, and another back to matrix. But that's unnecessarily verbose.

ludkinm commented 4 years ago

@ibell In this case you can do

 Eigen::VectorXdual x2 = x.array() + 0.5 * 0.001;

Although, be careful. Using the .array member function won't work inside expressions that take derivatives.

allanleal commented 4 years ago

Using the .array member function won't work inside expressions that take derivatives.

Hi @ludkinm , what do you mean above? Could you please give an example showing what should be the expected result and the actual one when using method .array() of Eigen on vectors/matrices?

ludkinm commented 4 years ago

I meant it would not compile. Say f(x) = sum( x*(x+1) )

return x.cwiseProduct(x.array() + 1.0).sum();

However, I realise if you construct the result into VectorXdual then it works:

return x.cwiseProduct(VectorXdual(x.array() + 1.0)).sum();

In my previous comment I guess the constructor for VectorXdual from VectorX is being called (explicitly), whereas VectorXdual::cwiseProduct(VectorXdual x) expects a VectorXdual. Since, the VectorXdual::array() returns a some kind of Eigen array object this is not converted to VectorXdual since there is no conversion operator (?)