[petsc-dev] Vertex assignment in DMPlexDistribute

[Michael Lange]

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?

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
>From 7a1aeca70d3c7454780576a35dc1ec09a4c2d479 Mon Sep 17 00:00:00 2001
From: Michael Lange <michael.la...@imperial.ac.uk>
Date: Sat, 16 Nov 2013 13:03:00 +0000
Subject: [PATCH] DMPlex: Fixing vertex assignment in DMPlexDistribute

When distributing a plex, only cells in the original partition are
assigned to the local proc. This commit extends this to vertices
by using the closure of the original partition.
---
 src/dm/impls/plex/plex.c |   21 ++++++++++++---------
 1 file changed, 12 insertions(+), 9 deletions(-)

diff --git a/src/dm/impls/plex/plex.c b/src/dm/impls/plex/plex.c
index eca507d..a8c8554 100644
--- a/src/dm/impls/plex/plex.c
+++ b/src/dm/impls/plex/plex.c
@@ -2909,8 +2909,8 @@ PetscErrorCode DMPlexDistribute(DM dm, const char partitioner[], PetscInt overla
   MPI_Comm               comm;
   const PetscInt         height = 0;
   PetscInt               dim, numRemoteRanks;
-  IS                     origCellPart,        cellPart,        part;
-  PetscSection           origCellPartSection, cellPartSection, partSection;
+  IS                     origCellPart,        origPart,        cellPart,        part;
+  PetscSection           origCellPartSection, origPartSection, cellPartSection, partSection;
   PetscSFNode           *remoteRanks;
   PetscSF                partSF, pointSF, coneSF;
   ISLocalToGlobalMapping renumbering;
@@ -3228,23 +3228,26 @@ PetscErrorCode DMPlexDistribute(DM dm, const char partitioner[], PetscInt overla
       rowners[p].index = -1;
     }
     if (origCellPart) {
-      /* Make sure cells in the original partition are not assigned to other procs */
-      const PetscInt *origCells;
+      /* Make sure points in the original partition are not assigned to other procs */
+      const PetscInt *origPoints;
 
-      ierr = ISGetIndices(origCellPart, &origCells);CHKERRQ(ierr);
+      ierr = DMPlexCreatePartitionClosure(dm, origCellPartSection, origCellPart, &origPartSection, &origPart);CHKERRQ(ierr);
+      ierr = ISGetIndices(origPart, &origPoints);CHKERRQ(ierr);
       for (p = 0; p < numProcs; ++p) {
         PetscInt dof, off, d;
 
-        ierr = PetscSectionGetDof(origCellPartSection, p, &dof);CHKERRQ(ierr);
-        ierr = PetscSectionGetOffset(origCellPartSection, p, &off);CHKERRQ(ierr);
+        ierr = PetscSectionGetDof(origPartSection, p, &dof);CHKERRQ(ierr);
+        ierr = PetscSectionGetOffset(origPartSection, p, &off);CHKERRQ(ierr);
         for (d = off; d < off+dof; ++d) {
-          rowners[origCells[d]].rank = p;
+          rowners[origPoints[d]].rank = p;
         }
       }
-      ierr = ISRestoreIndices(origCellPart, &origCells);CHKERRQ(ierr);
+      ierr = ISRestoreIndices(origPart, &origPoints);CHKERRQ(ierr);
     }
     ierr = ISDestroy(&origCellPart);CHKERRQ(ierr);
+    ierr = ISDestroy(&origPart);CHKERRQ(ierr);
     ierr = PetscSectionDestroy(&origCellPartSection);CHKERRQ(ierr);
+    ierr = PetscSectionDestroy(&origPartSection);CHKERRQ(ierr);
 
     ierr = PetscSFBcastBegin(pointSF, MPIU_2INT, rowners, lowners);CHKERRQ(ierr);
     ierr = PetscSFBcastEnd(pointSF, MPIU_2INT, rowners, lowners);CHKERRQ(ierr);
-- 
1.7.9.5