Actual source code: dmpleximpl.h
1: #if !defined(_PLEXIMPL_H)
2: #define _PLEXIMPL_H
4: #include <petscmat.h>
5: #include <petscdmplex.h>
6: #include <petscbt.h>
7: #include <petscsf.h>
8: #include <petsc/private/dmimpl.h>
10: #if defined(PETSC_HAVE_EXODUSII)
11: #include <exodusII.h>
12: #endif
14: PETSC_EXTERN PetscLogEvent DMPLEX_Interpolate;
15: PETSC_EXTERN PetscLogEvent DMPLEX_Partition;
16: PETSC_EXTERN PetscLogEvent DMPLEX_PartSelf;
17: PETSC_EXTERN PetscLogEvent DMPLEX_PartLabelInvert;
18: PETSC_EXTERN PetscLogEvent DMPLEX_PartLabelCreateSF;
19: PETSC_EXTERN PetscLogEvent DMPLEX_PartStratSF;
20: PETSC_EXTERN PetscLogEvent DMPLEX_CreatePointSF;
21: PETSC_EXTERN PetscLogEvent DMPLEX_Distribute;
22: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeCones;
23: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeLabels;
24: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeSF;
25: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeOverlap;
26: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeField;
27: PETSC_EXTERN PetscLogEvent DMPLEX_DistributeData;
28: PETSC_EXTERN PetscLogEvent DMPLEX_Migrate;
29: PETSC_EXTERN PetscLogEvent DMPLEX_InterpolateSF;
30: PETSC_EXTERN PetscLogEvent DMPLEX_GlobalToNaturalBegin;
31: PETSC_EXTERN PetscLogEvent DMPLEX_GlobalToNaturalEnd;
32: PETSC_EXTERN PetscLogEvent DMPLEX_NaturalToGlobalBegin;
33: PETSC_EXTERN PetscLogEvent DMPLEX_NaturalToGlobalEnd;
34: PETSC_EXTERN PetscLogEvent DMPLEX_Stratify;
35: PETSC_EXTERN PetscLogEvent DMPLEX_Symmetrize;
36: PETSC_EXTERN PetscLogEvent DMPLEX_Preallocate;
37: PETSC_EXTERN PetscLogEvent DMPLEX_ResidualFEM;
38: PETSC_EXTERN PetscLogEvent DMPLEX_JacobianFEM;
39: PETSC_EXTERN PetscLogEvent DMPLEX_InterpolatorFEM;
40: PETSC_EXTERN PetscLogEvent DMPLEX_InjectorFEM;
41: PETSC_EXTERN PetscLogEvent DMPLEX_IntegralFEM;
42: PETSC_EXTERN PetscLogEvent DMPLEX_CreateGmsh;
43: PETSC_EXTERN PetscLogEvent DMPLEX_RebalanceSharedPoints;
44: PETSC_EXTERN PetscLogEvent DMPLEX_CreateFromFile;
45: PETSC_EXTERN PetscLogEvent DMPLEX_BuildFromCellList;
46: PETSC_EXTERN PetscLogEvent DMPLEX_BuildCoordinatesFromCellList;
47: PETSC_EXTERN PetscLogEvent DMPLEX_LocatePoints;
48: PETSC_EXTERN PetscLogEvent DMPLEX_TopologyView;
49: PETSC_EXTERN PetscLogEvent DMPLEX_LabelsView;
50: PETSC_EXTERN PetscLogEvent DMPLEX_CoordinatesView;
51: PETSC_EXTERN PetscLogEvent DMPLEX_SectionView;
52: PETSC_EXTERN PetscLogEvent DMPLEX_GlobalVectorView;
53: PETSC_EXTERN PetscLogEvent DMPLEX_LocalVectorView;
54: PETSC_EXTERN PetscLogEvent DMPLEX_TopologyLoad;
55: PETSC_EXTERN PetscLogEvent DMPLEX_LabelsLoad;
56: PETSC_EXTERN PetscLogEvent DMPLEX_CoordinatesLoad;
57: PETSC_EXTERN PetscLogEvent DMPLEX_SectionLoad;
58: PETSC_EXTERN PetscLogEvent DMPLEX_GlobalVectorLoad;
59: PETSC_EXTERN PetscLogEvent DMPLEX_LocalVectorLoad;
60: PETSC_EXTERN PetscLogEvent DMPLEX_MetricEnforceSPD;
61: PETSC_EXTERN PetscLogEvent DMPLEX_MetricNormalize;
62: PETSC_EXTERN PetscLogEvent DMPLEX_MetricAverage;
63: PETSC_EXTERN PetscLogEvent DMPLEX_MetricIntersection;
65: /* Utility struct to store the contents of a Fluent file in memory */
66: typedef struct {
67: int index; /* Type of section */
68: unsigned int zoneID;
69: unsigned int first;
70: unsigned int last;
71: int type;
72: int nd; /* Either ND or element-type */
73: void *data;
74: } FluentSection;
76: struct _PetscGridHash {
77: PetscInt dim;
78: PetscReal lower[3]; /* The lower-left corner */
79: PetscReal upper[3]; /* The upper-right corner */
80: PetscReal extent[3]; /* The box size */
81: PetscReal h[3]; /* The subbox size */
82: PetscInt n[3]; /* The number of subboxes */
83: PetscSection cellSection; /* Offsets for cells in each subbox*/
84: IS cells; /* List of cells in each subbox */
85: DMLabel cellsSparse; /* Sparse storage for cell map */
86: };
88: typedef struct {
89: PetscBool isotropic; /* Is the metric isotropic? */
90: PetscBool uniform; /* Is the metric uniform? */
91: PetscBool restrictAnisotropyFirst; /* Should anisotropy or normalization come first? */
92: PetscBool noInsert; /* Should node insertion/deletion be turned off? */
93: PetscBool noSwap; /* Should facet swapping be turned off? */
94: PetscBool noMove; /* Should node movement be turned off? */
95: PetscBool noSurf; /* Should surface modification be turned off? */
96: PetscReal h_min, h_max; /* Minimum/maximum tolerated metric magnitudes */
97: PetscReal a_max; /* Maximum tolerated anisotropy */
98: PetscReal targetComplexity; /* Target metric complexity */
99: PetscReal p; /* Degree for L-p normalization methods */
100: PetscReal gradationFactor; /* Maximum tolerated length ratio for adjacent edges */
101: PetscReal hausdorffNumber; /* Max. distance between piecewise linear representation of boundary and reconstructed ideal boundary */
102: PetscInt numIter; /* Number of ParMmg mesh adaptation iterations */
103: PetscInt verbosity; /* Level of verbosity for remesher (-1 = no output, 10 = maximum) */
104: } DMPlexMetricCtx;
106: /* Point Numbering in Plex:
108: Points are numbered contiguously by stratum. Strate are organized as follows:
110: First Stratum: Cells [height 0]
111: Second Stratum: Vertices [depth 0]
112: Third Stratum: Faces [height 1]
113: Fourth Stratum: Edges [depth 1]
115: We do this so that the numbering of a cell-vertex mesh does not change after interpolation. Within a given stratum,
116: we allow additional segregation of by cell type.
117: */
118: typedef struct {
119: PetscInt refct;
121: PetscSection coneSection; /* Layout of cones (inedges for DAG) */
122: PetscInt *cones; /* Cone for each point */
123: 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 */
124: PetscSection supportSection; /* Layout of cones (inedges for DAG) */
125: PetscInt *supports; /* Cone for each point */
126: PetscBool refinementUniform; /* Flag for uniform cell refinement */
127: char *transformType; /* Type of transform for uniform cell refinement */
128: PetscReal refinementLimit; /* Maximum volume for refined cell */
129: PetscErrorCode (*refinementFunc)(const PetscReal [], PetscReal *); /* Function giving the maximum volume for refined cell */
130: PetscBool distDefault; /* Distribute the DM by default */
131: PetscInt overlap; /* Overlap of the partitions as passed to DMPlexDistribute() or DMPlexDistributeOverlap() */
132: DMPlexInterpolatedFlag interpolated;
133: DMPlexInterpolatedFlag interpolatedCollective;
135: PetscInt *facesTmp; /* Work space for faces operation */
137: /* Hierarchy */
138: PetscBool regularRefinement; /* This flag signals that we are a regular refinement of coarseMesh */
140: /* Generation */
141: char *tetgenOpts;
142: char *triangleOpts;
143: PetscPartitioner partitioner;
144: PetscBool partitionBalance; /* Evenly divide partition overlap when distributing */
145: PetscBool remeshBd;
147: /* Submesh */
148: DMLabel subpointMap; /* Label each original mesh point in the submesh with its depth, subpoint are the implicit numbering */
149: IS subpointIS; /* IS holding point number in the enclosing mesh of every point in the submesh chart */
150: PetscObjectState subpointState; /* The state of subpointMap when the subpointIS was last created */
152: /* Labels and numbering */
153: PetscObjectState depthState; /* State of depth label, so that we can determine if a user changes it */
154: PetscObjectState celltypeState; /* State of celltype label, so that we can determine if a user changes it */
155: IS globalVertexNumbers;
156: IS globalCellNumbers;
158: /* Constraints */
159: PetscSection anchorSection; /* maps constrained points to anchor points */
160: IS anchorIS; /* anchors indexed by the above section */
161: PetscErrorCode (*createanchors)(DM); /* automatically compute anchors (probably from tree constraints) */
162: PetscErrorCode (*computeanchormatrix)(DM,PetscSection,PetscSection,Mat);
164: /* Tree: automatically construct constraints for hierarchically non-conforming meshes */
165: PetscSection parentSection; /* dof == 1 if point has parent */
166: PetscInt *parents; /* point to parent */
167: PetscInt *childIDs; /* point to child ID */
168: PetscSection childSection; /* inverse of parent section */
169: PetscInt *children; /* point to children */
170: DM referenceTree; /* reference tree to which child ID's refer */
171: PetscErrorCode (*getchildsymmetry)(DM,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,PetscInt*);
173: /* MATIS support */
174: PetscSection subdomainSection;
176: /* Adjacency */
177: PetscBool useAnchors; /* Replace constrained points with their anchors in adjacency lists */
178: PetscErrorCode (*useradjacency)(DM,PetscInt,PetscInt*,PetscInt[],void*); /* User callback for adjacency */
179: void *useradjacencyctx; /* User context for callback */
181: /* Projection */
182: PetscInt maxProjectionHeight; /* maximum height of cells used in DMPlexProject functions */
183: PetscInt activePoint; /* current active point in iteration */
185: /* Output */
186: PetscInt vtkCellHeight; /* The height of cells for output, default is 0 */
187: PetscReal scale[NUM_PETSC_UNITS]; /* The scale for each SI unit */
189: /* Geometry */
190: PetscBool ignoreModel; /* Ignore the geometry model during refinement */
191: PetscReal minradius; /* Minimum distance from cell centroid to face */
192: PetscBool useHashLocation; /* Use grid hashing for point location */
193: PetscGridHash lbox; /* Local box for searching */
194: void (*coordFunc)(PetscInt, PetscInt, PetscInt, /* Function used to remap newly introduced vertices */
195: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
196: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
197: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
199: /* Neighbors */
200: PetscMPIInt* neighbors;
202: /* Metric */
203: DMPlexMetricCtx *metricCtx;
205: /* Debugging */
206: PetscBool printSetValues;
207: PetscInt printFEM;
208: PetscInt printL2;
209: PetscReal printTol;
210: } DM_Plex;
212: PETSC_INTERN PetscErrorCode DMPlexCopy_Internal(DM, PetscBool, DM);
214: PETSC_EXTERN PetscErrorCode DMPlexVTKWriteAll_VTU(DM,PetscViewer);
215: PETSC_EXTERN PetscErrorCode VecView_Plex_Local(Vec,PetscViewer);
216: PETSC_EXTERN PetscErrorCode VecView_Plex_Native(Vec,PetscViewer);
217: PETSC_EXTERN PetscErrorCode VecView_Plex(Vec,PetscViewer);
218: PETSC_EXTERN PetscErrorCode VecLoad_Plex_Local(Vec,PetscViewer);
219: PETSC_EXTERN PetscErrorCode VecLoad_Plex_Native(Vec,PetscViewer);
220: PETSC_EXTERN PetscErrorCode VecLoad_Plex(Vec,PetscViewer);
221: PETSC_INTERN PetscErrorCode DMPlexGetFieldType_Internal(DM, PetscSection, PetscInt, PetscInt *, PetscInt *, PetscViewerVTKFieldType *);
222: PETSC_INTERN PetscErrorCode DMPlexView_GLVis(DM,PetscViewer);
223: PETSC_INTERN PetscErrorCode DMSetUpGLVisViewer_Plex(PetscObject,PetscViewer);
224: #if defined(PETSC_HAVE_HDF5)
225: PETSC_EXTERN PetscErrorCode VecView_Plex_Local_HDF5(Vec, PetscViewer);
226: PETSC_EXTERN PetscErrorCode VecView_Plex_HDF5(Vec, PetscViewer);
227: PETSC_EXTERN PetscErrorCode VecLoad_Plex_HDF5(Vec, PetscViewer);
228: PETSC_EXTERN PetscErrorCode VecView_Plex_HDF5_Native(Vec, PetscViewer);
229: PETSC_EXTERN PetscErrorCode VecLoad_Plex_HDF5_Native(Vec, PetscViewer);
230: PETSC_EXTERN PetscErrorCode DMPlexView_HDF5(DM, PetscViewer);
231: PETSC_EXTERN PetscErrorCode DMPlexLoad_HDF5(DM, PetscViewer);
232: PETSC_INTERN PetscErrorCode DMPlexTopologyView_HDF5_Internal(DM, IS, PetscViewer);
233: PETSC_INTERN PetscErrorCode DMPlexCoordinatesView_HDF5_Internal(DM, PetscViewer);
234: PETSC_INTERN PetscErrorCode DMPlexLabelsView_HDF5_Internal(DM, IS, PetscViewer);
235: PETSC_INTERN PetscErrorCode DMPlexSectionView_HDF5_Internal(DM, PetscViewer, DM);
236: PETSC_INTERN PetscErrorCode DMPlexGlobalVectorView_HDF5_Internal(DM, PetscViewer, DM, Vec);
237: PETSC_INTERN PetscErrorCode DMPlexLocalVectorView_HDF5_Internal(DM, PetscViewer, DM, Vec);
238: PETSC_INTERN PetscErrorCode DMPlexTopologyLoad_HDF5_Internal(DM, PetscViewer, PetscSF*);
239: PETSC_INTERN PetscErrorCode DMPlexCoordinatesLoad_HDF5_Internal(DM, PetscViewer, PetscSF);
240: PETSC_INTERN PetscErrorCode DMPlexLabelsLoad_HDF5_Internal(DM, PetscViewer, PetscSF);
241: PETSC_INTERN PetscErrorCode DMPlexSectionLoad_HDF5_Internal(DM, PetscViewer, DM, PetscSF, PetscSF*, PetscSF*);
242: PETSC_INTERN PetscErrorCode DMPlexVecLoad_HDF5_Internal(DM, PetscViewer, DM, PetscSF, Vec);
243: PETSC_INTERN PetscErrorCode DMPlexView_HDF5_Internal(DM, PetscViewer);
244: PETSC_INTERN PetscErrorCode DMPlexLoad_HDF5_Internal(DM, PetscViewer);
245: PETSC_INTERN PetscErrorCode DMPlexLoad_HDF5_Xdmf_Internal(DM, PetscViewer);
246: PETSC_INTERN PetscErrorCode VecView_Plex_HDF5_Internal(Vec, PetscViewer);
247: PETSC_INTERN PetscErrorCode VecView_Plex_HDF5_Native_Internal(Vec, PetscViewer);
248: PETSC_INTERN PetscErrorCode VecView_Plex_Local_HDF5_Internal(Vec, PetscViewer);
249: PETSC_INTERN PetscErrorCode VecLoad_Plex_HDF5_Internal(Vec, PetscViewer);
250: PETSC_INTERN PetscErrorCode VecLoad_Plex_HDF5_Native_Internal(Vec, PetscViewer);
251: #endif
253: PETSC_INTERN PetscErrorCode DMPlexVecGetClosureAtDepth_Internal(DM, PetscSection, Vec, PetscInt, PetscInt, PetscInt *, PetscScalar *[]);
254: PETSC_INTERN PetscErrorCode DMPlexClosurePoints_Private(DM,PetscInt,const PetscInt[],IS*);
255: PETSC_INTERN PetscErrorCode DMSetFromOptions_NonRefinement_Plex(PetscOptionItems *, DM);
256: PETSC_INTERN PetscErrorCode DMCoarsen_Plex(DM, MPI_Comm, DM *);
257: PETSC_INTERN PetscErrorCode DMCoarsenHierarchy_Plex(DM, PetscInt, DM []);
258: PETSC_INTERN PetscErrorCode DMRefine_Plex(DM, MPI_Comm, DM *);
259: PETSC_INTERN PetscErrorCode DMRefineHierarchy_Plex(DM, PetscInt, DM []);
260: PETSC_INTERN PetscErrorCode DMAdaptLabel_Plex(DM, Vec, DMLabel, DMLabel, DM *);
261: PETSC_INTERN PetscErrorCode DMExtrude_Plex(DM, PetscInt, DM *);
262: PETSC_INTERN PetscErrorCode DMPlexInsertBoundaryValues_Plex(DM, PetscBool, Vec, PetscReal, Vec, Vec, Vec);
263: PETSC_INTERN PetscErrorCode DMPlexInsertTimeDerivativeBoundaryValues_Plex(DM, PetscBool, Vec, PetscReal, Vec, Vec, Vec);
264: PETSC_INTERN PetscErrorCode DMProjectFunctionLocal_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
265: 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);
266: 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);
267: 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);
268: 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);
269: PETSC_INTERN PetscErrorCode DMComputeL2Diff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
270: PETSC_INTERN PetscErrorCode DMComputeL2GradientDiff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[], const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,const PetscReal [],PetscReal *);
271: PETSC_INTERN PetscErrorCode DMComputeL2FieldDiff_Plex(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
272: PETSC_INTERN PetscErrorCode DMLocatePoints_Plex(DM, Vec, DMPointLocationType, PetscSF);
274: #if defined(PETSC_HAVE_EXODUSII)
275: PETSC_EXTERN PetscErrorCode DMView_PlexExodusII(DM, PetscViewer);
276: PETSC_INTERN PetscErrorCode VecView_PlexExodusII_Internal(Vec, PetscViewer);
277: PETSC_INTERN PetscErrorCode VecLoad_PlexExodusII_Internal(Vec, PetscViewer);
278: PETSC_INTERN PetscErrorCode VecViewPlex_ExodusII_Nodal_Internal(Vec, int, int, int);
279: PETSC_INTERN PetscErrorCode VecLoadPlex_ExodusII_Nodal_Internal(Vec, int, int, int);
280: PETSC_INTERN PetscErrorCode VecViewPlex_ExodusII_Zonal_Internal(Vec, int, int, int);
281: PETSC_INTERN PetscErrorCode VecLoadPlex_ExodusII_Zonal_Internal(Vec, int, int, int);
282: PETSC_INTERN PetscErrorCode EXOGetVarIndex_Internal(int, ex_entity_type, const char[], int*);
283: #endif
284: PETSC_INTERN PetscErrorCode DMPlexVTKGetCellType_Internal(DM,PetscInt,PetscInt,PetscInt*);
285: PETSC_INTERN PetscErrorCode DMPlexGetAdjacency_Internal(DM,PetscInt,PetscBool,PetscBool,PetscBool,PetscInt*,PetscInt*[]);
286: PETSC_INTERN PetscErrorCode DMPlexGetRawFaces_Internal(DM,DMPolytopeType,const PetscInt[],PetscInt*,const DMPolytopeType*[],const PetscInt*[],const PetscInt*[]);
287: PETSC_INTERN PetscErrorCode DMPlexRestoreRawFaces_Internal(DM,DMPolytopeType,const PetscInt[],PetscInt*,const DMPolytopeType*[],const PetscInt*[],const PetscInt*[]);
288: PETSC_INTERN PetscErrorCode DMPlexComputeCellType_Internal(DM, PetscInt, PetscInt, DMPolytopeType *);
289: PETSC_INTERN PetscErrorCode DMPlexVecSetFieldClosure_Internal(DM, PetscSection, Vec, PetscBool[], PetscInt, PetscInt, const PetscInt[], DMLabel, PetscInt, const PetscScalar[], InsertMode);
290: PETSC_INTERN PetscErrorCode DMPlexProjectConstraints_Internal(DM, Vec, Vec);
291: PETSC_EXTERN PetscErrorCode DMPlexCreateReferenceTree_SetTree(DM, PetscSection, PetscInt[], PetscInt[]);
292: PETSC_EXTERN PetscErrorCode DMPlexCreateReferenceTree_Union(DM,DM,const char *,DM*);
293: PETSC_EXTERN PetscErrorCode DMPlexComputeInterpolatorTree(DM,DM,PetscSF,PetscInt *,Mat);
294: PETSC_EXTERN PetscErrorCode DMPlexComputeInjectorTree(DM,DM,PetscSF,PetscInt *,Mat);
295: PETSC_EXTERN PetscErrorCode DMPlexAnchorsModifyMat(DM,PetscSection,PetscInt,PetscInt,const PetscInt[],const PetscInt ***,const PetscScalar[],PetscInt*,PetscInt*,PetscInt*[],PetscScalar*[],PetscInt[],PetscBool);
296: PETSC_EXTERN PetscErrorCode indicesPoint_private(PetscSection,PetscInt,PetscInt,PetscInt *,PetscBool,PetscInt,PetscInt []);
297: PETSC_EXTERN PetscErrorCode indicesPointFields_private(PetscSection,PetscInt,PetscInt,PetscInt [],PetscBool,PetscInt,PetscInt []);
298: PETSC_INTERN PetscErrorCode DMPlexLocatePoint_Internal(DM,PetscInt,const PetscScalar [],PetscInt,PetscInt *);
299: /* these two are PETSC_EXTERN just because of src/dm/impls/plex/tests/ex18.c */
300: PETSC_EXTERN PetscErrorCode DMPlexOrientInterface_Internal(DM);
302: /* Applications may use this function */
303: PETSC_EXTERN PetscErrorCode DMPlexCreateNumbering_Plex(DM, PetscInt, PetscInt, PetscInt, PetscInt *, PetscSF, IS *);
305: PETSC_INTERN PetscErrorCode DMPlexCreateCellNumbering_Internal(DM, PetscBool, IS *);
306: PETSC_INTERN PetscErrorCode DMPlexCreateVertexNumbering_Internal(DM, PetscBool, IS *);
307: PETSC_INTERN PetscErrorCode DMPlexRefine_Internal(DM, Vec, DMLabel, DMLabel, DM *);
308: PETSC_INTERN PetscErrorCode DMPlexCoarsen_Internal(DM, Vec, DMLabel, DMLabel, DM *);
309: PETSC_INTERN PetscErrorCode DMCreateMatrix_Plex(DM, Mat*);
311: PETSC_INTERN PetscErrorCode DMPlexGetOverlap_Plex(DM, PetscInt *);
312: PETSC_INTERN PetscErrorCode DMPlexDistributeGetDefault_Plex(DM, PetscBool *);
313: PETSC_INTERN PetscErrorCode DMPlexDistributeSetDefault_Plex(DM, PetscBool);
315: #if 1
316: static inline PetscInt DihedralInvert(PetscInt N, PetscInt a)
317: {
318: return (a <= 0) ? a : (N - a);
319: }
321: static inline PetscInt DihedralCompose(PetscInt N, PetscInt a, PetscInt b)
322: {
323: if (!N) return 0;
324: return (a >= 0) ?
325: ((b >= 0) ? ((a + b) % N) : -(((a - b - 1) % N) + 1)) :
326: ((b >= 0) ? -(((N - b - a - 1) % N) + 1) : ((N + b - a) % N));
327: }
329: static inline PetscInt DihedralSwap(PetscInt N, PetscInt a, PetscInt b)
330: {
331: return DihedralCompose(N,DihedralInvert(N,a),b);
332: }
333: #else
334: /* TODO
335: This is a reimplementation of the tensor dihedral symmetries using the new orientations.
336: These should be turned on when we convert to new-style orientations in p4est.
337: */
338: /* invert dihedral symmetry: return a^-1,
339: * using the representation described in
340: * DMPlexGetConeOrientation() */
341: static inline PetscInt DihedralInvert(PetscInt N, PetscInt a)
342: {
343: switch (N) {
344: case 0: return 0;
345: case 2: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_SEGMENT, 0, a);
346: case 4: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_QUADRILATERAL, 0, a);
347: case 8: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_HEXAHEDRON, 0, a);
348: default: SETERRABORT(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid celltype for DihedralInvert()");
349: }
350: return 0;
351: }
353: /* compose dihedral symmetry: return b * a,
354: * using the representation described in
355: * DMPlexGetConeOrientation() */
356: static inline PetscInt DihedralCompose(PetscInt N, PetscInt a, PetscInt b)
357: {
358: switch (N) {
359: case 0: return 0;
360: case 2: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_SEGMENT, b, a);
361: case 4: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_QUADRILATERAL, b, a);
362: case 8: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_HEXAHEDRON, b, a);
363: default: SETERRABORT(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid celltype for DihedralCompose()");
364: }
365: return 0;
366: }
368: /* swap dihedral symmetries: return b * a^-1,
369: * using the representation described in
370: * DMPlexGetConeOrientation() */
371: static inline PetscInt DihedralSwap(PetscInt N, PetscInt a, PetscInt b)
372: {
373: switch (N) {
374: case 0: return 0;
375: case 2: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_SEGMENT, b, a);
376: case 4: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_QUADRILATERAL, b, a);
377: case 8: return DMPolytopeTypeComposeOrientationInv(DM_POLYTOPE_HEXAHEDRON, b, a);
378: default: SETERRABORT(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid celltype for DihedralCompose()");
379: }
380: return 0;
381: }
382: #endif
383: PETSC_EXTERN PetscInt DMPolytopeConvertNewOrientation_Internal(DMPolytopeType, PetscInt);
384: PETSC_EXTERN PetscInt DMPolytopeConvertOldOrientation_Internal(DMPolytopeType, PetscInt);
385: PETSC_EXTERN PetscErrorCode DMPlexConvertOldOrientations_Internal(DM);
387: PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Internal(DM, PetscFormKey, IS, PetscReal, Vec, Vec, PetscReal, Vec, void *);
388: PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Hybrid_Internal(DM, PetscFormKey[], IS, PetscReal, Vec, Vec, PetscReal, Vec, void *);
389: PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Internal(DM, PetscFormKey, IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
390: PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Hybrid_Internal(DM, PetscFormKey[], IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
391: PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Action_Internal(DM, PetscFormKey, IS, PetscReal, PetscReal, Vec, Vec, Vec, Vec, void *);
392: PETSC_EXTERN PetscErrorCode DMPlexReconstructGradients_Internal(DM, PetscFV, PetscInt, PetscInt, Vec, Vec, Vec, Vec);
394: /* Matvec with A in row-major storage, x and y can be aliased */
395: static inline void DMPlex_Mult2D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
396: {
397: const PetscScalar z[2] = {x[0*ldx], x[1*ldx]};
398: y[0*ldx] = A[0]*z[0] + A[1]*z[1];
399: y[1*ldx] = A[2]*z[0] + A[3]*z[1];
400: (void)PetscLogFlops(6.0);
401: }
402: static inline void DMPlex_Mult3D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
403: {
404: const PetscScalar z[3] = {x[0*ldx], x[1*ldx], x[2*ldx]};
405: y[0*ldx] = A[0]*z[0] + A[1]*z[1] + A[2]*z[2];
406: y[1*ldx] = A[3]*z[0] + A[4]*z[1] + A[5]*z[2];
407: y[2*ldx] = A[6]*z[0] + A[7]*z[1] + A[8]*z[2];
408: (void)PetscLogFlops(15.0);
409: }
410: static inline void DMPlex_MultTranspose2D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
411: {
412: const PetscScalar z[2] = {x[0*ldx], x[1*ldx]};
413: y[0*ldx] = A[0]*z[0] + A[2]*z[1];
414: y[1*ldx] = A[1]*z[0] + A[3]*z[1];
415: (void)PetscLogFlops(6.0);
416: }
417: static inline void DMPlex_MultTranspose3D_Internal(const PetscScalar A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
418: {
419: const PetscScalar z[3] = {x[0*ldx], x[1*ldx], x[2*ldx]};
420: y[0*ldx] = A[0]*z[0] + A[3]*z[1] + A[6]*z[2];
421: y[1*ldx] = A[1]*z[0] + A[4]*z[1] + A[7]*z[2];
422: y[2*ldx] = A[2]*z[0] + A[5]*z[1] + A[8]*z[2];
423: (void)PetscLogFlops(15.0);
424: }
425: static inline void DMPlex_Mult2DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
426: {
427: const PetscScalar z[2] = {x[0*ldx], x[1*ldx]};
428: y[0*ldx] = A[0]*z[0] + A[1]*z[1];
429: y[1*ldx] = A[2]*z[0] + A[3]*z[1];
430: (void)PetscLogFlops(6.0);
431: }
432: static inline void DMPlex_Mult3DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
433: {
434: const PetscScalar z[3] = {x[0*ldx], x[1*ldx], x[2*ldx]};
435: y[0*ldx] = A[0]*z[0] + A[1]*z[1] + A[2]*z[2];
436: y[1*ldx] = A[3]*z[0] + A[4]*z[1] + A[5]*z[2];
437: y[2*ldx] = A[6]*z[0] + A[7]*z[1] + A[8]*z[2];
438: (void)PetscLogFlops(15.0);
439: }
440: static inline void DMPlex_MultAdd2DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
441: {
442: const PetscScalar z[2] = {x[0*ldx], x[1*ldx]};
443: y[0*ldx] += A[0]*z[0] + A[1]*z[1];
444: y[1*ldx] += A[2]*z[0] + A[3]*z[1];
445: (void)PetscLogFlops(6.0);
446: }
447: static inline void DMPlex_MultAdd3DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
448: {
449: const PetscScalar z[3] = {x[0*ldx], x[1*ldx], x[2*ldx]};
450: y[0*ldx] += A[0]*z[0] + A[1]*z[1] + A[2]*z[2];
451: y[1*ldx] += A[3]*z[0] + A[4]*z[1] + A[5]*z[2];
452: y[2*ldx] += A[6]*z[0] + A[7]*z[1] + A[8]*z[2];
453: (void)PetscLogFlops(15.0);
454: }
455: static inline void DMPlex_MultTransposeReal_Internal(const PetscReal A[], PetscInt m, PetscInt n, PetscInt ldx, const PetscScalar x[], PetscScalar y[])
456: {
457: PetscScalar z[3];
458: PetscInt i, j;
459: for (i = 0; i < m; ++i) z[i] = x[i*ldx];
460: for (j = 0; j < n; ++j) {
461: const PetscInt l = j*ldx;
462: y[l] = 0;
463: for (i = 0; i < m; ++i) {
464: y[l] += A[j*n+i]*z[i];
465: }
466: }
467: (void)PetscLogFlops(2*m*n);
468: }
469: static inline void DMPlex_MultTranspose2DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
470: {
471: const PetscScalar z[2] = {x[0*ldx], x[1*ldx]};
472: y[0*ldx] = A[0]*z[0] + A[2]*z[1];
473: y[1*ldx] = A[1]*z[0] + A[3]*z[1];
474: (void)PetscLogFlops(6.0);
475: }
476: static inline void DMPlex_MultTranspose3DReal_Internal(const PetscReal A[], PetscInt ldx, const PetscScalar x[], PetscScalar y[])
477: {
478: const PetscScalar z[3] = {x[0*ldx], x[1*ldx], x[2*ldx]};
479: y[0*ldx] = A[0]*z[0] + A[3]*z[1] + A[6]*z[2];
480: y[1*ldx] = A[1]*z[0] + A[4]*z[1] + A[7]*z[2];
481: y[2*ldx] = A[2]*z[0] + A[5]*z[1] + A[8]*z[2];
482: (void)PetscLogFlops(15.0);
483: }
485: static inline void DMPlex_MatMult2D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
486: {
487: #define PLEX_DIM__ 2
488: PetscScalar z[PLEX_DIM__];
489: for (PetscInt j = 0; j < n; ++j) {
490: for (int d = 0; d < PLEX_DIM__; ++d) z[d] = B[d*ldb+j];
491: DMPlex_Mult2D_Internal(A, 1, z, z);
492: for (int d = 0; d < PLEX_DIM__; ++d) C[d*ldb+j] = z[d];
493: }
494: (void)PetscLogFlops(8.0*n);
495: #undef PLEX_DIM__
496: }
497: static inline void DMPlex_MatMult3D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
498: {
499: #define PLEX_DIM__ 3
500: PetscScalar z[PLEX_DIM__];
501: for (PetscInt j = 0; j < n; ++j) {
502: for (int d = 0; d < PLEX_DIM__; ++d) z[d] = B[d*ldb+j];
503: DMPlex_Mult3D_Internal(A, 1, z, z);
504: for (int d = 0; d < PLEX_DIM__; ++d) C[d*ldb+j] = z[d];
505: }
506: (void)PetscLogFlops(8.0*n);
507: #undef PLEX_DIM__
508: }
509: static inline void DMPlex_MatMultTranspose2D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
510: {
511: #define PLEX_DIM__ 2
512: PetscScalar z[PLEX_DIM__];
513: for (PetscInt j = 0; j < n; ++j) {
514: for (int d = 0; d < PLEX_DIM__; ++d) z[d] = B[d*ldb+j];
515: DMPlex_MultTranspose2D_Internal(A, 1, z, z);
516: for (int d = 0; d < PLEX_DIM__; ++d) C[d*ldb+j] = z[d];
517: }
518: (void)PetscLogFlops(8.0*n);
519: #undef PLEX_DIM__
520: }
521: static inline void DMPlex_MatMultTranspose3D_Internal(const PetscScalar A[], PetscInt n, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
522: {
523: #define PLEX_DIM__ 3
524: PetscScalar z[PLEX_DIM__];
525: for (PetscInt j = 0; j < n; ++j) {
526: for (int d = 0; d < PLEX_DIM__; ++d) z[d] = B[d*ldb+j];
527: DMPlex_MultTranspose3D_Internal(A, 1, z, z);
528: for (int d = 0; d < PLEX_DIM__; ++d) C[d*ldb+j] = z[d];
529: }
530: (void)PetscLogFlops(8.0*n);
531: #undef PLEX_DIM__
532: }
534: static inline void DMPlex_MatMultLeft2D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
535: {
536: for (PetscInt j = 0; j < m; ++j) DMPlex_MultTranspose2D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
537: (void)PetscLogFlops(8.0*m);
538: }
539: static inline void DMPlex_MatMultLeft3D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
540: {
541: for (PetscInt j = 0; j < m; ++j) DMPlex_MultTranspose3D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
542: (void)PetscLogFlops(8.0*m);
543: }
544: static inline void DMPlex_MatMultTransposeLeft2D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
545: {
546: for (PetscInt j = 0; j < m; ++j) DMPlex_Mult2D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
547: (void)PetscLogFlops(8.0*m);
548: }
549: static inline void DMPlex_MatMultTransposeLeft3D_Internal(const PetscScalar A[], PetscInt m, PetscInt ldb, const PetscScalar B[], PetscScalar C[])
550: {
551: for (PetscInt j = 0; j < m; ++j) DMPlex_Mult3D_Internal(A, 1, &B[j*ldb], &C[j*ldb]);
552: (void)PetscLogFlops(8.0*m);
553: }
555: static inline void DMPlex_PTAP2DReal_Internal(const PetscReal P[], const PetscScalar A[], PetscScalar C[])
556: {
557: PetscScalar out[4];
558: PetscInt i, j, k, l;
559: for (i = 0; i < 2; ++i) {
560: for (j = 0; j < 2; ++j) {
561: out[i*2+j] = 0.;
562: for (k = 0; k < 2; ++k) {
563: for (l = 0; l < 2; ++l) {
564: out[i*2+j] += P[k*2+i]*A[k*2+l]*P[l*2+j];
565: }
566: }
567: }
568: }
569: for (i = 0; i < 2*2; ++i) C[i] = out[i];
570: (void)PetscLogFlops(48.0);
571: }
572: static inline void DMPlex_PTAP3DReal_Internal(const PetscReal P[], const PetscScalar A[], PetscScalar C[])
573: {
574: PetscScalar out[9];
575: PetscInt i, j, k, l;
576: for (i = 0; i < 3; ++i) {
577: for (j = 0; j < 3; ++j) {
578: out[i*3+j] = 0.;
579: for (k = 0; k < 3; ++k) {
580: for (l = 0; l < 3; ++l) {
581: out[i*3+j] += P[k*3+i]*A[k*3+l]*P[l*3+j];
582: }
583: }
584: }
585: }
586: for (i = 0; i < 2*2; ++i) C[i] = out[i];
587: (void)PetscLogFlops(243.0);
588: }
589: /* TODO Fix for aliasing of A and C */
590: static inline void DMPlex_PTAPReal_Internal(const PetscReal P[], PetscInt m, PetscInt n, const PetscScalar A[], PetscScalar C[])
591: {
592: PetscInt i, j, k, l;
593: for (i = 0; i < n; ++i) {
594: for (j = 0; j < n; ++j) {
595: C[i*n+j] = 0.;
596: for (k = 0; k < m; ++k) {
597: for (l = 0; l < m; ++l) {
598: C[i*n+j] += P[k*n+i]*A[k*m+l]*P[l*n+j];
599: }
600: }
601: }
602: }
603: (void)PetscLogFlops(243.0);
604: }
606: static inline void DMPlex_Transpose2D_Internal(PetscScalar A[])
607: {
608: PetscScalar tmp;
609: tmp = A[1]; A[1] = A[2]; A[2] = tmp;
610: }
611: static inline void DMPlex_Transpose3D_Internal(PetscScalar A[])
612: {
613: PetscScalar tmp;
614: tmp = A[1]; A[1] = A[3]; A[3] = tmp;
615: tmp = A[2]; A[2] = A[6]; A[6] = tmp;
616: tmp = A[5]; A[5] = A[7]; A[7] = tmp;
617: }
619: static inline void DMPlex_Invert2D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
620: {
621: const PetscReal invDet = 1.0/detJ;
623: invJ[0] = invDet*J[3];
624: invJ[1] = -invDet*J[1];
625: invJ[2] = -invDet*J[2];
626: invJ[3] = invDet*J[0];
627: (void)PetscLogFlops(5.0);
628: }
630: static inline void DMPlex_Invert3D_Internal(PetscReal invJ[], PetscReal J[], PetscReal detJ)
631: {
632: const PetscReal invDet = 1.0/detJ;
634: invJ[0*3+0] = invDet*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]);
635: invJ[0*3+1] = invDet*(J[0*3+2]*J[2*3+1] - J[0*3+1]*J[2*3+2]);
636: invJ[0*3+2] = invDet*(J[0*3+1]*J[1*3+2] - J[0*3+2]*J[1*3+1]);
637: invJ[1*3+0] = invDet*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]);
638: invJ[1*3+1] = invDet*(J[0*3+0]*J[2*3+2] - J[0*3+2]*J[2*3+0]);
639: invJ[1*3+2] = invDet*(J[0*3+2]*J[1*3+0] - J[0*3+0]*J[1*3+2]);
640: invJ[2*3+0] = invDet*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]);
641: invJ[2*3+1] = invDet*(J[0*3+1]*J[2*3+0] - J[0*3+0]*J[2*3+1]);
642: invJ[2*3+2] = invDet*(J[0*3+0]*J[1*3+1] - J[0*3+1]*J[1*3+0]);
643: (void)PetscLogFlops(37.0);
644: }
646: static inline void DMPlex_Det2D_Internal(PetscReal *detJ, const PetscReal J[])
647: {
648: *detJ = J[0]*J[3] - J[1]*J[2];
649: (void)PetscLogFlops(3.0);
650: }
652: static inline void DMPlex_Det3D_Internal(PetscReal *detJ, const PetscReal J[])
653: {
654: *detJ = (J[0*3+0]*(J[1*3+1]*J[2*3+2] - J[1*3+2]*J[2*3+1]) +
655: J[0*3+1]*(J[1*3+2]*J[2*3+0] - J[1*3+0]*J[2*3+2]) +
656: J[0*3+2]*(J[1*3+0]*J[2*3+1] - J[1*3+1]*J[2*3+0]));
657: (void)PetscLogFlops(12.0);
658: }
660: static inline void DMPlex_Det2D_Scalar_Internal(PetscReal *detJ, const PetscScalar J[])
661: {
662: *detJ = PetscRealPart(J[0])*PetscRealPart(J[3]) - PetscRealPart(J[1])*PetscRealPart(J[2]);
663: (void)PetscLogFlops(3.0);
664: }
666: static inline void DMPlex_Det3D_Scalar_Internal(PetscReal *detJ, const PetscScalar J[])
667: {
668: *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])) +
669: 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])) +
670: 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])));
671: (void)PetscLogFlops(12.0);
672: }
674: 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];}
676: 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;}
678: 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;}
680: 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);}
682: static inline PetscReal DMPlex_DistD_Internal(PetscInt dim, const PetscScalar *x, const PetscScalar *y) {PetscReal sum = 0.0; PetscInt d; for (d = 0; d < dim; ++d) sum += PetscRealPart(PetscConj(x[d] - y[d])*(x[d] - y[d])); return PetscSqrtReal(sum);}
684: PETSC_INTERN PetscErrorCode DMPlexGetPointDualSpaceFEM(DM,PetscInt,PetscInt,PetscDualSpace *);
685: PETSC_INTERN PetscErrorCode DMPlexGetIndicesPoint_Internal(PetscSection,PetscBool,PetscInt,PetscInt,PetscInt *,PetscBool,const PetscInt[],const PetscInt[],PetscInt[]);
686: PETSC_INTERN PetscErrorCode DMPlexGetIndicesPointFields_Internal(PetscSection,PetscBool,PetscInt,PetscInt,PetscInt[],PetscBool,const PetscInt***,PetscInt,const PetscInt[],PetscInt[]);
687: PETSC_INTERN PetscErrorCode DMPlexGetTransitiveClosure_Internal(DM, PetscInt, PetscInt, PetscBool, PetscInt *, PetscInt *[]);
689: PETSC_EXTERN PetscErrorCode DMPlexGetAllCells_Internal(DM, IS *);
690: PETSC_EXTERN PetscErrorCode DMSNESGetFEGeom(DMField, IS, PetscQuadrature, PetscBool, PetscFEGeom **);
691: PETSC_EXTERN PetscErrorCode DMSNESRestoreFEGeom(DMField, IS, PetscQuadrature, PetscBool, PetscFEGeom **);
692: PETSC_EXTERN PetscErrorCode DMPlexComputeResidual_Patch_Internal(DM, PetscSection, IS, PetscReal, Vec, Vec, Vec, void *);
693: PETSC_EXTERN PetscErrorCode DMPlexComputeJacobian_Patch_Internal(DM, PetscSection, PetscSection, IS, PetscReal, PetscReal, Vec, Vec, Mat, Mat, void *);
694: PETSC_INTERN PetscErrorCode DMCreateSubDomainDM_Plex(DM,DMLabel,PetscInt,IS*,DM*);
695: PETSC_INTERN PetscErrorCode DMPlexBasisTransformPoint_Internal(DM, DM, Vec, PetscInt, PetscBool[], PetscBool, PetscScalar *);
696: PETSC_EXTERN PetscErrorCode DMPlexBasisTransformPointTensor_Internal(DM, DM, Vec, PetscInt, PetscBool, PetscInt, PetscScalar *);
697: PETSC_INTERN PetscErrorCode DMPlexBasisTransformApplyReal_Internal(DM, const PetscReal[], PetscBool, PetscInt, const PetscReal *, PetscReal *, void *);
698: PETSC_INTERN PetscErrorCode DMPlexBasisTransformApply_Internal(DM, const PetscReal[], PetscBool, PetscInt, const PetscScalar *, PetscScalar *, void *);
699: PETSC_INTERN PetscErrorCode DMCreateNeumannOverlap_Plex(DM, IS*, Mat*, PetscErrorCode (**)(Mat, PetscReal, Vec, Vec, PetscReal, IS, void*), void **);
701: /* Functions in the vtable */
702: PETSC_INTERN PetscErrorCode DMCreateInterpolation_Plex(DM dmCoarse, DM dmFine, Mat *interpolation, Vec *scaling);
703: PETSC_INTERN PetscErrorCode DMCreateInjection_Plex(DM dmCoarse, DM dmFine, Mat *mat);
704: PETSC_INTERN PetscErrorCode DMCreateMassMatrix_Plex(DM dmCoarse, DM dmFine, Mat *mat);
705: PETSC_INTERN PetscErrorCode DMCreateMassMatrixLumped_Plex(DM dmCoarse, Vec *lm);
706: PETSC_INTERN PetscErrorCode DMCreateLocalSection_Plex(DM dm);
707: PETSC_INTERN PetscErrorCode DMCreateDefaultConstraints_Plex(DM dm);
708: PETSC_INTERN PetscErrorCode DMCreateMatrix_Plex(DM dm, Mat *J);
709: PETSC_INTERN PetscErrorCode DMCreateCoordinateDM_Plex(DM dm, DM *cdm);
710: PETSC_INTERN PetscErrorCode DMCreateCoordinateField_Plex(DM dm, DMField *field);
711: PETSC_INTERN PetscErrorCode DMClone_Plex(DM dm, DM *newdm);
712: PETSC_INTERN PetscErrorCode DMSetUp_Plex(DM dm);
713: PETSC_INTERN PetscErrorCode DMDestroy_Plex(DM dm);
714: PETSC_INTERN PetscErrorCode DMView_Plex(DM dm, PetscViewer viewer);
715: PETSC_INTERN PetscErrorCode DMLoad_Plex(DM dm, PetscViewer viewer);
716: PETSC_INTERN PetscErrorCode DMCreateSubDM_Plex(DM dm, PetscInt numFields, const PetscInt fields[], IS *is, DM *subdm);
717: PETSC_INTERN PetscErrorCode DMCreateSuperDM_Plex(DM dms[], PetscInt len, IS **is, DM *superdm);
719: #endif /* _PLEXIMPL_H */