Re: [eigen] Iterators with dense matrices: (was: Areas of eigen3 variation and using eigen2/3 in generic code)

[Benoit Jacob <jacob.benoit.1@xxxxxxxxx>]

2010/1/4 Hauke Heibel <hauke.heibel@xxxxxxxxxxxxxx>:
> On Mon, Jan 4, 2010 at 1:03 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
> wrote:
>>
>> 2010/1/3 Jesse Perla <jesseperla@xxxxxxxxx>:
>> > On Sun, Jan 3, 2010 at 8:59 PM, Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
>> > wrote:
>> >>
>> >> 2010/1/3 Paul C. Leopardi <paul.leopardi@xxxxxxxxxxxx>:
>> >>
>> >> So Eigen's approach has been to let the user write higher-level code in
>> >> terms of /expressions/ and /functors/. The user can provide his own
>> >> expressions and functors, which describe how coefficient(i,j) is
>> >> computed,
>> >> optionnally providing a SIMD implementation for which we have a
>> >> portable
>> >> infrastructure; and then the user uses that with Eigen's built-in
>> >> loops.
>> >
>> > Definitely the right approach.  The higher level and closer to linear
>> > algebra the better... But for cases where you have not implemented an
>> > algorithm and default generic algorithms exist, it is necessary to be
>> > able
>> > to use them without initially implementing extensions to Eigen.
>> > For example:  is there an Eigen equivalent to the following things I do
>> > all
>> > day long:
>> > std::vector<double> vec(4);
>> > //... fill it
>> > std::vector<double> vec2(4);
>> > std::transform(vec.begin(), vec.end(), vec2.begin(), [](double x){return
>> > my_func(x)}); //for some function
>>
>> For this, yes: first you would wrap my_func into a functor, and then you'd
>> call:
>>  vec2 = vec.unaryExpr(my_functor());
>>
>> Again, the advantages are vectorization and unrolling (we make it easy
>> to SIMD-enable your functor in a portable way).
>>
>> See:
>>
>> http://eigen.tuxfamily.org/dox/classEigen_1_1MatrixBase.html#a1d1835f182c0141dc073a1045c8a2e9e
>>
>> > //or do a binary search from the 2nd item for the number 1.0, with it
>> > sorted
>> > with some functor
>> > std::binary_search(vec.begin() + 1, vec.end(), 1.0,
>> > my_strict_weak_ordering_functor());
>>
>> OK, for this one, we don't provide anything, so indeed using STL
>> algorithms is useful.
>>
>> One remark: you could work by directly addressing the underlying array.
>>  vector.data()  returns a pointer to the vector's coeffs
>>  vector.size()  returns the number of coeffs
>>
>> STL algorithms can be applied directly to such an array, right?
>
> Totally right...
>
>>
>> OK, so, I'm not against it and it does seem useful. My next question
>> is: is it important that the methods returning begin and end iterators
>> are named begin(void) and end(void) ?
>>
>> Everybody: how about:
>>  start ---> head
>>  end ---> tail
>> Keep start and end(int) for compatibility in Eigen2Support
>> Remove end<int>() altogether
>
> That sounds like a reasonable change. To me the naming sounds even more
> clear as to what the underlying function does.
>
> Regarding iterators in general, I think we have to keep in mind that
> introducing them is not as trivial as renaming functions an returning the
> addresses of matrix/vector elements. If we want to support sparse and dense
> matrices/vectors in the same way, we really need an iterator class - correct
> me if I am wrong. Personally, I am already fine with coeffRef() or data()
> and size() functions. In the end we are dealing with linear algebra
> datatypes and not with containers.
>
> If we would take the iterator route, sooner or later people will want
> iterators for matrices as well. Or for blocks. Don't get me wrong - I am not
> completely against it but I think it has a lot of consequences and we will
> not be able to keep the performance guarantees we have right now.

Iterators would have to be indices-based, not address-based, in order
to work with non-DirectAccess expressions. So a Vector iterator would
store one m_index member; a general Matrix iterator would store m_row
and m_col; and one big weakness inherent to such an iterator API is
that we probably can't automatically take advantage of LinearAccess
when it is available, though we can still provide an alternate method
providing an index-based iterator making it the responsibility of the
user to only call it on a LinearAccess expression.

Let's be 100% clear that iterators for dense matrices/vectors would be
provided only for the user's convenience, and for the relatively rare
cases where Eigen's built-in loops don't provide what is needed.

I'm not proposing to rewrite all of Eigen to use iterator-based loops
everywhere possible!

My point is that, aside from the possible start()/end() API change,
this shouldn't affect at all people who don't need iterators.

Benoit