Thanks, sounds like something DMPLEX should be automatically optionally doing directly.
> On Oct 28, 2021, at 10:58 AM, Lawrence Mitchell <we...@gmx.li> wrote: > > Hi Barry et al., > >> On 28 Oct 2021, at 15:37, Barry Smith <bsm...@petsc.dev> wrote: >> >> >> >>> On Oct 28, 2021, at 10:31 AM, Matthew Knepley <knep...@gmail.com> wrote: >>> >>> On Thu, Oct 28, 2021 at 9:37 AM Barry Smith <bsm...@petsc.dev> wrote: >>> >>> Matt, >>> >>> How difficult would it be to rework DMPLEX to allow the use of VecGhost? >>> We have performance problems with GPUs with simple DMNETWORK models because >>> the code spends more time uselessly copying the local part of the vector to >>> another vector in global to local and local to global; more than 1/2 the >>> time of the total simulation. >>> >>> Firedrake already does this because they "vec ghost" their vectors by >>> default. Here is what you need: >>> >>> When you create the PetscSection, by default it orders the unknowns >>> according to the default point numbering. This >>> is what causes the ghost unknowns to be mixed in with the local unknowns. >>> However, PetscSection allows you to set >>> a point permutation >>> >>> >>> https://petsc.org/main/docs/manualpages/PetscSection/PetscSectionSetPermutation.html >>> >>> This determines the order of dogs by iterating through points in this >>> permutation, and you can put all shared points at the end. >> >> How do I know what are shared points to put at the end? Couldn't DMPLEX do >> this automatically with an option? Where is the Firedrake code that does >> this with DMPLEX so I can see it? > > The DM's "point sf" indicates shared points. To label them, do something like: > > DMCreateLabel(dm, "ghost_points"); > DMGetPointSF(dm, &pointsf); > PetscSFGetGraph(pointsf, NULL, &nleaves, &ilocal, NULL); > for (PetscInt p = 0; p < nleaves; p++) { > DMSetLabelValue(dm, "ghost_points", p, 1); > } > > Then you do something like (this is more pseudo-codey): > > PetscInt *permutation; > DMPlexGetChart(dm, &pStart, &pEnd); > PetscMalloc1(pEnd - pStart, &permutation); > PetscInt offsets[2] = {0, pEnd - pStart - nleaves}; > DMGetLabel(dm, &ghosts); > DMLabelCreateIndex(ghosts, pStart, pEnd); > for (PetscInt p = pStart, p < pEnd; p++) { > DMLabelHasPoint(ghosts, p, &has); > if (has) { > // this point is ghost point > permutation[offsets[1]++] = p; > } else { > permutation[offsets[0]++] = p; > } > } > DMLabelDestroyIndex(ghosts, pStart, pEnd); > ISCreateGeneral(..., permutation, &isperm); > > // Now whenever your do PetscSectionCreate, do > PetscSectionSetPermutation(..., isperm); > > And now ghost point dofs will appear after local ones. > > Note that this is probably more complicated for multi-field setups, depending > whether you are point major or field major. > > You can see what we actually do (if you like reading Cython) here: > https://github.com/firedrakeproject/firedrake/blob/master/firedrake/cython/dmcommon.pyx#L1734 > > It's more complicated because here we are doing some additional things: > > 1. We compute an RCM-order traversal for cells with DMPlexGetOrdering; > 2. Rather than ordering all plex points in the permutation, we walk the cells > in the RCM order and then greedily number points in the transitive closure > (so that in a section cell and vertex dofs from the same cell will be "close" > to each other in the final Vec). > > Thanks, > > Lawrence
