Hi there,

This has been puzzling me for some time and I'm guessing it's a bit of a newbie question.

I'm often writing functions that take as input several MatrixXf and return several MatrixXf. Until now I have been doing that the following way (simplified here for the example):

Contents of addsub.h

Code: Select all

#ifndef ADDSUB_H
#define ADDSUB_H

#include
#include
#include

// import most common Eigen types
USING_PART_OF_NAMESPACE_EIGEN

void addSub(MatrixXf a, MatrixXf b, MatrixXf* c, MatrixXf* d) {
  if(a.cols() != b.cols()) {
     // Throw a ValueError
  }

  if(a.rows() != b.rows()) {
     // Throw a ValueError
  }
 
  // Resize the outputs to fit the input
  (*c).resize(a.rows(), a.cols());
  (*d).resize(a.rows(), a.cols());
 
  (*c) = a - b;
  (*d) = a + b;

  return;
}

#endif  /* ADDSUB_H */


Contents of my_program.cpp:

Code: Select all

#include "addsub.h"

int main(int, char *[])
{
  MatrixXf a = MatrixXf::Constant(3,3, 19);
  MatrixXf b = MatrixXf::Constant(3,3, 14);
 
  MatrixXf c, d;
 
  addSub(a, b, &c, &d);
 
  std::cout << "cn" << c << std::endl;
  std::cout << "dn" << d << std::endl;

  // This is what I would like it to work with
  /*
  a.setConstant(12);
  b.setConstant(3);

  addSub(a.corner(Eigen::TopLeft, 2, 2),
         b.corner(Eigen::TopLeft, 2, 2),
         &(c.corner(Eigen::TopLeft, 2, 2)),
         &(d.corner(Eigen::TopLeft, 2, 2)));

  std::cout << "cn" << c << std::endl;
  std::cout << "dn" << d << std::endl;
  */
}


And here is the command to compile I use:

Code: Select all

g++ my_program.cpp -o my_program -I/path/to/eigen2


I would appreciate a lot any help to making work the last part (commented out) in the program. The reason I want to do this, is that I want my function to be able to work on both:
- uninitialized matrices (by resizing them and therefore allocating itself the memory)
- initialized matrices (where the memory has been previously allocated (like when working on buffers))

Thanks in advance.
ricard

rikrd wrote:

Code: Select all

void addSub(MatrixXf a, MatrixXf b, MatrixXf* c, MatrixXf* d) {


Here your addSub function takes a and b by value which means that when you call addSub, copies of a and b are made and addSub works on these copies. That's useless, inefficient. Moreover it would give you bad errors if instead of MatrixXd you had, say, Matrix4d. So you should rewrite this as:

Code: Select all

void addSub(const MatrixXf& a, const MatrixXf& b, MatrixXf* c, MatrixXf* d) {


and read this page: http://eigen.tuxfamily.org/dox/PassingByValue.html

rikrd wrote:

Code: Select all

  // Resize the outputs to fit the input
  (*c).resize(a.rows(), a.cols());
  (*d).resize(a.rows(), a.cols());
 
  (*c) = a - b;
  (*d) = a + b;


In order to make this work with both initialized and non-initialized matrices,you need to use set() like this:

Code: Select all

  c->set(a - b);
  d->set(a + b);


Notice set() takes care of resizing for you. Also notice that in C/C++ instead of (*c). you can just do c->.

Contents of my_program.cpp:

rikrd wrote:

Code: Select all

  // This is what I would like it to work with
  /*
  a.setConstant(12);
  b.setConstant(3);

  addSub(a.corner(Eigen::TopLeft, 2, 2),
         b.corner(Eigen::TopLeft, 2, 2),
         &(c.corner(Eigen::TopLeft, 2, 2)),
         &(d.corner(Eigen::TopLeft, 2, 2)));


Aaaah but then what you are passing to addSub aren't matrices, these are expressions. So you need to write addSub like that:

Code: Select all

template
void addSub(const MatrixBase& a, const MatrixBase& b,
 MatrixBase *c, MatrixBase *d)
{
  ...
}


However if addSub is declared like that, you can no longer call set() and resize() on c and d, because these methods are only of plain matrices and now c and d are expressions. So you'd have to initialize c and d with the correct size before calling addSub.

rikrd wrote:I would appreciate a lot any help to making work the last part (commented out) in the program. The reason I want to do this, is that I want my function to be able to work on both:
- uninitialized matrices (by resizing them and therefore allocating itself the memory)
- initialized matrices (where the memory has been previously allocated (like when working on buffers))

You have to make a compromise here:
- if you do like the first version i give you, then it works on both initialized and uninitialized matrices, but these are plain matrices not expressions so you can't let c and d be corner's.
- if you do like the 2nd version, c and d are expressions so they can be corner's but the thing is, even before you can call corner() on a matrix, this matrix has to be initialized.

Great! the explanation was fantastic and it helped me better understand how to use Eigen. Now I remade the example in order to fit my needs. But I get a compile warning which I don't understand, maybe you can help me get rid of it:

First of all the warning:

Code: Select all

my_program.cpp: In function ‘int main(int, char**)’:
my_program.cpp:24: warning: taking address of temporary
my_program.cpp:25: warning: taking address of temporary


Contents of my_program.cpp:

Code: Select all

#include "addsub.h"

using namespace std;

int main(int, char *[])
{
  MatrixXf a = MatrixXf::Constant(3,3, 19);
  MatrixXf b = MatrixXf::Constant(3,3, 14);
 
  MatrixXf c, d;
 
  addSub(a, b, &c, &d);
 
  cout
#include
#include

// import most common Eigen types
USING_PART_OF_NAMESPACE_EIGEN

template
void addSubBuffer(const MatrixBase& a, const MatrixBase& b,
                  MatrixBase *c, MatrixBase *d) {

  const int ncols = a.cols();
  const int nrows = a.rows();

  if(b.cols() != ncols || b.rows() != nrows) {
     // Throw a ValueError
  }

  if(c->cols() != ncols || c->rows() != nrows) {
     // Throw a ValueError
  }

  if(d->cols() != ncols || d->rows() != nrows) {
     // Throw a ValueError
  }
 
  (*c) = a - b;
  (*d) = a + b;

  return;
}


void addSub(const MatrixXf& a, const MatrixXf& b, MatrixXf* c, MatrixXf* d) {
  // Resize the outputs to fit the input
  c->resize(a.rows(), a.cols());
  d->resize(a.rows(), a.cols());
 
  addSubBuffer(a, b, c, d);

  return;
}

#endif  /* ADDSUB_H */


And to compile:

Code: Select all

g++ -o my_program my_program.cpp -Ipath/to/eigen


I think I've addressed all your comments Benoit.

PS: is there a better (as in more often used) name than *Buffer for these kind of methods that work on a preallocated matrix?

EDIT: Fix that I had forgotten the const MatrixXf& in the addSub() parameters.

Your compiler doesn't like these lines:

rikrd wrote:

Code: Select all

               &(c.corner(Eigen::TopLeft, 2, 2)),
               &(d.corner(Eigen::TopLeft, 2, 2)));


You can fix that by letting this function take references instead of pointers:

Code: Select all

template
void addSubBuffer(const MatrixBase& a, const MatrixBase& b,
                  MatrixBase& c, MatrixBase& d);


about the naming question, I would call that function matrixAddSub or something like that.

I hope not to be abusing from this forum, but just for the record...

With the help of ggael and bjacob I finally got this to work by turning into const references the output arguments of addSubBuffer() and using const_cast inside the function to actually modify them. This is not ideal, but it permits to do the calls I wanted to allow (calling a function on a part (Block) of a Matrix as well as on the whole Matrix).

Actually this is far from ideal, because to a user of the function it will be confusing that the const reference will not be treated as const and will actually be modified (it has to be clear in the doc!).

Contents of addsub.h:

Code: Select all

#ifndef ADDSUB_H
#define ADDSUB_H

#include
#include
#include

// import most common Eigen types
USING_PART_OF_NAMESPACE_EIGEN

/**
 *
 * The outputs a out_c and a out_d will be filled with the result
 * of the sum and the difference respectively of a a and a b.
 *
 */
template
void addSubBuffer(const MatrixBase& a, const MatrixBase& b,
                  const MatrixBase& out_c, const MatrixBase& out_d) {

  const int ncols = a.cols();
  const int nrows = a.rows();

  if(b.cols() != ncols || b.rows() != nrows) {
     // Throw a ValueError
  }

  if(out_c.cols() != ncols || out_c.rows() != nrows) {
     // Throw a ValueError
  }

  if(out_d.cols() != ncols || out_d.rows() != nrows) {
     // Throw a ValueError
  }
 
  MatrixBase& nonConst_out_c = const_cast&>(out_c);
  MatrixBase& nonConst_out_d = const_cast&>(out_d);

  nonConst_out_c = a - b;
  nonConst_out_d = a + b;

  return;
}


void addSub(const MatrixXf& a, const MatrixXf& b, MatrixXf& c, MatrixXf& d) {
  // Resize the outputs to fit the input
  c.resize(a.rows(), a.cols());
  d.resize(a.rows(), a.cols());
 
  addSubBuffer(a, b, c, d);

  return;
}

#endif  /* ADDSUB_H */


contents of my_program.cpp:

Code: Select all

#include "addsub.h"

using namespace std;

int main(int, char *[])
{
  MatrixXf a = MatrixXf::Constant(3, 3, 19);
  MatrixXf b = MatrixXf::Constant(3, 3, 14);
 
  MatrixXf c, d;
 
  addSub(a, b, c, d);
 
  cout << "cn" << c << endl;
  cout << "dn" << d << endl;

  // This is what I would like it to work with
 
  a.setConstant( 12 );
  b.setConstant( 3 );

  addSubBuffer(a.corner(Eigen::TopLeft, 2, 2),
                       b.corner(Eigen::TopLeft, 2, 2),
                       c.corner(Eigen::TopLeft, 2, 2),
                       d.corner(Eigen::TopLeft, 2, 2));

  cout << "cn" << c << endl;
  cout << "dn" << d << endl;
}



Hope this helps the next one running into this problem.

ricard