[petsc-dev] DMPlex overlap redistribution

[Matthew Knepley]

On Tue, Mar 12, 2024 at 6:49 PM Adrian Croucher <a.croucher at auckland.ac.nz>
wrote:

> hi Matt,
>
> Thanks for that. I think that is a useful change, even though it doesn't
> actually fix my problem (i.e. my code still gets an error). I think now
> that although the error message I'm getting is similar to the one that used
> to result from just a simple repeated stratification, the underlying cause
> is different.
>
> Looking more closely at the error message, it is saying that the depth 2
> stratum overlaps the depth 1 stratum, e.g.:
>
> [0]PETSC ERROR: New depth 2 range [0,141) overlaps with depth 1 range
> [0,141)
>
> That suggested to me that some of the strata in the re-distributed DMPlex
> are not contiguous. I checked that by looking at the cone and support sizes
> of the DMPlex points inside DMPlexStratify_Topological_Private(). On one of
> the processes they looked like this:
>
> p: 0, cone: 7, support: 0
> p: 1, cone: 7, support: 0
> p: 2, cone: 7, support: 0
> p: 3, cone: 6, support: 0
> p: 4, cone: 6, support: 0
> p: 5, cone: 6, support: 0
> p: 6, cone: 2, support: 0
> p: 7, cone: 1, support: 0
> p: 8, cone: 0, support: 1
> p: 9, cone: 0, support: 1
> p: 10, cone: 1, support: 2
> p: 11, cone: 1, support: 2
> p: 12, cone: 1, support: 1
> p: 13, cone: 1, support: 1
> p: 14, cone: 2, support: 0
> p: 15, cone: 2, support: 0
> p: 16, cone: 1, support: 0
> p: 17, cone: 1, support: 0
> p: 18, cone: 0, support: 1
> ...
>
> So here the cells (support 0) are not contiguous (and neither are the
> other strata) - some of the dual-porosity cells (14 - 17) are getting
> separated from the other cells (0 - 7). So I think that is messing up the
> identification of stratum bounds, because it's right, they do overlap.
>
> Before redistribution, the strata are contiguous and non-overlapping. So I
> think the redistribution is causing the trouble.
>
> When I first add the dual-porosity points to the DM I add them on the end
> of each stratum, after interior points and partition ghost points. That
> used to work, but perhaps something has changed so it doesn't any more? Are
> there any assumptions made about the order of storing interior points,
> partition ghost points and interior ghost points?
>
> Another possibility might be that it's getting confused by the
> dual-porosity faces, edges, and nodes I'm adding which have unusual
> cone/support sizes (e.g. faces with a cone of one edge, which in turn have
> a cone of one node). Again this used to work ok but perhaps something has
> changed so that it doesn't?
>
Thanks for finding this. This was indeed rewritten because the old
algorithm was too expensive at large sizes/parallelism. The function where
this breaks in DMPlexStratifyMigrationSF() in plexdistribute.c. The problem
is here:


https://urldefense.us/v3/__https://gitlab.com/petsc/petsc/-/blob/main/src/dm/impls/plex/plexdistribute.c?ref_type=heads*L915__;Iw!!G_uCfscf7eWS!Yk85ASiGW0qV4yJakUojWQTr2JxXr5K44nxh2tKUmC3i-hL5fa7GwIu28UGfXCrmRT91IEVostCpYPLUl3bR$ 

I treat DM_POLYTOPE_INTERIOR_GHOST as a cell no matter what. We need to
assign the interior  vertices, faces, and edges types that will give back
the right dimension. Vertices are easy, they should just be type POINT. I
can understand not wanting to use SEGMENT or QUAD for the interior edge and
face because they do not have full cones. However, I think nothing will
break (except explicit checks, which will not work anyway) if you give them
types SEGMENT and QUAD. Is this easy to try?

I could also make a new type for you at each dimension if faking it does
not work.

  Thanks,

     Matt

> - Adrian
> On 9/03/24 7:15 am, Matthew Knepley wrote:
>
> Okay, here are my first steps
>
> 1. I fixed the problem with repeated stratification
>
> knepley/fix-da-periodicity *$:/PETSc3/petsc/petsc-dev$ git diff
> diff --git a/src/dm/impls/plex/plex.c b/src/dm/impls/plex/plex.c
> index 0ffc372ea08..9d27d308538 100644
> --- a/src/dm/impls/plex/plex.c
> +++ b/src/dm/impls/plex/plex.c
> @@ -4439,6 +4439,7 @@ PetscErrorCode DMPlexStratify(DM dm)
>    PetscCall(PetscLogEventBegin(DMPLEX_Stratify, dm, 0, 0, 0));
>
>    // Create depth label
> +  PetscCall(DMRemoveLabel(dm, "depth", NULL));
>    PetscCall(DMCreateLabel(dm, "depth"));
>
>    PetscCall(DMPlexGetDepthLabel(dm, &label));
>
> I put it in an MR: https://urldefense.us/v3/__https://gitlab.com/petsc/petsc/-/merge_requests/7347__;!!G_uCfscf7eWS!Yk85ASiGW0qV4yJakUojWQTr2JxXr5K44nxh2tKUmC3i-hL5fa7GwIu28UGfXCrmRT91IEVostCpYFwOV8Q1$ 
>
> 2. I wrote a small program to load the mesh, and then restratify it. It is
> attached.
>
> 3. I ran stratification on this mesh, but I did not get an error:
>
> --
> Dr Adrian Croucher
> Senior Research Fellow
> Department of Engineering Science
> Waipapa Taumata Rau / University of Auckland, New Zealand
> email: a.croucher at auckland.ac.nz
> tel: +64 (0)9 923 4611
>
>

-- 
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

https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Yk85ASiGW0qV4yJakUojWQTr2JxXr5K44nxh2tKUmC3i-hL5fa7GwIu28UGfXCrmRT91IEVostCpYM8BJssE$  <https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Yk85ASiGW0qV4yJakUojWQTr2JxXr5K44nxh2tKUmC3i-hL5fa7GwIu28UGfXCrmRT91IEVostCpYC-76Ov-$ >
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20240312/bc72d19d/attachment-0001.html>