Re: [eigen] Eigen Types as Parameters for Functions

[ Thread Index | Date Index | More Archives ]

Nice work,

some questions/notes:
* Does A.row(3); change after calling foo1(A.row(3))? (i.e. is it
  copied back?) -- you could simply copy it back in the destructor, if
  match_helper<Derived>::type() is false_type.
  If you think back-copying causes too much performance penalty, you
  need to forbid this.
  Furthermore, it sometimes might be useful to mark Ref as a pure
  "out-ref" in this case, saving the copy from A.row() to m_object.
* If you pass, for example, a Ref<Vector4d>, don't you have a Vector4d
  (the m_object member) lying on the stack then? That might cause
  alignment issues if the stack is not aligned (I did not test it yet).


On 27.06.2012 10:25, Gael Guennebaud wrote:
Hi all,

I've implemented a solution fully described there:

feedback more than welcome!

For the lazy ones, here is a self-explanatory example of the proposed solution:

#include <iostream>
#include <Eigen/Dense>
using namespace Eigen;

// read-write access to 'a'
void foo1(Ref<VectorXf> a) { }

// read-only access to 'a'
void foo2(Ref<const VectorXf> a) { }

int main()
   VectorXf a(10);
   const VectorXf& ac(a);
   VectorBlock<VectorXf> ab(a,0,3);
   MatrixXf A(10,10);
   const VectorBlock<VectorXf> abc(a,0,3);

   //foo1(ac);       // does not compile because ac is const
   foo1(A.row(3));   // copied into a temp because innerstride!=1

   foo2(A*A.col(1)); // evaluated into a temp
   foo2(a+a);        // evaluated into a temp

   return 0;


On Mon, Feb 27, 2012 at 7:44 PM, Gael Guennebaud
<gael.guennebaud@xxxxxxxxx> wrote:

this more or less what Christian suggested but more restrictions because:
1 - the storage order really has to be known at compile time otherwise
we cannot construct an expression
2 - same for the alignment, especially for small fixed sizes.
Otherwise the alignment is already determined at runtime by checking
the pointers, so no need to specify it.
3 - a runtime inner-stride can really kill the performance so it is
important to let the user decide whether he can tolerate that or not.
4 - distinguishing at compile-time between vectors and matrices is
fundamental too

So if we go for a specific class (not Map) it could had 5 template parameters:
- typename scalar
- int rows
- int cols
- int storage_order
- int inner_stride = 1


On Mon, Feb 27, 2012 at 7:34 PM, Gael Guennebaud
<gael.guennebaud@xxxxxxxxx> wrote:
Alright, I thought I had already replied to this thread but I haven't
so here I go.

My general idea to offer the ability to write non template function
while still offering some flexibility on the input is really simple:
simply add constructors to Map<> that can take Eigen's objects such as
Matrix, Block, Map, etc. Then it is up to user to specify the level of
flexibility he wants for his arguments: dynamic/fixed sizes,
outer-stride? inner-stride? matrix/vector?

For instance:

void foo(const Map<MatrixXd,0,OuterStride<> >& arg);

should accept any MatrixXd, Matrix4d, Block, Map, etc.

There are some limitations though:
- No template is possible at all, in particular it is not possible to
templatize the scalar type, but this is fine to me since if the user
can tolerate some template parameter then why not allowing everything.
- The storage order has to be fixed too. I understand a runtime
storage order might be convenient, but this is really not the Eigen
philosophy, so this limitation will have to stay.

However, I have question regarding the handling of non compatible
expressions. For instance what should we do with:

foo(A+B);  // this is a pure expression without direct access


MatrixXcf C;
foo(C.real());  // C.real() has an inner-stride, but foo does not accept that

If we stick with Map<> I think that there is only one possible
solution: compilation error. Then it is up to the user to explicitly
evaluate it.

Another approach would be to add a variant of Map tailored for this
use case that would allow to silently evaluate the argument into a
temporary. This is what happen if the user declare his argument as a
const MatrixXd&.

Any opinion on which approach we should follow?


On Thu, Jan 26, 2012 at 3:12 PM, Christoph Hertzberg
<chtz@xxxxxxxxxxxxxxxxxxxxxxxx> wrote:
On 26.01.2012 14:50, Benoit Jacob wrote:

2012/1/26 Christoph Hertzberg<chtz@xxxxxxxxxxxxxxxxxxxxxxxx>:

On 26.01.2012 14:33, Benoit Jacob wrote:

2012/1/26 Christian Seiler<christian@xxxxxxxx>:

  - subclass of MatrixBase, i.e. can be used inside the function like a
   Matrix object
  - has private members containing the above-mentioned bits of
   information (pointer to first element, dimension, ...)
  - can be implicitly converetd to from any compatible MatrixBase (i.e.
   foo(mat) and foo(mat.block(...)) can be used directly)

Just note that this is only possible for those matrix expressions
whose coefficients are stored in memory with a simple enough layout
that can be described in terms of strides. That is what we call a
"Direct Access" matrix expression; the test is ExpressionType::Flags&

Yes, and I would say this would be the main application, if you have
complicated-enough functions (i.e. where you don't really bother about
evaluating expressions to a temporary if required).
And more importantly, when having parameters that are actually output
parameters, you are more or less stuck to direct access types anyway --
this is a very important point, where currently const-correctness is
(almost) entirely.

With that said, yes, we have been considering doing this, that's bug 58:

If someone wants to do this and needs mentoring, I should find time to
do at least that mentoring.

I'd be interested at least to join the "specification committee" ;)

Really, I haven't thought about this in a long time so your thoughts
on this subject are probably a lot fresher than mine. From the top of
my head, the starting point is to look at MapBase (which is inherited
both by Map and by Block<direct access type>) and ask: how much of the
templated stuff there can we replace by runtime variables (data
members), and try to make all direct-access types inherit the
resulting base class (in particular, Matrix should inherit it).

Yes, my thoughts were about the same. Actually, the only difference between
a Map/Block and a Matrix should be that the Matrix "owns" its memory, i.e.
has to deallocate it on destruction (with a special case of fixed size
matrices having static memory).

Furthermore, if one allows strides in the Matrix-class, some specialized
types, such as AlignedVector3 would become just a typedef. Currently, if
someone needs an aligned vector7, he'd have to to the same thing again --
also it might be interesting to store an aligned 3x3 matrix as the top 3
rows of a 4x3 matrix allowing a much more efficient matrix-vector product.


Dipl.-Inf. Christoph Hertzberg
Cartesium 0.051
Universität Bremen
Enrique-Schmidt-Straße 5
28359 Bremen

Tel: (+49) 421-218-64252

Dipl.-Inf. Christoph Hertzberg
Cartesium 0.049
Universität Bremen
Enrique-Schmidt-Straße 5
28359 Bremen

Tel: +49 (421) 218-64252

Mail converted by MHonArc 2.6.19+