1 #if !defined(_PLEXIMPL_H)
2 #define _PLEXIMPL_H
3 
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 
10 PETSC_EXTERN PetscLogEvent DMPLEX_Interpolate;
11 PETSC_EXTERN PetscLogEvent DMPLEX_Partition;
12 PETSC_EXTERN PetscLogEvent DMPLEX_PartSelf;
13 PETSC_EXTERN PetscLogEvent DMPLEX_PartLabelInvert;
14 PETSC_EXTERN PetscLogEvent DMPLEX_PartLabelCreateSF;
15 PETSC_EXTERN PetscLogEvent DMPLEX_PartStratSF;
16 PETSC_EXTERN PetscLogEvent DMPLEX_CreatePointSF;
17 PETSC_EXTERN PetscLogEvent DMPLEX_Distribute;
18 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeCones;
19 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeLabels;
20 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeSF;
21 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeOverlap;
22 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeField;
23 PETSC_EXTERN PetscLogEvent DMPLEX_DistributeData;
24 PETSC_EXTERN PetscLogEvent DMPLEX_Migrate;
25 PETSC_EXTERN PetscLogEvent DMPLEX_InterpolateSF;
26 PETSC_EXTERN PetscLogEvent DMPLEX_GlobalToNaturalBegin;
27 PETSC_EXTERN PetscLogEvent DMPLEX_GlobalToNaturalEnd;
28 PETSC_EXTERN PetscLogEvent DMPLEX_NaturalToGlobalBegin;
29 PETSC_EXTERN PetscLogEvent DMPLEX_NaturalToGlobalEnd;
30 PETSC_EXTERN PetscLogEvent DMPLEX_Stratify;
31 PETSC_EXTERN PetscLogEvent DMPLEX_Symmetrize;
32 PETSC_EXTERN PetscLogEvent DMPLEX_Preallocate;
33 PETSC_EXTERN PetscLogEvent DMPLEX_ResidualFEM;
34 PETSC_EXTERN PetscLogEvent DMPLEX_JacobianFEM;
35 PETSC_EXTERN PetscLogEvent DMPLEX_InterpolatorFEM;
36 PETSC_EXTERN PetscLogEvent DMPLEX_InjectorFEM;
37 PETSC_EXTERN PetscLogEvent DMPLEX_IntegralFEM;
38 PETSC_EXTERN PetscLogEvent DMPLEX_CreateGmsh;
39 PETSC_EXTERN PetscLogEvent DMPLEX_RebalanceSharedPoints;
40 PETSC_EXTERN PetscLogEvent DMPLEX_CreateFromFile;
41 PETSC_EXTERN PetscLogEvent DMPLEX_BuildFromCellList;
42 PETSC_EXTERN PetscLogEvent DMPLEX_BuildCoordinatesFromCellList;
43 PETSC_EXTERN PetscLogEvent DMPLEX_LocatePoints;
44 
45 typedef struct _DMPlexCellRefinerOps *DMPlexCellRefinerOps;
46 struct _DMPlexCellRefinerOps {
47   PetscErrorCode (*refine)(DMPlexCellRefiner, DMPolytopeType, PetscInt *, DMPolytopeType *[], PetscInt *[], PetscInt *[], PetscInt *[]);
48   PetscErrorCode (*mapsubcells)(DMPlexCellRefiner, DMPolytopeType, PetscInt, DMPolytopeType, PetscInt, PetscInt, PetscInt *, PetscInt *);
49   PetscErrorCode (*getaffinetransforms)(DMPlexCellRefiner, DMPolytopeType, PetscInt *, PetscReal *[], PetscReal *[], PetscReal *[]);
50   PetscErrorCode (*getaffinefacetransforms)(DMPlexCellRefiner, DMPolytopeType, PetscInt *, PetscReal *[], PetscReal *[], PetscReal *[], PetscReal *[]);
51   PetscErrorCode (*getcellvertices)(DMPlexCellRefiner, DMPolytopeType, PetscInt *, PetscReal *[]);
52   PetscErrorCode (*getsubcellvertices)(DMPlexCellRefiner, DMPolytopeType, DMPolytopeType, PetscInt, PetscInt *, PetscInt *[]);
53   PetscErrorCode (*mapcoords)(DMPlexCellRefiner, DMPolytopeType, DMPolytopeType, PetscInt, PetscInt, PetscInt, const PetscScalar[], PetscScalar[]);
54   PetscErrorCode (*setup)(DMPlexCellRefiner);
55   PetscErrorCode (*destroy)(DMPlexCellRefiner);
56 };
57 
58 struct _p_DMPlexCellRefiner {
59   PETSCHEADER(struct _DMPlexCellRefinerOps);
60   DM                    dm;          /* The original DM */
61   PetscBool             setupcalled;
62   DMPlexCellRefinerType type;
63   PetscInt              *ctOrder;    /* [i] = ct: An array with cell types in depth order */
64   PetscInt              *ctOrderInv; /* [ct] = i: An array with the ordinal numbers for each cell type */
65   PetscInt              *ctStart;    /* The number for the first cell of each polytope type in the original mesh, indexed by cell type */
66   PetscInt              *ctStartNew; /* The number for the first cell of each polytope type in the new mesh, indexed by cell type */
67   PetscInt              *offset;     /* [ct][ctNew]: The offset in the new point numbering of a point of type ctNew produced from an old point of type ct */
68   PetscFE               *coordFE;    /* Finite element for each cell type, used for localized coordinate interpolation */
69   PetscFEGeom           **refGeom;   /* Geometry of the reference cell for each cell type */
70   void                  *data;       /* refiner private data */
71 };
72 
73 /* Utility struct to store the contents of a Fluent file in memory */
74 typedef struct {
75   int          index;    /* Type of section */
76   unsigned int zoneID;
77   unsigned int first;
78   unsigned int last;
79   int          type;
80   int          nd;       /* Either ND or element-type */
81   void        *data;
82 } FluentSection;
83 
84 struct _PetscGridHash {
85   PetscInt     dim;
86   PetscReal    lower[3];    /* The lower-left corner */
87   PetscReal    upper[3];    /* The upper-right corner */
88   PetscReal    extent[3];   /* The box size */
89   PetscReal    h[3];        /* The subbox size */
90   PetscInt     n[3];        /* The number of subboxes */
91   PetscSection cellSection; /* Offsets for cells in each subbox*/
92   IS           cells;       /* List of cells in each subbox */
93   DMLabel      cellsSparse; /* Sparse storage for cell map */
94 };
95 
96 /* Point Numbering in Plex:
97 
98    Points are numbered contiguously by stratum. Strate are organized as follows:
99 
100    First Stratum:  Cells [height 0]
101    Second Stratum: Vertices [depth 0]
102    Third Stratum:  Faces [height 1]
103    Fourth Stratum: Edges [depth 1]
104 
105    We do this so that the numbering of a cell-vertex mesh does not change after interpolation. Within a given stratum,
106    we allow additional segregation of by cell type.
107 */
108 typedef struct {
109   PetscInt             refct;
110 
111   PetscSection         coneSection;       /* Layout of cones (inedges for DAG) */
112   PetscInt             maxConeSize;       /* Cached for fast lookup */
113   PetscInt            *cones;             /* Cone for each point */
114   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 */
115   PetscSection         supportSection;    /* Layout of cones (inedges for DAG) */
116   PetscInt             maxSupportSize;    /* Cached for fast lookup */
117   PetscInt            *supports;          /* Cone for each point */
118   PetscBool            refinementUniform; /* Flag for uniform cell refinement */
119   PetscReal            refinementLimit;   /* Maximum volume for refined cell */
120   PetscErrorCode     (*refinementFunc)(const PetscReal [], PetscReal *); /* Function giving the maximum volume for refined cell */
121   PetscInt             overlap;           /* Overlap of the partitions as passed to DMPlexDistribute() or DMPlexDistributeOverlap() */
122   DMPlexInterpolatedFlag interpolated;
123   DMPlexInterpolatedFlag interpolatedCollective;
124 
125   PetscInt            *facesTmp;          /* Work space for faces operation */
126 
127   /* Hierarchy */
128   DMPlexCellRefinerType cellRefiner;       /* Strategy for refining cells */
129   PetscBool             regularRefinement; /* This flag signals that we are a regular refinement of coarseMesh */
130 
131   /* Generation */
132   char                *tetgenOpts;
133   char                *triangleOpts;
134   PetscPartitioner     partitioner;
135   PetscBool            partitionBalance;  /* Evenly divide partition overlap when distributing */
136   PetscBool            remeshBd;
137 
138   /* Submesh */
139   DMLabel              subpointMap;       /* Label each original mesh point in the submesh with its depth, subpoint are the implicit numbering */
140   IS                   subpointIS;        /* IS holding point number in the enclosing mesh of every point in the submesh chart */
141   PetscObjectState     subpointState;     /* The state of subpointMap when the subpointIS was last created */
142 
143   /* Labels and numbering */
144   PetscObjectState     depthState;        /* State of depth label, so that we can determine if a user changes it */
145   PetscObjectState     celltypeState;     /* State of celltype label, so that we can determine if a user changes it */
146   IS                   globalVertexNumbers;
147   IS                   globalCellNumbers;
148 
149   /* Constraints */
150   PetscSection         anchorSection;      /* maps constrained points to anchor points */
151   IS                   anchorIS;           /* anchors indexed by the above section */
152   PetscErrorCode     (*createanchors)(DM); /* automatically compute anchors (probably from tree constraints) */
153   PetscErrorCode     (*computeanchormatrix)(DM,PetscSection,PetscSection,Mat);
154 
155   /* Tree: automatically construct constraints for hierarchically non-conforming meshes */
156   PetscSection         parentSection;     /* dof == 1 if point has parent */
157   PetscInt            *parents;           /* point to parent */
158   PetscInt            *childIDs;          /* point to child ID */
159   PetscSection         childSection;      /* inverse of parent section */
160   PetscInt            *children;          /* point to children */
161   DM                   referenceTree;     /* reference tree to which child ID's refer */
162   PetscErrorCode      (*getchildsymmetry)(DM,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*);
163 
164   /* MATIS support */
165   PetscSection         subdomainSection;
166 
167   /* Adjacency */
168   PetscBool            useAnchors;        /* Replace constrained points with their anchors in adjacency lists */
169   PetscErrorCode      (*useradjacency)(DM,PetscInt,PetscInt*,PetscInt[],void*); /* User callback for adjacency */
170   void                *useradjacencyctx;  /* User context for callback */
171 
172   /* Projection */
173   PetscInt             maxProjectionHeight; /* maximum height of cells used in DMPlexProject functions */
174   PetscInt             activePoint;         /* current active point in iteration */
175 
176   /* Output */
177   PetscInt             vtkCellHeight;            /* The height of cells for output, default is 0 */
178   PetscReal            scale[NUM_PETSC_UNITS];   /* The scale for each SI unit */
179 
180   /* Geometry */
181   PetscReal            minradius;         /* Minimum distance from cell centroid to face */
182   PetscBool            useHashLocation;   /* Use grid hashing for point location */
183   PetscGridHash        lbox;              /* Local box for searching */
184   void               (*coordFunc)(PetscInt, PetscInt, PetscInt, /* Function used to remap newly introduced vertices */
185                                   const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
186                                   const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
187                                   PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
188 
189   /* Neighbors */
190   PetscMPIInt*         neighbors;
191 
192   /* Debugging */
193   PetscBool            printSetValues;
194   PetscInt             printFEM;
195   PetscInt             printL2;
196   PetscReal            printTol;
197 } DM_Plex;
198 
199 PETSC_EXTERN PetscErrorCode DMPlexVTKWriteAll_VTU(DM,PetscViewer);
200 PETSC_EXTERN PetscErrorCode VecView_Plex_Local(Vec,PetscViewer);
201 PETSC_EXTERN PetscErrorCode VecView_Plex_Native(Vec,PetscViewer);
202 PETSC_EXTERN PetscErrorCode VecView_Plex(Vec,PetscViewer);
203 PETSC_EXTERN PetscErrorCode VecLoad_Plex_Local(Vec,PetscViewer);
204 PETSC_EXTERN PetscErrorCode VecLoad_Plex_Native(Vec,PetscViewer);
205 PETSC_EXTERN PetscErrorCode VecLoad_Plex(Vec,PetscViewer);
206 PETSC_INTERN PetscErrorCode DMPlexGetFieldType_Internal(DM, PetscSection, PetscInt, PetscInt *, PetscInt *, PetscViewerVTKFieldType *);
207 PETSC_INTERN PetscErrorCode DMPlexView_GLVis(DM,PetscViewer);
208 PETSC_INTERN PetscErrorCode DMSetUpGLVisViewer_Plex(PetscObject,PetscViewer);
209 #if defined(PETSC_HAVE_HDF5)
210 PETSC_EXTERN PetscErrorCode VecView_Plex_Local_HDF5(Vec, PetscViewer);
211 PETSC_EXTERN PetscErrorCode VecView_Plex_HDF5(Vec, PetscViewer);
212 PETSC_EXTERN PetscErrorCode VecLoad_Plex_HDF5(Vec, PetscViewer);
213 PETSC_EXTERN PetscErrorCode VecView_Plex_HDF5_Native(Vec, PetscViewer);
214 PETSC_EXTERN PetscErrorCode VecLoad_Plex_HDF5_Native(Vec, PetscViewer);
215 PETSC_EXTERN PetscErrorCode DMPlexView_HDF5(DM, PetscViewer);
216 PETSC_EXTERN PetscErrorCode DMPlexLoad_HDF5(DM, PetscViewer);
217 #endif
218 
219 PETSC_INTERN PetscErrorCode DMPlexVecGetClosureAtDepth_Internal(DM, PetscSection, Vec, PetscInt, PetscInt, PetscInt *, PetscScalar *[]);
220 PETSC_INTERN PetscErrorCode DMPlexClosurePoints_Private(DM,PetscInt,const PetscInt[],IS*);
221 PETSC_INTERN PetscErrorCode DMSetFromOptions_NonRefinement_Plex(PetscOptionItems *, DM);
222 PETSC_INTERN PetscErrorCode DMCoarsen_Plex(DM, MPI_Comm, DM *);
223 PETSC_INTERN PetscErrorCode DMCoarsenHierarchy_Plex(DM, PetscInt, DM []);
224 PETSC_INTERN PetscErrorCode DMRefine_Plex(DM, MPI_Comm, DM *);
225 PETSC_INTERN PetscErrorCode DMRefineHierarchy_Plex(DM, PetscInt, DM []);
226 PETSC_INTERN PetscErrorCode DMAdaptLabel_Plex(DM, DMLabel, DM *);
227 PETSC_INTERN PetscErrorCode DMAdaptMetric_Plex(DM, Vec, DMLabel, DM *);
228 PETSC_INTERN PetscErrorCode DMPlexInsertBoundaryValues_Plex(DM, PetscBool, Vec, PetscReal, Vec, Vec, Vec);
229 PETSC_INTERN PetscErrorCode DMPlexInsertTimeDerivativeBoundaryValues_Plex(DM, PetscBool, Vec, PetscReal, Vec, Vec, Vec);
230 PETSC_INTERN PetscErrorCode DMProjectFunctionLocal_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
231 PETSC_INTERN PetscErrorCode DMProjectFunctionLabelLocal_Plex(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
232 PETSC_INTERN PetscErrorCode DMProjectFieldLocal_Plex(DM,PetscReal,Vec,void (**)(PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],PetscReal,const PetscReal[],PetscInt,const PetscScalar[],PetscScalar[]),InsertMode,Vec);
233 PETSC_INTERN PetscErrorCode DMProjectFieldLabelLocal_Plex(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],Vec,void (**)(PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],PetscReal,const PetscReal[],PetscInt,const PetscScalar[],PetscScalar[]),InsertMode,Vec);
234 PETSC_INTERN PetscErrorCode DMProjectBdFieldLabelLocal_Plex(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],Vec,void (**)(PetscInt,PetscInt,PetscInt,const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],const PetscInt[],const PetscInt[],const PetscScalar[],const PetscScalar[],const PetscScalar[],PetscReal,const PetscReal[],const PetscReal[],PetscInt,const PetscScalar[],PetscScalar[]),InsertMode,Vec);
235 PETSC_INTERN PetscErrorCode DMComputeL2Diff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
236 PETSC_INTERN PetscErrorCode DMComputeL2GradientDiff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[], const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,const PetscReal [],PetscReal *);
237 PETSC_INTERN PetscErrorCode DMComputeL2FieldDiff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
238 PETSC_INTERN PetscErrorCode DMLocatePoints_Plex(DM, Vec, DMPointLocationType, PetscSF);
239 
240 PETSC_INTERN PetscErrorCode DMPlexLoadLabels_HDF5_Internal(DM, PetscViewer);
241 PETSC_INTERN PetscErrorCode DMPlexView_HDF5_Internal(DM, PetscViewer);
242 PETSC_INTERN PetscErrorCode DMPlexLoad_HDF5_Internal(DM, PetscViewer);
243 PETSC_INTERN PetscErrorCode DMPlexLoad_HDF5_Xdmf_Internal(DM, PetscViewer);
244 PETSC_INTERN PetscErrorCode VecView_Plex_HDF5_Internal(Vec, PetscViewer);
245 PETSC_INTERN PetscErrorCode VecView_Plex_HDF5_Native_Internal(Vec, PetscViewer);
246 PETSC_INTERN PetscErrorCode VecView_Plex_Local_HDF5_Internal(Vec, PetscViewer);
247 PETSC_INTERN PetscErrorCode VecLoad_Plex_HDF5_Internal(Vec, PetscViewer);
248 PETSC_INTERN PetscErrorCode VecLoad_Plex_HDF5_Native_Internal(Vec, PetscViewer);
249 /* TODO Make these INTERN */
250 PETSC_EXTERN PetscErrorCode DMPlexView_ExodusII_Internal(DM, int, PetscInt);
251 PETSC_EXTERN PetscErrorCode VecViewPlex_ExodusII_Nodal_Internal(Vec, int, int);
252 PETSC_EXTERN PetscErrorCode VecLoadPlex_ExodusII_Nodal_Internal(Vec, int, int);
253 PETSC_EXTERN PetscErrorCode VecViewPlex_ExodusII_Zonal_Internal(Vec, int, int);
254 PETSC_EXTERN PetscErrorCode VecLoadPlex_ExodusII_Zonal_Internal(Vec, int, int);
255 PETSC_INTERN PetscErrorCode DMPlexVTKGetCellType_Internal(DM,PetscInt,PetscInt,PetscInt*);
256 PETSC_INTERN PetscErrorCode DMPlexGetAdjacency_Internal(DM,PetscInt,PetscBool,PetscBool,PetscBool,PetscInt*,PetscInt*[]);
257 PETSC_INTERN PetscErrorCode DMPlexGetRawFaces_Internal(DM,DMPolytopeType,const PetscInt[],PetscInt*,const DMPolytopeType*[],const PetscInt*[],const PetscInt*[]);
258 PETSC_INTERN PetscErrorCode DMPlexRestoreRawFaces_Internal(DM,DMPolytopeType,const PetscInt[],PetscInt*,const DMPolytopeType*[],const PetscInt*[],const PetscInt*[]);
259 PETSC_INTERN PetscErrorCode CellRefinerInCellTest_Internal(DMPolytopeType, const PetscReal[], PetscBool *);
260 PETSC_INTERN PetscErrorCode DMPlexComputeCellType_Internal(DM, PetscInt, PetscInt, DMPolytopeType *);
261 PETSC_INTERN PetscErrorCode DMPlexCreateCellTypeOrder_Internal(DMPolytopeType, PetscInt *[], PetscInt *[]);
262 PETSC_INTERN PetscErrorCode DMPlexVecSetFieldClosure_Internal(DM, PetscSection, Vec, PetscBool[], PetscInt, PetscInt, const PetscInt[], DMLabel, PetscInt, const PetscScalar[], InsertMode);
263 PETSC_INTERN PetscErrorCode DMPlexProjectConstraints_Internal(DM, Vec, Vec);
264 PETSC_EXTERN PetscErrorCode DMPlexCreateReferenceTree_SetTree(DM, PetscSection, PetscInt[], PetscInt[]);
265 PETSC_EXTERN PetscErrorCode DMPlexCreateReferenceTree_Union(DM,DM,const char *,DM*);
266 PETSC_EXTERN PetscErrorCode DMPlexComputeInterpolatorTree(DM,DM,PetscSF,PetscInt *,Mat);
267 PETSC_EXTERN PetscErrorCode DMPlexComputeInjectorTree(DM,DM,PetscSF,PetscInt *,Mat);
268 PETSC_EXTERN PetscErrorCode DMPlexAnchorsModifyMat(DM,PetscSection,PetscInt,PetscInt,const PetscInt[],const PetscInt ***,const PetscScalar[],PetscInt*,PetscInt*,PetscInt*[],PetscScalar*[],PetscInt[],PetscBool);
269 PETSC_EXTERN PetscErrorCode indicesPoint_private(PetscSection,PetscInt,PetscInt,PetscInt *,PetscBool,PetscInt,PetscInt []);
270 PETSC_EXTERN PetscErrorCode indicesPointFields_private(PetscSection,PetscInt,PetscInt,PetscInt [],PetscBool,PetscInt,PetscInt []);
271 PETSC_INTERN PetscErrorCode DMPlexLocatePoint_Internal(DM,PetscInt,const PetscScalar [],PetscInt,PetscInt *);
272 /* these two are PETSC_EXTERN just because of src/dm/impls/plex/tests/ex18.c */
273 PETSC_EXTERN PetscErrorCode DMPlexOrientCell_Internal(DM,PetscInt,PetscInt,PetscBool);
274 PETSC_EXTERN PetscErrorCode DMPlexOrientInterface_Internal(DM);
275 
276 /* Applications may use this function */
277 PETSC_EXTERN PetscErrorCode DMPlexCreateNumbering_Plex(DM, PetscInt, PetscInt, PetscInt, PetscInt *, PetscSF, IS *);
278 
279 PETSC_INTERN PetscErrorCode DMPlexCreateCellNumbering_Internal(DM, PetscBool, IS *);
280 PETSC_INTERN PetscErrorCode DMPlexCreateVertexNumbering_Internal(DM, PetscBool, IS *);
281 PETSC_INTERN PetscErrorCode DMPlexRefine_Internal(DM, DMLabel, DM *);
282 PETSC_INTERN PetscErrorCode DMPlexCoarsen_Internal(DM, DMLabel, DM *);
283 PETSC_INTERN PetscErrorCode DMCreateMatrix_Plex(DM, Mat*);
284 
285 PETSC_INTERN PetscErrorCode DMPlexGetOverlap_Plex(DM, PetscInt *);
286 
287 /* invert dihedral symmetry: return a^-1,
288  * using the representation described in
289  * DMPlexGetConeOrientation() */
DihedralInvert(PetscInt N,PetscInt a)290 PETSC_STATIC_INLINE PetscInt DihedralInvert(PetscInt N, PetscInt a)
291 {
292   return (a <= 0) ? a : (N - a);
293 }
294 
295 /* invert dihedral symmetry: return b * a,
296  * using the representation described in
297  * DMPlexGetConeOrientation() */
DihedralCompose(PetscInt N,PetscInt a,PetscInt b)298 PETSC_STATIC_INLINE PetscInt DihedralCompose(PetscInt N, PetscInt a, PetscInt b)
299 {
300   if (!N) return 0;
301   return  (a >= 0) ?
302          ((b >= 0) ? ((a + b) % N) : -(((a - b - 1) % N) + 1)) :
303          ((b >= 0) ? -(((N - b - a - 1) % N) + 1) : ((N + b - a) % N));
304 }
305 
306 /* swap dihedral symmetries: return b * a^-1,
307  * using the representation described in
308  * DMPlexGetConeOrientation() */
DihedralSwap(PetscInt N,PetscInt a,PetscInt b)309 PETSC_STATIC_INLINE PetscInt DihedralSwap(PetscInt N, PetscInt a, PetscInt b)
310 {
311   return DihedralCompose(N,DihedralInvert(N,a),b);
312 }
313 
314 PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Internal(DM, IS , PetscReal, Vec, Vec, PetscReal, Vec, void *);
315 PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Hybrid_Internal(DM, IS , PetscReal, Vec, Vec, PetscReal, Vec, void *);
316 PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Internal(DM, IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
317 PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Hybrid_Internal(DM, IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
318 PETSC_EXTERN PetscErrorCode DMPlexReconstructGradients_Internal(DM, PetscFV, PetscInt, PetscInt, Vec, Vec, Vec, Vec);
319 
320 /* Matvec with A in row-major storage, x and y can be aliased */
DMPlex_Mult2D_Internal(const PetscScalar A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])321 PETSC_STATIC_INLINE void DMPlex_Mult2D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
322 {
323   PetscScalar z[2];
324   z[0] = x[0]; z[1] = x[ldx];
325   y[0]   = A[0]*z[0] + A[1]*z[1];
326   y[ldx] = A[2]*z[0] + A[3]*z[1];
327   (void)PetscLogFlops(6.0);
328 }
DMPlex_Mult3D_Internal(const PetscScalar A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])329 PETSC_STATIC_INLINE void DMPlex_Mult3D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
330 {
331   PetscScalar z[3];
332   z[0] = x[0]; z[1] = x[ldx]; z[2] = x[ldx*2];
333   y[0]     = A[0]*z[0] + A[1]*z[1] + A[2]*z[2];
334   y[ldx]   = A[3]*z[0] + A[4]*z[1] + A[5]*z[2];
335   y[ldx*2] = A[6]*z[0] + A[7]*z[1] + A[8]*z[2];
336   (void)PetscLogFlops(15.0);
337 }
DMPlex_MultTranspose2D_Internal(const PetscScalar A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])338 PETSC_STATIC_INLINE void DMPlex_MultTranspose2D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
339 {
340   PetscScalar z[2];
341   z[0] = x[0]; z[1] = x[ldx];
342   y[0]   = A[0]*z[0] + A[2]*z[1];
343   y[ldx] = A[1]*z[0] + A[3]*z[1];
344   (void)PetscLogFlops(6.0);
345 }
DMPlex_MultTranspose3D_Internal(const PetscScalar A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])346 PETSC_STATIC_INLINE void DMPlex_MultTranspose3D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
347 {
348   PetscScalar z[3];
349   z[0] = x[0]; z[1] = x[ldx]; z[2] = x[ldx*2];
350   y[0]     = A[0]*z[0] + A[3]*z[1] + A[6]*z[2];
351   y[ldx]   = A[1]*z[0] + A[4]*z[1] + A[7]*z[2];
352   y[ldx*2] = A[2]*z[0] + A[5]*z[1] + A[8]*z[2];
353   (void)PetscLogFlops(15.0);
354 }
DMPlex_Mult2DReal_Internal(const PetscReal A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])355 PETSC_STATIC_INLINE void DMPlex_Mult2DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
356 {
357   PetscScalar z[2];
358   z[0] = x[0]; z[1] = x[ldx];
359   y[0]   = A[0]*z[0] + A[1]*z[1];
360   y[ldx] = A[2]*z[0] + A[3]*z[1];
361   (void)PetscLogFlops(6.0);
362 }
DMPlex_Mult3DReal_Internal(const PetscReal A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])363 PETSC_STATIC_INLINE void DMPlex_Mult3DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
364 {
365   PetscScalar z[3];
366   z[0] = x[0]; z[1] = x[ldx]; z[2] = x[ldx*2];
367   y[0]     = A[0]*z[0] + A[1]*z[1] + A[2]*z[2];
368   y[ldx]   = A[3]*z[0] + A[4]*z[1] + A[5]*z[2];
369   y[ldx*2] = A[6]*z[0] + A[7]*z[1] + A[8]*z[2];
370   (void)PetscLogFlops(15.0);
371 }
DMPlex_MultTransposeReal_Internal(const PetscReal A[],PetscInt m,PetscInt n,PetscInt ldx,const PetscScalar x[],PetscScalar y[])372 PETSC_STATIC_INLINE void DMPlex_MultTransposeReal_Internal(const PetscReal A[], PetscInt m, PetscInt n, PetscInt ldx, const PetscScalar x[], PetscScalar y[])
373 {
374   PetscScalar z[3];
375   PetscInt    i, j;
376   for (i = 0; i < m; ++i) z[i] = x[i*ldx];
377   for (j = 0; j < n; ++j) {
378     const PetscInt l = j*ldx;
379     y[l] = 0;
380     for (i = 0; i < m; ++i) {
381       y[l] += A[j*n+i]*z[i];
382     }
383   }
384   (void)PetscLogFlops(2*m*n);
385 }
DMPlex_MultTranspose2DReal_Internal(const PetscReal A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])386 PETSC_STATIC_INLINE void DMPlex_MultTranspose2DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
387 {
388   PetscScalar z[2];
389   z[0] = x[0]; z[1] = x[ldx];
390   y[0]   = A[0]*z[0] + A[2]*z[1];
391   y[ldx] = A[1]*z[0] + A[3]*z[1];
392   (void)PetscLogFlops(6.0);
393 }
DMPlex_MultTranspose3DReal_Internal(const PetscReal A[],PetscInt ldx,const PetscScalar x[],PetscScalar y[])394 PETSC_STATIC_INLINE void DMPlex_MultTranspose3DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
395 {
396   PetscScalar z[3];
397   z[0] = x[0]; z[1] = x[ldx]; z[2] = x[ldx*2];
398   y[0]     = A[0]*z[0] + A[3]*z[1] + A[6]*z[2];
399   y[ldx]   = A[1]*z[0] + A[4]*z[1] + A[7]*z[2];
400   y[ldx*2] = A[2]*z[0] + A[5]*z[1] + A[8]*z[2];
401   (void)PetscLogFlops(15.0);
402 }
403 
DMPlex_MatMult2D_Internal(const PetscScalar A[],PetscInt n,PetscInt ldb,const PetscScalar B[],PetscScalar C[])404 PETSC_STATIC_INLINE void DMPlex_MatMult2D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
405 {
406   PetscInt j;
407   for (j = 0; j < n; ++j) {
408     PetscScalar z[2];
409     z[0] = B[0+j]; z[1] = B[1*ldb+j];
410     DMPlex_Mult2D_Internal(A, 1, z, z);
411     C[0+j] = z[0]; C[1*ldb+j] = z[1];
412   }
413   (void)PetscLogFlops(8.0*n);
414 }
DMPlex_MatMult3D_Internal(const PetscScalar A[],PetscInt n,PetscInt ldb,const PetscScalar B[],PetscScalar C[])415 PETSC_STATIC_INLINE void DMPlex_MatMult3D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
416 {
417   PetscInt j;
418   for (j = 0; j < n; ++j) {
419     PetscScalar z[3];
420     z[0] = B[0+j]; z[1] = B[1*ldb+j]; z[2] = B[2*ldb+j];
421     DMPlex_Mult3D_Internal(A, 1, z, z);
422     C[0+j] = z[0]; C[1*ldb+j] = z[1]; C[2*ldb+j] = z[2];
423   }
424   (void)PetscLogFlops(8.0*n);
425 }
DMPlex_MatMultTranspose2D_Internal(const PetscScalar A[],PetscInt n,PetscInt ldb,const PetscScalar B[],PetscScalar C[])426 PETSC_STATIC_INLINE void DMPlex_MatMultTranspose2D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
427 {
428   PetscInt j;
429   for (j = 0; j < n; ++j) {
430     PetscScalar z[2];
431     z[0] = B[0+j]; z[1] = B[1*ldb+j];
432     DMPlex_MultTranspose2D_Internal(A, 1, z, z);
433     C[0+j] = z[0]; C[1*ldb+j] = z[1];
434   }
435   (void)PetscLogFlops(8.0*n);
436 }
DMPlex_MatMultTranspose3D_Internal(const PetscScalar A[],PetscInt n,PetscInt ldb,const PetscScalar B[],PetscScalar C[])437 PETSC_STATIC_INLINE void DMPlex_MatMultTranspose3D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
438 {
439   PetscInt j;
440   for (j = 0; j < n; ++j) {
441     PetscScalar z[3];
442     z[0] = B[0+j]; z[1] = B[1*ldb+j]; z[2] = B[2*ldb+j];
443     DMPlex_MultTranspose3D_Internal(A, 1, z, z);
444     C[0+j] = z[0]; C[1*ldb+j] = z[1]; C[2*ldb+j] = z[2];
445   }
446   (void)PetscLogFlops(8.0*n);
447 }
448 
DMPlex_MatMultLeft2D_Internal(const PetscScalar A[],PetscInt m,PetscInt ldb,const PetscScalar B[],PetscScalar C[])449 PETSC_STATIC_INLINE void DMPlex_MatMultLeft2D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
450 {
451   PetscInt j;
452   for (j = 0; j < m; ++j) {
453     DMPlex_MultTranspose2D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
454   }
455   (void)PetscLogFlops(8.0*m);
456 }
DMPlex_MatMultLeft3D_Internal(const PetscScalar A[],PetscInt m,PetscInt ldb,const PetscScalar B[],PetscScalar C[])457 PETSC_STATIC_INLINE void DMPlex_MatMultLeft3D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
458 {
459   PetscInt j;
460   for (j = 0; j < m; ++j) {
461     DMPlex_MultTranspose3D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
462   }
463   (void)PetscLogFlops(8.0*m);
464 }
DMPlex_MatMultTransposeLeft2D_Internal(const PetscScalar A[],PetscInt m,PetscInt ldb,const PetscScalar B[],PetscScalar C[])465 PETSC_STATIC_INLINE void DMPlex_MatMultTransposeLeft2D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
466 {
467   PetscInt j;
468   for (j = 0; j < m; ++j) {
469     DMPlex_Mult2D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
470   }
471   (void)PetscLogFlops(8.0*m);
472 }
DMPlex_MatMultTransposeLeft3D_Internal(const PetscScalar A[],PetscInt m,PetscInt ldb,const PetscScalar B[],PetscScalar C[])473 PETSC_STATIC_INLINE void DMPlex_MatMultTransposeLeft3D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
474 {
475   PetscInt j;
476   for (j = 0; j < m; ++j) {
477     DMPlex_Mult3D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
478   }
479   (void)PetscLogFlops(8.0*m);
480 }
481 
DMPlex_Transpose2D_Internal(PetscScalar A[])482 PETSC_STATIC_INLINE void DMPlex_Transpose2D_Internal(PetscScalar A[])
483 {
484   PetscScalar tmp;
485   tmp = A[1]; A[1] = A[2]; A[2] = tmp;
486 }
DMPlex_Transpose3D_Internal(PetscScalar A[])487 PETSC_STATIC_INLINE void DMPlex_Transpose3D_Internal(PetscScalar A[])
488 {
489   PetscScalar tmp;
490   tmp = A[1]; A[1] = A[3]; A[3] = tmp;
491   tmp = A[2]; A[2] = A[6]; A[6] = tmp;
492   tmp = A[5]; A[5] = A[7]; A[7] = tmp;
493 }
494 
DMPlex_Invert2D_Internal(PetscReal invJ[],PetscReal J[],PetscReal detJ)495 PETSC_STATIC_INLINE void DMPlex_Invert2D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
496 {
497   const PetscReal invDet = 1.0/detJ;
498 
499   invJ[0] =  invDet*J[3];
500   invJ[1] = -invDet*J[1];
501   invJ[2] = -invDet*J[2];
502   invJ[3] =  invDet*J[0];
503   (void)PetscLogFlops(5.0);
504 }
505 
DMPlex_Invert3D_Internal(PetscReal invJ[],PetscReal J[],PetscReal detJ)506 PETSC_STATIC_INLINE void DMPlex_Invert3D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
507 {
508   const PetscReal invDet = 1.0/detJ;
509 
510   invJ[0*3+0] = invDet*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]);
511   invJ[0*3+1] = invDet*(J[0*3+2]*J[2*3+1] - J[0*3+1]*J[2*3+2]);
512   invJ[0*3+2] = invDet*(J[0*3+1]*J[1*3+2] - J[0*3+2]*J[1*3+1]);
513   invJ[1*3+0] = invDet*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]);
514   invJ[1*3+1] = invDet*(J[0*3+0]*J[2*3+2] - J[0*3+2]*J[2*3+0]);
515   invJ[1*3+2] = invDet*(J[0*3+2]*J[1*3+0] - J[0*3+0]*J[1*3+2]);
516   invJ[2*3+0] = invDet*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]);
517   invJ[2*3+1] = invDet*(J[0*3+1]*J[2*3+0] - J[0*3+0]*J[2*3+1]);
518   invJ[2*3+2] = invDet*(J[0*3+0]*J[1*3+1] - J[0*3+1]*J[1*3+0]);
519   (void)PetscLogFlops(37.0);
520 }
521 
DMPlex_Det2D_Internal(PetscReal * detJ,const PetscReal J[])522 PETSC_STATIC_INLINE void DMPlex_Det2D_Internal(PetscReal *detJ, const PetscReal J[])
523 {
524   *detJ = J[0]*J[3] - J[1]*J[2];
525   (void)PetscLogFlops(3.0);
526 }
527 
DMPlex_Det3D_Internal(PetscReal * detJ,const PetscReal J[])528 PETSC_STATIC_INLINE void DMPlex_Det3D_Internal(PetscReal *detJ, const PetscReal J[])
529 {
530   *detJ = (J[0*3+0]*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]) +
531            J[0*3+1]*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]) +
532            J[0*3+2]*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]));
533   (void)PetscLogFlops(12.0);
534 }
535 
DMPlex_Det2D_Scalar_Internal(PetscReal * detJ,const PetscScalar J[])536 PETSC_STATIC_INLINE void DMPlex_Det2D_Scalar_Internal(PetscReal *detJ, const PetscScalar J[])
537 {
538   *detJ = PetscRealPart(J[0])*PetscRealPart(J[3]) - PetscRealPart(J[1])*PetscRealPart(J[2]);
539   (void)PetscLogFlops(3.0);
540 }
541 
DMPlex_Det3D_Scalar_Internal(PetscReal * detJ,const PetscScalar J[])542 PETSC_STATIC_INLINE void DMPlex_Det3D_Scalar_Internal(PetscReal *detJ, const PetscScalar J[])
543 {
544   *detJ = (PetscRealPart(J[0*3+0])*(PetscRealPart(J[1*3+1])*PetscRealPart(J[2*3+2]) - PetscRealPart(J[1*3+2])*PetscRealPart(J[2*3+1])) +
545            PetscRealPart(J[0*3+1])*(PetscRealPart(J[1*3+2])*PetscRealPart(J[2*3+0]) - PetscRealPart(J[1*3+0])*PetscRealPart(J[2*3+2])) +
546            PetscRealPart(J[0*3+2])*(PetscRealPart(J[1*3+0])*PetscRealPart(J[2*3+1]) - PetscRealPart(J[1*3+1])*PetscRealPart(J[2*3+0])));
547   (void)PetscLogFlops(12.0);
548 }
549 
DMPlex_WaxpyD_Internal(PetscInt dim,PetscReal a,const PetscReal * x,const PetscReal * y,PetscReal * w)550 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];}
551 
DMPlex_DotD_Internal(PetscInt dim,const PetscScalar * x,const PetscReal * y)552 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;}
553 
DMPlex_DotRealD_Internal(PetscInt dim,const PetscReal * x,const PetscReal * y)554 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;}
555 
DMPlex_NormD_Internal(PetscInt dim,const PetscReal * x)556 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);}
557 
DMPlexFixFaceOrientations_Translate_Private(PetscInt ornt,PetscInt * start,PetscBool * reverse)558 PETSC_STATIC_INLINE PetscErrorCode DMPlexFixFaceOrientations_Translate_Private(PetscInt ornt, PetscInt *start, PetscBool *reverse)
559 {
560   PetscFunctionBegin;
561   *reverse = (ornt < 0) ? PETSC_TRUE : PETSC_FALSE;
562   *start = *reverse ? -(ornt+1) : ornt;
563   PetscFunctionReturn(0);
564 }
565 
DMPlexFixFaceOrientations_Combine_Private(PetscInt coneSize,PetscInt origStart,PetscBool origReverse,PetscInt rotateStart,PetscBool rotateReverse,PetscInt * newStart,PetscBool * newReverse)566 PETSC_STATIC_INLINE PetscErrorCode DMPlexFixFaceOrientations_Combine_Private(PetscInt coneSize, PetscInt origStart, PetscBool origReverse, PetscInt rotateStart, PetscBool rotateReverse, PetscInt *newStart, PetscBool *newReverse)
567 {
568   PetscFunctionBegin;
569   *newReverse = (origReverse == rotateReverse) ? PETSC_FALSE : PETSC_TRUE;
570   *newStart = rotateReverse ? (coneSize + rotateStart - origStart) : (coneSize + origStart - rotateStart);
571   *newStart %= coneSize;
572   PetscFunctionReturn(0);
573 }
574 
DMPlexFixFaceOrientations_TranslateBack_Private(PetscInt coneSize,PetscInt start,PetscBool reverse,PetscInt * ornt)575 PETSC_STATIC_INLINE PetscErrorCode DMPlexFixFaceOrientations_TranslateBack_Private(PetscInt coneSize, PetscInt start, PetscBool reverse, PetscInt *ornt)
576 {
577   PetscFunctionBegin;
578   if (coneSize < 3) {
579     /* edges just get flipped if start == 1 regardless direction */
580     *ornt = start ? -2 : 0;
581   } else {
582     *ornt = reverse ? -(start+1) : start;
583   }
584   PetscFunctionReturn(0);
585 }
586 
DMPlexFixFaceOrientations_Permute_Private(PetscInt n,const PetscInt arr[],PetscInt start,PetscBool reverse,PetscInt newarr[])587 PETSC_STATIC_INLINE PetscErrorCode DMPlexFixFaceOrientations_Permute_Private(PetscInt n, const PetscInt arr[], PetscInt start, PetscBool reverse, PetscInt newarr[])
588 {
589   PetscInt i;
590 
591   PetscFunctionBegin;
592   if (reverse) {for (i=0; i<n; i++) newarr[i] = arr[(n+start-i)%n];}
593   else         {for (i=0; i<n; i++) newarr[i] = arr[(start+i)%n];}
594   PetscFunctionReturn(0);
595 }
596 
597 PETSC_INTERN PetscErrorCode DMPlexGetPointDualSpaceFEM(DM,PetscInt,PetscInt,PetscDualSpace *);
598 PETSC_INTERN PetscErrorCode DMPlexGetIndicesPoint_Internal(PetscSection,PetscBool,PetscInt,PetscInt,PetscInt *,PetscBool,const PetscInt[],const PetscInt[],PetscInt[]);
599 PETSC_INTERN PetscErrorCode DMPlexGetIndicesPointFields_Internal(PetscSection,PetscBool,PetscInt,PetscInt,PetscInt[],PetscBool,const PetscInt***,PetscInt,const PetscInt[],PetscInt[]);
600 PETSC_INTERN PetscErrorCode DMPlexGetCompressedClosure(DM, PetscSection, PetscInt, PetscInt *, PetscInt **, PetscSection *, IS *, const PetscInt **);
601 PETSC_INTERN PetscErrorCode DMPlexRestoreCompressedClosure(DM, PetscSection, PetscInt, PetscInt *, PetscInt **, PetscSection *, IS *, const PetscInt **);
602 
603 PETSC_EXTERN PetscErrorCode DMSNESGetFEGeom(DMField, IS, PetscQuadrature, PetscBool, PetscFEGeom **);
604 PETSC_EXTERN PetscErrorCode DMSNESRestoreFEGeom(DMField, IS, PetscQuadrature, PetscBool, PetscFEGeom **);
605 PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Patch_Internal(DM, PetscSection, IS, PetscReal, Vec, Vec, Vec, void *);
606 PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Patch_Internal(DM, PetscSection, PetscSection, IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
607 PETSC_INTERN PetscErrorCode DMCreateSubDomainDM_Plex(DM,DMLabel,PetscInt,IS*,DM*);
608 PETSC_INTERN PetscErrorCode DMPlexBasisTransformPoint_Internal(DM, DM, Vec, PetscInt, PetscBool[], PetscBool, PetscScalar *);
609 PETSC_EXTERN PetscErrorCode DMPlexBasisTransformPointTensor_Internal(DM, DM, Vec, PetscInt, PetscBool, PetscInt, PetscScalar *);
610 PETSC_INTERN PetscErrorCode DMPlexBasisTransformApplyReal_Internal(DM, const PetscReal[], PetscBool, PetscInt, const PetscReal *, PetscReal *, void *);
611 PETSC_INTERN PetscErrorCode DMPlexBasisTransformApply_Internal(DM, const PetscReal[], PetscBool, PetscInt, const PetscScalar *, PetscScalar *, void *);
612 PETSC_INTERN PetscErrorCode DMCreateNeumannOverlap_Plex(DM, IS*, Mat*, PetscErrorCode (**)(Mat, PetscReal, Vec, Vec, PetscReal, IS, void*), void **);
613 
614 #endif /* _PLEXIMPL_H */
615