Commits

Jed Brown committed bb38c0a Merge

Merge branch 'jed/plex-f90-null-object'

* jed/plex-f90-null-object:
DMPlexMatSetClosure: allow NULL sections to use defaults
DMPlex Fortran: allow PETSC_NULL_OBJECT for sections
DMPlex F90: rename __ierr to ierr for use with macros

  • Participants
  • Parent commits 4098462, b4a4644

Comments (0)

Files changed (2)

File src/dm/impls/plex/f90-custom/zplexf90.c

 
 /* Definitions of Fortran Wrapper routines */
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetcone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetcone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   const PetscInt *v;
   PetscInt       n;
 
-  *__ierr = DMPlexGetConeSize(*dm, *p, &n);if (*__ierr) return;
-  *__ierr = DMPlexGetCone(*dm, *p, &v);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
+  *ierr = DMPlexGetConeSize(*dm, *p, &n);if (*ierr) return;
+  *ierr = DMPlexGetCone(*dm, *p, &v);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestorecone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestorecone_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
-  *__ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
+  *ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   const PetscInt *v;
   PetscInt       n;
 
-  *__ierr = DMPlexGetConeSize(*dm, *p, &n);if (*__ierr) return;
-  *__ierr = DMPlexGetConeOrientation(*dm, *p, &v);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
+  *ierr = DMPlexGetConeSize(*dm, *p, &n);if (*ierr) return;
+  *ierr = DMPlexGetConeOrientation(*dm, *p, &v);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestoreconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestoreconeorientation_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
-  *__ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
+  *ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetsupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetsupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   const PetscInt *v;
   PetscInt       n;
 
-  *__ierr = DMPlexGetSupportSize(*dm, *p, &n);if (*__ierr) return;
-  *__ierr = DMPlexGetSupport(*dm, *p, &v);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
+  *ierr = DMPlexGetSupportSize(*dm, *p, &n);if (*ierr) return;
+  *ierr = DMPlexGetSupport(*dm, *p, &v);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestoresupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestoresupport_(DM *dm, PetscInt *p, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
-  *__ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
+  *ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgettransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgettransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscInt *v = NULL;
   PetscInt n;
 
-  *__ierr = DMPlexGetTransitiveClosure(*dm, *p, *useCone, &n, &v);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n*2, ptr PETSC_F90_2PTR_PARAM(ptrd));
+  *ierr = DMPlexGetTransitiveClosure(*dm, *p, *useCone, &n, &v);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) v, PETSC_INT, 1, n*2, ptr PETSC_F90_2PTR_PARAM(ptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestoretransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestoretransitiveclosure_(DM *dm, PetscInt *p, PetscBool *useCone, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscInt *array;
 
-  *__ierr = F90Array1dAccess(ptr, PETSC_INT, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
-  *__ierr = DMPlexRestoreTransitiveClosure(*dm, *p, *useCone, NULL, &array);if (*__ierr) return;
-  *__ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
+  *ierr = F90Array1dAccess(ptr, PETSC_INT, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
+  *ierr = DMPlexRestoreTransitiveClosure(*dm, *p, *useCone, NULL, &array);if (*ierr) return;
+  *ierr = F90Array1dDestroy(ptr, PETSC_INT PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexvecgetclosure_(DM *dm, PetscSection *section, Vec *x, PetscInt *point, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexvecgetclosure_(DM *dm, PetscSection *section, Vec *x, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscScalar *v = NULL;
   PetscInt     n;
 
-  *__ierr = DMPlexVecGetClosure(*dm, *section, *x, *point, &n, &v);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) v, PETSC_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
+  CHKFORTRANNULLOBJECT(section);
+  *ierr = DMPlexVecGetClosure(*dm, *section, *x, *point, &n, &v);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) v, PETSC_SCALAR, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexvecrestoreclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexvecrestoreclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscScalar *array;
 
-  *__ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void **) &array PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
-  *__ierr = DMPlexVecRestoreClosure(*dm, *section, *v, *point, NULL, &array);if (*__ierr) return;
-  *__ierr = F90Array1dDestroy(ptr, PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
+  CHKFORTRANNULLOBJECT(section);
+  *ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void **) &array PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
+  *ierr = DMPlexVecRestoreClosure(*dm, *section, *v, *point, NULL, &array);if (*ierr) return;
+  *ierr = F90Array1dDestroy(ptr, PETSC_SCALAR PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexvecsetclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexvecsetclosure_(DM *dm, PetscSection *section, Vec *v, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscScalar *array;
 
-  *__ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
-  *__ierr = DMPlexVecSetClosure(*dm, *section, *v, *point, array, *mode);
+  CHKFORTRANNULLOBJECT(section);
+  *ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
+  *ierr = DMPlexVecSetClosure(*dm, *section, *v, *point, array, *mode);
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexmatsetclosure_(DM *dm, PetscSection *section, PetscSection *globalSection, Mat *A, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexmatsetclosure_(DM *dm, PetscSection *section, PetscSection *globalSection, Mat *A, PetscInt *point, F90Array1d *ptr, InsertMode *mode, int *ierr PETSC_F90_2PTR_PROTO(ptrd))
 {
   PetscScalar *array;
 
-  *__ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*__ierr) return;
-  *__ierr = DMPlexMatSetClosure(*dm, *section, *globalSection, *A, *point, array, *mode);
+  CHKFORTRANNULLOBJECT(section);
+  CHKFORTRANNULLOBJECT(globalSection);
+  *ierr = F90Array1dAccess(ptr, PETSC_SCALAR, (void**) &array PETSC_F90_2PTR_PARAM(ptrd));if (*ierr) return;
+  *ierr = DMPlexMatSetClosure(*dm, *section, *globalSection, *A, *point, array, *mode);
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetjoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetjoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt       *points;
   const PetscInt *coveredPoints;
   PetscInt       numCoveredPoints;
 
-  *__ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*__ierr) return;
-  *__ierr = DMPlexGetJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
+  *ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*ierr) return;
+  *ierr = DMPlexGetJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetfulljoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetfulljoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt       *points;
   const PetscInt *coveredPoints;
   PetscInt        numCoveredPoints;
 
-  *__ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*__ierr) return;
-  *__ierr = DMPlexGetFullJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
+  *ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*ierr) return;
+  *ierr = DMPlexGetFullJoin(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestorejoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestorejoin_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt *coveredPoints;
 
-  *__ierr = F90Array1dAccess(cptr, PETSC_INT, (void**) &coveredPoints PETSC_F90_2PTR_PARAM(cptrd));if (*__ierr) return;
-  *__ierr = DMPlexRestoreJoin(*dm, 0, NULL, NULL, (const PetscInt**) &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dDestroy(cptr, PETSC_INT PETSC_F90_2PTR_PARAM(cptrd));if (*__ierr) return;
+  *ierr = F90Array1dAccess(cptr, PETSC_INT, (void**) &coveredPoints PETSC_F90_2PTR_PARAM(cptrd));if (*ierr) return;
+  *ierr = DMPlexRestoreJoin(*dm, 0, NULL, NULL, (const PetscInt**) &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dDestroy(cptr, PETSC_INT PETSC_F90_2PTR_PARAM(cptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt       *points;
   const PetscInt *coveredPoints;
   PetscInt       numCoveredPoints;
 
-  *__ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*__ierr) return;
-  *__ierr = DMPlexGetMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
+  *ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*ierr) return;
+  *ierr = DMPlexGetMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexgetfullmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexgetfullmeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt       *points;
   const PetscInt *coveredPoints;
   PetscInt       numCoveredPoints;
 
-  *__ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*__ierr) return;
-  *__ierr = DMPlexGetFullMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
+  *ierr = F90Array1dAccess(pptr, PETSC_INT, (void**) &points PETSC_F90_2PTR_PARAM(pptrd));if (*ierr) return;
+  *ierr = DMPlexGetFullMeet(*dm, *numPoints, points, &numCoveredPoints, &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dCreate((void*) coveredPoints, PETSC_INT, 1, numCoveredPoints, cptr PETSC_F90_2PTR_PARAM(cptrd));
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexrestoremeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *__ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
+PETSC_EXTERN void PETSC_STDCALL dmplexrestoremeet_(DM *dm, PetscInt *numPoints, F90Array1d *pptr, F90Array1d *cptr, int *ierr PETSC_F90_2PTR_PROTO(pptrd) PETSC_F90_2PTR_PROTO(cptrd))
 {
   PetscInt *coveredPoints;
 
-  *__ierr = F90Array1dAccess(cptr, PETSC_INT, (void**) &coveredPoints PETSC_F90_2PTR_PARAM(cptrd));if (*__ierr) return;
-  *__ierr = DMPlexRestoreMeet(*dm, 0, NULL, NULL, (const PetscInt**) &coveredPoints);if (*__ierr) return;
-  *__ierr = F90Array1dDestroy(cptr, PETSC_INT PETSC_F90_2PTR_PARAM(cptrd));if (*__ierr) return;
+  *ierr = F90Array1dAccess(cptr, PETSC_INT, (void**) &coveredPoints PETSC_F90_2PTR_PARAM(cptrd));if (*ierr) return;
+  *ierr = DMPlexRestoreMeet(*dm, 0, NULL, NULL, (const PetscInt**) &coveredPoints);if (*ierr) return;
+  *ierr = F90Array1dDestroy(cptr, PETSC_INT PETSC_F90_2PTR_PARAM(cptrd));if (*ierr) return;
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexcreatesection_(DM *dm, PetscInt *dim, PetscInt *numFields, F90Array1d *ptrC, F90Array1d *ptrD, PetscInt *numBC, F90Array1d *ptrF, F90Array1d *ptrP, PetscSection *section, int *__ierr PETSC_F90_2PTR_PROTO(ptrCd) PETSC_F90_2PTR_PROTO(ptrDd) PETSC_F90_2PTR_PROTO(ptrFd) PETSC_F90_2PTR_PROTO(ptrPd))
+PETSC_EXTERN void PETSC_STDCALL dmplexcreatesection_(DM *dm, PetscInt *dim, PetscInt *numFields, F90Array1d *ptrC, F90Array1d *ptrD, PetscInt *numBC, F90Array1d *ptrF, F90Array1d *ptrP, PetscSection *section, int *ierr PETSC_F90_2PTR_PROTO(ptrCd) PETSC_F90_2PTR_PROTO(ptrDd) PETSC_F90_2PTR_PROTO(ptrFd) PETSC_F90_2PTR_PROTO(ptrPd))
 {
   PetscInt *numComp;
   PetscInt *numDof;
   PetscInt *bcField;
   IS       *bcPoints;
 
-  *__ierr = F90Array1dAccess(ptrC, PETSC_INT, (void**) &numComp PETSC_F90_2PTR_PARAM(ptrCd));if (*__ierr) return;
-  *__ierr = F90Array1dAccess(ptrD, PETSC_INT, (void**) &numDof PETSC_F90_2PTR_PARAM(ptrDd));if (*__ierr) return;
-  *__ierr = F90Array1dAccess(ptrF, PETSC_INT, (void**) &bcField PETSC_F90_2PTR_PARAM(ptrFd));if (*__ierr) return;
-  *__ierr = F90Array1dAccess(ptrP, PETSC_FORTRANADDR, (void**) &bcPoints PETSC_F90_2PTR_PARAM(ptrPd));if (*__ierr) return;
-  *__ierr = DMPlexCreateSection(*dm, *dim, *numFields, numComp, numDof, *numBC, bcField, bcPoints, section);
+  *ierr = F90Array1dAccess(ptrC, PETSC_INT, (void**) &numComp PETSC_F90_2PTR_PARAM(ptrCd));if (*ierr) return;
+  *ierr = F90Array1dAccess(ptrD, PETSC_INT, (void**) &numDof PETSC_F90_2PTR_PARAM(ptrDd));if (*ierr) return;
+  *ierr = F90Array1dAccess(ptrF, PETSC_INT, (void**) &bcField PETSC_F90_2PTR_PARAM(ptrFd));if (*ierr) return;
+  *ierr = F90Array1dAccess(ptrP, PETSC_FORTRANADDR, (void**) &bcPoints PETSC_F90_2PTR_PARAM(ptrPd));if (*ierr) return;
+  *ierr = DMPlexCreateSection(*dm, *dim, *numFields, numComp, numDof, *numBC, bcField, bcPoints, section);
 }
 
-PETSC_EXTERN void PETSC_STDCALL dmplexcomputecellgeometry_(DM *dm, PetscInt *cell, F90Array1d *ptrV, F90Array1d *ptrJ, F90Array1d *ptrIJ, PetscReal *detJ, int *__ierr PETSC_F90_2PTR_PROTO(ptrVd) PETSC_F90_2PTR_PROTO(ptrJd) PETSC_F90_2PTR_PROTO(ptrIJd))
+PETSC_EXTERN void PETSC_STDCALL dmplexcomputecellgeometry_(DM *dm, PetscInt *cell, F90Array1d *ptrV, F90Array1d *ptrJ, F90Array1d *ptrIJ, PetscReal *detJ, int *ierr PETSC_F90_2PTR_PROTO(ptrVd) PETSC_F90_2PTR_PROTO(ptrJd) PETSC_F90_2PTR_PROTO(ptrIJd))
 {
   PetscReal *v0;
   PetscReal *J;
   PetscReal *invJ;
 
-  *__ierr = F90Array1dAccess(ptrV,  PETSC_REAL, (void**) &v0 PETSC_F90_2PTR_PARAM(ptrVd));if (*__ierr) return;
-  *__ierr = F90Array1dAccess(ptrJ,  PETSC_REAL, (void**) &J PETSC_F90_2PTR_PARAM(ptrJd));if (*__ierr) return;
-  *__ierr = F90Array1dAccess(ptrIJ, PETSC_REAL, (void**) &invJ PETSC_F90_2PTR_PARAM(ptrIJd));if (*__ierr) return;
-  *__ierr = DMPlexComputeCellGeometry(*dm, *cell, v0, J, invJ, detJ);
+  *ierr = F90Array1dAccess(ptrV,  PETSC_REAL, (void**) &v0 PETSC_F90_2PTR_PARAM(ptrVd));if (*ierr) return;
+  *ierr = F90Array1dAccess(ptrJ,  PETSC_REAL, (void**) &J PETSC_F90_2PTR_PARAM(ptrJd));if (*ierr) return;
+  *ierr = F90Array1dAccess(ptrIJ, PETSC_REAL, (void**) &invJ PETSC_F90_2PTR_PARAM(ptrIJd));if (*ierr) return;
+  *ierr = DMPlexComputeCellGeometry(*dm, *cell, v0, J, invJ, detJ);
 }
 

File src/dm/impls/plex/plex.c

 
   Input Parameters:
 + dm - The DM
-. section - The section describing the layout in v
-. globalSection - The section describing the layout in v
+. section - The section describing the layout in v, or NULL to use the default section
+. globalSection - The section describing the layout in v, or NULL to use the default global section
 . A - The matrix
 . point - The sieve point in the DM
 . values - The array of values
 
   PetscFunctionBegin;
   PetscValidHeaderSpecific(dm, DM_CLASSID, 1);
+  if (!section) {ierr = DMGetDefaultSection(dm, &section);CHKERRQ(ierr);}
   PetscValidHeaderSpecific(section, PETSC_SECTION_CLASSID, 2);
+  if (!globalSection) {ierr = DMGetDefaultGlobalSection(dm, &globalSection);CHKERRQ(ierr);}
   PetscValidHeaderSpecific(globalSection, PETSC_SECTION_CLASSID, 3);
   PetscValidHeaderSpecific(A, MAT_CLASSID, 4);
   ierr = PetscSectionGetNumFields(section, &numFields);CHKERRQ(ierr);