Actual source code: dmpleximpl.h

petsc-3.6.3 2015-12-03
Report Typos and Errors
  1: #if !defined(_PLEXIMPL_H)
  2: #define _PLEXIMPL_H

  4: #include <petscmat.h>       /*I      "petscmat.h"          I*/
  5: #include <petscdmplex.h> /*I      "petscdmplex.h"    I*/
  6: #include <petscbt.h>
  7: #include <petscsf.h>
  8: #include <petsc/private/dmimpl.h>
  9: #include <petsc/private/isimpl.h>     /* for inline access to atlasOff */
 10: #include <../src/sys/utils/hash.h>

 12: PETSC_EXTERN PetscLogEvent DMPLEX_Interpolate, PETSCPARTITIONER_Partition, DMPLEX_Distribute, DMPLEX_DistributeCones, DMPLEX_DistributeLabels, DMPLEX_DistributeSF, DMPLEX_DistributeOverlap, DMPLEX_DistributeField, DMPLEX_DistributeData, DMPLEX_Migrate, DMPLEX_Stratify, DMPLEX_Preallocate, DMPLEX_ResidualFEM, DMPLEX_JacobianFEM, DMPLEX_InterpolatorFEM, DMPLEX_InjectorFEM, DMPLEX_IntegralFEM, DMPLEX_CreateGmsh;

 14: PETSC_EXTERN PetscBool      PetscPartitionerRegisterAllCalled;
 15: PETSC_EXTERN PetscErrorCode PetscPartitionerRegisterAll(void);
 16: PETSC_EXTERN PetscErrorCode PetscPartitionerSetTypeFromOptions_Internal(PetscPartitioner);

 18: typedef enum {REFINER_NOOP = 0,
 19:               REFINER_SIMPLEX_1D,
 20:               REFINER_SIMPLEX_2D,
 21:               REFINER_HYBRID_SIMPLEX_2D,
 22:               REFINER_HEX_2D,
 23:               REFINER_HYBRID_HEX_2D,
 24:               REFINER_SIMPLEX_3D,
 25:               REFINER_HYBRID_SIMPLEX_3D,
 26:               REFINER_HEX_3D,
 27:               REFINER_HYBRID_HEX_3D} CellRefiner;

 29: typedef struct _PetscPartitionerOps *PetscPartitionerOps;
 30: struct _PetscPartitionerOps {
 31:   PetscErrorCode (*setfromoptions)(PetscPartitioner);
 32:   PetscErrorCode (*setup)(PetscPartitioner);
 33:   PetscErrorCode (*view)(PetscPartitioner,PetscViewer);
 34:   PetscErrorCode (*destroy)(PetscPartitioner);
 35:   PetscErrorCode (*partition)(PetscPartitioner, DM, PetscInt, PetscInt, PetscInt[], PetscInt[], PetscSection, IS *);
 36: };

 38: struct _p_PetscPartitioner {
 39:   PETSCHEADER(struct _PetscPartitionerOps);
 40:   void           *data;             /* Implementation object */
 41:   PetscInt        height;           /* Height of points to partition into non-overlapping subsets */
 42: };

 44: typedef struct {
 45:   PetscInt dummy;
 46: } PetscPartitioner_Chaco;

 48: typedef struct {
 49:   PetscInt dummy;
 50: } PetscPartitioner_ParMetis;

 52: typedef struct {
 53:   PetscSection section;   /* Sizes for each partition */
 54:   IS           partition; /* Points in each partition */
 55: } PetscPartitioner_Shell;

 57: typedef struct {
 58:   PetscInt dummy;
 59: } PetscPartitioner_Simple;

 61: /* This is an integer map, in addition it is also a container class
 62:    Design points:
 63:      - Low storage is the most important design point
 64:      - We want flexible insertion and deletion
 65:      - We can live with O(log) query, but we need O(1) iteration over strata
 66: */
 67: struct _n_DMLabel {
 68:   PetscInt    refct;
 69:   char       *name;           /* Label name */
 70:   PetscInt    numStrata;      /* Number of integer values */
 71:   PetscInt   *stratumValues;  /* Value of each stratum */
 72:   /* Basic sorted array storage */
 73:   PetscBool  *arrayValid;     /* The array storage is valid (no additions need to be merged in) */
 74:   PetscInt   *stratumSizes;   /* Size of each stratum */
 75:   PetscInt  **points;         /* Points for each stratum, always sorted */
 76:   /* Hashtable for fast insertion */
 77:   PetscHashI *ht;             /* Hash table for fast insertion */
 78:   /* Index for fast search */
 79:   PetscInt    pStart, pEnd;   /* Bounds for index lookup */
 80:   PetscBT     bt;             /* A bit-wise index */
 81: };

 83: struct _n_PlexLabel {
 84:   DMLabel              label;
 85:   PetscBool            output;
 86:   struct _n_PlexLabel *next;
 87: };
 88: typedef struct _n_PlexLabel *PlexLabel;

 90: /* Utility struct to store the contents of a Gmsh file in memory */
 91: typedef struct {
 92:   PetscInt dim;      /* Entity dimension */
 93:   PetscInt id;       /* Element number */
 94:   PetscInt numNodes; /* Size of node array */
 95:   int nodes[8];      /* Node array */
 96:   PetscInt numTags;  /* Size of tag array */
 97:   int tags[4];       /* Tag array */
 98: } GmshElement;


101: /* Utility struct to store the contents of a Fluent file in memory */
102: typedef struct {
103:   int   index;    /* Type of section */
104:   int   zoneID;
105:   int   first;
106:   int   last;
107:   int   type;
108:   int   nd;       /* Either ND or element-type */
109:   void *data;
110: } FluentSection;

112: struct _n_Boundary {
113:   const char *name;
114:   const char *labelname;
115:   DMLabel     label;
116:   PetscBool   essential;
117:   PetscInt    field;
118:   PetscInt    numcomps;
119:   PetscInt   *comps;
120:   void      (*func)();
121:   PetscInt    numids;
122:   PetscInt   *ids;
123:   void       *ctx;
124:   DMBoundary  next;
125: };

127: typedef struct {
128:   PetscInt             refct;

130:   /* Sieve */
131:   PetscSection         coneSection;       /* Layout of cones (inedges for DAG) */
132:   PetscInt             maxConeSize;       /* Cached for fast lookup */
133:   PetscInt            *cones;             /* Cone for each point */
134:   PetscInt            *coneOrientations;  /* Orientation of each cone point, means cone traveral should start on point 'o', and if negative start on -(o+1) and go in reverse */
135:   PetscSection         supportSection;    /* Layout of cones (inedges for DAG) */
136:   PetscInt             maxSupportSize;    /* Cached for fast lookup */
137:   PetscInt            *supports;          /* Cone for each point */
138:   PetscBool            refinementUniform; /* Flag for uniform cell refinement */
139:   PetscReal            refinementLimit;   /* Maximum volume for refined cell */
140:   PetscInt             hybridPointMax[8]; /* Allow segregation of some points, each dimension has a divider (used in VTK output and refinement) */

142:   PetscInt            *facesTmp;          /* Work space for faces operation */

144:   /* Hierarchy */
145:   DM                   coarseMesh;        /* This mesh was obtained from coarse mesh using DMRefineHierarchy() */

147:   /* Generation */
148:   char                *tetgenOpts;
149:   char                *triangleOpts;
150:   PetscPartitioner     partitioner;

152:   /* Submesh */
153:   DMLabel              subpointMap;       /* Label each original mesh point in the submesh with its depth, subpoint are the implicit numbering */

155:   /* Labels and numbering */
156:   PlexLabel            labels;            /* Linked list of labels */
157:   DMLabel              depthLabel;
158:   IS                   globalVertexNumbers;
159:   IS                   globalCellNumbers;

161:   /* Constraints */
162:   PetscSection         anchorSection;      /* maps constrained points to anchor points */
163:   IS                   anchorIS;           /* anchors indexed by the above section */
164:   PetscErrorCode     (*createanchors)(DM); /* automatically compute anchors (probably from tree constraints) */
165:   PetscErrorCode     (*computeanchormatrix)(DM,PetscSection,PetscSection,Mat);

167:   /* Tree: automatically construct constraints for hierarchically non-conforming meshes */
168:   PetscSection         parentSection;     /* dof == 1 if point has parent */
169:   PetscInt            *parents;           /* point to parent */
170:   PetscInt            *childIDs;          /* point to child ID */
171:   PetscSection         childSection;      /* inverse of parent section */
172:   PetscInt            *children;          /* point to children */
173:   DM                   referenceTree;     /* reference tree to which child ID's refer */
174:   PetscErrorCode      (*getchildsymmetry)(DM,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*);

176:   /* Adjacency */
177:   PetscBool            useCone;           /* Use cone() first when defining adjacency */
178:   PetscBool            useClosure;        /* Use the transitive closure when defining adjacency */
179:   PetscBool            useAnchors;        /* Replace constrained points with their anchors in adjacency lists */

181:   /* Projection */
182:   PetscInt             maxProjectionHeight; /* maximum height of cells used in DMPlexProject functions */

184:   /* Output */
185:   PetscInt             vtkCellHeight;            /* The height of cells for output, default is 0 */
186:   PetscReal            scale[NUM_PETSC_UNITS];   /* The scale for each SI unit */

188:   /* Problem definition */
189:   DMBoundary           boundary;          /* List of boundary conditions */

191:   /* Geometry */
192:   PetscReal            minradius;         /* Minimum distance from cell centroid to face */


195:   /* Debugging */
196:   PetscBool            printSetValues;
197:   PetscInt             printFEM;
198:   PetscReal            printTol;
199: } DM_Plex;

201: PETSC_EXTERN PetscErrorCode DMPlexVTKWriteAll_VTU(DM,PetscViewer);
202: PETSC_EXTERN PetscErrorCode DMPlexVTKGetCellType(DM,PetscInt,PetscInt,PetscInt*);
203: PETSC_EXTERN PetscErrorCode VecView_Plex_Local(Vec,PetscViewer);
204: PETSC_EXTERN PetscErrorCode VecView_Plex(Vec,PetscViewer);
205: PETSC_EXTERN PetscErrorCode VecLoad_Plex_Local(Vec,PetscViewer);
206: PETSC_EXTERN PetscErrorCode VecLoad_Plex(Vec,PetscViewer);
207: PETSC_EXTERN PetscErrorCode DMPlexGetFieldType_Internal(DM, PetscSection, PetscInt, PetscInt *, PetscInt *, PetscViewerVTKFieldType *);
208: #if defined(PETSC_HAVE_HDF5)
209: PETSC_EXTERN PetscErrorCode VecView_Plex_Local_HDF5(Vec, PetscViewer);
210: PETSC_EXTERN PetscErrorCode VecView_Plex_HDF5(Vec, PetscViewer);
211: PETSC_EXTERN PetscErrorCode VecLoad_Plex_HDF5(Vec, PetscViewer);
212: PETSC_EXTERN PetscErrorCode DMPlexView_HDF5(DM, PetscViewer);
213: PETSC_EXTERN PetscErrorCode DMPlexLoad_HDF5(DM, PetscViewer);
214: #endif

216: PETSC_EXTERN PetscErrorCode DMPlexGetAdjacency_Internal(DM,PetscInt,PetscBool,PetscBool,PetscBool,PetscInt*,PetscInt*[]);
217: PETSC_EXTERN PetscErrorCode DMPlexGetFaces_Internal(DM,PetscInt,PetscInt,PetscInt*,PetscInt*,const PetscInt*[]);
218: PETSC_EXTERN PetscErrorCode DMPlexGetRawFaces_Internal(DM,PetscInt,PetscInt,const PetscInt[], PetscInt*,PetscInt*,const PetscInt*[]);
219: PETSC_EXTERN PetscErrorCode DMPlexRestoreFaces_Internal(DM,PetscInt,PetscInt,PetscInt*,PetscInt*,const PetscInt*[]);
220: PETSC_EXTERN PetscErrorCode DMPlexRefineUniform_Internal(DM,CellRefiner,DM*);
221: PETSC_EXTERN PetscErrorCode DMPlexGetCellRefiner_Internal(DM,CellRefiner*);
222: PETSC_EXTERN PetscErrorCode CellRefinerGetAffineTransforms_Internal(CellRefiner, PetscInt *, PetscReal *[], PetscReal *[], PetscReal *[]);
223: PETSC_EXTERN PetscErrorCode CellRefinerRestoreAffineTransforms_Internal(CellRefiner, PetscInt *, PetscReal *[], PetscReal *[], PetscReal *[]);
224: PETSC_EXTERN PetscErrorCode CellRefinerInCellTest_Internal(CellRefiner, const PetscReal[], PetscBool *);
225: PETSC_EXTERN PetscErrorCode DMPlexCreateGmsh_ReadElement(PetscViewer, PetscInt, PetscBool, PetscBool, GmshElement **);
226: PETSC_EXTERN PetscErrorCode DMPlexInvertCell_Internal(PetscInt, PetscInt, PetscInt[]);
227: PETSC_EXTERN PetscErrorCode DMPlexLocalizeCoordinate_Internal(DM, PetscInt, const PetscScalar[], const PetscScalar[], PetscScalar[]);
228: PETSC_EXTERN PetscErrorCode DMPlexLocalizeCoordinateReal_Internal(DM, PetscInt, const PetscReal[], const PetscReal[], PetscReal[]);
229: PETSC_EXTERN PetscErrorCode DMPlexLocalizeAddCoordinate_Internal(DM, PetscInt, const PetscScalar[], const PetscScalar[], PetscScalar[]);
230: PETSC_EXTERN PetscErrorCode DMPlexVecSetFieldClosure_Internal(DM, PetscSection, Vec, PetscBool[], PetscInt, const PetscScalar[], InsertMode);
231: PETSC_EXTERN PetscErrorCode DMPlexProjectConstraints_Internal(DM, Vec, Vec);

235: /* invert dihedral symmetry: return a^-1,
236:  * using the representation described in
237:  * DMPlexGetConeOrientation() */
238: PETSC_STATIC_INLINE PetscInt DihedralInvert(PetscInt N, PetscInt a)
239: {
240:   return (a <= 0) ? a : (N - a);
241: }

245: /* invert dihedral symmetry: return b * a,
246:  * using the representation described in
247:  * DMPlexGetConeOrientation() */
248: PETSC_STATIC_INLINE PetscInt DihedralCompose(PetscInt N, PetscInt a, PetscInt b)
249: {
250:   if (!N) return 0;
251:   return  (a >= 0) ?
252:          ((b >= 0) ? ((a + b) % N) : -(((a - b - 1) % N) + 1)) :
253:          ((b >= 0) ? -(((N - b - a - 1) % N) + 1) : ((N + b - a) % N));
254: }

258: /* swap dihedral symmetries: return b * a^-1,
259:  * using the representation described in
260:  * DMPlexGetConeOrientation() */
261: PETSC_STATIC_INLINE PetscInt DihedralSwap(PetscInt N, PetscInt a, PetscInt b)
262: {
263:   return DihedralCompose(N,DihedralInvert(N,a),b);
264: }

266: PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Internal(DM, PetscInt, PetscInt, PetscReal, Vec, Vec, Vec, void *);
267: PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Internal(DM, PetscInt, PetscInt, PetscReal, PetscReal, Vec, Vec, Mat, Mat,void *);

271: PETSC_STATIC_INLINE void DMPlex_Invert2D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
272: {
273:   const PetscReal invDet = 1.0/detJ;

275:   invJ[0] =  invDet*J[3];
276:   invJ[1] = -invDet*J[1];
277:   invJ[2] = -invDet*J[2];
278:   invJ[3] =  invDet*J[0];
279:   PetscLogFlops(5.0);
280: }

284: PETSC_STATIC_INLINE void DMPlex_Invert3D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
285: {
286:   const PetscReal invDet = 1.0/detJ;

288:   invJ[0*3+0] = invDet*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]);
289:   invJ[0*3+1] = invDet*(J[0*3+2]*J[2*3+1] - J[0*3+1]*J[2*3+2]);
290:   invJ[0*3+2] = invDet*(J[0*3+1]*J[1*3+2] - J[0*3+2]*J[1*3+1]);
291:   invJ[1*3+0] = invDet*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]);
292:   invJ[1*3+1] = invDet*(J[0*3+0]*J[2*3+2] - J[0*3+2]*J[2*3+0]);
293:   invJ[1*3+2] = invDet*(J[0*3+2]*J[1*3+0] - J[0*3+0]*J[1*3+2]);
294:   invJ[2*3+0] = invDet*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]);
295:   invJ[2*3+1] = invDet*(J[0*3+1]*J[2*3+0] - J[0*3+0]*J[2*3+1]);
296:   invJ[2*3+2] = invDet*(J[0*3+0]*J[1*3+1] - J[0*3+1]*J[1*3+0]);
297:   PetscLogFlops(37.0);
298: }

302: PETSC_STATIC_INLINE void DMPlex_Det2D_Internal(PetscReal *detJ, PetscReal J[])
303: {
304:   *detJ = J[0]*J[3] - J[1]*J[2];
305:   PetscLogFlops(3.0);
306: }

310: PETSC_STATIC_INLINE void DMPlex_Det3D_Internal(PetscReal *detJ, PetscReal J[])
311: {
312:   *detJ = (J[0*3+0]*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]) +
313:            J[0*3+1]*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]) +
314:            J[0*3+2]*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]));
315:   PetscLogFlops(12.0);
316: }

320: PETSC_STATIC_INLINE void DMPlex_WaxpyD_Internal(PetscInt dim, PetscReal a, const PetscReal *x, const PetscReal *y, PetscReal *w) {PetscInt d; for (d = 0; d < dim; ++d) w[d] = a*x[d] + y[d];}

324: PETSC_STATIC_INLINE PetscReal DMPlex_DotD_Internal(PetscInt dim, const PetscScalar *x, const PetscReal *y) {PetscReal sum = 0.0; PetscInt d; for (d = 0; d < dim; ++d) sum += PetscRealPart(x[d])*y[d]; return sum;}

328: PETSC_STATIC_INLINE PetscReal DMPlex_DotRealD_Internal(PetscInt dim, const PetscReal *x, const PetscReal *y) {PetscReal sum = 0.0; PetscInt d; for (d = 0; d < dim; ++d) sum += x[d]*y[d]; return sum;}

332: PETSC_STATIC_INLINE PetscReal DMPlex_NormD_Internal(PetscInt dim, const PetscReal *x) {PetscReal sum = 0.0; PetscInt d; for (d = 0; d < dim; ++d) sum += x[d]*x[d]; return PetscSqrtReal(sum);}


337: PETSC_STATIC_INLINE PetscErrorCode DMPlexGetLocalOffset_Private(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
338: {
340: #if defined(PETSC_USE_DEBUG)
341:   {
342:     PetscInt       dof;
344:     if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
345:     PetscSectionGetOffset(dm->defaultSection, point, start);
346:     PetscSectionGetDof(dm->defaultSection, point, &dof);
347:     *end = *start + dof;
348:   }
349: #else
350:   {
351:     const PetscSection s = dm->defaultSection;
352:     *start = s->atlasOff[point - s->pStart];
353:     *end   = *start + s->atlasDof[point - s->pStart];
354:   }
355: #endif
356:   return(0);
357: }

361: PETSC_STATIC_INLINE PetscErrorCode DMPlexGetLocalFieldOffset_Private(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
362: {
364: #if defined(PETSC_USE_DEBUG)
365:   {
366:     PetscInt       dof;
368:     if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
369:     PetscSectionGetFieldOffset(dm->defaultSection, point, field, start);
370:     PetscSectionGetFieldDof(dm->defaultSection, point, field, &dof);
371:     *end = *start + dof;
372:   }
373: #else
374:   {
375:     const PetscSection s = dm->defaultSection->field[field];
376:     *start = s->atlasOff[point - s->pStart];
377:     *end   = *start + s->atlasDof[point - s->pStart];
378:   }
379: #endif
380:   return(0);
381: }

385: PETSC_STATIC_INLINE PetscErrorCode DMPlexGetGlobalOffset_Private(DM dm, PetscInt point, PetscInt *start, PetscInt *end)
386: {
388: #if defined(PETSC_USE_DEBUG)
389:   {
391:     PetscInt       dof,cdof;
392:     if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
393:     if (!dm->defaultGlobalSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default global Section. It will be crated automatically by DMGetDefaultGlobalSection()");
394:     PetscSectionGetOffset(dm->defaultGlobalSection, point, start);
395:     PetscSectionGetDof(dm->defaultGlobalSection, point, &dof);
396:     PetscSectionGetConstraintDof(dm->defaultGlobalSection, point, &cdof);
397:     *end = *start + dof-cdof;
398:   }
399: #else
400:   {
401:     const PetscSection s    = dm->defaultGlobalSection;
402:     const PetscInt     dof  = s->atlasDof[point - s->pStart];
403:     const PetscInt     cdof = s->bc ? s->bc->atlasDof[point - s->bc->pStart] : 0;
404:     *start = s->atlasOff[point - s->pStart];
405:     *end   = *start + dof-cdof;
406:   }
407: #endif
408:   return(0);
409: }

413: PETSC_STATIC_INLINE PetscErrorCode DMPlexGetGlobalFieldOffset_Private(DM dm, PetscInt point, PetscInt field, PetscInt *start, PetscInt *end)
414: {
416: #if defined(PETSC_USE_DEBUG)
417:   {
418:     PetscInt       loff, lfoff, fdof, fcdof, ffcdof, f;

421:     if (!dm->defaultSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default Section, see DMSetDefaultSection()");
422:     if (!dm->defaultGlobalSection) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "DM must have a default global Section. It will be crated automatically by DMGetDefaultGlobalSection()");
423:     PetscSectionGetOffset(dm->defaultGlobalSection, point, start);
424:     PetscSectionGetOffset(dm->defaultSection, point, &loff);
425:     PetscSectionGetFieldOffset(dm->defaultSection, point, field, &lfoff);
426:     PetscSectionGetFieldDof(dm->defaultSection, point, field, &fdof);
427:     PetscSectionGetFieldConstraintDof(dm->defaultSection, point, field, &fcdof);
428:     *start = *start < 0 ? *start - (lfoff-loff) : *start + lfoff-loff;
429:     for (f = 0; f < field; ++f) {
430:       PetscSectionGetFieldConstraintDof(dm->defaultSection, point, f, &ffcdof);
431:       *start = *start < 0 ? *start + ffcdof : *start - ffcdof;
432:     }
433:     *end   = *start < 0 ? *start - (fdof-fcdof) : *start + fdof-fcdof;
434:   }
435: #else
436:   {
437:     const PetscSection s     = dm->defaultSection;
438:     const PetscSection fs    = dm->defaultSection->field[field];
439:     const PetscSection gs    = dm->defaultGlobalSection;
440:     const PetscInt     loff  = s->atlasOff[point - s->pStart];
441:     const PetscInt     goff  = gs->atlasOff[point - s->pStart];
442:     const PetscInt     lfoff = fs->atlasOff[point - s->pStart];
443:     const PetscInt     fdof  = fs->atlasDof[point - s->pStart];
444:     const PetscInt     fcdof = fs->bc ? fs->bc->atlasDof[point - fs->bc->pStart] : 0;
445:     PetscInt           ffcdof = 0, f;

447:     for (f = 0; f < field; ++f) {
448:       const PetscSection ffs = dm->defaultSection->field[f];
449:       ffcdof += ffs->bc ? ffs->bc->atlasDof[point - ffs->bc->pStart] : 0;
450:     }
451:     *start = goff + (goff < 0 ? loff-lfoff + ffcdof : lfoff-loff - ffcdof);
452:     *end   = *start < 0 ? *start - (fdof-fcdof) : *start + fdof-fcdof;
453:   }
454: #endif
455:   return(0);
456: }

458: #endif /* _PLEXIMPL_H */