Hi Matt,
I think there is a misunderstanding here. I am referring to the case
where DMPlexDistribute() is run with overlap=1 (which is not the case in
SNES ex12) and vertices in the overlap/halo region are assigned to the
wrong rank. This can lead to a case where a proc may own a vertex that
is not in its original (non-overlapping) partition, although the
attached cell is not owned and will be marked as "ghost" by
DMPlexConstructGhostCells().
To illustrate this, I have attached an example consisting of a unit
square with 3 faces in each dimension and a section with only vertex
dofs. If run with two ranks, rank 1 will own all its vertices (13
roots), whereas rank 0 only owns vertices not in the overlap/halo of
rank 1 (3 roots). My understanding is that, since the original partition
splits the square along its diagonal, the vertex distribution should be
10 to 6 with the 4 diagonal vertices assigned to rank 1 and all other
vertices assigned according to the original partition. Is this correct,
or am I missing something here?
Many thanks for all your help
Michael
On 16/11/13 13:54, Matthew Knepley wrote:
On Sat, Nov 16, 2013 at 7:22 AM, Michael Lange
<[email protected] <mailto:[email protected]>>
wrote:
Hi,
I notice that, when creating the point SF for the parallel
partition in DMPlexDistribute, cells are assigned to procs
according to the original partition but vertices aren't. Was this
done by design or is this a bug?
If this were true, there would be no communication for the P1 test of
SNES ex12. Here is running it with
-interpolate 1 and -dm_view ::ascii_info_detail
PetscSF Object: 2 MPI processes
type: basic
sort=rank-order
[0] Number of roots=19, leaves=5, remote ranks=1
[0] 4 <- (1,6)
[0] 5 <- (1,8)
[0] 7 <- (1,9)
[0] 10 <- (1,13)
[0] 11 <- (1,17)
[1] Number of roots=19, leaves=0, remote ranks=0
[0] Roots referenced by my leaves, by rank
[0] 1: 5 edges
[0] 4 <- 6
[0] 5 <- 8
[0] 7 <- 9
[0] 10 <- 13
[0] 11 <- 17
[1] Roots referenced by my leaves, by rank
So there are 3 vertices and 2 edges in the point SF.
Matt
In case it is a bug, I have attached a patch that fixes this by
using the closure of the original partition instead.
Thanks and kind regards
Michael
--
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
static char help[] = "Small test problem. \n";
#include <petscdmplex.h>
#include <petscsnes.h>
#undef __FUNCT__
#define __FUNCT__ "main"
int main(int argc, char **argv)
{
DM dm,boundary,pardm,ghostdm;
PetscSF sf;
PetscSection section;
MPI_Comm comm=MPI_COMM_WORLD;
PetscInt comm_size,comm_rank,numGhostCells;
DMLabel lbl;
PetscErrorCode ierr;
ierr = PetscInitialize(&argc, &argv, NULL, help);CHKERRQ(ierr);
ierr = MPI_Comm_size(comm, &comm_size);CHKERRQ(ierr);
ierr = MPI_Comm_rank(comm, &comm_rank);CHKERRQ(ierr);
PetscInt dim = 2, overlap = 1;
PetscReal lower[2] = {0.0, 0.0};
PetscReal upper[2] = {1.0, 1.0};
PetscInt faces[2] = {3, 3};
PetscInt comp[1] = { 1 };
PetscInt dof[4] = {1, 0, 0, 0};
/* Build a unit mesh (3x3 square) */
ierr = DMPlexCreate(comm, &boundary);CHKERRQ(ierr);
ierr = DMSetType(boundary, DMPLEX);CHKERRQ(ierr);
ierr = DMPlexSetDimension(boundary, dim-1);CHKERRQ(ierr);
ierr = DMPlexCreateSquareBoundary(boundary, lower, upper, faces);CHKERRQ(ierr);
ierr = DMPlexGenerate(boundary, NULL, PETSC_TRUE, &dm);CHKERRQ(ierr);
ierr = DMDestroy(&boundary);CHKERRQ(ierr);
/* Label exterior facets */
ierr = DMPlexCreateLabel(dm, "exterior_facets");CHKERRQ(ierr);
ierr = DMPlexGetLabel(dm, "exterior_facets", &lbl);CHKERRQ(ierr);
ierr = DMPlexMarkBoundaryFaces(dm, lbl);CHKERRQ(ierr);
if (comm_size > 1) {
ierr = DMPlexDistribute(dm,NULL,overlap,&sf,&pardm);CHKERRQ(ierr);
if (pardm) {
ierr = DMDestroy(&dm);CHKERRQ(ierr);
dm = pardm;
}
/* Mark ghost cells */
ierr = DMPlexConstructGhostCells(dm, "exterior_facets", &numGhostCells, &ghostdm);CHKERRQ(ierr);
if (ghostdm) {
ierr = DMDestroy(&dm);CHKERRQ(ierr);
dm = ghostdm;
}
ierr = DMView(dm, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
}
ierr = DMPlexCreateSection(dm, dim, 1, comp, dof, 0, NULL, NULL, §ion);CHKERRQ(ierr);
ierr = DMSetDefaultSection(dm, section);CHKERRQ(ierr);
ierr = DMGetDefaultSF(dm, &sf);CHKERRQ(ierr);CHKERRQ(ierr);
ierr = PetscSFView(sf, PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr);
ierr = PetscFinalize();
return 0;
}
