[petsc-dev] DM field interfaces

[Barry Smith]

  See bottom.

On Feb 26, 2012, at 9:33 PM, Jed Brown wrote:

> On Sat, Feb 25, 2012 at 15:27, Barry Smith <bsmith at mcs.anl.gov> wrote:
> 
>   Let's back up a step before everyone gets totally confused (and worry less about exact syntax to use instead focus on concepts)
> 
> Proposal 1: Jed
> 
>    // Conservative gas dynamics with fields ordered as [rho, rho*u, rho*v, rho*w, E]
>   DMGetFieldSplitting(dm,"conservative-vector",&fscons); // creates if it doesn't exist
>   DMFieldSplittingAddField(dm,fscons,"rho",{0});
>   DMFieldSplittingAddField(dm,fscons,"momentum",{1,2,3});
> 
>    This seems to be premised on having an underlying "canonical" set of fields that one can define unions of to make new fields. For example with DMDA the underlying
>    fields are simply the 0th, 1st, 2nd etc DOF at each node. For a simple staggered grid 0 canonical may be the pressure centers, 1 the x velocities, 2 the y velocities
> 
> Alright, lets abandon this "canonical field" idea, it was broken to begin with and less flexible, as you say.
>  
> 
> Proposal 2: Barry
> 
>     Use IS's to define the fields, not unions of canonical fields. Reason: more general than proposal 1
> 
> Great, my concern is that we would like DM to be able to declare multiple named partitions. For example, so that we can have a conservative-scalar splitting, a conservative vector splitting (keeps vector components together), primitive variable splittings, local-characteristic splittings (analytically diagonalize the flux Jacobian for hyperbolic problems), etc.
> 
> We can think of these transformations as symmetric preconditioning. Suppose that we want to solve with J, operating on conservative variables. Suppose P converts conservative to primitive variables and C = P^{-1} converts primitive to conservative. I'm using matrix notation for because this would often be applied to perturbations. We need to solve a linear system with J_c, the Jacobian acting on conservative variables, but the preconditioner wants to operate in primitive variables, so the user assembles
> 
> J_p = P J_c C.
> 
> (Usually they do this only approximately and it's done without looking at J_c.) The conservative system is
> 
> J_c = C J_p P
> 
> At the algebraic level, I think this logic belongs in a new preconditioner that applies a user-defined pre- and post-transformation, the concern is how to forward DM-based information into PCFieldSplit and PCMG operating on J_p.
>  
>  
>  Because Jed's proposal 1 violates the use of IS for indexing paradigm and is less general than Proposal 2 I'd like to drop that Proposal 1; but it does bring up an issue. I think we need another IS implementation that is more "general" than stride for handling sets of DOFs of freedom associated with those canonical sets of fields. So it is easy to access canonical fields that are not simply stride. For example ISCreateCanonical(MPI_Comm,PetscInt ncanonical,PetscInt *fields), we can even have "predefined" ones available for use such as IS_CANONICAL_(MPI_Comm,3,{1,2,3}) Thus we can preserve some of the simplicity of proposal 1.  Jed's (again not worrying about syntax)
> 
> DMFieldSplittingAddField(dm,fscons,"momentum",{1,2,3});  becomes
> 
> DMFieldSplittingAddField(dm,fscons,"momentum",IS_CANONICAL_(comm,3,{1,2,3})
> 
> Note that for a Vec with a given block size the canonical fields are simply the stride fields so IS_CANONICAL_(comm,1,{start}) is pretty much ISStride().
> 
> Ignoring syntax, I think this is fine. 
> 
> Proposal 3: Jed
> 
>    Observation: Using IS does not give all desired generality since the new field may be a linear combination of DOFs. Or, yikes, a nonlinear transformation of DOF.
>    Thus we may need something like DMTransformVec() to do these transformations.
> 
> 
> So this moves us beyond the traditional block linear algebra methods but if we could solve it we'd have some very nice new powerful capabilities on our hands.
> 
> Question 1:  Do we need to solve this for the current generation of PCFieldsplit?
> Question 2:  Do we want to solve it in PCFieldSplit or should that come later in a different PC, PCBasisChangeThingy or something?
> 
> Yes, but what DM should PCBasisChangeThingy forward to the inner PC?

   Presumably a "new" DM that represents the new set of equations ... and to start we will need to do very specific examples before coming up with a general paradigm.

   But it seems best to me to make this orthogonal to PCFieldSplit.

   Barry