[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