xref: /petsc/include/petscdm.h (revision ac09b9214d23ea9ad238aa607de9fa447fd4e91b)
1 /*
2       Objects to manage the interactions between the mesh data structures and the algebraic objects
3 */
4 #if !defined(PETSCDM_H)
5 #define PETSCDM_H
6 #include <petscmat.h>
7 #include <petscdmtypes.h>
8 #include <petscfetypes.h>
9 #include <petscdstypes.h>
10 #include <petscdmlabel.h>
11 
12 /* SUBMANSEC = DM */
13 
14 PETSC_EXTERN PetscErrorCode DMInitializePackage(void);
15 
16 PETSC_EXTERN PetscClassId DM_CLASSID;
17 PETSC_EXTERN PetscClassId DMLABEL_CLASSID;
18 
19 #define DMLOCATEPOINT_POINT_NOT_FOUND -367
20 
21 /*J
22     DMType - String with the name of a PETSc DM
23 
24    Level: beginner
25 
26 .seealso: `DMSetType()`, `DM`
27 J*/
28 typedef const char* DMType;
29 #define DMDA        "da"
30 #define DMCOMPOSITE "composite"
31 #define DMSLICED    "sliced"
32 #define DMSHELL     "shell"
33 #define DMPLEX      "plex"
34 #define DMREDUNDANT "redundant"
35 #define DMPATCH     "patch"
36 #define DMMOAB      "moab"
37 #define DMNETWORK   "network"
38 #define DMFOREST    "forest"
39 #define DMP4EST     "p4est"
40 #define DMP8EST     "p8est"
41 #define DMSWARM     "swarm"
42 #define DMPRODUCT   "product"
43 #define DMSTAG      "stag"
44 
45 PETSC_EXTERN const char *const DMBoundaryTypes[];
46 PETSC_EXTERN const char *const DMBoundaryConditionTypes[];
47 PETSC_EXTERN PetscFunctionList DMList;
48 PETSC_EXTERN DMGeneratorFunctionList DMGenerateList;
49 PETSC_EXTERN PetscErrorCode DMCreate(MPI_Comm,DM*);
50 PETSC_EXTERN PetscErrorCode DMClone(DM,DM*);
51 PETSC_EXTERN PetscErrorCode DMSetType(DM, DMType);
52 PETSC_EXTERN PetscErrorCode DMGetType(DM, DMType *);
53 PETSC_EXTERN PetscErrorCode DMRegister(const char[],PetscErrorCode (*)(DM));
54 PETSC_EXTERN PetscErrorCode DMRegisterDestroy(void);
55 
56 PETSC_EXTERN PetscErrorCode DMView(DM,PetscViewer);
57 PETSC_EXTERN PetscErrorCode DMLoad(DM,PetscViewer);
58 PETSC_EXTERN PetscErrorCode DMDestroy(DM*);
59 PETSC_EXTERN PetscErrorCode DMCreateGlobalVector(DM,Vec*);
60 PETSC_EXTERN PetscErrorCode DMCreateLocalVector(DM,Vec*);
61 PETSC_EXTERN PetscErrorCode DMGetLocalVector(DM,Vec *);
62 PETSC_EXTERN PetscErrorCode DMRestoreLocalVector(DM,Vec *);
63 PETSC_EXTERN PetscErrorCode DMGetGlobalVector(DM,Vec *);
64 PETSC_EXTERN PetscErrorCode DMRestoreGlobalVector(DM,Vec *);
65 PETSC_EXTERN PetscErrorCode DMClearGlobalVectors(DM);
66 PETSC_EXTERN PetscErrorCode DMClearLocalVectors(DM);
67 PETSC_EXTERN PetscErrorCode DMHasNamedGlobalVector(DM,const char*,PetscBool*);
68 PETSC_EXTERN PetscErrorCode DMGetNamedGlobalVector(DM,const char*,Vec*);
69 PETSC_EXTERN PetscErrorCode DMRestoreNamedGlobalVector(DM,const char*,Vec*);
70 PETSC_EXTERN PetscErrorCode DMHasNamedLocalVector(DM,const char*,PetscBool*);
71 PETSC_EXTERN PetscErrorCode DMGetNamedLocalVector(DM,const char*,Vec*);
72 PETSC_EXTERN PetscErrorCode DMRestoreNamedLocalVector(DM,const char*,Vec*);
73 PETSC_EXTERN PetscErrorCode DMGetLocalToGlobalMapping(DM,ISLocalToGlobalMapping*);
74 PETSC_EXTERN PetscErrorCode DMCreateFieldIS(DM,PetscInt*,char***,IS**);
75 PETSC_EXTERN PetscErrorCode DMGetBlockSize(DM,PetscInt*);
76 PETSC_EXTERN PetscErrorCode DMCreateColoring(DM,ISColoringType,ISColoring*);
77 PETSC_EXTERN PetscErrorCode DMCreateMatrix(DM,Mat*);
78 PETSC_EXTERN PetscErrorCode DMSetMatrixPreallocateSkip(DM,PetscBool);
79 PETSC_EXTERN PetscErrorCode DMSetMatrixPreallocateOnly(DM,PetscBool);
80 PETSC_EXTERN PetscErrorCode DMSetMatrixStructureOnly(DM,PetscBool);
81 PETSC_EXTERN PetscErrorCode DMCreateInterpolation(DM,DM,Mat*,Vec*);
82 PETSC_EXTERN PetscErrorCode DMCreateRestriction(DM,DM,Mat*);
83 PETSC_EXTERN PetscErrorCode DMCreateInjection(DM,DM,Mat*);
84 PETSC_EXTERN PetscErrorCode DMCreateMassMatrix(DM,DM,Mat*);
85 PETSC_EXTERN PetscErrorCode DMCreateMassMatrixLumped(DM,Vec*);
86 PETSC_EXTERN PetscErrorCode DMGetWorkArray(DM,PetscInt,MPI_Datatype,void*);
87 PETSC_EXTERN PetscErrorCode DMRestoreWorkArray(DM,PetscInt,MPI_Datatype,void*);
88 PETSC_EXTERN PetscErrorCode DMRefine(DM,MPI_Comm,DM*);
89 PETSC_EXTERN PetscErrorCode DMCoarsen(DM,MPI_Comm,DM*);
90 PETSC_EXTERN PetscErrorCode DMGetCoarseDM(DM,DM*);
91 PETSC_EXTERN PetscErrorCode DMSetCoarseDM(DM,DM);
92 PETSC_EXTERN PetscErrorCode DMGetFineDM(DM,DM*);
93 PETSC_EXTERN PetscErrorCode DMSetFineDM(DM,DM);
94 PETSC_EXTERN PetscErrorCode DMRefineHierarchy(DM,PetscInt,DM[]);
95 PETSC_EXTERN PetscErrorCode DMCoarsenHierarchy(DM,PetscInt,DM[]);
96 PETSC_EXTERN PetscErrorCode DMCoarsenHookAdd(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,Mat,Vec,Mat,DM,void*),void*);
97 PETSC_EXTERN PetscErrorCode DMCoarsenHookRemove(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,Mat,Vec,Mat,DM,void*),void*);
98 PETSC_EXTERN PetscErrorCode DMRefineHookAdd(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,Mat,DM,void*),void*);
99 PETSC_EXTERN PetscErrorCode DMRefineHookRemove(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,Mat,DM,void*),void*);
100 PETSC_EXTERN PetscErrorCode DMRestrict(DM,Mat,Vec,Mat,DM);
101 PETSC_EXTERN PetscErrorCode DMInterpolate(DM,Mat,DM);
102 PETSC_EXTERN PetscErrorCode DMInterpolateSolution(DM,DM,Mat,Vec,Vec);
103 PETSC_EXTERN PetscErrorCode DMExtrude(DM,PetscInt,DM*);
104 PETSC_EXTERN PetscErrorCode DMSetFromOptions(DM);
105 PETSC_EXTERN PetscErrorCode DMViewFromOptions(DM,PetscObject,const char[]);
106 
107 PETSC_EXTERN PetscErrorCode DMGenerate(DM, const char [], PetscBool , DM *);
108 PETSC_EXTERN PetscErrorCode DMGenerateRegister(const char[],PetscErrorCode (*)(DM,PetscBool,DM*),PetscErrorCode (*)(DM,PetscReal*,DM*),PetscErrorCode (*)(DM,Vec,DMLabel,DMLabel,DM*),PetscInt);
109 PETSC_EXTERN PetscErrorCode DMGenerateRegisterAll(void);
110 PETSC_EXTERN PetscErrorCode DMGenerateRegisterDestroy(void);
111 PETSC_EXTERN PetscErrorCode DMAdaptLabel(DM,DMLabel,DM*);
112 PETSC_EXTERN PetscErrorCode DMAdaptMetric(DM, Vec, DMLabel, DMLabel, DM *);
113 
114 PETSC_EXTERN PetscErrorCode DMSetUp(DM);
115 PETSC_EXTERN PetscErrorCode DMCreateInterpolationScale(DM,DM,Mat,Vec*);
116 PETSC_EXTERN PETSC_DEPRECATED_FUNCTION("Use DMDACreateAggregates() or DMCreateRestriction() (since version 3.12)") PetscErrorCode DMCreateAggregates(DM,DM,Mat*);
117 PETSC_EXTERN PetscErrorCode DMGlobalToLocalHookAdd(DM,PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),void*);
118 PETSC_EXTERN PetscErrorCode DMLocalToGlobalHookAdd(DM,PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),PetscErrorCode (*)(DM,Vec,InsertMode,Vec,void*),void*);
119 PETSC_EXTERN PetscErrorCode DMGlobalToLocal(DM,Vec,InsertMode,Vec);
120 PETSC_EXTERN PetscErrorCode DMGlobalToLocalBegin(DM,Vec,InsertMode,Vec);
121 PETSC_EXTERN PetscErrorCode DMGlobalToLocalEnd(DM,Vec,InsertMode,Vec);
122 PETSC_EXTERN PetscErrorCode DMLocalToGlobal(DM,Vec,InsertMode,Vec);
123 PETSC_EXTERN PetscErrorCode DMLocalToGlobalBegin(DM,Vec,InsertMode,Vec);
124 PETSC_EXTERN PetscErrorCode DMLocalToGlobalEnd(DM,Vec,InsertMode,Vec);
125 PETSC_EXTERN PetscErrorCode DMLocalToLocalBegin(DM,Vec,InsertMode,Vec);
126 PETSC_EXTERN PetscErrorCode DMLocalToLocalEnd(DM,Vec,InsertMode,Vec);
127 PETSC_EXTERN PetscErrorCode DMConvert(DM,DMType,DM*);
128 
129 /* Topology support */
130 PETSC_EXTERN PetscErrorCode DMGetDimension(DM,PetscInt*);
131 PETSC_EXTERN PetscErrorCode DMSetDimension(DM,PetscInt);
132 PETSC_EXTERN PetscErrorCode DMGetDimPoints(DM,PetscInt,PetscInt*,PetscInt*);
133 PETSC_EXTERN PetscErrorCode DMGetUseNatural(DM,PetscBool*);
134 PETSC_EXTERN PetscErrorCode DMSetUseNatural(DM,PetscBool);
135 
136 /* Coordinate support */
137 PETSC_EXTERN PetscErrorCode DMGetCoordinateDM(DM,DM*);
138 PETSC_EXTERN PetscErrorCode DMSetCoordinateDM(DM,DM);
139 PETSC_EXTERN PetscErrorCode DMGetCoordinateDim(DM,PetscInt*);
140 PETSC_EXTERN PetscErrorCode DMSetCoordinateDim(DM,PetscInt);
141 PETSC_EXTERN PetscErrorCode DMGetCoordinateSection(DM,PetscSection*);
142 PETSC_EXTERN PetscErrorCode DMSetCoordinateSection(DM,PetscInt,PetscSection);
143 PETSC_EXTERN PetscErrorCode DMGetCoordinates(DM,Vec*);
144 PETSC_EXTERN PetscErrorCode DMSetCoordinates(DM,Vec);
145 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocal(DM,Vec*);
146 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocalSetUp(DM);
147 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocalNoncollective(DM,Vec*);
148 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocalTuple(DM,IS,PetscSection*,Vec*);
149 PETSC_EXTERN PetscErrorCode DMSetCoordinatesLocal(DM,Vec);
150 PETSC_EXTERN PetscErrorCode DMLocatePoints(DM,Vec,DMPointLocationType,PetscSF*);
151 PETSC_EXTERN PetscErrorCode DMGetPeriodicity(DM,PetscBool*,const PetscReal**,const PetscReal**,const DMBoundaryType**);
152 PETSC_EXTERN PetscErrorCode DMSetPeriodicity(DM,PetscBool,const PetscReal[],const PetscReal[],const DMBoundaryType[]);
153 PETSC_EXTERN PetscErrorCode DMLocalizeCoordinate(DM, const PetscScalar[], PetscBool, PetscScalar[]);
154 PETSC_EXTERN PetscErrorCode DMLocalizeCoordinates(DM);
155 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocalized(DM,PetscBool*);
156 PETSC_EXTERN PetscErrorCode DMGetCoordinatesLocalizedLocal(DM,PetscBool*);
157 PETSC_EXTERN PetscErrorCode DMGetNeighbors(DM,PetscInt*,const PetscMPIInt**);
158 PETSC_EXTERN PetscErrorCode DMGetCoordinateField(DM,DMField*);
159 PETSC_EXTERN PetscErrorCode DMSetCoordinateField(DM,DMField);
160 PETSC_EXTERN PetscErrorCode DMGetBoundingBox(DM,PetscReal[],PetscReal[]);
161 PETSC_EXTERN PetscErrorCode DMGetLocalBoundingBox(DM,PetscReal[],PetscReal[]);
162 PETSC_EXTERN PetscErrorCode DMProjectCoordinates(DM,PetscFE);
163 
164 /* block hook interface */
165 PETSC_EXTERN PetscErrorCode DMSubDomainHookAdd(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,VecScatter,VecScatter,DM,void*),void*);
166 PETSC_EXTERN PetscErrorCode DMSubDomainHookRemove(DM,PetscErrorCode (*)(DM,DM,void*),PetscErrorCode (*)(DM,VecScatter,VecScatter,DM,void*),void*);
167 PETSC_EXTERN PetscErrorCode DMSubDomainRestrict(DM,VecScatter,VecScatter,DM);
168 
169 PETSC_EXTERN PetscErrorCode DMSetOptionsPrefix(DM,const char []);
170 PETSC_EXTERN PetscErrorCode DMAppendOptionsPrefix(DM,const char []);
171 PETSC_EXTERN PetscErrorCode DMGetOptionsPrefix(DM,const char*[]);
172 PETSC_EXTERN PetscErrorCode DMSetVecType(DM,VecType);
173 PETSC_EXTERN PetscErrorCode DMGetVecType(DM,VecType*);
174 PETSC_EXTERN PetscErrorCode DMSetMatType(DM,MatType);
175 PETSC_EXTERN PetscErrorCode DMGetMatType(DM,MatType*);
176 PETSC_EXTERN PetscErrorCode DMSetISColoringType(DM,ISColoringType);
177 PETSC_EXTERN PetscErrorCode DMGetISColoringType(DM,ISColoringType*);
178 PETSC_EXTERN PetscErrorCode DMSetApplicationContext(DM,void*);
179 PETSC_EXTERN PetscErrorCode DMSetApplicationContextDestroy(DM,PetscErrorCode (*)(void**));
180 PETSC_EXTERN PetscErrorCode DMGetApplicationContext(DM,void*);
181 PETSC_EXTERN PetscErrorCode DMSetVariableBounds(DM,PetscErrorCode (*)(DM,Vec,Vec));
182 PETSC_EXTERN PetscErrorCode DMHasVariableBounds(DM,PetscBool *);
183 PETSC_EXTERN PetscErrorCode DMHasColoring(DM,PetscBool *);
184 PETSC_EXTERN PetscErrorCode DMHasCreateRestriction(DM,PetscBool *);
185 PETSC_EXTERN PetscErrorCode DMHasCreateInjection(DM,PetscBool *);
186 PETSC_EXTERN PetscErrorCode DMComputeVariableBounds(DM,Vec,Vec);
187 
188 PETSC_EXTERN PetscErrorCode DMCreateSubDM(DM, PetscInt, const PetscInt[], IS *, DM *);
189 PETSC_EXTERN PetscErrorCode DMCreateSuperDM(DM[], PetscInt, IS **, DM *);
190 PETSC_EXTERN PetscErrorCode DMCreateSectionSubDM(DM,PetscInt,const PetscInt[],IS*,DM*);
191 PETSC_EXTERN PetscErrorCode DMCreateSectionSuperDM(DM[],PetscInt,IS**,DM*);
192 PETSC_EXTERN PetscErrorCode DMCreateFieldDecomposition(DM,PetscInt*,char***,IS**,DM**);
193 PETSC_EXTERN PetscErrorCode DMCreateDomainDecomposition(DM,PetscInt*,char***,IS**,IS**,DM**);
194 PETSC_EXTERN PetscErrorCode DMCreateDomainDecompositionScatters(DM,PetscInt,DM*,VecScatter**,VecScatter**,VecScatter**);
195 
196 PETSC_EXTERN PetscErrorCode DMGetRefineLevel(DM,PetscInt*);
197 PETSC_EXTERN PetscErrorCode DMSetRefineLevel(DM,PetscInt);
198 PETSC_EXTERN PetscErrorCode DMGetCoarsenLevel(DM,PetscInt*);
199 PETSC_EXTERN PetscErrorCode DMSetCoarsenLevel(DM,PetscInt);
200 PETSC_EXTERN PetscErrorCode DMFinalizePackage(void);
201 
202 PETSC_EXTERN PetscErrorCode VecGetDM(Vec, DM*);
203 PETSC_EXTERN PetscErrorCode VecSetDM(Vec, DM);
204 PETSC_EXTERN PetscErrorCode MatGetDM(Mat, DM*);
205 PETSC_EXTERN PetscErrorCode MatSetDM(Mat, DM);
206 PETSC_EXTERN PetscErrorCode MatFDColoringUseDM(Mat,MatFDColoring);
207 
208 typedef struct NLF_DAAD* NLF;
209 
210 #define DM_FILE_CLASSID 1211221
211 
212 /* FEM support */
213 PETSC_EXTERN PetscErrorCode DMPrintCellVector(PetscInt, const char [], PetscInt, const PetscScalar []);
214 PETSC_EXTERN PetscErrorCode DMPrintCellMatrix(PetscInt, const char [], PetscInt, PetscInt, const PetscScalar []);
215 PETSC_EXTERN PetscErrorCode DMPrintLocalVec(DM, const char [], PetscReal, Vec);
216 
217 PETSC_EXTERN PetscErrorCode DMSetNullSpaceConstructor(DM, PetscInt, PetscErrorCode (*)(DM, PetscInt, PetscInt, MatNullSpace *));
218 PETSC_EXTERN PetscErrorCode DMGetNullSpaceConstructor(DM, PetscInt, PetscErrorCode (**)(DM, PetscInt, PetscInt, MatNullSpace *));
219 PETSC_EXTERN PetscErrorCode DMSetNearNullSpaceConstructor(DM, PetscInt, PetscErrorCode (*)(DM, PetscInt, PetscInt, MatNullSpace *));
220 PETSC_EXTERN PetscErrorCode DMGetNearNullSpaceConstructor(DM, PetscInt, PetscErrorCode (**)(DM, PetscInt, PetscInt, MatNullSpace *));
221 
222 PETSC_EXTERN PetscErrorCode DMGetSection(DM, PetscSection *); /* Use DMGetLocalSection() in new code (since v3.12) */
223 PETSC_EXTERN PetscErrorCode DMSetSection(DM, PetscSection);   /* Use DMSetLocalSection() in new code (since v3.12) */
224 PETSC_EXTERN PetscErrorCode DMGetLocalSection(DM, PetscSection *);
225 PETSC_EXTERN PetscErrorCode DMSetLocalSection(DM, PetscSection);
226 PETSC_EXTERN PetscErrorCode DMGetGlobalSection(DM, PetscSection *);
227 PETSC_EXTERN PetscErrorCode DMSetGlobalSection(DM, PetscSection);
228 static inline PETSC_DEPRECATED_FUNCTION("Use DMGetSection() (since v3.9)") PetscErrorCode DMGetDefaultSection(DM dm, PetscSection *s) {return DMGetSection(dm,s);}
229 static inline PETSC_DEPRECATED_FUNCTION("Use DMSetSection() (since v3.9)") PetscErrorCode DMSetDefaultSection(DM dm, PetscSection s) {return DMSetSection(dm,s);}
230 static inline PETSC_DEPRECATED_FUNCTION("Use DMGetGlobalSection() (since v3.9)") PetscErrorCode DMGetDefaultGlobalSection(DM dm, PetscSection *s) {return DMGetGlobalSection(dm,s);}
231 static inline PETSC_DEPRECATED_FUNCTION("Use DMSetGlobalSection() (since v3.9)") PetscErrorCode DMSetDefaultGlobalSection(DM dm, PetscSection s) {return DMSetGlobalSection(dm,s);}
232 
233 PETSC_EXTERN PetscErrorCode DMGetSectionSF(DM, PetscSF*);
234 PETSC_EXTERN PetscErrorCode DMSetSectionSF(DM, PetscSF);
235 PETSC_EXTERN PetscErrorCode DMCreateSectionSF(DM, PetscSection, PetscSection);
236 static inline PETSC_DEPRECATED_FUNCTION("Use DMGetSectionSF() (since v3.12)") PetscErrorCode DMGetDefaultSF(DM dm, PetscSF *s) {return DMGetSectionSF(dm,s);}
237 static inline PETSC_DEPRECATED_FUNCTION("Use DMSetSectionSF() (since v3.12)") PetscErrorCode DMSetDefaultSF(DM dm, PetscSF s) {return DMSetSectionSF(dm,s);}
238 static inline PETSC_DEPRECATED_FUNCTION("Use DMCreateSectionSF() (since v3.12)") PetscErrorCode DMCreateDefaultSF(DM dm, PetscSection l, PetscSection g) {return DMCreateSectionSF(dm,l,g);}
239 PETSC_EXTERN PetscErrorCode DMGetPointSF(DM, PetscSF *);
240 PETSC_EXTERN PetscErrorCode DMSetPointSF(DM, PetscSF);
241 PETSC_EXTERN PetscErrorCode DMGetNaturalSF(DM, PetscSF *);
242 PETSC_EXTERN PetscErrorCode DMSetNaturalSF(DM, PetscSF);
243 
244 PETSC_EXTERN PetscErrorCode DMGetDefaultConstraints(DM, PetscSection *, Mat *, Vec *);
245 PETSC_EXTERN PetscErrorCode DMSetDefaultConstraints(DM, PetscSection, Mat, Vec);
246 
247 PETSC_EXTERN PetscErrorCode DMGetOutputDM(DM, DM *);
248 PETSC_EXTERN PetscErrorCode DMGetOutputSequenceNumber(DM, PetscInt *, PetscReal *);
249 PETSC_EXTERN PetscErrorCode DMSetOutputSequenceNumber(DM, PetscInt, PetscReal);
250 PETSC_EXTERN PetscErrorCode DMOutputSequenceLoad(DM, PetscViewer, const char *, PetscInt, PetscReal *);
251 
252 PETSC_EXTERN PetscErrorCode DMGetNumFields(DM, PetscInt *);
253 PETSC_EXTERN PetscErrorCode DMSetNumFields(DM, PetscInt);
254 PETSC_EXTERN PetscErrorCode DMGetField(DM, PetscInt, DMLabel *, PetscObject *);
255 PETSC_EXTERN PetscErrorCode DMSetField(DM, PetscInt, DMLabel, PetscObject);
256 PETSC_EXTERN PetscErrorCode DMAddField(DM, DMLabel, PetscObject);
257 PETSC_EXTERN PetscErrorCode DMSetFieldAvoidTensor(DM, PetscInt, PetscBool);
258 PETSC_EXTERN PetscErrorCode DMGetFieldAvoidTensor(DM, PetscInt, PetscBool *);
259 PETSC_EXTERN PetscErrorCode DMClearFields(DM);
260 PETSC_EXTERN PetscErrorCode DMCopyFields(DM, DM);
261 PETSC_EXTERN PetscErrorCode DMGetAdjacency(DM, PetscInt, PetscBool *, PetscBool *);
262 PETSC_EXTERN PetscErrorCode DMSetAdjacency(DM, PetscInt, PetscBool, PetscBool);
263 PETSC_EXTERN PetscErrorCode DMGetBasicAdjacency(DM, PetscBool *, PetscBool *);
264 PETSC_EXTERN PetscErrorCode DMSetBasicAdjacency(DM, PetscBool, PetscBool);
265 
266 PETSC_EXTERN PetscErrorCode DMGetNumDS(DM, PetscInt *);
267 PETSC_EXTERN PetscErrorCode DMGetDS(DM, PetscDS *);
268 PETSC_EXTERN PetscErrorCode DMGetCellDS(DM, PetscInt, PetscDS *);
269 PETSC_EXTERN PetscErrorCode DMGetRegionDS(DM, DMLabel, IS *, PetscDS *);
270 PETSC_EXTERN PetscErrorCode DMSetRegionDS(DM, DMLabel, IS, PetscDS);
271 PETSC_EXTERN PetscErrorCode DMGetRegionNumDS(DM, PetscInt, DMLabel *, IS *, PetscDS *);
272 PETSC_EXTERN PetscErrorCode DMSetRegionNumDS(DM, PetscInt, DMLabel, IS, PetscDS);
273 PETSC_EXTERN PetscErrorCode DMFindRegionNum(DM, PetscDS, PetscInt *);
274 PETSC_EXTERN PetscErrorCode DMCreateFEDefault(DM, PetscInt, const char[], PetscInt, PetscFE *);
275 PETSC_EXTERN PetscErrorCode DMCreateDS(DM);
276 PETSC_EXTERN PetscErrorCode DMClearDS(DM);
277 PETSC_EXTERN PetscErrorCode DMCopyDS(DM, DM);
278 PETSC_EXTERN PetscErrorCode DMCopyDisc(DM, DM);
279 PETSC_EXTERN PetscErrorCode DMComputeExactSolution(DM, PetscReal, Vec, Vec);
280 PETSC_EXTERN PetscErrorCode DMGetNumAuxiliaryVec(DM, PetscInt *);
281 PETSC_EXTERN PetscErrorCode DMGetAuxiliaryVec(DM, DMLabel, PetscInt, PetscInt, Vec *);
282 PETSC_EXTERN PetscErrorCode DMSetAuxiliaryVec(DM, DMLabel, PetscInt, PetscInt, Vec);
283 PETSC_EXTERN PetscErrorCode DMGetAuxiliaryLabels(DM, DMLabel[], PetscInt[], PetscInt[]);
284 PETSC_EXTERN PetscErrorCode DMCopyAuxiliaryVec(DM, DM);
285 
286 /*MC
287   DMInterpolationInfo - Structure for holding information about interpolation on a mesh
288 
289   Level: intermediate
290 
291   Synopsis:
292     comm   - The communicator
293     dim    - The spatial dimension of points
294     nInput - The number of input points
295     points - The input point coordinates
296     cells  - The cell containing each point
297     n      - The number of local points
298     coords - The point coordinates
299     dof    - The number of components to interpolate
300 
301 .seealso: `DMInterpolationCreate()`, `DMInterpolationEvaluate()`, `DMInterpolationAddPoints()`
302 M*/
303 struct _DMInterpolationInfo {
304   MPI_Comm   comm;
305   PetscInt   dim;    /*1 The spatial dimension of points */
306   PetscInt   nInput; /* The number of input points */
307   PetscReal *points; /* The input point coordinates */
308   PetscInt  *cells;  /* The cell containing each point */
309   PetscInt   n;      /* The number of local points */
310   Vec        coords; /* The point coordinates */
311   PetscInt   dof;    /* The number of components to interpolate */
312 };
313 typedef struct _DMInterpolationInfo *DMInterpolationInfo;
314 
315 PETSC_EXTERN PetscErrorCode DMInterpolationCreate(MPI_Comm, DMInterpolationInfo *);
316 PETSC_EXTERN PetscErrorCode DMInterpolationSetDim(DMInterpolationInfo, PetscInt);
317 PETSC_EXTERN PetscErrorCode DMInterpolationGetDim(DMInterpolationInfo, PetscInt *);
318 PETSC_EXTERN PetscErrorCode DMInterpolationSetDof(DMInterpolationInfo, PetscInt);
319 PETSC_EXTERN PetscErrorCode DMInterpolationGetDof(DMInterpolationInfo, PetscInt *);
320 PETSC_EXTERN PetscErrorCode DMInterpolationAddPoints(DMInterpolationInfo, PetscInt, PetscReal[]);
321 PETSC_EXTERN PetscErrorCode DMInterpolationSetUp(DMInterpolationInfo, DM, PetscBool, PetscBool);
322 PETSC_EXTERN PetscErrorCode DMInterpolationGetCoordinates(DMInterpolationInfo, Vec *);
323 PETSC_EXTERN PetscErrorCode DMInterpolationGetVector(DMInterpolationInfo, Vec *);
324 PETSC_EXTERN PetscErrorCode DMInterpolationRestoreVector(DMInterpolationInfo, Vec *);
325 PETSC_EXTERN PetscErrorCode DMInterpolationEvaluate(DMInterpolationInfo, DM, Vec, Vec);
326 PETSC_EXTERN PetscErrorCode DMInterpolationDestroy(DMInterpolationInfo *);
327 
328 PETSC_EXTERN PetscErrorCode DMCreateLabel(DM, const char []);
329 PETSC_EXTERN PetscErrorCode DMGetLabelValue(DM, const char[], PetscInt, PetscInt *);
330 PETSC_EXTERN PetscErrorCode DMSetLabelValue(DM, const char[], PetscInt, PetscInt);
331 PETSC_EXTERN PetscErrorCode DMClearLabelValue(DM, const char[], PetscInt, PetscInt);
332 PETSC_EXTERN PetscErrorCode DMGetLabelSize(DM, const char[], PetscInt *);
333 PETSC_EXTERN PetscErrorCode DMGetLabelIdIS(DM, const char[], IS *);
334 PETSC_EXTERN PetscErrorCode DMGetStratumSize(DM, const char [], PetscInt, PetscInt *);
335 PETSC_EXTERN PetscErrorCode DMGetStratumIS(DM, const char [], PetscInt, IS *);
336 PETSC_EXTERN PetscErrorCode DMSetStratumIS(DM, const char [], PetscInt, IS);
337 PETSC_EXTERN PetscErrorCode DMClearLabelStratum(DM, const char[], PetscInt);
338 PETSC_EXTERN PetscErrorCode DMGetLabelOutput(DM, const char[], PetscBool *);
339 PETSC_EXTERN PetscErrorCode DMSetLabelOutput(DM, const char[], PetscBool);
340 PETSC_EXTERN PetscErrorCode DMGetFirstLabeledPoint(DM, DM, DMLabel, PetscInt, const PetscInt *, PetscInt, PetscInt *, PetscDS *);
341 
342 /*E
343    DMCopyLabelsMode - Determines how DMCopyLabels() behaves when there is a DMLabel in the source and destination DMs with the same name
344 
345    Level: advanced
346 
347 $ DM_COPY_LABELS_REPLACE  - replace label in destination by label from source
348 $ DM_COPY_LABELS_KEEP     - keep destination label
349 $ DM_COPY_LABELS_FAIL     - throw error
350 
351 E*/
352 typedef enum {DM_COPY_LABELS_REPLACE, DM_COPY_LABELS_KEEP, DM_COPY_LABELS_FAIL} DMCopyLabelsMode;
353 PETSC_EXTERN const char *const DMCopyLabelsModes[];
354 
355 PETSC_EXTERN PetscErrorCode DMGetNumLabels(DM, PetscInt *);
356 PETSC_EXTERN PetscErrorCode DMGetLabelName(DM, PetscInt, const char **);
357 PETSC_EXTERN PetscErrorCode DMHasLabel(DM, const char [], PetscBool *);
358 PETSC_EXTERN PetscErrorCode DMGetLabel(DM, const char *, DMLabel *);
359 PETSC_EXTERN PetscErrorCode DMSetLabel(DM, DMLabel);
360 PETSC_EXTERN PetscErrorCode DMGetLabelByNum(DM, PetscInt, DMLabel *);
361 PETSC_EXTERN PetscErrorCode DMAddLabel(DM, DMLabel);
362 PETSC_EXTERN PetscErrorCode DMRemoveLabel(DM, const char [], DMLabel *);
363 PETSC_EXTERN PetscErrorCode DMRemoveLabelBySelf(DM, DMLabel *, PetscBool);
364 PETSC_EXTERN PetscErrorCode DMCopyLabels(DM, DM, PetscCopyMode, PetscBool, DMCopyLabelsMode emode);
365 PETSC_EXTERN PetscErrorCode DMCompareLabels(DM, DM, PetscBool *, char **);
366 
367 PETSC_EXTERN PetscErrorCode DMAddBoundary(DM, DMBoundaryConditionType, const char[], DMLabel, PetscInt, const PetscInt[], PetscInt, PetscInt, const PetscInt[], void (*)(void), void (*)(void), void *, PetscInt *);
368 PETSC_EXTERN PetscErrorCode DMIsBoundaryPoint(DM, PetscInt, PetscBool *);
369 
370 PETSC_EXTERN PetscErrorCode DMProjectFunction(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void**,InsertMode,Vec);
371 PETSC_EXTERN PetscErrorCode DMProjectFunctionLocal(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void**,InsertMode,Vec);
372 PETSC_EXTERN PetscErrorCode DMProjectFunctionLabel(DM, PetscReal, DMLabel, PetscInt, const PetscInt[], PetscInt, const PetscInt[], PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **, InsertMode, Vec);
373 PETSC_EXTERN PetscErrorCode DMProjectFunctionLabelLocal(DM,PetscReal,DMLabel,PetscInt,const PetscInt[],PetscInt,const PetscInt[],PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,InsertMode,Vec);
374 PETSC_EXTERN PetscErrorCode DMProjectFieldLocal(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);
375 PETSC_EXTERN PetscErrorCode DMProjectFieldLabelLocal(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);
376 PETSC_EXTERN PetscErrorCode DMProjectBdFieldLabelLocal(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);
377 PETSC_EXTERN PetscErrorCode DMComputeL2Diff(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
378 PETSC_EXTERN PetscErrorCode DMComputeL2GradientDiff(DM, PetscReal, PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], const PetscReal [], PetscInt, PetscScalar *, void *), void **, Vec, const PetscReal [], PetscReal *);
379 PETSC_EXTERN PetscErrorCode DMComputeL2FieldDiff(DM,PetscReal,PetscErrorCode(**)(PetscInt,PetscReal,const PetscReal[],PetscInt,PetscScalar *,void *),void **,Vec,PetscReal *);
380 PETSC_EXTERN PetscErrorCode DMComputeError(DM, Vec, PetscReal[], Vec *);
381 PETSC_EXTERN PetscErrorCode DMHasBasisTransform(DM,PetscBool*);
382 PETSC_EXTERN PetscErrorCode DMCopyTransform(DM, DM);
383 
384 PETSC_EXTERN PetscErrorCode DMGetCompatibility(DM,DM,PetscBool*,PetscBool*);
385 
386 PETSC_EXTERN PetscErrorCode DMMonitorSet(DM, PetscErrorCode (*)(DM, void *), void *, PetscErrorCode (*)(void**));
387 PETSC_EXTERN PetscErrorCode DMMonitorCancel(DM);
388 PETSC_EXTERN PetscErrorCode DMMonitorSetFromOptions(DM, const char[], const char[], const char[], PetscErrorCode (*)(DM, void *), PetscErrorCode (*)(DM, PetscViewerAndFormat *), PetscBool *);
389 PETSC_EXTERN PetscErrorCode DMMonitor(DM);
390 
391 static inline PetscInt DMPolytopeTypeGetDim(DMPolytopeType ct)
392 {
393   switch (ct) {
394     case DM_POLYTOPE_POINT:
395       return 0;
396     case DM_POLYTOPE_SEGMENT:
397     case DM_POLYTOPE_POINT_PRISM_TENSOR:
398       return 1;
399     case DM_POLYTOPE_TRIANGLE:
400     case DM_POLYTOPE_QUADRILATERAL:
401     case DM_POLYTOPE_SEG_PRISM_TENSOR:
402       return 2;
403     case DM_POLYTOPE_TETRAHEDRON:
404     case DM_POLYTOPE_HEXAHEDRON:
405     case DM_POLYTOPE_TRI_PRISM:
406     case DM_POLYTOPE_TRI_PRISM_TENSOR:
407     case DM_POLYTOPE_QUAD_PRISM_TENSOR:
408     case DM_POLYTOPE_PYRAMID:
409       return 3;
410     default: return -1;
411   }
412 }
413 
414 static inline PetscInt DMPolytopeTypeGetConeSize(DMPolytopeType ct)
415 {
416   switch (ct) {
417     case DM_POLYTOPE_POINT:              return 0;
418     case DM_POLYTOPE_SEGMENT:            return 2;
419     case DM_POLYTOPE_POINT_PRISM_TENSOR: return 2;
420     case DM_POLYTOPE_TRIANGLE:           return 3;
421     case DM_POLYTOPE_QUADRILATERAL:      return 4;
422     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return 4;
423     case DM_POLYTOPE_TETRAHEDRON:        return 4;
424     case DM_POLYTOPE_HEXAHEDRON:         return 6;
425     case DM_POLYTOPE_TRI_PRISM:          return 5;
426     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return 5;
427     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return 6;
428     case DM_POLYTOPE_PYRAMID:            return 5;
429     default: return -1;
430   }
431 }
432 
433 static inline PetscInt DMPolytopeTypeGetNumVertices(DMPolytopeType ct)
434 {
435   switch (ct) {
436     case DM_POLYTOPE_POINT:              return 1;
437     case DM_POLYTOPE_SEGMENT:            return 2;
438     case DM_POLYTOPE_POINT_PRISM_TENSOR: return 2;
439     case DM_POLYTOPE_TRIANGLE:           return 3;
440     case DM_POLYTOPE_QUADRILATERAL:      return 4;
441     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return 4;
442     case DM_POLYTOPE_TETRAHEDRON:        return 4;
443     case DM_POLYTOPE_HEXAHEDRON:         return 8;
444     case DM_POLYTOPE_TRI_PRISM:          return 6;
445     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return 6;
446     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return 8;
447     case DM_POLYTOPE_PYRAMID:            return 5;
448     default: return -1;
449   }
450 }
451 
452 static inline DMPolytopeType DMPolytopeTypeSimpleShape(PetscInt dim, PetscBool simplex)
453 {
454   return dim == 0 ? DM_POLYTOPE_POINT :
455         (dim == 1 ? DM_POLYTOPE_SEGMENT :
456         (dim == 2 ? (simplex ? DM_POLYTOPE_TRIANGLE : DM_POLYTOPE_QUADRILATERAL) :
457         (dim == 3 ? (simplex ? DM_POLYTOPE_TETRAHEDRON : DM_POLYTOPE_HEXAHEDRON) : DM_POLYTOPE_UNKNOWN)));
458 }
459 
460 static inline PetscInt DMPolytopeTypeGetNumArrangments(DMPolytopeType ct)
461 {
462   switch (ct) {
463     case DM_POLYTOPE_POINT:              return 1;
464     case DM_POLYTOPE_SEGMENT:            return 2;
465     case DM_POLYTOPE_POINT_PRISM_TENSOR: return 2;
466     case DM_POLYTOPE_TRIANGLE:           return 6;
467     case DM_POLYTOPE_QUADRILATERAL:      return 8;
468     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return 4;
469     case DM_POLYTOPE_TETRAHEDRON:        return 24;
470     case DM_POLYTOPE_HEXAHEDRON:         return 48;
471     case DM_POLYTOPE_TRI_PRISM:          return 12;
472     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return 12;
473     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return 16;
474     case DM_POLYTOPE_PYRAMID:            return 8;
475     default: return -1;
476   }
477 }
478 
479 /* An arrangement is a face order combined with an orientation for each face */
480 static inline const PetscInt *DMPolytopeTypeGetArrangment(DMPolytopeType ct, PetscInt o)
481 {
482   static const PetscInt pntArr[1*2] = {0, 0};
483   /* a: swap */
484   static const PetscInt segArr[2*2*2] = {
485     1, 0,  0, 0, /* -1: a */
486     0, 0,  1, 0, /*  0: e */};
487   /* a: swap first two
488      b: swap last two */
489   static const PetscInt triArr[6*3*2] = {
490     0, -1,  2, -1,  1, -1, /* -3: b */
491     2, -1,  1, -1,  0, -1, /* -2: aba */
492     1, -1,  0, -1,  2, -1, /* -1: a */
493     0,  0,  1,  0,  2,  0, /*  0: identity */
494     1,  0,  2,  0,  0,  0, /*  1: ba */
495     2,  0,  0,  0,  1,  0, /*  2: ab */};
496   /* a: forward cyclic permutation
497      b: swap first and last pairs */
498   static const PetscInt quadArr[8*4*2] = {
499     1, -1,  0, -1,  3, -1,  2, -1, /* -4: b */
500     0, -1,  3, -1,  2, -1,  1, -1, /* -3: b a^3 = a b */
501     3, -1,  2, -1,  1, -1,  0, -1, /* -2: b a^2 = a^2 b */
502     2, -1,  1, -1,  0, -1,  3, -1, /* -1: b a   = a^3 b */
503     0,  0,  1,  0,  2,  0,  3,  0, /*  0: identity */
504     1,  0,  2,  0,  3,  0,  0,  0, /*  1: a */
505     2,  0,  3,  0,  0,  0,  1,  0, /*  2: a^2 */
506     3,  0,  0,  0,  1,  0,  2,  0, /*  3: a^3 */};
507   /* r: rotate 180
508      b: swap top and bottom segments */
509   static const PetscInt tsegArr[4*4*2] = {
510     1, -1,  0, -1,  3, -1,  2, -1, /* -2: r b */
511     0, -1,  1, -1,  3,  0,  2,  0, /* -1: r */
512     0,  0,  1,  0,  2,  0,  3,  0, /*  0: identity */
513     1,  0,  0,  0,  2, -1,  3, -1, /*  1: b */};
514   /* https://en.wikiversity.org/wiki/Symmetric_group_S4 */
515   static const PetscInt tetArr[24*4*2] = {
516     3, -2,  2, -3,  0, -1,  1, -1, /* -12: (1324)   p22 */
517     3, -1,  1, -3,  2, -1,  0, -1, /* -11: (14)     p21 */
518     3, -3,  0, -3,  1, -1,  2, -1, /* -10: (1234)   p18 */
519     2, -1,  3, -1,  1, -3,  0, -2, /*  -9: (1423)   p17 */
520     2, -3,  0, -1,  3, -2,  1, -3, /*  -8: (1342)   p13 */
521     2, -2,  1, -2,  0, -2,  3, -2, /*  -7: (24)     p14 */
522     1, -2,  0, -2,  2, -2,  3, -1, /*  -6: (34)     p6  */
523     1, -1,  3, -3,  0, -3,  2, -2, /*  -5: (1243)   p10 */
524     1, -3,  2, -1,  3, -1,  0, -3, /*  -4: (1432)   p9  */
525     0, -3,  1, -1,  3, -3,  2, -3, /*  -3: (12)     p1  */
526     0, -2,  2, -2,  1, -2,  3, -3, /*  -2: (23)     p2  */
527     0, -1,  3, -2,  2, -3,  1, -2, /*  -1: (13)     p5  */
528     0,  0,  1,  0,  2,  0,  3,  0, /*   0: ()       p0  */
529     0,  1,  3,  1,  1,  2,  2,  0, /*   1: (123)    p4  */
530     0,  2,  2,  1,  3,  0,  1,  2, /*   2: (132)    p3  */
531     1,  2,  0,  1,  3,  1,  2,  2, /*   3: (12)(34) p7  */
532     1,  0,  2,  0,  0,  0,  3,  1, /*   4: (243)    p8  */
533     1,  1,  3,  2,  2,  2,  0,  0, /*   5: (143)    p11 */
534     2,  1,  3,  0,  0,  2,  1,  0, /*   6: (13)(24) p16 */
535     2,  2,  1,  1,  3,  2,  0,  2, /*   7: (142)    p15 */
536     2,  0,  0,  0,  1,  0,  3,  2, /*   8: (234)    p12 */
537     3,  2,  2,  2,  1,  1,  0,  1, /*   9: (14)(23) p23 */
538     3,  0,  0,  2,  2,  1,  1,  1, /*  10: (134)    p19 */
539     3,  1,  1,  2,  0,  1,  2,  1  /*  11: (124)    p20 */};
540   /* Each rotation determines a permutation of the four diagonals, and this defines the isomorphism with S_4 */
541   static const PetscInt hexArr[48*6*2] = {
542     2, -3,  3, -2,  4, -2,  5, -3,  1, -3,  0, -1, /* -24: reflect bottom and use -3 on top */
543     4, -2,  5, -2,  0, -1,  1, -4,  3, -2,  2, -3, /* -23: reflect bottom and use -3 on top */
544     5, -3,  4, -1,  1, -2,  0, -3,  3, -4,  2, -1, /* -22: reflect bottom and use -3 on top */
545     3, -1,  2, -4,  4, -4,  5, -1,  0, -4,  1, -4, /* -21: reflect bottom and use -3 on top */
546     3, -3,  2, -2,  5, -1,  4, -4,  1, -1,  0, -3, /* -20: reflect bottom and use -3 on top */
547     4, -4,  5, -4,  1, -4,  0, -1,  2, -4,  3, -1, /* -19: reflect bottom and use -3 on top */
548     2, -1,  3, -4,  5, -3,  4, -2,  0, -2,  1, -2, /* -18: reflect bottom and use -3 on top */
549     5, -1,  4, -3,  0, -3,  1, -2,  2, -2,  3, -3, /* -17: reflect bottom and use -3 on top */
550     4, -3,  5, -1,  3, -2,  2, -4,  1, -4,  0, -4, /* -16: reflect bottom and use -3 on top */
551     5, -4,  4, -4,  3, -4,  2, -2,  0, -3,  1, -1, /* -15: reflect bottom and use -3 on top */
552     3, -4,  2, -1,  1, -1,  0, -4,  4, -4,  5, -4, /* -14: reflect bottom and use -3 on top */
553     2, -2,  3, -3,  0, -2,  1, -3,  4, -2,  5, -2, /* -13: reflect bottom and use -3 on top */
554     1, -3,  0, -1,  4, -1,  5, -4,  3, -1,  2, -4, /* -12: reflect bottom and use -3 on top */
555     1, -1,  0, -3,  5, -4,  4, -1,  2, -1,  3, -4, /* -11: reflect bottom and use -3 on top */
556     5, -2,  4, -2,  2, -2,  3, -4,  1, -2,  0, -2, /* -10: reflect bottom and use -3 on top */
557     1, -2,  0, -2,  2, -1,  3, -1,  4, -1,  5, -3, /*  -9: reflect bottom and use -3 on top */
558     4, -1,  5, -3,  2, -4,  3, -2,  0, -1,  1, -3, /*  -8: reflect bottom and use -3 on top */
559     3, -2,  2, -3,  0, -4,  1, -1,  5, -1,  4, -3, /*  -7: reflect bottom and use -3 on top */
560     1, -4,  0, -4,  3, -1,  2, -1,  5, -4,  4, -4, /*  -6: reflect bottom and use -3 on top */
561     2, -4,  3, -1,  1, -3,  0, -2,  5, -3,  4, -1, /*  -5: reflect bottom and use -3 on top */
562     0, -4,  1, -4,  4, -3,  5, -2,  2, -3,  3, -2, /*  -4: reflect bottom and use -3 on top */
563     0, -3,  1, -1,  3, -3,  2, -3,  4, -3,  5, -1, /*  -3: reflect bottom and use -3 on top */
564     0, -2,  1, -2,  5, -2,  4, -3,  3, -3,  2, -2, /*  -2: reflect bottom and use -3 on top */
565     0, -1,  1, -3,  2, -3,  3, -3,  5, -2,  4, -2, /*  -1: reflect bottom and use -3 on top */
566     0,  0,  1,  0,  2,  0,  3,  0,  4,  0,  5,  0, /*   0: identity */
567     0,  1,  1,  3,  5,  3,  4,  0,  2,  0,  3,  1, /*   1: 90  rotation about z */
568     0,  2,  1,  2,  3,  0,  2,  0,  5,  3,  4,  1, /*   2: 180 rotation about z */
569     0,  3,  1,  1,  4,  0,  5,  3,  3,  0,  2,  1, /*   3: 270 rotation about z */
570     2,  3,  3,  2,  1,  0,  0,  3,  4,  3,  5,  1, /*   4: 90  rotation about x */
571     1,  3,  0,  1,  3,  2,  2,  2,  4,  2,  5,  2, /*   5: 180 rotation about x */
572     3,  1,  2,  0,  0,  1,  1,  2,  4,  1,  5,  3, /*   6: 270 rotation about x */
573     4,  0,  5,  0,  2,  1,  3,  3,  1,  1,  0,  3, /*   7: 90  rotation about y */
574     1,  1,  0,  3,  2,  2,  3,  2,  5,  1,  4,  3, /*   8: 180 rotation about y */
575     5,  1,  4,  3,  2,  3,  3,  1,  0,  0,  1,  0, /*   9: 270 rotation about y */
576     1,  0,  0,  0,  5,  1,  4,  2,  3,  2,  2,  3, /*  10: 180 rotation about x+y */
577     1,  2,  0,  2,  4,  2,  5,  1,  2,  2,  3,  3, /*  11: 180 rotation about x-y */
578     2,  1,  3,  0,  0,  3,  1,  0,  5,  0,  4,  0, /*  12: 180 rotation about y+z */
579     3,  3,  2,  2,  1,  2,  0,  1,  5,  2,  4,  2, /*  13: 180 rotation about y-z */
580     5,  3,  4,  1,  3,  1,  2,  3,  1,  3,  0,  1, /*  14: 180 rotation about z+x */
581     4,  2,  5,  2,  3,  3,  2,  1,  0,  2,  1,  2, /*  15: 180 rotation about z-x */
582     5,  0,  4,  0,  0,  0,  1,  3,  3,  1,  2,  0, /*  16: 120 rotation about x+y+z (v0v6) */
583     2,  0,  3,  1,  5,  0,  4,  3,  1,  0,  0,  0, /*  17: 240 rotation about x+y+z (v0v6) */
584     4,  3,  5,  1,  1,  1,  0,  2,  3,  3,  2,  2, /*  18: 120 rotation about x+y-z (v4v2) */
585     3,  2,  2,  3,  5,  2,  4,  1,  0,  1,  1,  3, /*  19: 240 rotation about x+y-z (v4v2) */
586     3,  0,  2,  1,  4,  1,  5,  2,  1,  2,  0,  2, /*  20: 120 rotation about x-y+z (v1v5) */
587     5,  2,  4,  2,  1,  3,  0,  0,  2,  3,  3,  2, /*  21: 240 rotation about x-y+z (v1v5) */
588     4,  1,  5,  3,  0,  2,  1,  1,  2,  1,  3,  0, /*  22: 120 rotation about x-y-z (v7v3) */
589     2,  2,  3,  3,  4,  3,  5,  0,  0,  3,  1,  1, /*  23: 240 rotation about x-y-z (v7v3) */
590   };
591   static const PetscInt tripArr[12*5*2] = {
592     1, -3,  0, -1,  3, -1,  4, -1,  2, -1, /* -6: reflect bottom and top */
593     1, -1,  0, -3,  4, -1,  2, -1,  3, -1, /* -5: reflect bottom and top */
594     1, -2,  0, -2,  2, -1,  3, -1,  4, -1, /* -4: reflect bottom and top */
595     0, -3,  1, -1,  3, -3,  2, -3,  4, -3, /* -3: reflect bottom and top */
596     0, -2,  1, -2,  4, -3,  3, -3,  2, -3, /* -2: reflect bottom and top */
597     0, -1,  1, -3,  2, -3,  4, -3,  3, -3, /* -1: reflect bottom and top */
598     0,  0,  1,  0,  2,  0,  3,  0,  4,  0, /*  0: identity */
599     0,  1,  1,  2,  4,  0,  2,  0,  3,  0, /*  1: 120 rotation about z */
600     0,  2,  1,  1,  3,  0,  4,  0,  2,  0, /*  2: 240 rotation about z */
601     1,  1,  0,  2,  2,  2,  4,  2,  3,  2, /*  3: 180 rotation about y of 0 */
602     1,  0,  0,  0,  4,  2,  3,  2,  2,  2, /*  4: 180 rotation about y of 1 */
603     1,  2,  0,  1,  3,  2,  2,  2,  4,  2, /*  5: 180 rotation about y of 2 */
604   };
605   /* a: rotate 120 about z
606      b: swap top and bottom segments
607      r: reflect */
608   static const PetscInt ttriArr[12*5*2] = {
609     1, -3,  0, -3,  2, -2,  4, -2,  3, -2, /* -6: r b a^2 */
610     1, -2,  0, -2,  4, -2,  3, -2,  2, -2, /* -5: r b a */
611     1, -1,  0, -1,  3, -2,  2, -2,  4, -2, /* -4: r b */
612     0, -3,  1, -3,  2, -1,  4, -1,  3, -1, /* -3: r a^2 */
613     0, -2,  1, -2,  4, -1,  3, -1,  2, -1, /* -2: r a */
614     0, -1,  1, -1,  3, -1,  2, -1,  4, -1, /* -1: r */
615     0,  0,  1,  0,  2,  0,  3,  0,  4,  0, /*  0: identity */
616     0,  1,  1,  1,  3,  0,  4,  0,  2,  0, /*  1: a */
617     0,  2,  1,  2,  4,  0,  2,  0,  3,  0, /*  2: a^2 */
618     1,  0,  0,  0,  2,  1,  3,  1,  4,  1, /*  3: b */
619     1,  1,  0,  1,  3,  1,  4,  1,  2,  1, /*  4: b a */
620     1,  2,  0,  2,  4,  1,  2,  1,  3,  1, /*  5: b a^2 */
621   };
622   /* a: rotate 90 about z
623      b: swap top and bottom segments
624      r: reflect */
625   static const PetscInt tquadArr[16*6*2] = {
626     1, -4,  0, -4,  3, -2,  2, -2,  5, -2,  4, -2, /* -8: r b a^3 */
627     1, -3,  0, -3,  2, -2,  5, -2,  4, -2,  3, -2, /* -7: r b a^2 */
628     1, -2,  0, -2,  5, -2,  4, -2,  3, -2,  2, -2, /* -6: r b a */
629     1, -1,  0, -1,  4, -2,  3, -2,  2, -2,  5, -2, /* -5: r b */
630     0, -4,  1, -4,  3, -1,  2, -1,  5, -1,  4, -1, /* -4: r a^3 */
631     0, -3,  1, -3,  2, -1,  5, -1,  4, -1,  3, -1, /* -3: r a^2 */
632     0, -2,  1, -2,  5, -1,  4, -1,  3, -1,  2, -1, /* -2: r a */
633     0, -1,  1, -1,  4, -1,  3, -1,  2, -1,  5, -1, /* -1: r */
634     0,  0,  1,  0,  2,  0,  3,  0,  4,  0,  5,  0, /*  0: identity */
635     0,  1,  1,  1,  3,  0,  4,  0,  5,  0,  2,  0, /*  1: a */
636     0,  2,  1,  2,  4,  0,  5,  0,  2,  0,  3,  0, /*  2: a^2 */
637     0,  3,  1,  3,  5,  0,  2,  0,  3,  0,  4,  0, /*  3: a^3 */
638     1,  0,  0,  0,  2,  1,  3,  1,  4,  1,  5,  1, /*  4: b */
639     1,  1,  0,  1,  3,  1,  4,  1,  5,  1,  2,  1, /*  5: b a */
640     1,  2,  0,  2,  4,  1,  5,  1,  2,  1,  3,  1, /*  6: b a^2 */
641     1,  3,  0,  3,  5,  1,  2,  1,  3,  1,  4,  1, /*  7: b a^3 */
642   };
643   static const PetscInt pyrArr[8*5*2] = {
644     0, -4,  2, -3,  1, -3,  4, -3,  3, -3, /* -4: Reflect bottom face */
645     0, -3,  3, -3,  2, -3,  1, -3,  4, -3, /* -3: Reflect bottom face */
646     0, -2,  4, -3,  3, -3,  2, -3,  1, -3, /* -2: Reflect bottom face */
647     0, -1,  1, -3,  4, -3,  3, -3,  2, -3, /* -1: Reflect bottom face */
648     0,  0,  1,  0,  2,  0,  3,  0,  4,  0, /*  0: identity */
649     0,  1,  4,  0,  1,  0,  2,  0,  3,  0, /*  1:  90 rotation about z */
650     0,  2,  3,  0,  4,  0,  1,  0,  2,  0, /*  2: 180 rotation about z */
651     0,  3,  2,  0,  3,  0,  4,  0,  1,  0, /*  3: 270 rotation about z */
652   };
653   switch (ct) {
654     case DM_POLYTOPE_POINT:              return pntArr;
655     case DM_POLYTOPE_SEGMENT:            return &segArr[(o+1)*2*2];
656     case DM_POLYTOPE_POINT_PRISM_TENSOR: return &segArr[(o+1)*2*2];
657     case DM_POLYTOPE_TRIANGLE:           return &triArr[(o+3)*3*2];
658     case DM_POLYTOPE_QUADRILATERAL:      return &quadArr[(o+4)*4*2];
659     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return &tsegArr[(o+2)*4*2];
660     case DM_POLYTOPE_TETRAHEDRON:        return &tetArr[(o+12)*4*2];
661     case DM_POLYTOPE_HEXAHEDRON:         return &hexArr[(o+24)*6*2];
662     case DM_POLYTOPE_TRI_PRISM:          return &tripArr[(o+6)*5*2];
663     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return &ttriArr[(o+6)*5*2];
664     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return &tquadArr[(o+8)*6*2];
665     case DM_POLYTOPE_PYRAMID:            return &pyrArr[(o+4)*5*2];
666     default: return NULL;
667   }
668 }
669 
670 /* A vertex arrangment is a vertex order */
671 static inline const PetscInt *DMPolytopeTypeGetVertexArrangment(DMPolytopeType ct, PetscInt o)
672 {
673   static const PetscInt pntVerts[1]    = {0};
674   static const PetscInt segVerts[2*2]  = {
675     1, 0,
676     0, 1};
677   static const PetscInt triVerts[6*3]  = {
678     1, 0, 2,
679     0, 2, 1,
680     2, 1, 0,
681     0, 1, 2,
682     1, 2, 0,
683     2, 0, 1};
684   static const PetscInt quadVerts[8*4]  = {
685     2, 1, 0, 3,
686     1, 0, 3, 2,
687     0, 3, 2, 1,
688     3, 2, 1, 0,
689     0, 1, 2, 3,
690     1, 2, 3, 0,
691     2, 3, 0, 1,
692     3, 0, 1, 2};
693   static const PetscInt tsegVerts[4*4]  = {
694     3, 2, 1, 0,
695     1, 0, 3, 2,
696     0, 1, 2, 3,
697     2, 3, 0, 1};
698   static const PetscInt tetVerts[24*4] = {
699     2, 3, 1, 0, /* -12: (1324)   p22 */
700     3, 1, 2, 0, /* -11: (14)     p21 */
701     1, 2, 3, 0, /* -10: (1234)   p18 */
702     3, 2, 0, 1, /*  -9: (1423)   p17 */
703     2, 0, 3, 1, /*  -8: (1342)   p13 */
704     0, 3, 2, 1, /*  -7: (24)     p14 */
705     0, 1, 3, 2, /*  -6: (34)     p6  */
706     1, 3, 0, 2, /*  -5: (1243)   p10 */
707     3, 0, 1, 2, /*  -4: (1432    p9  */
708     1, 0, 2, 3, /*  -3: (12)     p1  */
709     0, 2, 1, 3, /*  -2: (23)     p2  */
710     2, 1, 0, 3, /*  -1: (13)     p5  */
711     0, 1, 2, 3, /*   0: ()       p0  */
712     1, 2, 0, 3, /*   1: (123)    p4  */
713     2, 0, 1, 3, /*   2: (132)    p3  */
714     1, 0, 3, 2, /*   3: (12)(34) p7  */
715     0, 3, 1, 2, /*   4: (243)    p8  */
716     3, 1, 0, 2, /*   5: (143)    p11 */
717     2, 3, 0, 1, /*   6: (13)(24) p16 */
718     3, 0, 2, 1, /*   7: (142)    p15 */
719     0, 2, 3, 1, /*   8: (234)    p12 */
720     3, 2, 1, 0, /*   9: (14)(23) p23 */
721     2, 1, 3, 0, /*  10: (134)    p19 */
722     1, 3, 2, 0  /*  11: (124)    p20 */};
723   static const PetscInt hexVerts[48*8] = {
724     3,  0,  4,  5,  2,  6,  7,  1, /* -24: reflected 23 */
725     3,  5,  6,  2,  0,  1,  7,  4, /* -23: reflected 22 */
726     4,  0,  1,  7,  5,  6,  2,  3, /* -22: reflected 21 */
727     6,  7,  1,  2,  5,  3,  0,  4, /* -21: reflected 20 */
728     1,  2,  6,  7,  0,  4,  5,  3, /* -20: reflected 19 */
729     6,  2,  3,  5,  7,  4,  0,  1, /* -19: reflected 18 */
730     4,  5,  3,  0,  7,  1,  2,  6, /* -18: reflected 17 */
731     1,  7,  4,  0,  2,  3,  5,  6, /* -17: reflected 16 */
732     2,  3,  5,  6,  1,  7,  4,  0, /* -16: reflected 15 */
733     7,  4,  0,  1,  6,  2,  3,  5, /* -15: reflected 14 */
734     7,  1,  2,  6,  4,  5,  3,  0, /* -14: reflected 13 */
735     0,  4,  5,  3,  1,  2,  6,  7, /* -13: reflected 12 */
736     5,  4,  7,  6,  3,  2,  1,  0, /* -12: reflected 11 */
737     7,  6,  5,  4,  1,  0,  3,  2, /* -11: reflected 10 */
738     0,  1,  7,  4,  3,  5,  6,  2, /* -10: reflected  9 */
739     4,  7,  6,  5,  0,  3,  2,  1, /*  -9: reflected  8 */
740     5,  6,  2,  3,  4,  0,  1,  7, /*  -8: reflected  7 */
741     2,  6,  7,  1,  3,  0,  4,  5, /*  -7: reflected  6 */
742     6,  5,  4,  7,  2,  1,  0,  3, /*  -6: reflected  5 */
743     5,  3,  0,  4,  6,  7,  1,  2, /*  -5: reflected  4 */
744     2,  1,  0,  3,  6,  5,  4,  7, /*  -4: reflected  3 */
745     1,  0,  3,  2,  7,  6,  5,  4, /*  -3: reflected  2 */
746     0,  3,  2,  1,  4,  7,  6,  5, /*  -2: reflected  1 */
747     3,  2,  1,  0,  5,  4,  7,  6, /*  -1: reflected  0 */
748     0,  1,  2,  3,  4,  5,  6,  7, /*   0: identity */
749     1,  2,  3,  0,  7,  4,  5,  6, /*   1: 90  rotation about z */
750     2,  3,  0,  1,  6,  7,  4,  5, /*   2: 180 rotation about z */
751     3,  0,  1,  2,  5,  6,  7,  4, /*   3: 270 rotation about z */
752     4,  0,  3,  5,  7,  6,  2,  1, /*   4: 90  rotation about x */
753     7,  4,  5,  6,  1,  2,  3,  0, /*   5: 180 rotation about x */
754     1,  7,  6,  2,  0,  3,  5,  4, /*   6: 270 rotation about x */
755     3,  2,  6,  5,  0,  4,  7,  1, /*   7: 90  rotation about y */
756     5,  6,  7,  4,  3,  0,  1,  2, /*   8: 180 rotation about y */
757     4,  7,  1,  0,  5,  3,  2,  6, /*   9: 270 rotation about y */
758     4,  5,  6,  7,  0,  1,  2,  3, /*  10: 180 rotation about x+y */
759     6,  7,  4,  5,  2,  3,  0,  1, /*  11: 180 rotation about x-y */
760     3,  5,  4,  0,  2,  1,  7,  6, /*  12: 180 rotation about y+z */
761     6,  2,  1,  7,  5,  4,  0,  3, /*  13: 180 rotation about y-z */
762     1,  0,  4,  7,  2,  6,  5,  3, /*  14: 180 rotation about z+x */
763     6,  5,  3,  2,  7,  1,  0,  4, /*  15: 180 rotation about z-x */
764     0,  4,  7,  1,  3,  2,  6,  5, /*  16: 120 rotation about x+y+z (v0v6) */
765     0,  3,  5,  4,  1,  7,  6,  2, /*  17: 240 rotation about x+y+z (v0v6) */
766     5,  3,  2,  6,  4,  7,  1,  0, /*  18: 120 rotation about x+y-z (v4v2) */
767     7,  6,  2,  1,  4,  0,  3,  5, /*  19: 240 rotation about x+y-z (v4v2) */
768     2,  1,  7,  6,  3,  5,  4,  0, /*  20: 120 rotation about x-y+z (v1v5) */
769     7,  1,  0,  4,  6,  5,  3,  2, /*  21: 240 rotation about x-y+z (v1v5) */
770     2,  6,  5,  3,  1,  0,  4,  7, /*  22: 120 rotation about x-y-z (v7v3) */
771     5,  4,  0,  3,  6,  2,  1,  7, /*  23: 240 rotation about x-y-z (v7v3) */
772   };
773   static const PetscInt tripVerts[12*6] = {
774     4,  3,  5,  2,  1,  0, /* -6: reflect bottom and top */
775     5,  4,  3,  1,  0,  2, /* -5: reflect bottom and top */
776     3,  5,  4,  0,  2,  1, /* -4: reflect bottom and top */
777     1,  0,  2,  5,  4,  3, /* -3: reflect bottom and top */
778     0,  2,  1,  3,  5,  4, /* -2: reflect bottom and top */
779     2,  1,  0,  4,  3,  5, /* -1: reflect bottom and top */
780     0,  1,  2,  3,  4,  5, /*  0: identity */
781     1,  2,  0,  5,  3,  4, /*  1: 120 rotation about z */
782     2,  0,  1,  4,  5,  3, /*  2: 240 rotation about z */
783     4,  5,  3,  2,  0,  1, /*  3: 180 rotation about y of 0 */
784     3,  4,  5,  0,  1,  2, /*  4: 180 rotation about y of 1 */
785     5,  3,  4,  1,  2,  0, /*  5: 180 rotation about y of 2 */
786   };
787   static const PetscInt ttriVerts[12*6] = {
788     4,  3,  5,  1,  0,  2, /* -6: r b a^2 */
789     3,  5,  4,  0,  2,  1, /* -5: r b a */
790     5,  4,  3,  2,  1,  0, /* -4: r b */
791     1,  0,  2,  4,  3,  5, /* -3: r a^2 */
792     0,  2,  1,  3,  5,  4, /* -2: r a */
793     2,  1,  0,  5,  4,  3, /* -1: r */
794     0,  1,  2,  3,  4,  5, /*  0: identity */
795     1,  2,  0,  4,  5,  3, /*  1: a */
796     2,  0,  1,  5,  3,  4, /*  2: a^2 */
797     3,  4,  5,  0,  1,  2, /*  3: b */
798     4,  5,  3,  1,  2,  0, /*  4: b a */
799     5,  3,  4,  2,  0,  1, /*  5: b a^2 */
800   };
801   /* a: rotate 90 about z
802      b: swap top and bottom segments
803      r: reflect */
804   static const PetscInt tquadVerts[16*8] = {
805     6,  5,  4,  7,  2,  1,  0,  3, /* -8: r b a^3 */
806     5,  4,  7,  6,  1,  0,  3,  2, /* -7: r b a^2 */
807     4,  7,  6,  5,  0,  3,  2,  1, /* -6: r b a */
808     7,  6,  5,  4,  3,  2,  1,  0, /* -5: r b */
809     2,  1,  0,  3,  6,  5,  4,  7, /* -4: r a^3 */
810     1,  0,  3,  2,  5,  4,  7,  6, /* -3: r a^2 */
811     0,  3,  2,  1,  4,  7,  6,  5, /* -2: r a */
812     3,  2,  1,  0,  7,  6,  5,  4, /* -1: r */
813     0,  1,  2,  3,  4,  5,  6,  7, /*  0: identity */
814     1,  2,  3,  0,  5,  6,  7,  4, /*  1: a */
815     2,  3,  0,  1,  6,  7,  4,  5, /*  2: a^2 */
816     3,  0,  1,  2,  7,  4,  5,  6, /*  3: a^3 */
817     4,  5,  6,  7,  0,  1,  2,  3, /*  4: b */
818     5,  6,  7,  4,  1,  2,  3,  0, /*  5: b a */
819     6,  7,  4,  5,  2,  3,  0,  1, /*  6: b a^2 */
820     7,  4,  5,  6,  3,  0,  1,  2, /*  7: b a^3 */
821   };
822   static const PetscInt pyrVerts[8*5] = {
823     2,  1,  0,  3,  4, /* -4: Reflect bottom face */
824     1,  0,  3,  2,  4, /* -3: Reflect bottom face */
825     0,  3,  2,  1,  4, /* -2: Reflect bottom face */
826     3,  2,  1,  0,  4, /* -1: Reflect bottom face */
827     0,  1,  2,  3,  4, /*  0: identity */
828     1,  2,  3,  0,  4, /*  1:  90 rotation about z */
829     2,  3,  0,  1,  4, /*  2: 180 rotation about z */
830     3,  0,  1,  2,  4, /*  3: 270 rotation about z */
831   };
832   switch (ct) {
833     case DM_POLYTOPE_POINT:              return pntVerts;
834     case DM_POLYTOPE_SEGMENT:            return &segVerts[(o+1)*2];
835     case DM_POLYTOPE_POINT_PRISM_TENSOR: return &segVerts[(o+1)*2];
836     case DM_POLYTOPE_TRIANGLE:           return &triVerts[(o+3)*3];
837     case DM_POLYTOPE_QUADRILATERAL:      return &quadVerts[(o+4)*4];
838     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return &tsegVerts[(o+2)*4];
839     case DM_POLYTOPE_TETRAHEDRON:        return &tetVerts[(o+12)*4];
840     case DM_POLYTOPE_HEXAHEDRON:         return &hexVerts[(o+24)*8];
841     case DM_POLYTOPE_TRI_PRISM:          return &tripVerts[(o+6)*6];
842     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return &ttriVerts[(o+6)*6];
843     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return &tquadVerts[(o+8)*8];
844     case DM_POLYTOPE_PYRAMID:            return &pyrVerts[(o+4)*5];
845     default: return NULL;
846   }
847 }
848 
849 /* This is orientation o1 acting on orientation o2 */
850 static inline PetscInt DMPolytopeTypeComposeOrientation(DMPolytopeType ct, PetscInt o1, PetscInt o2)
851 {
852   static const PetscInt segMult[2*2] = {
853      0, -1,
854     -1,  0};
855   static const PetscInt triMult[6*6] = {
856      0,  2,  1, -3, -1, -2,
857      1,  0,  2, -2, -3, -1,
858      2,  1,  0, -1, -2, -3,
859     -3, -2, -1,  0,  1,  2,
860     -2, -1, -3,  1,  2,  0,
861     -1, -3, -2,  2,  0,  1};
862   static const PetscInt quadMult[8*8] = {
863      0,  3,  2,  1, -4, -1, -2, -3,
864      1,  0,  3,  2, -3, -4, -1, -2,
865      2,  1,  0,  3, -2, -3, -4, -1,
866      3,  2,  1,  0, -1, -2, -3, -4,
867     -4, -3, -2, -1,  0,  1,  2,  3,
868     -3, -2, -1, -4,  1,  2,  3,  0,
869     -2, -1, -4, -3,  2,  3,  0,  1,
870     -1, -4, -3, -2,  3,  0,  1,  2};
871   static const PetscInt tsegMult[4*4] = {
872      0,  1, -2, -1,
873      1,  0, -1, -2,
874     -2, -1,  0,  1,
875     -1, -2,  1,  0};
876   static const PetscInt tetMult[24*24] = {
877     3, 2, 7, 0, 5, 10, 9, 8, 1, 6, 11, 4, -12, -7, -5, -9, -10, -2, -6, -1, -11, -3, -4, -8,
878     4, 0, 8, 1, 3, 11, 10, 6, 2, 7, 9, 5, -11, -9, -4, -8, -12, -1, -5, -3, -10, -2, -6, -7,
879     5, 1, 6, 2, 4, 9, 11, 7, 0, 8, 10, 3, -10, -8, -6, -7, -11, -3, -4, -2, -12, -1, -5, -9,
880     0, 8, 4, 3, 11, 1, 6, 2, 10, 9, 5, 7, -9, -4, -11, -12, -1, -8, -3, -10, -5, -6, -7, -2,
881     1, 6, 5, 4, 9, 2, 7, 0, 11, 10, 3, 8, -8, -6, -10, -11, -3, -7, -2, -12, -4, -5, -9, -1,
882     2, 7, 3, 5, 10, 0, 8, 1, 9, 11, 4, 6, -7, -5, -12, -10, -2, -9, -1, -11, -6, -4, -8, -3,
883     6, 5, 1, 9, 2, 4, 0, 11, 7, 3, 8, 10, -6, -10, -8, -3, -7, -11, -12, -4, -2, -9, -1, -5,
884     7, 3, 2, 10, 0, 5, 1, 9, 8, 4, 6, 11, -5, -12, -7, -2, -9, -10, -11, -6, -1, -8, -3, -4,
885     8, 4, 0, 11, 1, 3, 2, 10, 6, 5, 7, 9, -4, -11, -9, -1, -8, -12, -10, -5, -3, -7, -2, -6,
886     9, 11, 10, 6, 8, 7, 3, 5, 4, 0, 2, 1, -3, -1, -2, -6, -4, -5, -9, -7, -8, -12, -10, -11,
887     10, 9, 11, 7, 6, 8, 4, 3, 5, 1, 0, 2, -2, -3, -1, -5, -6, -4, -8, -9, -7, -11, -12, -10,
888     11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5, -6, -7, -8, -9, -10, -11, -12,
889     -12, -11, -10, -9, -8, -7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11,
890     -11, -10, -12, -8, -7, -9, -5, -4, -6, -2, -1, -3, 1, 2, 0, 4, 5, 3, 7, 8, 6, 10, 11, 9,
891     -10, -12, -11, -7, -9, -8, -4, -6, -5, -1, -3, -2, 2, 0, 1, 5, 3, 4, 8, 6, 7, 11, 9, 10,
892     -9, -5, -1, -12, -2, -4, -3, -11, -7, -6, -8, -10, 3, 10, 8, 0, 7, 11, 9, 4, 2, 6, 1, 5,
893     -8, -4, -3, -11, -1, -6, -2, -10, -9, -5, -7, -12, 4, 11, 6, 1, 8, 9, 10, 5, 0, 7, 2, 3,
894     -7, -6, -2, -10, -3, -5, -1, -12, -8, -4, -9, -11, 5, 9, 7, 2, 6, 10, 11, 3, 1, 8, 0, 4,
895     -3, -8, -4, -6, -11, -1, -9, -2, -10, -12, -5, -7, 6, 4, 11, 9, 1, 8, 0, 10, 5, 3, 7, 2,
896     -2, -7, -6, -5, -10, -3, -8, -1, -12, -11, -4, -9, 7, 5, 9, 10, 2, 6, 1, 11, 3, 4, 8, 0,
897     -1, -9, -5, -4, -12, -2, -7, -3, -11, -10, -6, -8, 8, 3, 10, 11, 0, 7, 2, 9, 4, 5, 6, 1,
898     -6, -2, -7, -3, -5, -10, -12, -8, -1, -9, -11, -4, 9, 7, 5, 6, 10, 2, 3, 1, 11, 0, 4, 8,
899     -5, -1, -9, -2, -4, -12, -11, -7, -3, -8, -10, -6, 10, 8, 3, 7, 11, 0, 4, 2, 9, 1, 5, 6,
900     -4, -3, -8, -1, -6, -11, -10, -9, -2, -7, -12, -5, 11, 6, 4, 8, 9, 1, 5, 0, 10, 2, 3, 7,
901     };
902   static const PetscInt hexMult[48*48] = {
903     18, 2, 5, 22, 21, 8, 16, 0, 13, 6, 11, 3, 15, 9, 4, 23, 12, 1, 19, 10, 7, 20, 14, 17, -24, -10, -20, -16, -12, -21, -4, -5, -18, -13, -15, -8, -2, -11, -14, -7, -3, -22, -6, -17, -19, -9, -1, -23,
904     8, 20, 19, 2, 5, 23, 0, 17, 11, 1, 15, 7, 13, 4, 10, 18, 3, 14, 21, 9, 12, 22, 6, 16, -23, -13, -17, -7, -8, -19, -16, -12, -22, -2, -14, -5, -10, -15, -11, -4, -20, -9, -21, -3, -6, -18, -24, -1,
905     2, 17, 23, 8, 0, 19, 5, 20, 1, 11, 9, 14, 12, 6, 3, 16, 10, 7, 22, 15, 13, 21, 4, 18, -22, -14, -19, -5, -15, -17, -10, -2, -23, -12, -13, -7, -16, -8, -4, -11, -24, -3, -18, -9, -1, -21, -20, -6,
906     21, 5, 2, 16, 18, 0, 22, 8, 4, 12, 3, 11, 14, 7, 13, 20, 6, 10, 17, 1, 9, 23, 15, 19, -21, -8, -18, -15, -4, -24, -12, -14, -20, -7, -16, -10, -11, -2, -5, -13, -6, -19, -3, -23, -22, -1, -9, -17,
907     16, 8, 0, 21, 22, 2, 18, 5, 12, 4, 1, 10, 9, 15, 6, 19, 13, 11, 23, 3, 14, 17, 7, 20, -20, -16, -24, -10, -2, -18, -11, -7, -21, -14, -8, -15, -12, -4, -13, -5, -9, -23, -1, -19, -17, -3, -6, -22,
908     5, 19, 20, 0, 8, 17, 2, 23, 10, 3, 7, 15, 6, 12, 11, 22, 1, 9, 16, 14, 4, 18, 13, 21, -19, -5, -22, -14, -16, -23, -8, -11, -17, -4, -7, -13, -15, -10, -12, -2, -21, -6, -20, -1, -9, -24, -18, -3,
909     22, 0, 8, 18, 16, 5, 21, 2, 6, 13, 10, 1, 7, 14, 12, 17, 4, 3, 20, 11, 15, 19, 9, 23, -18, -15, -21, -8, -11, -20, -2, -13, -24, -5, -10, -16, -4, -12, -7, -14, -1, -17, -9, -22, -23, -6, -3, -19,
910     0, 23, 17, 5, 2, 20, 8, 19, 3, 10, 14, 9, 4, 13, 1, 21, 11, 15, 18, 7, 6, 16, 12, 22, -17, -7, -23, -13, -10, -22, -15, -4, -19, -11, -5, -14, -8, -16, -2, -12, -18, -1, -24, -6, -3, -20, -21, -9,
911     10, 13, 6, 1, 11, 12, 3, 4, 8, 0, 22, 18, 19, 23, 5, 15, 2, 21, 9, 16, 17, 7, 20, 14, -16, -24, -10, -20, -23, -8, -19, -6, -15, -3, -21, -18, -22, -17, -9, -1, -14, -12, -7, -4, -11, -13, -5, -2,
912     1, 4, 12, 10, 3, 6, 11, 13, 0, 8, 16, 21, 17, 20, 2, 14, 5, 18, 7, 22, 19, 9, 23, 15, -15, -21, -8, -18, -17, -10, -22, -3, -16, -6, -24, -20, -19, -23, -1, -9, -5, -4, -13, -12, -2, -7, -14, -11,
913     14, 10, 3, 9, 7, 1, 15, 11, 17, 23, 0, 5, 16, 22, 20, 6, 19, 8, 12, 2, 21, 4, 18, 13, -14, -19, -5, -22, -3, -13, -9, -20, -7, -21, -23, -17, -6, -1, -24, -18, -12, -16, -2, -8, -10, -4, -11, -15,
914     7, 3, 10, 15, 14, 11, 9, 1, 20, 19, 5, 0, 18, 21, 17, 4, 23, 2, 13, 8, 22, 6, 16, 12, -13, -17, -7, -23, -9, -14, -3, -24, -5, -18, -22, -19, -1, -6, -20, -21, -2, -10, -12, -15, -16, -11, -4, -8,
915     13, 14, 15, 12, 4, 9, 6, 7, 21, 22, 23, 20, 2, 0, 18, 3, 16, 17, 1, 19, 8, 11, 5, 10, -12, -9, -11, -6, -21, -4, -24, -22, -2, -23, -3, -1, -20, -18, -19, -17, -16, -14, -15, -13, -5, -8, -10, -7,
916     6, 9, 7, 4, 12, 14, 13, 15, 16, 18, 17, 19, 0, 2, 22, 1, 21, 23, 3, 20, 5, 10, 8, 11, -11, -6, -12, -9, -20, -2, -18, -17, -4, -19, -1, -3, -21, -24, -23, -22, -8, -7, -10, -5, -13, -16, -15, -14,
917     3, 12, 4, 11, 1, 13, 10, 6, 2, 5, 21, 16, 23, 19, 0, 9, 8, 22, 15, 18, 20, 14, 17, 7, -10, -20, -16, -24, -22, -15, -17, -1, -8, -9, -18, -21, -23, -19, -3, -6, -13, -2, -5, -11, -4, -14, -7, -12,
918     20, 16, 18, 23, 17, 21, 19, 22, 14, 15, 4, 6, 3, 1, 7, 0, 9, 12, 2, 13, 11, 5, 10, 8, -9, -11, -6, -12, -14, -3, -13, -10, -1, -8, -2, -4, -7, -5, -16, -15, -23, -20, -22, -18, -24, -19, -17, -21,
919     11, 6, 13, 3, 10, 4, 1, 12, 5, 2, 18, 22, 20, 17, 8, 7, 0, 16, 14, 21, 23, 15, 19, 9, -8, -18, -15, -21, -19, -16, -23, -9, -10, -1, -20, -24, -17, -22, -6, -3, -7, -11, -14, -2, -12, -5, -13, -4,
920     9, 11, 1, 14, 15, 3, 7, 10, 23, 17, 2, 8, 21, 18, 19, 13, 20, 5, 4, 0, 16, 12, 22, 6, -7, -23, -13, -17, -1, -5, -6, -21, -14, -20, -19, -22, -9, -3, -18, -24, -11, -8, -4, -16, -15, -2, -12, -10,
921     19, 21, 22, 17, 23, 16, 20, 18, 9, 7, 12, 13, 1, 3, 15, 2, 14, 4, 0, 6, 10, 8, 11, 5, -6, -12, -9, -11, -7, -1, -5, -15, -3, -16, -4, -2, -14, -13, -8, -10, -19, -21, -17, -24, -18, -23, -22, -20,
922     15, 1, 11, 7, 9, 10, 14, 3, 19, 20, 8, 2, 22, 16, 23, 12, 17, 0, 6, 5, 18, 13, 21, 4, -5, -22, -14, -19, -6, -7, -1, -18, -13, -24, -17, -23, -3, -9, -21, -20, -4, -15, -11, -10, -8, -12, -2, -16,
923     4, 15, 14, 6, 13, 7, 12, 9, 18, 16, 20, 23, 5, 8, 21, 11, 22, 19, 10, 17, 0, 3, 2, 1, -4, -1, -2, -3, -24, -12, -21, -19, -11, -17, -6, -9, -18, -20, -22, -23, -15, -5, -16, -7, -14, -10, -8, -13,
924     17, 18, 16, 19, 20, 22, 23, 21, 7, 9, 6, 4, 10, 11, 14, 5, 15, 13, 8, 12, 1, 0, 3, 2, -3, -4, -1, -2, -13, -9, -14, -16, -6, -15, -12, -11, -5, -7, -10, -8, -22, -24, -23, -21, -20, -17, -19, -18,
925     12, 7, 9, 13, 6, 15, 4, 14, 22, 21, 19, 17, 8, 5, 16, 10, 18, 20, 11, 23, 2, 1, 0, 3, -2, -3, -4, -1, -18, -11, -20, -23, -12, -22, -9, -6, -24, -21, -17, -19, -10, -13, -8, -14, -7, -15, -16, -5,
926     23, 22, 21, 20, 19, 18, 17, 16, 15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5, -6, -7, -8, -9, -10, -11, -12, -13, -14, -15, -16, -17, -18, -19, -20, -21, -22, -23, -24,
927     -24, -23, -22, -21, -20, -19, -18, -17, -16, -15, -14, -13, -12, -11, -10, -9, -8, -7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23,
928     -13, -8, -10, -14, -7, -16, -5, -15, -23, -22, -20, -18, -9, -6, -17, -11, -19, -21, -12, -24, -3, -2, -1, -4, 1, 2, 3, 0, 17, 10, 19, 22, 11, 21, 8, 5, 23, 20, 16, 18, 9, 12, 7, 13, 6, 14, 15, 4,
929     -18, -19, -17, -20, -21, -23, -24, -22, -8, -10, -7, -5, -11, -12, -15, -6, -16, -14, -9, -13, -2, -1, -4, -3, 2, 3, 0, 1, 12, 8, 13, 15, 5, 14, 11, 10, 4, 6, 9, 7, 21, 23, 22, 20, 19, 16, 18, 17,
930     -5, -16, -15, -7, -14, -8, -13, -10, -19, -17, -21, -24, -6, -9, -22, -12, -23, -20, -11, -18, -1, -4, -3, -2, 3, 0, 1, 2, 23, 11, 20, 18, 10, 16, 5, 8, 17, 19, 21, 22, 14, 4, 15, 6, 13, 9, 7, 12,
931     -16, -2, -12, -8, -10, -11, -15, -4, -20, -21, -9, -3, -23, -17, -24, -13, -18, -1, -7, -6, -19, -14, -22, -5, 4, 21, 13, 18, 5, 6, 0, 17, 12, 23, 16, 22, 2, 8, 20, 19, 3, 14, 10, 9, 7, 11, 1, 15,
932     -20, -22, -23, -18, -24, -17, -21, -19, -10, -8, -13, -14, -2, -4, -16, -3, -15, -5, -1, -7, -11, -9, -12, -6, 5, 11, 8, 10, 6, 0, 4, 14, 2, 15, 3, 1, 13, 12, 7, 9, 18, 20, 16, 23, 17, 22, 21, 19,
933     -10, -12, -2, -15, -16, -4, -8, -11, -24, -18, -3, -9, -22, -19, -20, -14, -21, -6, -5, -1, -17, -13, -23, -7, 6, 22, 12, 16, 0, 4, 5, 20, 13, 19, 18, 21, 8, 2, 17, 23, 10, 7, 3, 15, 14, 1, 11, 9,
934     -12, -7, -14, -4, -11, -5, -2, -13, -6, -3, -19, -23, -21, -18, -9, -8, -1, -17, -15, -22, -24, -16, -20, -10, 7, 17, 14, 20, 18, 15, 22, 8, 9, 0, 19, 23, 16, 21, 5, 2, 6, 10, 13, 1, 11, 4, 12, 3,
935     -21, -17, -19, -24, -18, -22, -20, -23, -15, -16, -5, -7, -4, -2, -8, -1, -10, -13, -3, -14, -12, -6, -11, -9, 8, 10, 5, 11, 13, 2, 12, 9, 0, 7, 1, 3, 6, 4, 15, 14, 22, 19, 21, 17, 23, 18, 16, 20,
936     -4, -13, -5, -12, -2, -14, -11, -7, -3, -6, -22, -17, -24, -20, -1, -10, -9, -23, -16, -19, -21, -15, -18, -8, 9, 19, 15, 23, 21, 14, 16, 0, 7, 8, 17, 20, 22, 18, 2, 5, 12, 1, 4, 10, 3, 13, 6, 11,
937     -7, -10, -8, -5, -13, -15, -14, -16, -17, -19, -18, -20, -1, -3, -23, -2, -22, -24, -4, -21, -6, -11, -9, -12, 10, 5, 11, 8, 19, 1, 17, 16, 3, 18, 0, 2, 20, 23, 22, 21, 7, 6, 9, 4, 12, 15, 14, 13,
938     -14, -15, -16, -13, -5, -10, -7, -8, -22, -23, -24, -21, -3, -1, -19, -4, -17, -18, -2, -20, -9, -12, -6, -11, 11, 8, 10, 5, 20, 3, 23, 21, 1, 22, 2, 0, 19, 17, 18, 16, 15, 13, 14, 12, 4, 7, 9, 6,
939     -8, -4, -11, -16, -15, -12, -10, -2, -21, -20, -6, -1, -19, -22, -18, -5, -24, -3, -14, -9, -23, -7, -17, -13, 12, 16, 6, 22, 8, 13, 2, 23, 4, 17, 21, 18, 0, 5, 19, 20, 1, 9, 11, 14, 15, 10, 3, 7,
940     -15, -11, -4, -10, -8, -2, -16, -12, -18, -24, -1, -6, -17, -23, -21, -7, -20, -9, -13, -3, -22, -5, -19, -14, 13, 18, 4, 21, 2, 12, 8, 19, 6, 20, 22, 16, 5, 0, 23, 17, 11, 15, 1, 7, 9, 3, 10, 14,
941     -2, -5, -13, -11, -4, -7, -12, -14, -1, -9, -17, -22, -18, -21, -3, -15, -6, -19, -8, -23, -20, -10, -24, -16, 14, 20, 7, 17, 16, 9, 21, 2, 15, 5, 23, 19, 18, 22, 0, 8, 4, 3, 12, 11, 1, 6, 13, 10,
942     -11, -14, -7, -2, -12, -13, -4, -5, -9, -1, -23, -19, -20, -24, -6, -16, -3, -22, -10, -17, -18, -8, -21, -15, 15, 23, 9, 19, 22, 7, 18, 5, 14, 2, 20, 17, 21, 16, 8, 0, 13, 11, 6, 3, 10, 12, 4, 1,
943     -1, -24, -18, -6, -3, -21, -9, -20, -4, -11, -15, -10, -5, -14, -2, -22, -12, -16, -19, -8, -7, -17, -13, -23, 16, 6, 22, 12, 9, 21, 14, 3, 18, 10, 4, 13, 7, 15, 1, 11, 17, 0, 23, 5, 2, 19, 20, 8,
944     -23, -1, -9, -19, -17, -6, -22, -3, -7, -14, -11, -2, -8, -15, -13, -18, -5, -4, -21, -12, -16, -20, -10, -24, 17, 14, 20, 7, 10, 19, 1, 12, 23, 4, 9, 15, 3, 11, 6, 13, 0, 16, 8, 21, 22, 5, 2, 18,
945     -6, -20, -21, -1, -9, -18, -3, -24, -11, -4, -8, -16, -7, -13, -12, -23, -2, -10, -17, -15, -5, -19, -14, -22, 18, 4, 21, 13, 15, 22, 7, 10, 16, 3, 6, 12, 14, 9, 11, 1, 20, 5, 19, 0, 8, 23, 17, 2,
946     -17, -9, -1, -22, -23, -3, -19, -6, -13, -5, -2, -11, -10, -16, -7, -20, -14, -12, -24, -4, -15, -18, -8, -21, 19, 15, 23, 9, 1, 17, 10, 6, 20, 13, 7, 14, 11, 3, 12, 4, 8, 22, 0, 18, 16, 2, 5, 21,
947     -22, -6, -3, -17, -19, -1, -23, -9, -5, -13, -4, -12, -15, -8, -14, -21, -7, -11, -18, -2, -10, -24, -16, -20, 20, 7, 17, 14, 3, 23, 11, 13, 19, 6, 15, 9, 10, 1, 4, 12, 5, 18, 2, 22, 21, 0, 8, 16,
948     -3, -18, -24, -9, -1, -20, -6, -21, -2, -12, -10, -15, -13, -7, -4, -17, -11, -8, -23, -16, -14, -22, -5, -19, 21, 13, 18, 4, 14, 16, 9, 1, 22, 11, 12, 6, 15, 7, 3, 10, 23, 2, 17, 8, 0, 20, 19, 5,
949     -9, -21, -20, -3, -6, -24, -1, -18, -12, -2, -16, -8, -14, -5, -11, -19, -4, -15, -22, -10, -13, -23, -7, -17, 22, 12, 16, 6, 7, 18, 15, 11, 21, 1, 13, 4, 9, 14, 10, 3, 19, 8, 20, 2, 5, 17, 23, 0,
950     -19, -3, -6, -23, -22, -9, -17, -1, -14, -7, -12, -4, -16, -10, -5, -24, -13, -2, -20, -11, -8, -21, -15, -18, 23, 9, 19, 15, 11, 20, 3, 4, 17, 12, 14, 7, 1, 10, 13, 6, 2, 21, 5, 16, 18, 8, 0, 22,
951     };
952   static const PetscInt tripMult[12*12] = {
953     1, 0, 2, 3, 5, 4, -6, -4, -5, -2, -3, -1,
954     0, 2, 1, 4, 3, 5, -5, -6, -4, -3, -1, -2,
955     2, 1, 0, 5, 4, 3, -4, -5, -6, -1, -2, -3,
956     4, 3, 5, 0, 2, 1, -3, -1, -2, -5, -6, -4,
957     3, 5, 4, 1, 0, 2, -2, -3, -1, -6, -4, -5,
958     5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5, -6,
959     -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5,
960     -4, -6, -5, -2, -1, -3, 1, 2, 0, 5, 3, 4,
961     -5, -4, -6, -1, -3, -2, 2, 0, 1, 4, 5, 3,
962     -3, -2, -1, -6, -5, -4, 3, 4, 5, 0, 1, 2,
963     -1, -3, -2, -5, -4, -6, 4, 5, 3, 2, 0, 1,
964     -2, -1, -3, -4, -6, -5, 5, 3, 4, 1, 2, 0,
965   };
966   static const PetscInt ttriMult[12*12] = {
967     0, 2, 1, 3, 5, 4, -6, -4, -5, -3, -1, -2,
968     1, 0, 2, 4, 3, 5, -5, -6, -4, -2, -3, -1,
969     2, 1, 0, 5, 4, 3, -4, -5, -6, -1, -2, -3,
970     3, 5, 4, 0, 2, 1, -3, -1, -2, -6, -4, -5,
971     4, 3, 5, 1, 0, 2, -2, -3, -1, -5, -6, -4,
972     5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5, -6,
973     -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5,
974     -5, -4, -6, -2, -1, -3, 1, 2, 0, 4, 5, 3,
975     -4, -6, -5, -1, -3, -2, 2, 0, 1, 5, 3, 4,
976     -3, -2, -1, -6, -5, -4, 3, 4, 5, 0, 1, 2,
977     -2, -1, -3, -5, -4, -6, 4, 5, 3, 1, 2, 0,
978     -1, -3, -2, -4, -6, -5, 5, 3, 4, 2, 0, 1,
979   };
980   static const PetscInt tquadMult[16*16] = {
981     0, 3, 2, 1, 4, 7, 6, 5, -8, -5, -6, -7, -4, -1, -2, -3,
982     1, 0, 3, 2, 5, 4, 7, 6, -7, -8, -5, -6, -3, -4, -1, -2,
983     2, 1, 0, 3, 6, 5, 4, 7, -6, -7, -8, -5, -2, -3, -4, -1,
984     3, 2, 1, 0, 7, 6, 5, 4, -5, -6, -7, -8, -1, -2, -3, -4,
985     4, 7, 6, 5, 0, 3, 2, 1, -4, -1, -2, -3, -8, -5, -6, -7,
986     5, 4, 7, 6, 1, 0, 3, 2, -3, -4, -1, -2, -7, -8, -5, -6,
987     6, 5, 4, 7, 2, 1, 0, 3, -2, -3, -4, -1, -6, -7, -8, -5,
988     7, 6, 5, 4, 3, 2, 1, 0, -1, -2, -3, -4, -5, -6, -7, -8,
989     -8, -7, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 7,
990     -7, -6, -5, -8, -3, -2, -1, -4, 1, 2, 3, 0, 5, 6, 7, 4,
991     -6, -5, -8, -7, -2, -1, -4, -3, 2, 3, 0, 1, 6, 7, 4, 5,
992     -5, -8, -7, -6, -1, -4, -3, -2, 3, 0, 1, 2, 7, 4, 5, 6,
993     -4, -3, -2, -1, -8, -7, -6, -5, 4, 5, 6, 7, 0, 1, 2, 3,
994     -3, -2, -1, -4, -7, -6, -5, -8, 5, 6, 7, 4, 1, 2, 3, 0,
995     -2, -1, -4, -3, -6, -5, -8, -7, 6, 7, 4, 5, 2, 3, 0, 1,
996     -1, -4, -3, -2, -5, -8, -7, -6, 7, 4, 5, 6, 3, 0, 1, 2,
997   };
998   static const PetscInt pyrMult[8*8] = {
999     0, 3, 2, 1, -4, -1, -2, -3,
1000     1, 0, 3, 2, -3, -4, -1, -2,
1001     2, 1, 0, 3, -2, -3, -4, -1,
1002     3, 2, 1, 0, -1, -2, -3, -4,
1003     -4, -3, -2, -1, 0, 1, 2, 3,
1004     -3, -2, -1, -4, 1, 2, 3, 0,
1005     -2, -1, -4, -3, 2, 3, 0, 1,
1006     -1, -4, -3, -2, 3, 0, 1, 2,
1007   };
1008   switch (ct) {
1009     case DM_POLYTOPE_POINT:              return 0;
1010     case DM_POLYTOPE_SEGMENT:
1011     case DM_POLYTOPE_POINT_PRISM_TENSOR: return segMult[(o1+1)*2+o2+1];
1012     case DM_POLYTOPE_TRIANGLE:           return triMult[(o1+3)*6+o2+3];
1013     case DM_POLYTOPE_QUADRILATERAL:      return quadMult[(o1+4)*8+o2+4];
1014     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return tsegMult[(o1+2)*4+o2+2];
1015     case DM_POLYTOPE_TETRAHEDRON:        return tetMult[(o1+12)*24+o2+12];
1016     case DM_POLYTOPE_HEXAHEDRON:         return hexMult[(o1+24)*48+o2+24];
1017     case DM_POLYTOPE_TRI_PRISM:          return tripMult[(o1+6)*12+o2+6];
1018     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return ttriMult[(o1+6)*12+o2+6];
1019     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return tquadMult[(o1+8)*16+o2+8];
1020     case DM_POLYTOPE_PYRAMID:            return pyrMult[(o1+4)*8+o2+4];
1021     default: return 0;
1022   }
1023 }
1024 
1025 /* This is orientation o1 acting on orientation o2^{-1} */
1026 static inline PetscInt DMPolytopeTypeComposeOrientationInv(DMPolytopeType ct, PetscInt o1, PetscInt o2)
1027 {
1028   static const PetscInt triInv[6]    = {-3, -2, -1, 0, 2, 1};
1029   static const PetscInt quadInv[8]   = {-4, -3, -2, -1, 0, 3, 2, 1};
1030   static const PetscInt tetInv[24]   = {-9, -11, -4, -12, -5, -7, -6, -8, -10, -3, -2, -1, 0, 2, 1, 3, 8, 10, 6, 11, 4, 9, 5, 7};
1031   static const PetscInt hexInv[48]   = {-17, -18, -20, -19, -22, -21, -23, -24, -15, -16, -14, -13, -11, -12, -10, -9, -8, -5, -6, -7, -4, -3, -2, -1,
1032                                           0,   3,   2,   1,   6,   5,   4,   9,   8,   7,  10,  11,  12,  13,  14, 15, 17, 16, 19, 18, 21, 20, 23, 22};
1033   static const PetscInt tripInv[12]  = {-5, -6, -4, -3, -2, -1, 0, 2, 1, 3, 4, 5};
1034   static const PetscInt ttriInv[12]  = {-6, -5, -4, -3, -2, -1, 0, 2, 1, 3, 5, 4};
1035   static const PetscInt tquadInv[16] = {-8, -7, -6, -5, -4, -3, -2, -1, 0, 3, 2, 1, 4, 7, 6, 5};
1036   static const PetscInt pyrInv[8]    = {-4, -3, -2, -1, 0, 3, 2, 1};
1037   switch (ct) {
1038     case DM_POLYTOPE_POINT:              return 0;
1039     case DM_POLYTOPE_SEGMENT:
1040     case DM_POLYTOPE_POINT_PRISM_TENSOR: return DMPolytopeTypeComposeOrientation(ct, o1, o2);
1041     case DM_POLYTOPE_TRIANGLE:           return DMPolytopeTypeComposeOrientation(ct, o1, triInv[o2+3]);
1042     case DM_POLYTOPE_QUADRILATERAL:      return DMPolytopeTypeComposeOrientation(ct, o1, quadInv[o2+4]);
1043     case DM_POLYTOPE_SEG_PRISM_TENSOR:   return DMPolytopeTypeComposeOrientation(ct, o1, o2);
1044     case DM_POLYTOPE_TETRAHEDRON:        return DMPolytopeTypeComposeOrientation(ct, o1, tetInv[o2+12]);
1045     case DM_POLYTOPE_HEXAHEDRON:         return DMPolytopeTypeComposeOrientation(ct, o1, hexInv[o2+24]);
1046     case DM_POLYTOPE_TRI_PRISM:          return DMPolytopeTypeComposeOrientation(ct, o1, tripInv[o2+6]);
1047     case DM_POLYTOPE_TRI_PRISM_TENSOR:   return DMPolytopeTypeComposeOrientation(ct, o1, ttriInv[o2+6]);
1048     case DM_POLYTOPE_QUAD_PRISM_TENSOR:  return DMPolytopeTypeComposeOrientation(ct, o1, tquadInv[o2+8]);
1049     case DM_POLYTOPE_PYRAMID:            return DMPolytopeTypeComposeOrientation(ct, o1, pyrInv[o2+4]);
1050     default: return 0;
1051   }
1052 }
1053 
1054 PETSC_EXTERN PetscErrorCode DMPolytopeMatchOrientation(DMPolytopeType, const PetscInt[], const PetscInt[], PetscInt *, PetscBool *);
1055 PETSC_EXTERN PetscErrorCode DMPolytopeMatchVertexOrientation(DMPolytopeType, const PetscInt[], const PetscInt[], PetscInt *, PetscBool *);
1056 PETSC_EXTERN PetscErrorCode DMPolytopeGetOrientation(DMPolytopeType, const PetscInt[], const PetscInt[], PetscInt *);
1057 PETSC_EXTERN PetscErrorCode DMPolytopeGetVertexOrientation(DMPolytopeType, const PetscInt[], const PetscInt[], PetscInt *);
1058 PETSC_EXTERN PetscErrorCode DMPolytopeInCellTest(DMPolytopeType, const PetscReal[], PetscBool *);
1059 
1060 #endif
1061