On Thu, Nov 21, 2013 at 3:23 PM, Geoffrey Irving <[email protected]> wrote:
> On Wed, Nov 20, 2013 at 3:52 PM, Matthew Knepley <[email protected]> > wrote: > > On Tue, Nov 19, 2013 at 5:27 PM, Geoffrey Irving <[email protected]> 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? Definitely we should guard the feBd loop to only look at co-dimension 1 things. I think we need to allow all points in boundary labels. Matt > > Geoffrey > -- 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
