On Thu, 16 Feb 2023 at 16:43, Matthew Knepley <knep...@gmail.com> wrote: > > On Thu, Feb 16, 2023 at 10:54 AM Lawrence Mitchell <we...@gmx.li> wrote: >> >> Hi Blaise, >> >> On Thu, 16 Feb 2023 at 15:17, Blaise Bourdin <bour...@mcmaster.ca> wrote: >> > >> > Hi, >> > >> > I am trying to implement a non-local finite elements reconstruction >> > operator in parallel. >> > >> > Given a dmplex distributed with an overlap, is there a way to figure out >> > which cells are in the overlap and which are not? >> >> Yes. Get the pointSF of the DM, and the cell chart >> >> DMPlexGetPointSF(dm, &sf); >> DMPlexGetHeightStratum(dm, 0, &cstart, &cend); >> >> Now get the graph (specifically ilocal of the sf): >> >> PetscSFGetGraph(sf, NULL, &nleaves, &ilocal, NULL); >> >> Now any value of ilocal that lies in [cstart, cend) is a cell which is >> not owned by this process (i.e. in the overlap). Note that ilocal can >> be NULL which just means it is the identity map [0, ..., nleaves), so >> you just intersect [cstart, cend) with [0, nleaves) in that case to >> find the overlap cells. >> >> But that is very unlikely to be true, so: > > > Note that you can use > > PetscFindInt(nleaves, ilocal, cell, &loc);
Modulo argument order, is it not PetscFindInt(cell, nleaves, ilocal, &loc)? I guess one should probably document that PetscSFSetGraph ensures that ilocal is sorted. Lawrence