[petsc-dev] Is DMPlexLabelComplete necessary in snes ex12?

Geoffrey Irving irving at naml.us
Thu Nov 21 15:23:39 CST 2013


On Wed, Nov 20, 2013 at 3:52 PM, Matthew Knepley <knepley at gmail.com> wrote:
> On Tue, Nov 19, 2013 at 5:27 PM, Geoffrey Irving <irving at naml.us> wrote:
>>
>> This question isn't about changing ex12, I'm just making sure I
>> understand what's going on:
>>
>> SetupSection in snes ex12 calls DMPlexLabelComplete on the "marker"
>> label if Dirichlet conditions are chosen (not if Neumann is chosen).
>> Specifically,
>>
>>     if (user->bcType == DIRICHLET) {ierr = DMPlexLabelComplete(dm,
>> label);CHKERRQ(ierr);}
>>
>> It then uses the same label ("marker") to create an index set:
>>
>>     if (user->bcType == DIRICHLET) {ierr  = DMPlexGetStratumIS(dm,
>> bdLabel, 1, &bcPoints[0]);CHKERRQ(ierr);}
>>
>> I believe 1 is the codimension.  If the index set is restricted to
>> codimension 1, and DMPlexLabelComplete labels vertices based on edge
>> labels, isn't the DMPlexLabelComplete call unnecessary?
>
>
> 1) The value 1 is just a convention that I use for marking the boundary. In a real simulation,
>      there are probably many makers and they are user specified. This is what happens in
>      TS ex11.
>
> 2) I should probably complete the label in the Neumann case, but I just forgot. The completion
>      is not necessary the way my defaults mesh generators work, but is sometimes necessary
>      when reading stuff in.

Actually, it looks like calling DMPlexLabelComplete in the Neumann
case breaks with an error later.  In plexfem.c, if feBd is nonzero,
DMPlexComputeCellGeometry is called on all "points", where points are
actually whatever the DMLabel "boundary" spits out.
DMPlexComputeGeometry bails with an exception if called on anything
with dimension 0:

[0]PETSC ERROR: DMPlexComputeCellGeometry() line 783 in
src/dm/impls/plex/plexgeometry.c Unsupported dimension 0 in cell 2 for
element geometry computation

Should the "points" variable be "faces"?  Should I restrict the loops
in plexfem.c to only touch codimension 1 elements?

Geoffrey



More information about the petsc-dev mailing list