[petsc-dev] attaching field information to PetscLayout?

[Barry Smith]

   In preparation for adding  field information to PetscLayout (as discussed below in the previous email) I have fixed up PETSc-dev so that VecDuplicate(), MatDuplicate() and MatGetVecs() results in new objects that reference the previous PetscLayouts (rather than produce new ones).  Also, inspired by a comment of Dmitry I have put the ISLocalToGlobalMapping mapping and bmapping objects into the PetscLayout so that each individual Vec and Mat doesn't need to mess with them. Please report any problems and if you find places where new vectors are created directly rather than by duplication etc report this as well. Note that these changes are all good things to do independent of putting field information into the PetscLayout.

   At this point I could add a PetscLayoutSetFields() or similar beasty. Some possibilities

      PetscLayoutAddField(PetscLayout,name?,IS) to allow incrementally putting them in or 
  
      PetscLayoutSetFields(PetscLayout,n,names[]?,IS[]) to set them all at once

     and how would we get them into the PetscLayout that we usually don't access directly. Have a 

      PetscVecSetFields(Vec,n,names[]?,IS[])  and PetscMatSetFields(Mat,n,names[]?,IS[]).  

      One natural place to have set this information up is with the DMs.
          For example DMDA with dof  > 1 would have be default dof stride fields; to allow other possibilities there could be DMDASetFields(). Also Vecs and Mats with bs > 1 could be default get similar bs stride fields by default.
          DMComposite would by default have a field for each of the composed DMs perhaps. 
           Likely each DM object should keep the root PetscLayouts that it uses to create its vectors and matrices so all of its creations share a common layout rather than the current model where each newly created vec or mat from a DM gets a new PetscLayout (hmm, maybe I should fix the first).


     Another issue is that there are other places where field information should eventually be propagated. For example VecGetSubVector() and MatGetSubMatrix/Matrices() new vectors matrices should get this information from their parents (this way nesting of fieldsplit/bjacobi/asm/etc  solvers will become much easier and more automatic). And redoing all the VecStrideXXX() routines to be more general for fields.

   Comments, suggestions, corrections?

   Barry







On May 14, 2011, at 11:51 AM, Barry Smith wrote:

> 
>    I have a proposal for a moderate/major change for how we handle knowledge of subfields/splits in PETSc.
> 
>    Currently one can call PCFieldSplitSetIS() or PCFieldSplitSetFields() to indicate the splitting of vectors(fields) into subvectors(subfields).  But, in fact, the knowledge of what splits make sense really comes in the generation of the vectors and matrices and it is bad that the user needs to stash that information somewhere and then bring it out and attach it to the PCFieldSplit. It is especially annoying if one is using a fieldsplit inside, say a multigrid preconditioner, and one has to somehow get that information down into the inner PC. 
> 
>   I propose an alternative. All Vecs and Mats currently have associated PetscLayouts and when a Vec or Mat is duplicated the new Vec or Mat has that same PetscLayout. I propose each PetscLayout have associated with it an optional set of IS indicating the subfields/subvectors (note that this is sort of already true when you set a blocksize larger than 1 the PCFieldSplit will by default use the strides for subfields). Thus PCFieldSplit can get directly from its matrix the appropriate default splits.   Similarly the VecStrideXXX() operations can be made more general becoming VecSubVecXXX(). We also add VecGetSubVec() to get appropriated sized subvectors easily.
> 
>   The end result will be a relatively small change in the PETSc API and for most users little or absolutely no change, but the power to compose easily will really increase.
> 
>   Comments, objections, improvements?
> 
>    Do this before or after the next release? (I guess the conventional wisdom answer would be after?).
> 
> 
>    Barry
> 
> 
> BTW: It is true that parallel vectors and matrices share the PetscLayout from their parent object (via VecDuplicate/MatDuplicate) but this is not true to sequential Vec and Mat. Before adding Field information to PetscLayouts I would need to change all the sequential duplicates to get a reference to their parents PetscLayout instead of creating a new one.
>