[petsc-dev] Recent error checking in DMCreateDefaultSF breaks parallel DMs

Lawrence Mitchell lawrence.mitchell at imperial.ac.uk
Wed May 28 10:16:01 CDT 2014


1cc89ce (DM: Added check for erroneous input when creating default SF)

Checks that the constrained point does not constrain more dofs than the
unconstrained point has.  However, this breaks in parallel, since
PetscSectionCreateGlobalSection says:

"  Note: This gives negative sizes and offsets to points not owned by
this process"

Hence the following example breakage:

$ cat > foo.py <<\EOF
import petsc4py
import sys
petsc4py.init(sys.argv)
from petsc4py import PETSc


bdy = PETSc.DMPlex().create()
bdy.setDimension(1)
bdy.createSquareBoundary([0, 0], [1, 1], [3, 3])

dm = PETSc.DMPlex().generate(bdy)

dm.distribute()
sec = dm.createSection(1, [1], [1, 0, 0])

dm.setDefaultSection(sec)
gsec = dm.getDefaultGlobalSection()

dm.createDefaultSF(sec, gsec)
EOF

$ mpiexec -n 2 python foo.py
Traceback (most recent call last):
  File "foo.py", line 19, in <module>
    dm.createDefaultSF(sec, gsec)
  File "DM.pyx", line 268, in petsc4py.PETSc.DM.createDefaultSF
(src/petsc4py.PETSc.c:179320)
petsc4py.PETSc.Error: error code 63
[0] DMCreateDefaultSF() line 3065 in
/data/lmitche1/src/petsc/src/dm/interface/dm.c
[0] Argument out of range
[0] Point 9 has 0 constraints > -2 dof


Maybe this patch?

diff --git a/src/dm/interface/dm.c b/src/dm/interface/dm.c
index d830c0c..71bf7ad 100644
--- a/src/dm/interface/dm.c
+++ b/src/dm/interface/dm.c
@@ -3062,7 +3062,7 @@ PetscErrorCode DMCreateDefaultSF(DM dm,
PetscSection localSection, PetscSection

     ierr     = PetscSectionGetDof(globalSection, p, &gdof);CHKERRQ(ierr);
     ierr     = PetscSectionGetConstraintDof(globalSection, p,
&gcdof);CHKERRQ(ierr);
-    if (gcdof > gdof) SETERRQ3(PETSC_COMM_SELF,
PETSC_ERR_ARG_OUTOFRANGE, "Point %d has %d constraints > %d dof", p,
gcdof, gdof);
+    if (gcdof && gcdof > gdof) SETERRQ3(PETSC_COMM_SELF,
PETSC_ERR_ARG_OUTOFRANGE, "Point %d has %d constraints > %d dof", p,
gcdof, gdof);
     nleaves += gdof < 0 ? -(gdof+1)-gcdof : gdof-gcdof;
   }
   ierr = PetscMalloc1(nleaves, &local);CHKERRQ(ierr);

Cheers,

Lawrence



More information about the petsc-dev mailing list