Re: [boost] Re: multidimensional arrays (was Boost Mathematicians)

From: Matthias Schabel (boost_at_[hidden])
Date: 2004-07-12 12:35:22

> Maybe so, but I wouldn't start there.  The way any particular class
> template is defined should be one of the last concerns in the design
> of any generic framework.  It's all about concepts and algorithms.

I agree in principle; on the other hand, it is often quite difficult to  
tease out appropriately
complete and orthogonal concepts until one has had significant  
experience in practice
with the problem domain - that requires an actual working  
implementation against which
concepts can be tested.

> If by that you mean driving towards something based on the STL
> (reversible) container requirements in tables 65 and 66, then I
> strongly disagree.  Those are bloated concepts and little or no
> generic code has been written to take advantage of them.

I'm not familiar with tables 65 and 66; I just think that the interface  
should support as
much of the existing standard for one-dimensional containers as  
possible as there is
existing code which would presumably work with little modification if  
that is done.

>> so that algorithms which are unconcerned with the order of traversal
>> can use an efficient 1D iterator interface
>
> That, I agree, is important.
>
>>  - for example, if I want
>> to replace every element of a matrix with its square, I should be
>> able to do this :
>>
>> matrix_type mat;
>>
>> for (matrix_type::iterator it=mat.begin();it!=mat.end();++it) *it =
>> (*it)*(*it);
>
> That's potentially rather inefficient with all those .end()
> invocations.  Anyway, I don't think low-level iteration over elements
> should be required by the concepts to be exposed as the primary
> begin() and end() iterators of a matrix, but through some kind of
> indirection.  Since I am predisposed to using a free-function
> interface:
>
>    std::for_each(begin(elements(mat)), end(elements(mat)), _1 = _1*_1);

I'm curious - what's the rationale for preferring begin(elements(mat))  
to mat.begin()  ?
While the syntax you describe looks fine to me, I think it's also  
important to remember that
many numericists for whom this sort of library would be interesting are  
not interested in
advanced coding techniques, per se, and will likely be put off by the  
need to learn lambda
expressions, etc... and will generally prefer to not have to completely  
rework to the logic
of all their code to port it effectively.  There is a significant  
contingent of programmers who are
clinging to FORTRAN because 1) there is a lot of existing code written  
in it, 2) they understand
the syntax, 3) there is a perception that the added language machinery  
in C++ leads to
performance compromises.  If we're hoping to gain broad acceptance, it  
seems to me that
the interface specification may need to be a bit fatter rather than  
lean to soften the transition.

> I have no problem with any of these details being part of any specific
> array class.  Whether or not they should be part of the library's
> concept requirements is another question (on which I haven't formed
> an opinion yet).

That's fine - I'm primarily interested in starting a dialogue.  There  
continues to be a proliferation
of different vector/matrix/tensor/array libraries, each of which is  
tailored to some particular domain
or another.  The existence of a standard interface (or several levels  
of interfaces supporting
various functionality) would go a long way toward facilitating the  
growth of a coherent and
interoperable set of numerical libraries in C++.

>> of arbitrary indices with no need to reorder the underlying storage
>> (that is mat.transpose(0,1) causes all subsequent calls to mat(i,j) to
>> provide the transposed element without moving the elements
>> themselves).
>
> I'm not sure that's a nice property, since it means that the
> traversal order of the matrix is not a compile-time property.

This is actually precisely the kind of issue that I'm hoping to  
side-step by allowing flexible implementation -
I imagine different requirements will call for different concepts to be  
supported.  In some applications
having static traversal order may not be as important as being able to  
dynamically reconfigure array access,
whereas in others it may be critical.  However, for another large class  
of algorithms it may well be irrelevant
and, for those applications, it should be an unspecified implementation  
detail.

> Your matrix implementation might be very nice, but from my point of
> view, using any specific implementation as the basis of a new linear
> algebra library design is starting from the wrong end of things... so
> for me, at least, looking at it in detail would be premature.

I understand your point here and agree that algorithms should only  
require that their arguments support the
relevant concepts.  I'm not arguing that my particular implementation  
should be considered as the basis for a
new MTL, just that it would be hugely beneficial if a new linear  
algebra library was based on a well considered
and reasonably mature multidimensional array interface/concept set.  I  
think that it is critically important to
decouple the 2D array aspects of a linear algebra library from the  
actual linear algebra aspects - there are many
algorithms which one might want to apply to a matrix which have nothing  
to do with linear algebra.  It would be
very nice to see Boost put forth at least the start of a solution to  
the ongoing problem of proliferating array and
matrix libraries.

------------------------------------------------------------------------ 
---------------------------
Matthias Schabel, Ph.D.                           	
Utah Center for Advanced Imaging Research
729 Arapeen Drive
Salt Lake City, UT  84108
801-587-9413 (work)
801-585-3592 (fax)
801-706-5760 (cell)
801-484-0811 (home)
mschabel at ucair med utah edu