Hi Barry et al.,
> On 28 Oct 2021, at 15:37, Barry Smith <[email protected]> wrote:
>
>
>
>> On Oct 28, 2021, at 10:31 AM, Matthew Knepley <[email protected]> wrote:
>>
>> On Thu, Oct 28, 2021 at 9:37 AM Barry Smith <[email protected]> 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
