> On 19 Mar 2021, at 14:21, Jed Brown <[email protected]> wrote:
>
> Notice how the permutations are contained within the vertices {0, ..., 8},
> edges {9, ..., 24}, and cells {25, ..., 32}. I would like to get rid of that
> restriction, but you've said it would have significant non-local consequences
> so I haven't tried.
What we do, and I'm not sure if this is exactly what you want, is the following.
Use DMPlexGetOrdering to compute permutatations of points with RCM. As you say,
this is stratified.
Then, we traverse cells in this RCM ordering and compute a permutation of the
full plex chart greedily.
Pseudo-code:
permutation = empty(pstart - pend)
i = 0
for cell in reordered_cells:
points = plex.getclosure(cell)
for point in points:
if not seen(point):
permutation[i] = point
i += 1
Now when we create a section, we say
PetscSectionSetPermutation(sec, permutation)
And then when the section offsets are calculated, the points are visited and
numbered in the order provided by the permutation.
This seems to work, here are the sparsity patterns for a Q3 mass matrix on an
unstructured quad mesh (from gmsh). One with RCM ordering applied, the other
without).
I guess the only issue is now one needs to always remember to traverse cells in
the reordered_cells order (although if you call DMPlexPermute that is
automatic?).
Cheers,
Lawrence
