<div dir="ltr"><div dir="ltr">On Wed, Mar 31, 2021 at 3:36 PM Nicolas Barral <<a href="mailto:nicolas.barral@math.u-bordeaux.fr">nicolas.barral@math.u-bordeaux.fr</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">On 31/03/2021 21:24, Matthew Knepley wrote:<br>
>     Ok so that's worth clarifying, but my problem is that dm has a<br>
>     0-overlap<br>
>     before DMPlexDistributeOverlap and odm has a 1-overlap even though I<br>
>     passed "0". If I understand what you just wrote, adding "0" more levels<br>
>     of overlap to 0-overlap, that should still be 0 overlap ?<br>
> <br>
>     Yet when I look at the code of DMPlexDistributeOverlap, you're flagging<br>
>     points to be added to the overlap whatever "k" is.<br>
> <br>
> <br>
> Crap. There is a bug handling 0. Evidently, no one ever asks for overlap 0.<br>
> Will fix.<br>
> <br>
ok now I'm sure I understand what you explained :) Thanks Matt.<br>
<br>
Now we can look at the other question: I need to be able to pass the <br>
non-overlapping mesh to the remesher. I can either maintain 2 plexes, or <br>
trim the overlapping plex when I create the arrays I pass to the <br>
remesher. I'm not sure which is the best/worst ?<br></blockquote><div><br></div><div>I would start with two plexes since it is very easy.</div><div><br></div><div>Trimming the plex would essentially make another plex anyway, but it is not hard using DMPlexFilter().</div><div><br></div><div>  Thanks,</div><div><br></div><div>     Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
Thanks<br>
<br>
-- <br>
Nicolas<br>
<br>
<br>
>    Thanks,<br>
> <br>
>       Matt<br>
> <br>
>     Sample code:<br>
>     static char help[] = "Tests plex distribution and overlaps.\n";<br>
> <br>
>     #include <petsc/private/dmpleximpl.h><br>
> <br>
>     int main (int argc, char * argv[]) {<br>
> <br>
>         DM             dm, odm;<br>
>         MPI_Comm       comm;<br>
>         PetscErrorCode ierr;<br>
> <br>
>         ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr)<br>
>     return ierr;<br>
>         comm = PETSC_COMM_WORLD;<br>
> <br>
>         ierr = DMPlexCreateBoxMesh(comm, 2, PETSC_TRUE, NULL, NULL, NULL,<br>
>     NULL, PETSC_TRUE, &dm);CHKERRQ(ierr);<br>
>         ierr = DMSetFromOptions(dm);CHKERRQ(ierr);<br>
> <br>
> <br>
> <br>
>         ierr = PetscObjectSetName((PetscObject) dm, "DM<br>
>     before");CHKERRQ(ierr);<br>
>         ierr = DMViewFromOptions(dm, NULL, "-before_dm_view");CHKERRQ(ierr);<br>
>         DMPlexDistributeOverlap(dm, 0, NULL, &odm);<br>
>         ierr = PetscObjectSetName((PetscObject) odm, "DM<br>
>     after");CHKERRQ(ierr);<br>
>         ierr = DMViewFromOptions(odm, NULL, "-after_dm_view");CHKERRQ(ierr);<br>
> <br>
> <br>
>         ierr = DMDestroy(&dm);CHKERRQ(ierr);<br>
>         ierr = DMDestroy(&odm);CHKERRQ(ierr);<br>
>         ierr = PetscFinalize();<br>
>         return ierr;<br>
>     }<br>
> <br>
>     % mpiexec -n 2 ./test_overlapV3 -dm_plex_box_faces 5,5 -dm_distribute<br>
>     -before_dm_view -after_dm_view<br>
>     DM Object: DM before 2 MPI processes<br>
>         type: plex<br>
>     DM before in 2 dimensions:<br>
>         0-cells: 21 21<br>
>         1-cells: 45 45<br>
>         2-cells: 25 25<br>
>     Labels:<br>
>         depth: 3 strata with value/size (0 (21), 1 (45), 2 (25))<br>
>         celltype: 3 strata with value/size (0 (21), 1 (45), 3 (25))<br>
>         marker: 1 strata with value/size (1 (21))<br>
>         Face Sets: 1 strata with value/size (1 (10))<br>
>     DM Object: DM after 2 MPI processes<br>
>         type: plex<br>
>     DM after in 2 dimensions:<br>
>         0-cells: 29 29<br>
>         1-cells: 65 65<br>
>         2-cells: 37 37<br>
>     Labels:<br>
>         depth: 3 strata with value/size (0 (29), 1 (65), 2 (37))<br>
>         celltype: 3 strata with value/size (0 (29), 1 (65), 3 (37))<br>
>         marker: 1 strata with value/size (1 (27))<br>
>         Face Sets: 1 strata with value/size (1 (13))<br>
> <br>
> <br>
> <br>
>      ><br>
>      >    Thanks,<br>
>      ><br>
>      >       Matt<br>
>      ><br>
>      >     Thanks,<br>
>      ><br>
>      >     --<br>
>      >     Nicolas<br>
>      ><br>
>      >      ><br>
>      >      >    Thanks,<br>
>      >      ><br>
>      >      >       Matt<br>
>      >      ><br>
>      >      >         if (!dm) {printf("Big problem\n"); dm = odm;}<br>
>      >      >         else     {DMDestroy(&odm);}<br>
>      >      >         ierr = PetscObjectSetName((PetscObject) dm, "Initial<br>
>      >      >     DM");CHKERRQ(ierr);<br>
>      >      >         ierr = DMViewFromOptions(dm, NULL,<br>
>      >      >     "-initial_dm_view");CHKERRQ(ierr);<br>
>      >      ><br>
>      >      ><br>
>      >      >         ierr = DMDestroy(&dm);CHKERRQ(ierr);<br>
>      >      >         ierr = PetscFinalize();<br>
>      >      >         return ierr;<br>
>      >      >     }<br>
>      >      ><br>
>      >      >     called with mpiexec -n 2 ./test_overlapV3 -initial_dm_view<br>
>      >      >     -dm_plex_box_faces 5,5 -dm_distribute<br>
>      >      ><br>
>      >      >     gives:<br>
>      >      >     DM Object: Initial DM 2 MPI processes<br>
>      >      >         type: plex<br>
>      >      >     Initial DM in 2 dimensions:<br>
>      >      >         0-cells: 29 29<br>
>      >      >         1-cells: 65 65<br>
>      >      >         2-cells: 37 37<br>
>      >      >     Labels:<br>
>      >      >         depth: 3 strata with value/size (0 (29), 1 (65), 2<br>
>     (37))<br>
>      >      >         celltype: 3 strata with value/size (0 (29), 1<br>
>     (65), 3 (37))<br>
>      >      >         marker: 1 strata with value/size (1 (27))<br>
>      >      >         Face Sets: 1 strata with value/size (1 (13))<br>
>      >      ><br>
>      >      >     which is not what I expect ?<br>
>      >      ><br>
>      >      >     Thanks,<br>
>      >      ><br>
>      >      >     --<br>
>      >      >     Nicolas<br>
>      >      ><br>
>      >      >     On 31/03/2021 19:02, Matthew Knepley wrote:<br>
>      >      >      > Alright, I think the problem had to do with keeping<br>
>     track<br>
>      >     of what<br>
>      >      >     DM you<br>
>      >      >      > were looking at. This code increases the overlap of an<br>
>      >     initial DM:<br>
>      >      >      ><br>
>      >      >      > static char help[] = "Tests plex distribution and<br>
>      >     overlaps.\n";<br>
>      >      >      ><br>
>      >      >      > #include <petsc/private/dmpleximpl.h><br>
>      >      >      ><br>
>      >      >      > int main (int argc, char * argv[]) {<br>
>      >      >      ><br>
>      >      >      >    DM             dm, dm2;<br>
>      >      >      >    PetscInt       overlap;<br>
>      >      >      >    MPI_Comm       comm;<br>
>      >      >      >    PetscErrorCode ierr;<br>
>      >      >      ><br>
>      >      >      >    ierr = PetscInitialize(&argc, &argv, NULL,<br>
>     help);if (ierr)<br>
>      >      >     return ierr;<br>
>      >      >      >    comm = PETSC_COMM_WORLD;<br>
>      >      >      ><br>
>      >      >      >    ierr = DMPlexCreateBoxMesh(comm, 2, PETSC_TRUE,<br>
>     NULL,<br>
>      >     NULL, NULL,<br>
>      >      >      > NULL, PETSC_TRUE, &dm);CHKERRQ(ierr);<br>
>      >      >      >    ierr = DMSetFromOptions(dm);CHKERRQ(ierr);<br>
>      >      >      >    ierr = PetscObjectSetName((PetscObject) dm, "Initial<br>
>      >      >     DM");CHKERRQ(ierr);<br>
>      >      >      >    ierr = DMViewFromOptions(dm, NULL,<br>
>      >      >     "-initial_dm_view");CHKERRQ(ierr);<br>
>      >      >      ><br>
>      >      >      >    ierr = DMPlexGetOverlap(dm, &overlap);CHKERRQ(ierr);<br>
>      >      >      >    ierr = DMPlexDistributeOverlap(dm, overlap+1, NULL,<br>
>      >      >     &dm2);CHKERRQ(ierr);<br>
>      >      >      >    ierr = PetscObjectSetName((PetscObject) dm2,<br>
>     "More Overlap<br>
>      >      >      > DM");CHKERRQ(ierr);<br>
>      >      >      >    ierr = DMViewFromOptions(dm2, NULL,<br>
>      >      >     "-over_dm_view");CHKERRQ(ierr);<br>
>      >      >      ><br>
>      >      >      >    ierr = DMDestroy(&dm2);CHKERRQ(ierr);<br>
>      >      >      >    ierr = DMDestroy(&dm);CHKERRQ(ierr);<br>
>      >      >      >    ierr = PetscFinalize();<br>
>      >      >      >    return ierr;<br>
>      >      >      > }<br>
>      >      >      ><br>
>      >      >      > and when we run it we get the expected result<br>
>      >      >      ><br>
>      >      >      > master *:~/Downloads/tmp/Nicolas$<br>
>      >     /PETSc3/petsc/apple/bin/mpiexec<br>
>      >      >     -n 2<br>
>      >      >      > ./test_overlap -initial_dm_view -dm_plex_box_faces 5,5<br>
>      >      >     -dm_distribute<br>
>      >      >      > -dm_distribute_overlap 1 -over_dm_view<br>
>      >      >      > DM Object: Initial DM 2 MPI processes<br>
>      >      >      >    type: plex<br>
>      >      >      > Initial DM in 2 dimensions:<br>
>      >      >      >    0-cells: 29 29<br>
>      >      >      >    1-cells: 65 65<br>
>      >      >      >    2-cells: 37 37<br>
>      >      >      > Labels:<br>
>      >      >      >    depth: 3 strata with value/size (0 (29), 1 (65),<br>
>     2 (37))<br>
>      >      >      >    celltype: 3 strata with value/size (0 (29), 1<br>
>     (65), 3 (37))<br>
>      >      >      >    marker: 1 strata with value/size (1 (27))<br>
>      >      >      >    Face Sets: 1 strata with value/size (1 (13))<br>
>      >      >      > DM Object: More Overlap DM 2 MPI processes<br>
>      >      >      >    type: plex<br>
>      >      >      > More Overlap DM in 2 dimensions:<br>
>      >      >      >    0-cells: 36 36<br>
>      >      >      >    1-cells: 85 85<br>
>      >      >      >    2-cells: 50 50<br>
>      >      >      > Labels:<br>
>      >      >      >    depth: 3 strata with value/size (0 (36), 1 (85),<br>
>     2 (50))<br>
>      >      >      >    celltype: 3 strata with value/size (0 (36), 1<br>
>     (85), 3 (50))<br>
>      >      >      >    marker: 1 strata with value/size (1 (40))<br>
>      >      >      >    Face Sets: 1 strata with value/size (1 (20))<br>
>      >      >      ><br>
>      >      >      >    Thanks,<br>
>      >      >      ><br>
>      >      >      >       Matt<br>
>      >      >      ><br>
>      >      >      > On Wed, Mar 31, 2021 at 12:57 PM Matthew Knepley<br>
>      >      >     <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>><br>
>     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>>><br>
>      >     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>><br>
>     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>>>><br>
>      >      >      > <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a><br>
>     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a><br>
>     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>>><br>
>      >     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>><br>
>     <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a> <mailto:<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>>>>>> wrote:<br>
>      >      >      ><br>
>      >      >      >     Okay, let me show a really simple example that<br>
>     gives<br>
>      >     the expected<br>
>      >      >      >     result before I figure out what is going wrong for<br>
>      >     you. This code<br>
>      >      >      ><br>
>      >      >      >     static char help[] = "Tests plex distribution and<br>
>      >     overlaps.\n";<br>
>      >      >      ><br>
>      >      >      >     #include <petsc/private/dmpleximpl.h><br>
>      >      >      ><br>
>      >      >      >     int main (int argc, char * argv[]) {<br>
>      >      >      >        DM                    dm;<br>
>      >      >      >        MPI_Comm       comm;<br>
>      >      >      >        PetscErrorCode ierr;<br>
>      >      >      ><br>
>      >      >      >        ierr = PetscInitialize(&argc, &argv, NULL,<br>
>     help);if<br>
>      >     (ierr)<br>
>      >      >     return<br>
>      >      >      >     ierr;<br>
>      >      >      >        comm = PETSC_COMM_WORLD;<br>
>      >      >      ><br>
>      >      >      >        ierr = DMPlexCreateBoxMesh(comm, 2,<br>
>     PETSC_TRUE, NULL,<br>
>      >      >     NULL, NULL,<br>
>      >      >      >     NULL, PETSC_TRUE, &dm);CHKERRQ(ierr);<br>
>      >      >      >        ierr = DMSetFromOptions(dm);CHKERRQ(ierr);<br>
>      >      >      >        ierr = PetscObjectSetName((PetscObject) dm,<br>
>     "Initial<br>
>      >      >      >     DM");CHKERRQ(ierr);<br>
>      >      >      >        ierr = DMViewFromOptions(dm, NULL,<br>
>      >      >     "-initial_dm_view");CHKERRQ(ierr);<br>
>      >      >      >        ierr = DMDestroy(&dm);CHKERRQ(ierr);<br>
>      >      >      >        ierr = PetscFinalize();<br>
>      >      >      >        return ierr;<br>
>      >      >      >     }<br>
>      >      >      ><br>
>      >      >      >     can do all the overlap tests. For example, you<br>
>     can run<br>
>      >     it naively<br>
>      >      >      >     and get a serial mesh<br>
>      >      >      ><br>
>      >      >      >     master *:~/Downloads/tmp/Nicolas$<br>
>      >      >     /PETSc3/petsc/apple/bin/mpiexec -n<br>
>      >      >      >     2 ./test_overlap -initial_dm_view<br>
>     -dm_plex_box_faces 5,5<br>
>      >      >      >     DM Object: Initial DM 2 MPI processes<br>
>      >      >      >        type: plex<br>
>      >      >      >     Initial DM in 2 dimensions:<br>
>      >      >      >        0-cells: 36 0<br>
>      >      >      >        1-cells: 85 0<br>
>      >      >      >        2-cells: 50 0<br>
>      >      >      >     Labels:<br>
>      >      >      >        celltype: 3 strata with value/size (0 (36),<br>
>     3 (50),<br>
>      >     1 (85))<br>
>      >      >      >        depth: 3 strata with value/size (0 (36), 1<br>
>     (85), 2<br>
>      >     (50))<br>
>      >      >      >        marker: 1 strata with value/size (1 (40))<br>
>      >      >      >        Face Sets: 1 strata with value/size (1 (20))<br>
>      >      >      ><br>
>      >      >      >     Then run it telling Plex to distribute after<br>
>     creating<br>
>      >     the mesh<br>
>      >      >      ><br>
>      >      >      >     master *:~/Downloads/tmp/Nicolas$<br>
>      >      >     /PETSc3/petsc/apple/bin/mpiexec -n<br>
>      >      >      >     2 ./test_overlap -initial_dm_view<br>
>     -dm_plex_box_faces 5,5<br>
>      >      >     -dm_distribute<br>
>      >      >      >     DM Object: Initial DM 2 MPI processes<br>
>      >      >      >        type: plex<br>
>      >      >      >     Initial DM in 2 dimensions:<br>
>      >      >      >        0-cells: 21 21<br>
>      >      >      >        1-cells: 45 45<br>
>      >      >      >        2-cells: 25 25<br>
>      >      >      >     Labels:<br>
>      >      >      >        depth: 3 strata with value/size (0 (21), 1<br>
>     (45), 2<br>
>      >     (25))<br>
>      >      >      >        celltype: 3 strata with value/size (0 (21),<br>
>     1 (45),<br>
>      >     3 (25))<br>
>      >      >      >        marker: 1 strata with value/size (1 (21))<br>
>      >      >      >        Face Sets: 1 strata with value/size (1 (10))<br>
>      >      >      ><br>
>      >      >      >     The get the same thing back with overlap = 0<br>
>      >      >      ><br>
>      >      >      >     master *:~/Downloads/tmp/Nicolas$<br>
>      >      >     /PETSc3/petsc/apple/bin/mpiexec -n<br>
>      >      >      >     2 ./test_overlap -initial_dm_view<br>
>     -dm_plex_box_faces 5,5<br>
>      >      >      >     -dm_distribute -dm_distribute_overlap 0<br>
>      >      >      >     DM Object: Initial DM 2 MPI processes<br>
>      >      >      >        type: plex<br>
>      >      >      >     Initial DM in 2 dimensions:<br>
>      >      >      >        0-cells: 21 21<br>
>      >      >      >        1-cells: 45 45<br>
>      >      >      >        2-cells: 25 25<br>
>      >      >      >     Labels:<br>
>      >      >      >        depth: 3 strata with value/size (0 (21), 1<br>
>     (45), 2<br>
>      >     (25))<br>
>      >      >      >        celltype: 3 strata with value/size (0 (21),<br>
>     1 (45),<br>
>      >     3 (25))<br>
>      >      >      >        marker: 1 strata with value/size (1 (21))<br>
>      >      >      >        Face Sets: 1 strata with value/size (1 (10))<br>
>      >      >      ><br>
>      >      >      >     and get larger local meshes with overlap = 1<br>
>      >      >      ><br>
>      >      >      >     master *:~/Downloads/tmp/Nicolas$<br>
>      >      >     /PETSc3/petsc/apple/bin/mpiexec -n<br>
>      >      >      >     2 ./test_overlap -initial_dm_view<br>
>     -dm_plex_box_faces 5,5<br>
>      >      >      >     -dm_distribute -dm_distribute_overlap 1<br>
>      >      >      >     DM Object: Initial DM 2 MPI processes<br>
>      >      >      >        type: plex<br>
>      >      >      >     Initial DM in 2 dimensions:<br>
>      >      >      >        0-cells: 29 29<br>
>      >      >      >        1-cells: 65 65<br>
>      >      >      >        2-cells: 37 37<br>
>      >      >      >     Labels:<br>
>      >      >      >        depth: 3 strata with value/size (0 (29), 1<br>
>     (65), 2<br>
>      >     (37))<br>
>      >      >      >        celltype: 3 strata with value/size (0 (29),<br>
>     1 (65),<br>
>      >     3 (37))<br>
>      >      >      >        marker: 1 strata with value/size (1 (27))<br>
>      >      >      >        Face Sets: 1 strata with value/size (1 (13))<br>
>      >      >      ><br>
>      >      >      >        Thanks,<br>
>      >      >      ><br>
>      >      >      >           Matt<br>
>      >      >      ><br>
>      >      >      >     On Wed, Mar 31, 2021 at 12:22 PM Nicolas Barral<br>
>      >      >      >     <<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>><br>
>      >      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>>>> wrote:<br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      >         @+<br>
>      >      >      ><br>
>      >      >      >         --<br>
>      >      >      >         Nicolas<br>
>      >      >      ><br>
>      >      >      >         On 31/03/2021 17:51, Matthew Knepley wrote:<br>
>      >      >      >          > On Sat, Mar 27, 2021 at 9:27 AM Nicolas<br>
>     Barral<br>
>      >      >      >          > <<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>><br>
>      >      >      >         <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>>><br>
>      >      >      >          ><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>><br>
>      >      >      >         <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>><br>
>      >      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>><br>
>      >     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a><br>
>     <mailto:<a href="mailto:nicolas.barral@math.u-bordeaux.fr" target="_blank">nicolas.barral@math.u-bordeaux.fr</a>>>>>>> wrote:<br>
>      >      >      >          ><br>
>      >      >      >          >     Hi all,<br>
>      >      >      >          ><br>
>      >      >      >          >     First, I'm not sure I understand<br>
>     what the<br>
>      >     overlap<br>
>      >      >      >         parameter in<br>
>      >      >      >          >     DMPlexDistributeOverlap does. I<br>
>     tried the<br>
>      >     following:<br>
>      >      >      >         generate a small<br>
>      >      >      >          >     mesh on 1 rank with<br>
>     DMPlexCreateBoxMesh, then<br>
>      >      >     distribute<br>
>      >      >      >         it with<br>
>      >      >      >          >     DMPlexDistribute. At this point I<br>
>     have two nice<br>
>      >      >      >         partitions, with shared<br>
>      >      >      >          >     vertices and no overlapping cells.<br>
>     Then I call<br>
>      >      >      >         DMPlexDistributeOverlap<br>
>      >      >      >          >     with the overlap parameter set to 0<br>
>     or 1,<br>
>      >     and get the<br>
>      >      >      >         same resulting<br>
>      >      >      >          >     plex in both cases. Why is that ?<br>
>      >      >      >          ><br>
>      >      >      >          ><br>
>      >      >      >          > The overlap parameter says how many cell<br>
>      >     adjacencies to go<br>
>      >      >      >         out. You<br>
>      >      >      >          > should not get the same<br>
>      >      >      >          > mesh out. We have lots of examples that use<br>
>      >     this. If<br>
>      >      >     you send<br>
>      >      >      >         your small<br>
>      >      >      >          > example, I can probably<br>
>      >      >      >          > tell you what is happening.<br>
>      >      >      >          ><br>
>      >      >      ><br>
>      >      >      >         Ok so I do have a small example on that and the<br>
>      >     DMClone<br>
>      >      >     thing I<br>
>      >      >      >         set up<br>
>      >      >      >         to understand! I attach it to the email.<br>
>      >      >      ><br>
>      >      >      >         For the overlap, you can change the overlap<br>
>      >     constant at<br>
>      >      >     the top<br>
>      >      >      >         of the<br>
>      >      >      >         file. With OVERLAP=0 or 1, the distributed<br>
>      >     overlapping mesh<br>
>      >      >      >         (shown using<br>
>      >      >      >         -over_dm_view, it's DMover) are the same, and<br>
>      >     different<br>
>      >      >     from the<br>
>      >      >      >         mesh<br>
>      >      >      >         before distributing the overlap (shown using<br>
>      >      >     -distrib_dm_view). For<br>
>      >      >      >         larger overlap values they're different.<br>
>      >      >      ><br>
>      >      >      >         The process is:<br>
>      >      >      >         1/ create a DM dm on 1 rank<br>
>      >      >      >         2/ clone dm into dm2<br>
>      >      >      >         3/ distribute dm<br>
>      >      >      >         4/ clone dm into dm3<br>
>      >      >      >         5/ distribute dm overlap<br>
>      >      >      ><br>
>      >      >      >         I print all the DMs after each step. dm has a<br>
>      >     distributed<br>
>      >      >      >         overlap, dm2<br>
>      >      >      >         is not distributed, dm3 is distributed but<br>
>     without<br>
>      >      >     overlap. Since<br>
>      >      >      >         distribute and distributeOverlap create new<br>
>     DMs, I<br>
>      >     don't seem<br>
>      >      >      >         have a<br>
>      >      >      >         problem with the shallow copies.<br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      >          >     Second, I'm wondering what would be<br>
>     a good<br>
>      >     way to<br>
>      >      >     handle<br>
>      >      >      >         two overlaps<br>
>      >      >      >          >     and associated local vectors. In my<br>
>     adaptation<br>
>      >      >     code, the<br>
>      >      >      >         remeshing<br>
>      >      >      >          >     library requires a non-overlapping mesh,<br>
>      >     while the<br>
>      >      >      >         refinement criterion<br>
>      >      >      >          >     computation is based on hessian<br>
>      >     computations, which<br>
>      >      >      >         require a layer of<br>
>      >      >      >          >     overlap. What I can do is clone the<br>
>     dm before<br>
>      >      >      >         distributing the overlap,<br>
>      >      >      >          >     then manage two independent plex<br>
>     objects with<br>
>      >      >     their own<br>
>      >      >      >         local sections<br>
>      >      >      >          >     etc. and copy/trim local vectors<br>
>     manually.<br>
>      >     Is there a<br>
>      >      >      >         more automatic<br>
>      >      >      >          >     way<br>
>      >      >      >          >     to do this ?<br>
>      >      >      >          ><br>
>      >      >      >          ><br>
>      >      >      >          > DMClone() is a shallow copy, so that<br>
>     will not work.<br>
>      >      >     You would<br>
>      >      >      >         maintain<br>
>      >      >      >          > two different Plexes, overlapping<br>
>      >      >      >          > and non-overlapping, with their own sections<br>
>      >     and vecs. Are<br>
>      >      >      >         you sure you<br>
>      >      >      >          > need to keep around the non-overlapping one?<br>
>      >      >      >          > Maybe if I understood what operations<br>
>     you want<br>
>      >     to work, I<br>
>      >      >      >         could say<br>
>      >      >      >          > something more definitive.<br>
>      >      >      >          ><br>
>      >      >      >         I need to be able to pass the<br>
>     non-overlapping mesh<br>
>      >     to the<br>
>      >      >      >         remesher. I<br>
>      >      >      >         can either maintain 2 plexes, or trim the<br>
>     overlapping<br>
>      >      >     plex when<br>
>      >      >      >         I create<br>
>      >      >      >         the arrays I give to the remesher. I'm not sure<br>
>      >     which is the<br>
>      >      >      >         best/worst ?<br>
>      >      >      ><br>
>      >      >      >         Thanks<br>
>      >      >      ><br>
>      >      >      >         --<br>
>      >      >      >         Nicolas<br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      >          >    Thanks,<br>
>      >      >      >          ><br>
>      >      >      >          >       Matt<br>
>      >      >      >          ><br>
>      >      >      >          >     Thanks<br>
>      >      >      >          ><br>
>      >      >      >          >     --<br>
>      >      >      >          >     Nicolas<br>
>      >      >      >          ><br>
>      >      >      >          ><br>
>      >      >      >          ><br>
>      >      >      >          > --<br>
>      >      >      >          > What most experimenters take for granted<br>
>     before<br>
>      >     they<br>
>      >      >     begin their<br>
>      >      >      >          > experiments is infinitely more<br>
>     interesting than any<br>
>      >      >     results<br>
>      >      >      >         to which<br>
>      >      >      >          > their experiments lead.<br>
>      >      >      >          > -- Norbert Wiener<br>
>      >      >      >          ><br>
>      >      >      >          > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
>      >      >      >         <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      >     --<br>
>      >      >      >     What most experimenters take for granted before<br>
>     they<br>
>      >     begin their<br>
>      >      >      >     experiments is infinitely more interesting than any<br>
>      >     results<br>
>      >      >     to which<br>
>      >      >      >     their experiments lead.<br>
>      >      >      >     -- Norbert Wiener<br>
>      >      >      ><br>
>      >      >      > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
>      >      >      >     <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      ><br>
>      >      >      > --<br>
>      >      >      > What most experimenters take for granted before<br>
>     they begin<br>
>      >     their<br>
>      >      >      > experiments is infinitely more interesting than any<br>
>      >     results to which<br>
>      >      >      > their experiments lead.<br>
>      >      >      > -- Norbert Wiener<br>
>      >      >      ><br>
>      >      >      > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
>      >      >     <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
>      >      ><br>
>      >      ><br>
>      >      ><br>
>      >      > --<br>
>      >      > What most experimenters take for granted before they begin<br>
>     their<br>
>      >      > experiments is infinitely more interesting than any<br>
>     results to which<br>
>      >      > their experiments lead.<br>
>      >      > -- Norbert Wiener<br>
>      >      ><br>
>      >      > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
>      >     <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
>      ><br>
>      ><br>
>      ><br>
>      > --<br>
>      > What most experimenters take for granted before they begin their<br>
>      > experiments is infinitely more interesting than any results to which<br>
>      > their experiments lead.<br>
>      > -- Norbert Wiener<br>
>      ><br>
>      > <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
>     <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
> <br>
> <br>
> <br>
> -- <br>
> What most experimenters take for granted before they begin their <br>
> experiments is infinitely more interesting than any results to which <br>
> their experiments lead.<br>
> -- Norbert Wiener<br>
> <br>
> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a> <<a href="http://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div><div><br></div><div><a href="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>