Re: slicing and ranges of matrices

From: jhrwalter (walter_at_[hidden])
Date: 2002-01-19 10:20:55

--- In boost_at_y..., Toon Knapen <toon.knapen_at_s...> wrote:
> Looking at the slicing and ranging possibilities of both 
multi_array and 
> ublas, I find the concepts as they are currently defined 
overlapping and 
> limited in respect to e.g. Blitz++.
> 
> The most important difference with Blitz++ is that slicing in 
Blitz++ 
> speak means also reducing the order of the array. In ublas & 
> multi_array, slicing means a 'strided range'. I think it is 
important to 
> define concepts in conformance with Blitz++ as it's still one of 
the 
> default array libraries out there.

A slice in ublas is conforming to the valarray slice ;-). A range is 
the specialization of a slice with stride = 1. As we intentionally 
chose to distinguish vectors and matrices (and *not* to model vectors 
as m x 1 or 1 x n matrices) we needed to introduce different proxy 
classes for order reduction.
 
To cite the Blitz++ documentation w.r.t. it's current 
status: ,Blitz++ currently provides dense vectors and 
multidimensional arrays. Some other features (matrices, random number 
generators) are partially implemented.'
 
> But another concept I'm missing, and which is in neither of the 
three 
> libraries, is creating a sub-array (of the same order as the array) 
that 
> contains a specific selection of indices for every dimension of the 
> array. e.g. in the 2D case I would like to do
> matrix< double > A(8,8);
> boost::array<3> i1 = { 0,1,6 };
> boost::array<5> i2 = { 5,6,7 };
> sub_matrix sub_A = submatrix(A, i1,i2 );
> 
> where sub_A contains all element that are in columns 5,6 and 7 and 
in 
> rows 0,1 and 6. The difference with a range is thus that this sub-
array 
> is not a continuous block in the original array. The diff with the 
> current definition of slicing is that the selection does not need 
to be 
> strided.
> 
> Is there a reason this kind of subarray does not exist ? 
 
Interesting question. Should this sub matrix support dense matrices 
only or sparse matrices also (the old question of genericity ;-)?
 
> With this kind 
> of functionality though, I would never have to write loops over my 
> matrices again as I currently rely on ublas to do this for me. 
> Currently, my only loop is for assembling my FEM matrix from all 
local 
> element matrices which comes down basically to assigning a sub-
array 
> with non-strided indices.
 
Sounds rather cool. Would you like to post some sample code?
 
Regards
 
Joerg