[petsc-users] Agglomeration for Multigrid on Unstructured Meshes

[Lawrence Mitchell]

> On 2 Jun 2020, at 09:35, Matthew Knepley <knepley at gmail.com> wrote:
> 
> On Tue, Jun 2, 2020 at 4:25 AM Lawrence Mitchell <wencel at gmail.com> wrote:
> Hi Dave,
> 
> > On 2 Jun 2020, at 05:43, Dave May <dave.mayhem23 at gmail.com> wrote:
> > 
> > 
> > 
> > On Tue 2. Jun 2020 at 03:30, Matthew Knepley <knepley at gmail.com> wrote:
> > On Mon, Jun 1, 2020 at 7:03 PM Danyang Su <danyang.su at gmail.com> wrote:
> > Thanks Jed for the quick response. Yes I am asking about the repartitioning of coarse grids in geometric multigrid for unstructured mesh. I am happy with AMG. Thanks for letting me know.
> > 
> > All the pieces are there, we just have not had users asking for this, and it will take some work to put together.
> > 
> > Matt - I created a branch for you and Lawrence last year which added full support for PLEX within Telescope. This implementation was not a fully automated algmoeration strategy - it utilized the partition associated with the DM returned from DMGetCoarseDM. Hence the job of building the distributed coarse hierarchy was let to the user.
> > 
> > I’m pretty sure that code got merged into master as the branch also contained several bug mixes for Telescope. Or am I mistaken?
> 
> I think you're right. I didn't manage to get the redistribution of the DMPlex object done last summer (it's bubbling up again).
> 
> As I see it, for redistributed geometric multigrid on plexes, the missing piece is a function:
> 
> DMPlexRedistributeOntoComm(DM old, MPI_Comm comm, DM *new)
> 
> I went down a rabbit hole of trying to do this, since I actually think this replaced the current interface to DMPlexDistribute, which is
> 
> DMPlexDistribute(DM old, PetscInt overlap, PetscSF *pointDistSF, DM *new)
> 
> Where the new DM comes out on the same communicator as the old DM, just with a different partition.
> 
> This has lots of follow-on consequences, for example, one can't easily load on P processes and then compute on Q.
> 
> Unfortunately, collectiveness over MPI_Comm(old) is baked into the redistribution routines everywhere, and I didn't manage to finish things.
> 
> Yes, I remember thinking this out. I believe the conclusion was that redistribution should happen on the large comm, which
> some fraction of processes getting no cells. Then at the end we call one new function which copies that DM onto the smaller comm.

I think this kind of works, but I worry that it makes the following usecase really painful.

1. Parallel load on P processes
2. Redistribute onto Q > P processes
... Compute

3. (Occasionally), checkpoint by mapping onto P processes and dumping

Perhaps I don't understand your proposal properly, but I think you're saying that the _input_ DM in redistribute defines the processes that participate in the redistribution (and the _output_ DM just controls the size).

This would then require me in my example to load on Q processes (not P) or else do some futzing around.

Or do I have the wrong end of the stick?

I think I would prefer an interface where I go from

commOld -> commNew

with the number of partitions determined by size(commNew)

Then the restriction you have on the two communicators is one of:

i. commOld \subseteq commNew

ii. commOld \supseteq commNew

I guess if you wanted to be really fancy you could make an intercommunicator of commOld \cap commNew, but now you have to carry another communicator around with you for data transfer etc...

Perhaps internally you do this by making a DM on whichever the smaller of the two communicators is that is empty on a bunch of ranks and then copying in (or copying out) appropriately.

I think a first step here is to break the assumption everywhere in the redistribution interface that the size of commOld determines the number of partitions you want, and instead move to using size(commNew). Then for now we just pass commNew=commOld.

Lawrence