Actual source code: ex11.c
1: static char help[] = "Mesh Orientation Tutorial\n\n";
3: #include <petscdmplex.h>
4: #include <petscdmplextransform.h>
6: typedef struct {
7: PetscBool genArr; /* Generate all possible cell arrangements */
8: PetscBool refArr; /* Refine all possible cell arrangements */
9: PetscBool printTable; /* Print the CAyley table */
10: PetscInt orntBounds[2]; /* Bounds for the orientation check */
11: PetscInt numOrnt; /* Number of specific orientations specified, or -1 for all orientations */
12: PetscInt ornts[48]; /* Specific orientations if specified */
13: PetscInt initOrnt; /* Initial orientation for starting mesh */
14: } AppCtx;
16: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
17: {
18: PetscInt n = 2;
19: PetscBool flg;
21: PetscFunctionBeginUser;
22: options->genArr = PETSC_FALSE;
23: options->refArr = PETSC_FALSE;
24: options->printTable = PETSC_FALSE;
25: options->orntBounds[0] = PETSC_INT_MIN;
26: options->orntBounds[1] = PETSC_INT_MAX;
27: options->numOrnt = -1;
28: options->initOrnt = 0;
30: PetscOptionsBegin(comm, "", "Mesh Orientation Tutorials Options", "DMPLEX");
31: PetscCall(PetscOptionsBool("-gen_arrangements", "Flag for generating all arrangements of the cell", "ex11.c", options->genArr, &options->genArr, NULL));
32: PetscCall(PetscOptionsBool("-ref_arrangements", "Flag for refining all arrangements of the cell", "ex11.c", options->refArr, &options->refArr, NULL));
33: PetscCall(PetscOptionsBool("-print_table", "Print the Cayley table", "ex11.c", options->printTable, &options->printTable, NULL));
34: PetscCall(PetscOptionsIntArray("-ornt_bounds", "Bounds for orientation checks", "ex11.c", options->orntBounds, &n, NULL));
35: n = 48;
36: PetscCall(PetscOptionsIntArray("-ornts", "Specific orientations for checks", "ex11.c", options->ornts, &n, &flg));
37: if (flg) {
38: options->numOrnt = n;
39: PetscCall(PetscSortInt(n, options->ornts));
40: }
41: PetscCall(PetscOptionsInt("-init_ornt", "Initial orientation for starting mesh", "ex11.c", options->initOrnt, &options->initOrnt, NULL));
42: PetscOptionsEnd();
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: static PetscBool ignoreOrnt(AppCtx *user, PetscInt o)
47: {
48: PetscInt loc;
49: PetscErrorCode ierr;
51: if (user->numOrnt < 0) return PETSC_FALSE;
52: ierr = PetscFindInt(o, user->numOrnt, user->ornts, &loc);
53: if (loc < 0 || ierr) return PETSC_TRUE;
54: return PETSC_FALSE;
55: }
57: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
58: {
59: PetscFunctionBeginUser;
60: PetscCall(DMCreate(comm, dm));
61: PetscCall(DMSetType(*dm, DMPLEX));
62: PetscCall(DMSetFromOptions(*dm));
63: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: static PetscErrorCode CheckCellVertices(DM dm, PetscInt cell, PetscInt o)
68: {
69: DMPolytopeType ct;
70: const PetscInt *arrVerts;
71: PetscInt *closure = NULL;
72: PetscInt Ncl, cl, Nv, vStart, vEnd, v;
73: MPI_Comm comm;
75: PetscFunctionBeginUser;
76: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
77: PetscCall(DMPlexGetCellType(dm, cell, &ct));
78: PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
79: PetscCall(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
80: for (cl = 0, Nv = 0; cl < Ncl * 2; cl += 2) {
81: const PetscInt vertex = closure[cl];
83: if (vertex < vStart || vertex >= vEnd) continue;
84: closure[Nv++] = vertex;
85: }
86: PetscCheck(Nv == DMPolytopeTypeGetNumVertices(ct), comm, PETSC_ERR_ARG_WRONG, "Cell %" PetscInt_FMT " has %" PetscInt_FMT " vertices != %" PetscInt_FMT " vertices in a %s", cell, Nv, DMPolytopeTypeGetNumVertices(ct), DMPolytopeTypes[ct]);
87: arrVerts = DMPolytopeTypeGetVertexArrangement(ct, o);
88: for (v = 0; v < Nv; ++v) {
89: PetscCheck(closure[v] == arrVerts[v] + vStart, comm, PETSC_ERR_ARG_WRONG, "Cell %" PetscInt_FMT " vertex[%" PetscInt_FMT "]: %" PetscInt_FMT " should be %" PetscInt_FMT " for arrangement %" PetscInt_FMT, cell, v, closure[v], arrVerts[v] + vStart, o);
90: }
91: PetscCall(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
92: PetscFunctionReturn(PETSC_SUCCESS);
93: }
95: /* Transform cell with group operation o */
96: static PetscErrorCode ReorientCell(DM dm, PetscInt cell, PetscInt o, PetscBool swapCoords)
97: {
98: DM cdm;
99: Vec coordinates;
100: PetscScalar *coords, *ccoords = NULL;
101: PetscInt *closure = NULL;
102: PetscInt cdim, d, Nc, Ncl, cl, vStart, vEnd, Nv;
104: PetscFunctionBegin;
105: /* Change vertex coordinates so that it plots as we expect */
106: PetscCall(DMGetCoordinateDM(dm, &cdm));
107: PetscCall(DMGetCoordinateDim(dm, &cdim));
108: PetscCall(DMGetCoordinatesLocal(dm, &coordinates));
109: PetscCall(DMPlexVecGetClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords));
110: /* Reorient cone */
111: PetscCall(DMPlexOrientPoint(dm, cell, o));
112: /* Finish resetting coordinates */
113: if (swapCoords) {
114: PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
115: PetscCall(VecGetArrayWrite(coordinates, &coords));
116: PetscCall(DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
117: for (cl = 0, Nv = 0; cl < Ncl * 2; cl += 2) {
118: const PetscInt vertex = closure[cl];
119: PetscScalar *vcoords;
121: if (vertex < vStart || vertex >= vEnd) continue;
122: PetscCall(DMPlexPointLocalRef(cdm, vertex, coords, &vcoords));
123: for (d = 0; d < cdim; ++d) vcoords[d] = ccoords[Nv * cdim + d];
124: ++Nv;
125: }
126: PetscCall(DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure));
127: PetscCall(VecRestoreArrayWrite(coordinates, &coords));
128: }
129: PetscCall(DMPlexVecRestoreClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords));
130: PetscFunctionReturn(PETSC_SUCCESS);
131: }
133: static PetscErrorCode GenerateArrangements(DM dm, AppCtx *user)
134: {
135: DM odm;
136: DMPolytopeType ct;
137: PetscInt No, o;
138: const char *name;
140: PetscFunctionBeginUser;
141: if (!user->genArr) PetscFunctionReturn(PETSC_SUCCESS);
142: PetscCall(PetscObjectGetName((PetscObject)dm, &name));
143: PetscCall(DMPlexGetCellType(dm, 0, &ct));
144: No = DMPolytopeTypeGetNumArrangements(ct) / 2;
145: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
146: if (ignoreOrnt(user, o)) continue;
147: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm));
148: PetscCall(ReorientCell(odm, 0, o, PETSC_TRUE));
149: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o));
150: PetscCall(DMViewFromOptions(odm, NULL, "-gen_dm_view"));
151: PetscCall(CheckCellVertices(odm, 0, o));
152: PetscCall(DMDestroy(&odm));
153: }
154: PetscFunctionReturn(PETSC_SUCCESS);
155: }
157: static PetscErrorCode VerifyCayleyTable(DM dm, AppCtx *user)
158: {
159: DM dm1, dm2;
160: DMPolytopeType ct;
161: const PetscInt *refcone, *cone;
162: PetscInt No, o1, o2, o3, o4;
163: PetscBool equal;
164: const char *name;
166: PetscFunctionBeginUser;
167: if (!user->genArr) PetscFunctionReturn(PETSC_SUCCESS);
168: PetscCall(PetscObjectGetName((PetscObject)dm, &name));
169: PetscCall(DMPlexGetCellType(dm, 0, &ct));
170: PetscCall(DMPlexGetCone(dm, 0, &refcone));
171: No = DMPolytopeTypeGetNumArrangements(ct) / 2;
172: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Cayley Table for %s\n", DMPolytopeTypes[ct]));
173: for (o1 = PetscMax(-No, user->orntBounds[0]); o1 < PetscMin(No, user->orntBounds[1]); ++o1) {
174: for (o2 = PetscMax(-No, user->orntBounds[0]); o2 < PetscMin(No, user->orntBounds[1]); ++o2) {
175: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1));
176: PetscCall(DMPlexOrientPoint(dm1, 0, o2));
177: PetscCall(DMPlexCheckFaces(dm1, 0));
178: PetscCall(DMPlexOrientPoint(dm1, 0, o1));
179: PetscCall(DMPlexCheckFaces(dm1, 0));
180: o3 = DMPolytopeTypeComposeOrientation(ct, o1, o2);
181: /* First verification */
182: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2));
183: PetscCall(DMPlexOrientPoint(dm2, 0, o3));
184: PetscCall(DMPlexCheckFaces(dm2, 0));
185: PetscCall(DMPlexEqual(dm1, dm2, &equal));
186: if (!equal) {
187: PetscCall(DMViewFromOptions(dm1, NULL, "-error_dm_view"));
188: PetscCall(DMViewFromOptions(dm2, NULL, "-error_dm_view"));
189: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cayley table error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != %" PetscInt_FMT, DMPolytopeTypes[ct], o1, o2, o3);
190: }
191: /* Second verification */
192: PetscCall(DMPlexGetCone(dm1, 0, &cone));
193: PetscCall(DMPolytopeGetOrientation(ct, refcone, cone, &o4));
194: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", o4));
195: PetscCheck(o3 == o4, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cayley table error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " = %" PetscInt_FMT " != %" PetscInt_FMT, DMPolytopeTypes[ct], o1, o2, o3, o4);
196: PetscCall(DMDestroy(&dm1));
197: PetscCall(DMDestroy(&dm2));
198: }
199: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
200: }
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: static PetscErrorCode VerifyInverse(DM dm, AppCtx *user)
205: {
206: DM dm1, dm2;
207: DMPolytopeType ct;
208: const PetscInt *refcone, *cone;
209: PetscInt No, o, oi, o2;
210: PetscBool equal;
211: const char *name;
213: PetscFunctionBeginUser;
214: if (!user->genArr) PetscFunctionReturn(PETSC_SUCCESS);
215: PetscCall(PetscObjectGetName((PetscObject)dm, &name));
216: PetscCall(DMPlexGetCellType(dm, 0, &ct));
217: PetscCall(DMPlexGetCone(dm, 0, &refcone));
218: No = DMPolytopeTypeGetNumArrangements(ct) / 2;
219: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Inverse table for %s\n", DMPolytopeTypes[ct]));
220: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
221: if (ignoreOrnt(user, o)) continue;
222: oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o);
223: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1));
224: PetscCall(DMPlexOrientPoint(dm1, 0, o));
225: PetscCall(DMPlexCheckFaces(dm1, 0));
226: PetscCall(DMPlexOrientPoint(dm1, 0, oi));
227: PetscCall(DMPlexCheckFaces(dm1, 0));
228: /* First verification */
229: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2));
230: PetscCall(DMPlexEqual(dm1, dm2, &equal));
231: if (!equal) {
232: PetscCall(DMViewFromOptions(dm1, NULL, "-error_dm_view"));
233: PetscCall(DMViewFromOptions(dm2, NULL, "-error_dm_view"));
234: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inverse error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != 0", DMPolytopeTypes[ct], o, oi);
235: }
236: /* Second verification */
237: PetscCall(DMPlexGetCone(dm1, 0, &cone));
238: PetscCall(DMPolytopeGetOrientation(ct, refcone, cone, &o2));
239: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", oi));
240: PetscCheck(o2 == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inverse error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " = %" PetscInt_FMT " != 0", DMPolytopeTypes[ct], o, oi, o2);
241: PetscCall(DMDestroy(&dm1));
242: PetscCall(DMDestroy(&dm2));
243: }
244: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
245: PetscFunctionReturn(PETSC_SUCCESS);
246: }
248: /* Suppose that point p has the same arrangement as o from canonical, compare the subcells to canonical subcells */
249: static PetscErrorCode CheckSubcells(DM dm, DM odm, PetscInt p, PetscInt o, AppCtx *user)
250: {
251: DMPlexTransform tr, otr;
252: DMPolytopeType ct;
253: DMPolytopeType *rct;
254: const PetscInt *cone, *ornt, *ocone, *oornt;
255: PetscInt *rsize, *rcone, *rornt;
256: PetscInt Nct, n, oi, debug = 0;
258: PetscFunctionBeginUser;
259: PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr));
260: PetscCall(DMPlexTransformSetDM(tr, dm));
261: PetscCall(DMPlexTransformSetFromOptions(tr));
262: PetscCall(DMPlexTransformSetUp(tr));
264: PetscCall(DMPlexTransformCreate(PetscObjectComm((PetscObject)odm), &otr));
265: PetscCall(DMPlexTransformSetDM(otr, odm));
266: PetscCall(DMPlexTransformSetFromOptions(otr));
267: PetscCall(DMPlexTransformSetUp(otr));
269: PetscCall(DMPlexGetCellType(dm, p, &ct));
270: PetscCall(DMPlexGetCone(dm, p, &cone));
271: PetscCall(DMPlexGetConeOrientation(dm, p, &ornt));
272: PetscCall(DMPlexGetCone(odm, p, &ocone));
273: PetscCall(DMPlexGetConeOrientation(odm, p, &oornt));
274: oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o);
275: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "Orientation %" PetscInt_FMT "\n", oi));
277: PetscCall(DMPlexTransformCellTransform(tr, ct, p, NULL, &Nct, &rct, &rsize, &rcone, &rornt));
278: for (n = 0; n < Nct; ++n) {
279: DMPolytopeType ctNew = rct[n];
280: PetscInt ro;
282: if (debug) PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking type %s\n", DMPolytopeTypes[ctNew]));
283: for (PetscInt r = 0; r < rsize[n]; ++r) {
284: const PetscInt *qcone, *qornt, *oqcone, *oqornt;
285: PetscInt pNew, opNew, oo, pr, fo;
286: PetscBool restore = PETSC_TRUE;
288: PetscCall(DMPlexTransformGetTargetPoint(tr, ct, ctNew, p, r, &pNew));
289: PetscCall(DMPlexTransformGetCone(tr, pNew, &qcone, &qornt));
290: if (debug) {
291: PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking replica %" PetscInt_FMT " (%" PetscInt_FMT ")\n Original Cone", r, pNew));
292: for (PetscInt c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", qcone[c], qornt[c]));
293: PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
294: }
295: for (ro = 0; ro < rsize[n]; ++ro) {
296: PetscBool found;
298: PetscCall(DMPlexTransformGetTargetPoint(otr, ct, ctNew, p, ro, &opNew));
299: PetscCall(DMPlexTransformGetConeOriented(otr, opNew, o, &oqcone, &oqornt));
300: PetscCall(DMPolytopeMatchOrientation(ctNew, oqcone, qcone, &oo, &found));
301: if (found) break;
302: PetscCall(DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt));
303: }
304: if (debug) {
305: PetscCall(PetscPrintf(PETSC_COMM_SELF, " Checking transform replica %" PetscInt_FMT " (%" PetscInt_FMT ") (%" PetscInt_FMT ")\n Transform Cone", ro, opNew, o));
306: for (PetscInt c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", oqcone[c], oqornt[c]));
307: PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
308: PetscCall(PetscPrintf(PETSC_COMM_SELF, " Matched %" PetscInt_FMT "\n", oo));
309: }
310: if (ro == rsize[n]) {
311: /* The tetrahedron has 3 pairs of opposing edges, and any pair can be connected by the interior segment */
312: if (ct == DM_POLYTOPE_TETRAHEDRON) {
313: /* The segment in a tetrahedron does not map into itself under the group action */
314: if (ctNew == DM_POLYTOPE_SEGMENT) {
315: restore = PETSC_FALSE;
316: ro = r;
317: oo = 0;
318: }
319: /* The last four interior faces do not map into themselves under the group action */
320: if (r > 3 && ctNew == DM_POLYTOPE_TRIANGLE) {
321: restore = PETSC_FALSE;
322: ro = r;
323: oo = 0;
324: }
325: /* The last four interior faces do not map into themselves under the group action */
326: if (r > 3 && ctNew == DM_POLYTOPE_TETRAHEDRON) {
327: restore = PETSC_FALSE;
328: ro = r;
329: oo = 0;
330: }
331: }
332: PetscCheck(!restore, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Unable to find matching %s %" PetscInt_FMT " orientation for cell orientation %" PetscInt_FMT, DMPolytopeTypes[ctNew], r, o);
333: }
334: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", %" PetscInt_FMT ", ", ro, oo));
335: PetscCall(DMPlexTransformGetSubcellOrientation(tr, ct, p, oi, ctNew, r, 0, &pr, &fo));
336: if (!user->printTable) {
337: PetscCheck(pr == ro, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Choose wrong replica %" PetscInt_FMT " != %" PetscInt_FMT, pr, ro);
338: PetscCheck(fo == oo, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Choose wrong orientation %" PetscInt_FMT " != %" PetscInt_FMT, fo, oo);
339: }
340: PetscCall(DMPlexTransformRestoreCone(tr, pNew, &qcone, &qornt));
341: if (restore) PetscCall(DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt));
342: }
343: if (user->printTable) PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
344: }
345: PetscCall(DMPlexTransformDestroy(&tr));
346: PetscCall(DMPlexTransformDestroy(&otr));
347: PetscFunctionReturn(PETSC_SUCCESS);
348: }
350: static PetscErrorCode RefineArrangements(DM dm, AppCtx *user)
351: {
352: DM odm, rdm;
353: DMPolytopeType ct;
354: PetscInt No, o;
355: const char *name;
357: PetscFunctionBeginUser;
358: if (!user->refArr) PetscFunctionReturn(PETSC_SUCCESS);
359: PetscCall(PetscObjectGetName((PetscObject)dm, &name));
360: PetscCall(DMPlexGetCellType(dm, 0, &ct));
361: No = DMPolytopeTypeGetNumArrangements(ct) / 2;
362: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
363: if (ignoreOrnt(user, o)) continue;
364: PetscCall(CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm));
365: if (user->initOrnt) PetscCall(ReorientCell(odm, 0, user->initOrnt, PETSC_FALSE));
366: PetscCall(ReorientCell(odm, 0, o, PETSC_TRUE));
367: PetscCall(DMViewFromOptions(odm, NULL, "-orig_dm_view"));
368: PetscCall(DMRefine(odm, MPI_COMM_NULL, &rdm));
369: PetscCall(DMSetFromOptions(rdm));
370: PetscCall(PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o));
371: PetscCall(DMViewFromOptions(rdm, NULL, "-ref_dm_view"));
372: PetscCall(CheckSubcells(dm, odm, 0, o, user));
373: PetscCall(DMDestroy(&odm));
374: PetscCall(DMDestroy(&rdm));
375: }
376: PetscFunctionReturn(PETSC_SUCCESS);
377: }
379: int main(int argc, char **argv)
380: {
381: DM dm;
382: AppCtx user;
384: PetscFunctionBeginUser;
385: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
386: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
387: PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
388: if (user.initOrnt) {
389: PetscCall(ReorientCell(dm, 0, user.initOrnt, PETSC_FALSE));
390: PetscCall(DMViewFromOptions(dm, NULL, "-ornt_dm_view"));
391: }
392: PetscCall(GenerateArrangements(dm, &user));
393: PetscCall(VerifyCayleyTable(dm, &user));
394: PetscCall(VerifyInverse(dm, &user));
395: PetscCall(RefineArrangements(dm, &user));
396: PetscCall(DMDestroy(&dm));
397: PetscCall(PetscFinalize());
398: return 0;
399: }
401: /*TEST
403: testset:
404: args: -dm_coord_space 0 -dm_plex_reference_cell_domain -gen_arrangements \
405: -gen_dm_view ::ascii_latex -dm_plex_view_tikzscale 0.5
407: test:
408: suffix: segment
409: args: -dm_plex_cell segment \
410: -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0
412: test:
413: suffix: triangle
414: args: -dm_plex_cell triangle \
415: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
417: test:
418: suffix: quadrilateral
419: args: -dm_plex_cell quadrilateral \
420: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
422: test:
423: suffix: tensor_segment
424: args: -dm_plex_cell tensor_quad \
425: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
427: test:
428: suffix: tetrahedron
429: args: -dm_plex_cell tetrahedron \
430: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
432: test:
433: suffix: hexahedron
434: args: -dm_plex_cell hexahedron \
435: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 0.3
437: test:
438: suffix: triangular_prism
439: args: -dm_plex_cell triangular_prism \
440: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
442: test:
443: suffix: tensor_triangular_prism
444: args: -dm_plex_cell tensor_triangular_prism \
445: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
447: test:
448: suffix: tensor_quadrilateral_prism
449: args: -dm_plex_cell tensor_quadrilateral_prism \
450: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
452: test:
453: suffix: pyramid
454: args: -dm_plex_cell pyramid \
455: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
457: testset:
458: args: -dm_coord_space 0 -dm_plex_reference_cell_domain -ref_arrangements -dm_plex_check_all \
459: -ref_dm_view ::ascii_latex -dm_plex_view_tikzscale 1.0
461: test:
462: suffix: ref_segment
463: args: -dm_plex_cell segment \
464: -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0
466: test:
467: suffix: ref_triangle
468: args: -dm_plex_cell triangle \
469: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
471: test:
472: suffix: ref_quadrilateral
473: args: -dm_plex_cell quadrilateral \
474: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
476: test:
477: suffix: ref_tensor_segment
478: args: -dm_plex_cell tensor_quad \
479: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
481: test:
482: suffix: ref_tetrahedron
483: args: -dm_plex_cell tetrahedron \
484: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
486: test:
487: suffix: ref_hexahedron
488: args: -dm_plex_cell hexahedron \
489: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
491: test:
492: suffix: ref_triangular_prism
493: args: -dm_plex_cell triangular_prism \
494: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
496: test:
497: suffix: ref_tensor_triangular_prism
498: args: -dm_plex_cell tensor_triangular_prism \
499: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
501: test:
502: suffix: ref_tensor_quadrilateral_prism
503: args: -dm_plex_cell tensor_quadrilateral_prism \
504: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
506: # ToBox should recreate the coordinate space since the cell shape changes
507: testset:
508: args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -ref_arrangements -dm_plex_check_all
510: test:
511: suffix: tobox_triangle
512: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \
513: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
515: test:
516: suffix: tobox_tensor_segment
517: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \
518: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
520: test:
521: suffix: tobox_tetrahedron
522: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \
523: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
525: test:
526: suffix: tobox_triangular_prism
527: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism \
528: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
530: test:
531: suffix: tobox_tensor_triangular_prism
532: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \
533: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
535: test:
536: suffix: tobox_tensor_quadrilateral_prism
537: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism \
538: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
540: testset:
541: args: -dm_coord_space 0 -dm_plex_transform_type refine_alfeld -ref_arrangements -dm_plex_check_all
543: test:
544: suffix: alfeld_triangle
545: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \
546: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
548: test:
549: suffix: alfeld_tetrahedron
550: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \
551: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
553: testset:
554: args: -dm_plex_transform_type refine_sbr -ref_arrangements -dm_plex_check_all
556: # This splits edge 1 of the triangle, and reflects about the added edge
557: test:
558: suffix: sbr_triangle_0
559: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -ornts -2,0 \
560: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
562: # This splits edge 0 of the triangle, and reflects about the added edge
563: test:
564: suffix: sbr_triangle_1
565: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 1 -ornts -3,0 \
566: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
568: # This splits edge 2 of the triangle, and reflects about the added edge
569: test:
570: suffix: sbr_triangle_2
571: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 2 -ornts -1,0 \
572: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
574: # This splits edges 1 and 2 of the triangle
575: test:
576: suffix: sbr_triangle_3
577: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -ornts 0 \
578: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
580: # This splits edges 1 and 0 of the triangle
581: test:
582: suffix: sbr_triangle_4
583: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -ornts 0 \
584: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
586: # This splits edges 0 and 1 of the triangle
587: test:
588: suffix: sbr_triangle_5
589: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 1 -ornts 0 \
590: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
592: # This splits edges 0 and 2 of the triangle
593: test:
594: suffix: sbr_triangle_6
595: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 1 -ornts 0 \
596: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
598: # This splits edges 2 and 0 of the triangle
599: test:
600: suffix: sbr_triangle_7
601: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 2 -ornts 0 \
602: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
604: # This splits edges 2 and 1 of the triangle
605: test:
606: suffix: sbr_triangle_8
607: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 2 -ornts 0 \
608: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
610: testset:
611: args: -dm_plex_transform_type refine_boundary_layer -dm_plex_transform_bl_splits 2 -ref_arrangements -dm_plex_check_all
613: test:
614: suffix: bl_tensor_segment
615: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \
616: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0 -dm_plex_view_tikzscale 1.0
618: # The subcell check is broken because at orientation 3, the internal triangles do not get properly permuted for the check
619: test:
620: suffix: bl_tensor_triangular_prism
621: requires: TODO
622: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \
623: -ref_dm_view ::ascii_latex -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 1.0
625: TEST*/