dolfin team mailing list archive

[Johan Hake <johan.hake@xxxxxxxxx>]

Hello!

As pointed out by Garth in a couple of commit messages lately we need to 
consider the slicing interface for vectors and matrices in PyDOLFIN.

Now that the parallel interface are starting to get more and more robust 
(great work Garth!) we need to consider what:

  v[6], v[4:], v[index_array] # returning a scalar and a vector
  A[5,6], A[5:,6], A[5:,6:] # returning a scalar, a vector and a matrix

should return for serial and parallel structures.

I guess Garth comments that the Matrix slicing "is just a toy" has some merits 
as it uses a dense representation of the Matrix (even though it is erasing 
zero entries). Regardless if we decide to discard the Matrix slicing or not we 
need to decide how off process indices should be treated.

As it is now, a user need to check whether an index is owned by the local 
processor or not, by using GenericVector::owns_index, before he access the 
data. 

An other scenario would be to always allow a user to write:

  v[5], v[4:] and v[index_array]

and then return for example None (for single value access) and a reduced 
vector (only including local values) for the two other alternatives. Then the 
user need to deal with the result knowing how we treat this in parallel.

Matrix slicing is more tricky especially now when GenericMatrix::resize is 
removed. There are probably ways of getting around this by building a 
SparsityPattern based on the slice. But I am not sure it makes any senses to 
do a matrix slice for a parallel structure, as it is only the rows that are 
local. Then we need to decide if the Matrix slicing should be:

  1) removed 
  2) reimplemented but disabled in parallel 
  3) reimplemented for both serial and parallel cases using for example 
     a similar approach as sketched above for vectors, that is always
     returning something.

Best regards,

Johan