Re: [eigen] The "nested" concept

[Thomas Capricelli <orzel@xxxxxxxxxxxxxxx>]

Guys! This was very interesting. Typically the kind of thing I would see fit in the 'developers' documentation'!

Thomas
-- 
Thomas Capricelli <orzel@xxxxxxxxxxxxxxx>
http://www.freehackers.org/thomas

In data martedì 29 giugno 2010 23:27:45, Gael Guennebaud ha scritto:
> Hi,
> 
> let's pick a simple example:
> 
> MatrixXf A, B, C;
> C = A.transpose() + B;
> 
> here A.transpose() returns an object of type:
> 
> Transpose<MatrixXf>
> 
> where the MatrixXf A is nested into the Transpose expression. The
> nested type tells how to store the nested object. Here
> MatrixXf::Nested boils down to a MatrixXf&, and thus "A.transpose()"
> stores a reference to A.
> 
> Now, you might wonder why we introduced such a nesting type mechanism
> and not always use a reference. There are two main reasons:
> 
> 1 - Unlike Matrix or Array, other expressions are lightweight and
> better nested by value. In the previous example, "A.transpose() + B"
> returns an object of type CwiseBinaryOp<ei_scalar_sum_op<float>,
> Transpose<MatrixXf>, MatrixXf>  storing both sides of the addition as
> follow:
> 
> Transpose<MatrixXf>::Nested lhs; // left hand side
> MatrixXf::Nested rhs; // right hand side
> 
> which boils down to:
> 
> const Transpose<MatrixXf> lhs; // nesting by value
> const MatrixXf& rhs; // nesting by reference
> 
> Nesting by value small object avoids temporary headache when a
> function has to return complex expressions, e.g.:
> 
> template<typename A, typename B>
> CwiseBinaryOp<ei_scalar_sum_op<float>, Transpose<A>, B> adjoint(const
> A& a, const B& b)
> {
>   return a.transpose() + b;
> }
> 
> If the temporary "a.transpose()" was stored by reference by the
> CwiseBinaryOp expression, then you would end up with a segfault
> because the "a.transpose()" temporary is destroyed just before the
> function returns, and so the CwiseBinaryOp expression would store a
> reference to dead object...
> 
> 2 - Another reason is to allow to handle expressions which has to be
> evaluated into a temporary before being used. For instance, in the
> following example:
> 
> d = a * b + c;
> 
> for performance reason, the matrix product a * b has to be evaluated
> into a temporary before evaluating the addition. This is achieved as
> follow. Here we build the expression of type:
> 
> CwiseBinaryOp<ei_scalar_sum_op<float>, Product<MatrixXf,MatrixXf>, MatrixXf>
> 
> which stores a Product<MatrixXf,MatrixXf>::Nested for the lhs, and
> Product<MatrixXf,MatrixXf>::Nested is ... MatrixXf !!
> 
> Something more complicated:
> 
> (a + b) * c
> 
> Here, if c is not too small, it is better to evaluate (a+b) into a
> temporary before doing the matrix product, otherwise, a+b would be
> computed c.cols() times... To this end we have a ei_nested<> helper
> class to determine the ideal nesting type. In Product, we have
> something like:
> 
> ei_nested<CwiseBinaryOp<ei_scalar_sum_op<float>, type_of_a,
> type_of_b>, type_of_c::ColsAtCompileTime>::type
> 
> giving us the nesting type of the left hand side of the product (here
> a MatrixXf if a, b, and c are MatrixXf). For the right hand side here
> we have:
> 
> ei_nested<CwiseBinaryOp<ei_scalar_sum_op<float>, type_of_c,
> type_of_a_plus_b::RowsAtCompileTime>::type
> 
> which gives us a MatrixXf&.
> 
> ei_nested<> determines whether the nested expression has to be
> evaluated or not in function of an estimation of the evaluation cost
> of one coefficient. This cost is automatically computed by the
> expressions in the ei_traits<> specializations.
> 
> Yes, this is a bit complicated, but what is very important here, and
> that we have to better document, is that when you write a generic
> function taking, e.g., a MatrixBase<Derived> object you should really
> honor the nesting type of the Derived class:
> 
> template<typename Derived>
> void foo(const MatrixBase<Derived>& _x)
> {
>    typename Derived::Nested x(_x.derived());
> 
>   // use x safely
> }
> 
> Actually, is you use the argument more than once, you should even use
> the ei_nested<> helper:
> 
> template<typename Derived>
> typename Derived::Scalar foo(const MatrixBase<Derived>& _x)
> {
>    typename ei_nested<Derived,2>::type x(_x.derived());
> 
>    return (x + x.adjoint()).maxCoeff();
> }
> 
> if you don't do so, and call:
> 
> foo(a*b);
> 
> then the expensive product a*b will be computed twice !!
> 
> gael
> 
> On Tue, Jun 29, 2010 at 9:29 PM, Manoj Rajagopalan <rmanoj@xxxxxxxxx> wrote:
> > Hi,
> >
> >  I see typedefs with names like "MatrixTypeNested". What is this nesting
> > concept?
> >
> > thanks,
> > Manoj
> >
> >
> >
> 
> 
>