Re: [eigen] about .lazy()

[ Thread Index | Date Index | More Archives ]

On Fri, 14 Aug 2009, Gael Guennebaud wrote:

I know this is not a very nice situation, and my suggestion to solve this mess is to remove the .lazy() function.

That looks like a good idea. Not sure whether the rest of my email is useful ....

c - [ .lazy() ] covers two different features at once:

 c1 - it means that the result does not alias with the operands of
the product, but for that purpose it makes more sense to control that
via a special operator=, like res.noalias() = ...

 c2 - it also means that the product should not be evaluated
immediately, but evaluated as a standard expression. However, in
practice it is (almost) never a good idea to do so, and when it is not
the case the speed difference is negligible.

I agree that these need to be de-coupled. Apart from that they're separate things, there is an important difference. The first case (no aliasing) is something that the compiler cannot find out by itself, so we need to have a function for that. In constrast, the second case is something that the library / compiler writers probably know better then the user, so ideally the compiler should decide here.

Incidentally, I prefer Boost's syntax "noalias(res) = ..."

For large matrices, my last statement is obviously true.

For large *dense* matrices, yes. Probably not for large sparse matrices, though this probably does not make a difference for this discussion.

So if you wonder what happens for small matrices, here is a benchmark for small fixed and dynamic sizes matrices which evaluates D = C + A*B; using three different strategies:

The obvious thing to note is that Eigen 2.1 is doing so much better than Eigen 2.0, a testament to the work that you put into it. Thank you!

Finally, there is also the question whether operator+= and -= should
be "no-alias" by default ? This because I think that in 99% of the
case, when you write:

m += <product>

it's very unlikely that m is one of the operand of the product. The
drawback is that might be confusing for the user, (because operator=
and operator+= would behave differently wrt aliasing).

If you make operator+= no-alias by default, then A = A + A*A would do what you expect but A += A*A would not, correct? That would be very confusing. I don't like it; I'd always expect A = A + stuff to have the same result (though perhaps arrived at less efficiently) as A += stuff.

Just to echo Benoit who wrote that predictability is important here ...

My main issue with using C++ libraries for linear algebra is that a lot is going on under the hood. Things like how expressions are evaluated and which expressions lead to the introduction of temporaries are important for performance and it's difficult to find out exactly what's going on (partially because I'm still learning about the advanced usage of templates that's so important for eigen). I understand that this is part of the price one pays for using a high-level library, but for me it's important that this is well-documented.


Mail converted by MHonArc 2.6.19+