6

Writing a function implementing a Jacobi Matrix for Newton's method I have noticed a really nasty error.

Calling the function 

auto DF = [T](VectorXd y){
    return PhiAndW(y(0),y(1),T).second - MatrixXd::Identity(2,2);
 };

only returns the value of PhiAndW(y(0),y(1),T).second and omits the subtraction of MatrixXd::Identity(2,2). But if I change the code to 

auto DF = [T](VectorXd y){
    MatrixXd mat = PhiAndW(y(0),y(1),T).second - MatrixXd::Identity(2,2);
    return mat;
 };

everything works seamlessly.

I have tried to reproduce it and this is not the exact same behavior but it also doesn't behave as wanted:

Consider the following code:

MatrixXd FF(MatrixXd y){
  return y;
}

int other(){

  auto DF = [](MatrixXd y){
    MatrixXd test = FF(y)  - MatrixXd::Identity(2,2);
    return test;
  };

  std::cout << DF(MatrixXd::Ones(2,2)) <<std::endl;
  std::cout << std::endl;
  std::cout << (MatrixXd::Ones(2,2) - MatrixXd::Identity(2,2))<< std::endl;

}

It will print

>  1 0
>  0 1 
>
>  1 0 
>  0 1

to the console.

However if I change to function DF to 

  auto DF = [](MatrixXd y){
    return FF(y)  - MatrixXd::Identity(2,2);
  };

The console will print 

>  2.22045e-15 1.63042e-322
>  2.37152e-322    -0.999998

for the second matrix. (Which is just some uninitialized junk from memory).

Could someone explain what is happening with my code and my example problem. I have truly no idea what is going on here. I am especially interested why saving the result of the calculation in a temporary variable fixes the problem.

David
  • 143
  • 7
  • 4
    Please read the section about using `auto` on this page: https://eigen.tuxfamily.org/dox-devel/TopicPitfalls.html – chtz Jan 03 '20 at 23:43
  • (in this case `auto` is implicitly the return type of the lambda) – Marc Glisse Jan 03 '20 at 23:46
  • Did you try `auto DF = []() -> MatrixXd {...}` ? This should force the evaluation of return value. – cplusplusrat Jan 03 '20 at 23:46
  • That's another nice example that `auto` might not be that cool like everyone else (except me) seems to believe. ;-) I assume an explicit cast would fix the first sample as well: `auto DF = [](MatrixXd y){ return MatrixXd(FF(y) - MatrixXd::Identity(2,2)); };` – Scheff's Cat Jan 04 '20 at 08:27
  • 1
    @Scheff [You're not the only one.](https://stackoverflow.com/a/32758566/2899559) – Avi Ginsburg Jan 05 '20 at 08:55

1 Answers1

4

Since the comments have pretty much solved my issue (thank you very much) I thought I'd go ahead and answer my question so that other people see that this problem has been solved.

Why does the problem occur?

The problem is that, for example the result type of an Eigen multiplication of two matrices is not an Eigen matrix, but some internal object which represents the multiplication and references the two matrices we are trying to multiply.

Hence, if we use the auto keyword the compiler will very likely not give the variable we are setting the type MatrixXd but the type of some internal object.

For more information refer to the Eigen documentation which explicitly states:

In short: do not use the auto keywords with Eigen's expressions, unless you are 100% sure about what you are doing. In particular, do not use the auto keyword as a replacement for a Matrix<> type

How do I prevent it from happening?

  • do not use the auto keyword, use explicit types.
  • for a lambda function always specify the return type: Write auto DF = []() -> MatrixXd {...} instead of auto DF = []() {...}
David
  • 143
  • 7