February 2012 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] DM field interfaces
by Barry Smith 25 Feb '12

25 Feb '12

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 Proposal 2: Barry Use IS's to define the fields, not unions of canonical fields. Reason: more general than proposal 1 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. Comment: Barry I very much like the universality in PETSc of using IS __ALWAYS__ to grab a set of indexed things and not using another indexing technique like 0 to indicate the "first field" and 1 to indicate the "second field". 1) Matt seems to like this too since he very recently lobbied for simplifying PCFIELDSPLIT to ONLY use the IS for indicating parts and to remove the (somewhat redundant) PetscInt nfields; PetscInt *fields; which we will do. 2) We currently have VecStrideXXXXX(vec,start,....) that index using the 0th field, 1st field etc Does not match the paradigm 3) VecGetSubVec(Vec,IS,Vec*) uses IS to pull out pieces so matches the paradigm. Way to go Jed! 4) There may be other places that do not match the paradigm 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(). We could push this even further and assume that the underlying object being indexed in has "named" fields, we could still use those named fields to index via the IS if we had and ISCreateNamed(comm,"Name") and helper IS_NAMED_(comm,"Name") so one could do VecGetSubVec(vec,IS_NAMED_(comm,"T"),&subvec); etc. Ok, so this solves the original issue of indicating sets of indices for DM, useful for standard use of block methods for matrices :-). It allows general IS, canonical fields and named fields all under the same interface. ---------------- Now to the hard stuff Jed dumped on us. 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? Question 3: Can we go ahead with deciding on exact syntax etc for Proposal 2 and implementing before we resolve the Proposal 3 issue? When I see DMTransformVec() I get worried, we already have DMGetInterpolation/Restriction() that maps between two DMs, cannot we not somehow extend that instead of putting in more methods? (I know they return Mat which is a linear operator and if I had listened to Matt many years ago everything in PETSc would be a nonlinear operator so everything would be fine). If at all possibly I'd to reuse/fix up those old interfaces to handle this need. Barry On Feb 25, 2012, at 1:33 PM, Jed Brown wrote: > On Sat, Feb 25, 2012 at 13:10, Dmitry Karpeev <karpeev(a)mcs.anl.gov> wrote: > I wrote this before the thread moved to petsc-dev, so there is some overlap with the below. > A given problem with the same "topology" -- a mesh and a total index layout -- might define several useful splits -- collections of fields. Each field might have further useful splits. I have useful examples from my interaction with Moose/libMesh, although these are rather general. A field is now represented by an IS. We could try to put netsted split definitions into ISs (e.g., ISNest), but I think it's much more natural to keep all nestedness inside DMs -- we already have examples of this in the coarsening, refinement for MG and a pretty well understood interaction of this sort of DM hierarchy. with a hierarchy of PCs. > > Thus, I think it would be useful for a DM to return a split as a list of DMs, which can be further split, if necessary: > DMGetFieldDMs(DM dm, PetscInt *fieldcount, DM **fielddms); > The field names can simply be the corresponding DM names, etc. The fielddms can be lightweight, sharing the underlying > topology. That's up to the implementation to arrange, though. This addresses nested splits. > > There are two concepts, defining what the space for a given split _is_ and how to get the part of the current state is needed to get a vector living in that space. The IS is a rudimentary way to define what part is needed. (You can use VecGetSubVector(X,isfield,&Xfield) or other methods to reference that part of the state.) Unfortunately, it's not sufficient even for a linear change of basis. An IS combined with a PF is sufficient for non-mixed spaces. Mixed/staggered spaces need even more because the transformation is defined through quadrature and projection. Maybe the transformation should be entirely handled by the DM instead of giving the user objects that apply the transformation. > > But is the interface going to be DMTransformVec(dm,X,Y,dmt,Yt) that transforms the perturbation Y to the state X into transformed variables according to dmt? Or should we be getting back an object that does this sort of transformation, so that we can amortize some setup coming from matching dmt relative to dm? (Note that we also need an inverse transformation.) > > > We could use the same mehanism to implement alternative splits on the same topology -- these are different DMs sharing the topology internally. Instead of pushing and popping splits, they can be retrieved by name as > DMGetSplitDM(DM dm, const char *splitname, DM *splitdm); > Then calling DMGetFieldDMs on the splitdm retrieves the corresponding fields. > Again, splitdm and dm can share common data structures, as arranged by the particular DM implementation. > The advantage of this approach over push/pop is that we can simultaneously use different splits of the same topology. > > It is rather easy to apply the same scheme and have DMs serve up (overlapping) domain decompositions useful for (G)ASM. > I'm now implementing both of these schemes (fieldsplit and subdomains) on top of libMesh, but the code resides in the libMesh source tree to avoid a circular dependence: as designed right now libMesh depends on PETSc. The API has to live in DM, though. >

1 0

Re: [petsc-dev] DM field interfaces
by Jed Brown 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 13:10, Dmitry Karpeev <karpeev(a)mcs.anl.gov> wrote: > I wrote this before the thread moved to petsc-dev, so there is some > overlap with the below. > A given problem with the same "topology" -- a mesh and a total index > layout -- might define several useful splits -- collections of fields. Each > field might have further useful splits. I have useful examples from my > interaction with Moose/libMesh, although these are rather general. A field > is now represented by an IS. We could try to put netsted split definitions > into ISs (e.g., ISNest), but I think it's much more natural to keep all > nestedness inside DMs -- we already have examples of this in the > coarsening, refinement for MG and a pretty well understood interaction of > this sort of DM hierarchy. with a hierarchy of PCs. > > Thus, I think it would be useful for a DM to return a split as a list of > DMs, which can be further split, if necessary: > DMGetFieldDMs(DM dm, PetscInt *fieldcount, DM **fielddms); > The field names can simply be the corresponding DM names, etc. The > fielddms can be lightweight, sharing the underlying > topology. That's up to the implementation to arrange, though. This > addresses nested splits. > There are two concepts, defining what the space for a given split _is_ and how to get the part of the current state is needed to get a vector living in that space. The IS is a rudimentary way to define what part is needed. (You can use VecGetSubVector(X,isfield,&Xfield) or other methods to reference that part of the state.) Unfortunately, it's not sufficient even for a linear change of basis. An IS combined with a PF is sufficient for non-mixed spaces. Mixed/staggered spaces need even more because the transformation is defined through quadrature and projection. Maybe the transformation should be entirely handled by the DM instead of giving the user objects that apply the transformation. But is the interface going to be DMTransformVec(dm,X,Y,dmt,Yt) that transforms the perturbation Y to the state X into transformed variables according to dmt? Or should we be getting back an object that does this sort of transformation, so that we can amortize some setup coming from matching dmt relative to dm? (Note that we also need an inverse transformation.) > > We could use the same mehanism to implement alternative splits on the same > topology -- these are different DMs sharing the topology internally. > Instead of pushing and popping splits, they can be retrieved by name as > DMGetSplitDM(DM dm, const char *splitname, DM *splitdm); > Then calling DMGetFieldDMs on the splitdm retrieves the corresponding > fields. > Again, splitdm and dm can share common data structures, as arranged by > the particular DM implementation. > The advantage of this approach over push/pop is that we can simultaneously > use different splits of the same topology. > > It is rather easy to apply the same scheme and have DMs serve up > (overlapping) domain decompositions useful for (G)ASM. > I'm now implementing both of these schemes (fieldsplit and subdomains) on > top of libMesh, but the code resides in the libMesh source tree to avoid a > circular dependence: as designed right now libMesh depends on PETSc. The > API has to live in DM, though. >

1 0

Re: [petsc-dev] DM field interfaces
by Jed Brown 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 12:36, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > Jed, It also makes it impossible to handle bizare things where for one > spitting you want a subset of one thing actually in another thing. So why > not use IS to define the names, this gives complete flexibility. Also this > above is not adding fields: you are just giving names to subsets of fields > indicated by the original numbering of fields. Which is fine, but has > confusing names (you are not AddField just naming it). > I'm defining collections of fields, then adding/naming field in those collections. You want these collections defined for visualization anyway. > My suggestion is just a more or less a concise way of doing what you are > doing without the need for YET ANOTHER CLASS (fcons). > That is just a token, it could be an int. There is a very real concept of multiple (possibly nonlinear) bases for model state. How do you want to refer to these collections? Just strings? Instead of making these collections, we could just define an arbitrarily large (thus redundant) list of "observables" (quantities that can be computed from model state). Primitive variables are one such collection of observables. Then when asking for splits, the caller could provide the subset that they wanted. This isn't good for automatically ensuring that the bases span and are non-redundant, but it would be the best interface for visualization. I'm not sure how you would want to define the equation of state (which would normally provide all the primitive variables from the state) or its inverse. > > > > Then PCFieldSplit and eventual nonlinear and time-stepping splitting > methods would call > > > > DMGetFieldSplitting(dm,&nsplits,&splits); > > // choose split number to use > > split = splits[splitidx]; > > DMFieldSplittingGetName(dm,split,&splitname); > > > DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); > > I do not like this approach based on requesting by name, just push your > name and then pop it. This way something "above" in the call stack can push > to the right one then something further down can use it without having to > know its name (sort of like your "default" split) but the push allows > changing the default. You kind of assume that the guts is suppose to know > what name to use? // choose split number to use how is some PC deep inside > suppose to know this? > > I realize my proposal was not complete on detail. I'm suggesting changes > to your detail. > Push/pop is fine with me. But note that the GUI should be able to get a list of options at any level in which a decision makes sense. It's possible that we would want to split in terms of different variables in that PCComposite+2*PCFieldSplit example.

1 0

Re: [petsc-dev] DM field interfaces
by Dmitry Karpeev 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 12:18 PM, Matthew Knepley <knepley(a)gmail.com> wrote: > On Sat, Feb 25, 2012 at 12:13 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > >> [migrating from petsc-maint] >> >> On Sat, Feb 25, 2012 at 11:45, Barry Smith <bsmith(a)mcs.anl.gov> wrote: >> >>> > I agree, but it would be quite helpful to refactor DM_DA have a section >>> > that holds topology and geometry, plus the (very small) extra >>> information >>> > that depends on the number of dofs per point. Then the common part >>> could be >>> > shared, making all the DMDAs coming out of the split easier to work >>> with. >>> > >>> >>> this is true but completely irrelevent here >>> >> >> My understanding is that much of the reason Matt wants to use >> PetscSection in places is because he feels like the DM implementation is >> quite heavy, so if you want to address a subset, making a new DM would >> involve heavy copies. I don't want PetscSection because I think the >> algebraic components should have access to the richer information in a DM, >> but I also don't want full copies of heavy DMs. If we can implement >> topology/geometry sharing without much complexity, we shouldn't need to use >> PetscSection. >> > > This is incorrect. I want to use PetscSection because DM just cannot do > arbitrary data layout. > I would argue that generating arbitrary layouts belongs in particular DM implementations, but the generic interface has to allow for split retrieval so that they can be used by various PCs (PCFieldSplit, PC(G)ASM) and even SNES (splitting primal and dual variables in versions of SNESVI). Dmitry. > Can we please > stop making examples with vertex unknowns on structures meshes? They have > nothing left to tell us. > > Matt > > >> >> >>> >>> > What interface are we going to use for splitting out a DM for a vector >>> > field in a larger system? How do we get the bs=2 velocity DM out of the >>> > lid-driven cavity's bs=4 velocity, vorticity, and temperature DM? >>> > >>> >>> again this is irrelevent here because the same issue is true for the >>> ISs. >>> >>> But >>> >>> I explained this in my previous email: DMPushSubDMType(DM,"name of a >>> type of splitting like u,v,p or velocity,p"); DMPopSubDMType(); it is not >>> perfect but a start. See that email for more details. >> >> >> // 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}); >> DMFieldSplittingAddField(dm,fscons,"Energy",{4}); >> >> DMGetFieldSplitting(dm,"primitive-vector",&fsprim); >> >> DMFieldSplittingSetPF(dm,fsprim,UserEquationOfStatePF,UserEquationOfStateInversePF); >> DMFieldSplittingAddField(dm,fsprim,"velocity",{0,1,2}); >> DMFieldSplittingAddField(dm,fsprim,"pressure",{3}); >> DMFieldSplittingAddField(dm,fsprim,"Temperature",{4}); >> >> If a transform is not specified, the indices would address the original >> state variables. >> >> Then PCFieldSplit and eventual nonlinear and time-stepping splitting >> methods would call >> >> DMGetFieldSplitting(dm,&nsplits,&splits); >> // choose split number to use >> split = splits[splitidx]; >> DMFieldSplittingGetName(dm,split,&splitname); >> >> DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); >> >> If pf is not NULL, apply PF to the state variables, then isfields[i] >> denotes what to extract in order to live in dmfields[i], pfinv is needed to >> put a contribution back into the coupled dm. >> > > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener >

1 0

Re: [petsc-dev] DM field interfaces
by Dmitry Karpeev 25 Feb '12

25 Feb '12

I wrote this before the thread moved to petsc-dev, so there is some overlap with the below. A given problem with the same "topology" -- a mesh and a total index layout -- might define several useful splits -- collections of fields. Each field might have further useful splits. I have useful examples from my interaction with Moose/libMesh, although these are rather general. A field is now represented by an IS. We could try to put netsted split definitions into ISs (e.g., ISNest), but I think it's much more natural to keep all nestedness inside DMs -- we already have examples of this in the coarsening, refinement for MG and a pretty well understood interaction of this sort of DM hierarchy. with a hierarchy of PCs. Thus, I think it would be useful for a DM to return a split as a list of DMs, which can be further split, if necessary: DMGetFieldDMs(DM dm, PetscInt *fieldcount, DM **fielddms); The field names can simply be the corresponding DM names, etc. The fielddms can be lightweight, sharing the underlying topology. That's up to the implementation to arrange, though. This addresses nested splits. We could use the same mehanism to implement alternative splits on the same topology -- these are different DMs sharing the topology internally. Instead of pushing and popping splits, they can be retrieved by name as DMGetSplitDM(DM dm, const char *splitname, DM *splitdm); Then calling DMGetFieldDMs on the splitdm retrieves the corresponding fields. Again, splitdm and dm can share common data structures, as arranged by the particular DM implementation. The advantage of this approach over push/pop is that we can simultaneously use different splits of the same topology. It is rather easy to apply the same scheme and have DMs serve up (overlapping) domain decompositions useful for (G)ASM. I'm now implementing both of these schemes (fieldsplit and subdomains) on top of libMesh, but the code resides in the libMesh source tree to avoid a circular dependence: as designed right now libMesh depends on PETSc. The API has to live in DM, though. On Sat, Feb 25, 2012 at 12:36 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > On Feb 25, 2012, at 12:13 PM, Jed Brown wrote: > > > [migrating from petsc-maint] > > > > > What interface are we going to use for splitting out a DM for a vector > > > field in a larger system? How do we get the bs=2 velocity DM out of the > > > lid-driven cavity's bs=4 velocity, vorticity, and temperature DM? > > > > > > > again this is irrelevent here because the same issue is true for the > ISs. > > > > But > > > > I explained this in my previous email: DMPushSubDMType(DM,"name of a > type of splitting like u,v,p or velocity,p"); DMPopSubDMType(); it is not > perfect but a start. See that email for more details. > > > > // 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}); > > DMFieldSplittingAddField(dm,fscons,"Energy",{4}); > > > > DMGetFieldSplitting(dm,"primitive-vector",&fsprim); > > > DMFieldSplittingSetPF(dm,fsprim,UserEquationOfStatePF,UserEquationOfStateInversePF); > > DMFieldSplittingAddField(dm,fsprim,"velocity",{0,1,2}); > > DMFieldSplittingAddField(dm,fsprim,"pressure",{3}); > > DMFieldSplittingAddField(dm,fsprim,"Temperature",{4}); > > > > If a transform is not specified, the indices would address the original > state variables. > This API is for constructing splits, rather than pulling them out. The two APIs do not conflict, but certain implementations might already have "naturally defined" splits (e.g., in libMesh by variable name). Dmitry. > > Jed, It also makes it impossible to handle bizare things where for one > spitting you want a subset of one thing actually in another thing. So why > not use IS to define the names, this gives complete flexibility. Also this > above is not adding fields: you are just giving names to subsets of fields > indicated by the original numbering of fields. Which is fine, but has > confusing names (you are not AddField just naming it). My suggestion is > just a more or less a concise way of doing what you are doing without the > need for YET ANOTHER CLASS (fcons). > > > > > Then PCFieldSplit and eventual nonlinear and time-stepping splitting > methods would call > > > > DMGetFieldSplitting(dm,&nsplits,&splits); > > // choose split number to use > > split = splits[splitidx]; > > DMFieldSplittingGetName(dm,split,&splitname); > > > DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); > > I do not like this approach based on requesting by name, just push your > name and then pop it. This way something "above" in the call stack can push > to the right one then something further down can use it without having to > know its name (sort of like your "default" split) but the push allows > changing the default. You kind of assume that the guts is suppose to know > what name to use? // choose split number to use how is some PC deep inside > suppose to know this? > > I realize my proposal was not complete on detail. I'm suggesting changes > to your detail. > > Barry > > > > > If pf is not NULL, apply PF to the state variables, then isfields[i] > denotes what to extract in order to live in dmfields[i], pfinv is needed to > put a contribution back into the coupled dm. > >

1 0

Re: [petsc-dev] DM field interfaces
by Barry Smith 25 Feb '12

25 Feb '12

On Feb 25, 2012, at 12:13 PM, Jed Brown wrote: > [migrating from petsc-maint] > > > What interface are we going to use for splitting out a DM for a vector > > field in a larger system? How do we get the bs=2 velocity DM out of the > > lid-driven cavity's bs=4 velocity, vorticity, and temperature DM? > > > > again this is irrelevent here because the same issue is true for the ISs. > > But > > I explained this in my previous email: DMPushSubDMType(DM,"name of a type of splitting like u,v,p or velocity,p"); DMPopSubDMType(); it is not perfect but a start. See that email for more details. > > // 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}); > DMFieldSplittingAddField(dm,fscons,"Energy",{4}); > > DMGetFieldSplitting(dm,"primitive-vector",&fsprim); > DMFieldSplittingSetPF(dm,fsprim,UserEquationOfStatePF,UserEquationOfStateInversePF); > DMFieldSplittingAddField(dm,fsprim,"velocity",{0,1,2}); > DMFieldSplittingAddField(dm,fsprim,"pressure",{3}); > DMFieldSplittingAddField(dm,fsprim,"Temperature",{4}); > > If a transform is not specified, the indices would address the original state variables. Jed, It also makes it impossible to handle bizare things where for one spitting you want a subset of one thing actually in another thing. So why not use IS to define the names, this gives complete flexibility. Also this above is not adding fields: you are just giving names to subsets of fields indicated by the original numbering of fields. Which is fine, but has confusing names (you are not AddField just naming it). My suggestion is just a more or less a concise way of doing what you are doing without the need for YET ANOTHER CLASS (fcons). > > Then PCFieldSplit and eventual nonlinear and time-stepping splitting methods would call > > DMGetFieldSplitting(dm,&nsplits,&splits); > // choose split number to use > split = splits[splitidx]; > DMFieldSplittingGetName(dm,split,&splitname); > DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); I do not like this approach based on requesting by name, just push your name and then pop it. This way something "above" in the call stack can push to the right one then something further down can use it without having to know its name (sort of like your "default" split) but the push allows changing the default. You kind of assume that the guts is suppose to know what name to use? // choose split number to use how is some PC deep inside suppose to know this? I realize my proposal was not complete on detail. I'm suggesting changes to your detail. Barry > > If pf is not NULL, apply PF to the state variables, then isfields[i] denotes what to extract in order to live in dmfields[i], pfinv is needed to put a contribution back into the coupled dm.

1 0

Re: [petsc-dev] DM field interfaces
by Jed Brown 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 12:23, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > I now realize a weakness in my PF transformation suggestion: it doesn't > work with vector elements (Raviart-Thomas or Nedelec) or staggered grids. I > don't know a way to include that with the ability to transform from > conservative to primitive variables. This is unfortunate because many > important preconditioning methods operate in different variables than the > nonlinear system is defined in. This is to say, that sort of transformation can't reasonably be _implemented_ using only DM base, it has to be implemented using DMXXXSet* methods.

1 0

Re: [petsc-dev] DM field interfaces
by Jed Brown 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 12:18, Matthew Knepley <knepley(a)gmail.com> wrote: > This is incorrect. I want to use PetscSection because DM just cannot do > arbitrary data layout. Can we please > stop making examples with vertex unknowns on structures meshes? They have > nothing left to tell us. > This has nothing to do with DM, that's DMDA. I now realize a weakness in my PF transformation suggestion: it doesn't work with vector elements (Raviart-Thomas or Nedelec) or staggered grids. I don't know a way to include that with the ability to transform from conservative to primitive variables. This is unfortunate because many important preconditioning methods operate in different variables than the nonlinear system is defined in.

1 0

Re: [petsc-dev] DM field interfaces
by Matthew Knepley 25 Feb '12

25 Feb '12

On Sat, Feb 25, 2012 at 12:13 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > [migrating from petsc-maint] > > On Sat, Feb 25, 2012 at 11:45, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > >> > I agree, but it would be quite helpful to refactor DM_DA have a section >> > that holds topology and geometry, plus the (very small) extra >> information >> > that depends on the number of dofs per point. Then the common part >> could be >> > shared, making all the DMDAs coming out of the split easier to work >> with. >> > >> >> this is true but completely irrelevent here >> > > My understanding is that much of the reason Matt wants to use PetscSection > in places is because he feels like the DM implementation is quite heavy, so > if you want to address a subset, making a new DM would involve heavy > copies. I don't want PetscSection because I think the algebraic components > should have access to the richer information in a DM, but I also don't want > full copies of heavy DMs. If we can implement topology/geometry sharing > without much complexity, we shouldn't need to use PetscSection. > This is incorrect. I want to use PetscSection because DM just cannot do arbitrary data layout. Can we please stop making examples with vertex unknowns on structures meshes? They have nothing left to tell us. Matt > > >> >> > What interface are we going to use for splitting out a DM for a vector >> > field in a larger system? How do we get the bs=2 velocity DM out of the >> > lid-driven cavity's bs=4 velocity, vorticity, and temperature DM? >> > >> >> again this is irrelevent here because the same issue is true for the >> ISs. >> >> But >> >> I explained this in my previous email: DMPushSubDMType(DM,"name of a >> type of splitting like u,v,p or velocity,p"); DMPopSubDMType(); it is not >> perfect but a start. See that email for more details. > > > // 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}); > DMFieldSplittingAddField(dm,fscons,"Energy",{4}); > > DMGetFieldSplitting(dm,"primitive-vector",&fsprim); > > DMFieldSplittingSetPF(dm,fsprim,UserEquationOfStatePF,UserEquationOfStateInversePF); > DMFieldSplittingAddField(dm,fsprim,"velocity",{0,1,2}); > DMFieldSplittingAddField(dm,fsprim,"pressure",{3}); > DMFieldSplittingAddField(dm,fsprim,"Temperature",{4}); > > If a transform is not specified, the indices would address the original > state variables. > > Then PCFieldSplit and eventual nonlinear and time-stepping splitting > methods would call > > DMGetFieldSplitting(dm,&nsplits,&splits); > // choose split number to use > split = splits[splitidx]; > DMFieldSplittingGetName(dm,split,&splitname); > > DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); > > If pf is not NULL, apply PF to the state variables, then isfields[i] > denotes what to extract in order to live in dmfields[i], pfinv is needed to > put a contribution back into the coupled dm. > -- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener

1 0

DM field interfaces
by Jed Brown 25 Feb '12

25 Feb '12

[migrating from petsc-maint] On Sat, Feb 25, 2012 at 11:45, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > I agree, but it would be quite helpful to refactor DM_DA have a section > > that holds topology and geometry, plus the (very small) extra information > > that depends on the number of dofs per point. Then the common part could > be > > shared, making all the DMDAs coming out of the split easier to work with. > > > > this is true but completely irrelevent here > My understanding is that much of the reason Matt wants to use PetscSection in places is because he feels like the DM implementation is quite heavy, so if you want to address a subset, making a new DM would involve heavy copies. I don't want PetscSection because I think the algebraic components should have access to the richer information in a DM, but I also don't want full copies of heavy DMs. If we can implement topology/geometry sharing without much complexity, we shouldn't need to use PetscSection. > > > What interface are we going to use for splitting out a DM for a vector > > field in a larger system? How do we get the bs=2 velocity DM out of the > > lid-driven cavity's bs=4 velocity, vorticity, and temperature DM? > > > > again this is irrelevent here because the same issue is true for the ISs. > > But > > I explained this in my previous email: DMPushSubDMType(DM,"name of a > type of splitting like u,v,p or velocity,p"); DMPopSubDMType(); it is not > perfect but a start. See that email for more details. // 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}); DMFieldSplittingAddField(dm,fscons,"Energy",{4}); DMGetFieldSplitting(dm,"primitive-vector",&fsprim); DMFieldSplittingSetPF(dm,fsprim,UserEquationOfStatePF,UserEquationOfStateInversePF); DMFieldSplittingAddField(dm,fsprim,"velocity",{0,1,2}); DMFieldSplittingAddField(dm,fsprim,"pressure",{3}); DMFieldSplittingAddField(dm,fsprim,"Temperature",{4}); If a transform is not specified, the indices would address the original state variables. Then PCFieldSplit and eventual nonlinear and time-stepping splitting methods would call DMGetFieldSplitting(dm,&nsplits,&splits); // choose split number to use split = splits[splitidx]; DMFieldSplittingGetName(dm,split,&splitname); DMFieldSplittingGetFields(dm,split,&pf,&pfinv,&nfields,&isfields,&dmfields); If pf is not NULL, apply PF to the state variables, then isfields[i] denotes what to extract in order to live in dmfields[i], pfinv is needed to put a contribution back into the coupled dm.

1 0