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