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;
22: options->genArr = PETSC_FALSE;
23: options->refArr = PETSC_FALSE;
24: options->printTable = PETSC_FALSE;
25: options->orntBounds[0] = PETSC_MIN_INT;
26: options->orntBounds[1] = PETSC_MAX_INT;
27: options->numOrnt = -1;
28: options->initOrnt = 0;
30: PetscOptionsBegin(comm, "", "Mesh Orientation Tutorials Options", "DMPLEX");
31: PetscOptionsBool("-gen_arrangements", "Flag for generating all arrangements of the cell", "ex11.c", options->genArr, &options->genArr, NULL);
32: PetscOptionsBool("-ref_arrangements", "Flag for refining all arrangements of the cell", "ex11.c", options->refArr, &options->refArr, NULL);
33: PetscOptionsBool("-print_table", "Print the Cayley table", "ex11.c", options->printTable, &options->printTable, NULL);
34: PetscOptionsIntArray("-ornt_bounds", "Bounds for orientation checks", "ex11.c", options->orntBounds, &n, NULL);
35: n = 48;
36: PetscOptionsIntArray("-ornts", "Specific orientations for checks", "ex11.c", options->ornts, &n, &flg);
37: if (flg) {
38: options->numOrnt = n;
39: PetscSortInt(n, options->ornts);
40: }
41: PetscOptionsInt("-init_ornt", "Initial orientation for starting mesh", "ex11.c", options->initOrnt, &options->initOrnt, NULL);
42: PetscOptionsEnd();
43: return 0;
44: }
46: static PetscBool ignoreOrnt(AppCtx *user, PetscInt o)
47: {
48: PetscInt loc;
51: if (user->numOrnt < 0) return PETSC_FALSE;
52: 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: {
60: DMCreate(comm, dm);
61: DMSetType(*dm, DMPLEX);
62: DMSetFromOptions(*dm);
63: DMViewFromOptions(*dm, NULL, "-dm_view");
64: return 0;
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;
76: PetscObjectGetComm((PetscObject)dm, &comm);
77: DMPlexGetCellType(dm, cell, &ct);
78: DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
79: 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: }
87: arrVerts = DMPolytopeTypeGetVertexArrangment(ct, o);
88: for (v = 0; v < Nv; ++v) {
90: }
91: DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure);
92: return 0;
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: /* Change vertex coordinates so that it plots as we expect */
105: DMGetCoordinateDM(dm, &cdm);
106: DMGetCoordinateDim(dm, &cdim);
107: DMGetCoordinatesLocal(dm, &coordinates);
108: DMPlexVecGetClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords);
109: /* Reorient cone */
110: DMPlexOrientPoint(dm, cell, o);
111: /* Finish resetting coordinates */
112: if (swapCoords) {
113: DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);
114: VecGetArrayWrite(coordinates, &coords);
115: DMPlexGetTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure);
116: for (cl = 0, Nv = 0; cl < Ncl * 2; cl += 2) {
117: const PetscInt vertex = closure[cl];
118: PetscScalar *vcoords;
120: if (vertex < vStart || vertex >= vEnd) continue;
121: DMPlexPointLocalRef(cdm, vertex, coords, &vcoords);
122: for (d = 0; d < cdim; ++d) vcoords[d] = ccoords[Nv * cdim + d];
123: ++Nv;
124: }
125: DMPlexRestoreTransitiveClosure(dm, cell, PETSC_TRUE, &Ncl, &closure);
126: VecRestoreArrayWrite(coordinates, &coords);
127: }
128: DMPlexVecRestoreClosure(cdm, NULL, coordinates, cell, &Nc, &ccoords);
129: return 0;
130: }
132: static PetscErrorCode GenerateArrangments(DM dm, AppCtx *user)
133: {
134: DM odm;
135: DMPolytopeType ct;
136: PetscInt No, o;
137: const char *name;
140: if (!user->genArr) return 0;
141: PetscObjectGetName((PetscObject)dm, &name);
142: DMPlexGetCellType(dm, 0, &ct);
143: No = DMPolytopeTypeGetNumArrangments(ct) / 2;
144: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
145: if (ignoreOrnt(user, o)) continue;
146: CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm);
147: ReorientCell(odm, 0, o, PETSC_TRUE);
148: PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o);
149: DMViewFromOptions(odm, NULL, "-gen_dm_view");
150: CheckCellVertices(odm, 0, o);
151: DMDestroy(&odm);
152: }
153: return 0;
154: }
156: static PetscErrorCode VerifyCayleyTable(DM dm, AppCtx *user)
157: {
158: DM dm1, dm2;
159: DMPolytopeType ct;
160: const PetscInt *refcone, *cone;
161: PetscInt No, o1, o2, o3, o4;
162: PetscBool equal;
163: const char *name;
166: if (!user->genArr) return 0;
167: PetscObjectGetName((PetscObject)dm, &name);
168: DMPlexGetCellType(dm, 0, &ct);
169: DMPlexGetCone(dm, 0, &refcone);
170: No = DMPolytopeTypeGetNumArrangments(ct) / 2;
171: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "Cayley Table for %s\n", DMPolytopeTypes[ct]);
172: for (o1 = PetscMax(-No, user->orntBounds[0]); o1 < PetscMin(No, user->orntBounds[1]); ++o1) {
173: for (o2 = PetscMax(-No, user->orntBounds[0]); o2 < PetscMin(No, user->orntBounds[1]); ++o2) {
174: CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1);
175: DMPlexOrientPoint(dm1, 0, o2);
176: DMPlexCheckFaces(dm1, 0);
177: DMPlexOrientPoint(dm1, 0, o1);
178: DMPlexCheckFaces(dm1, 0);
179: o3 = DMPolytopeTypeComposeOrientation(ct, o1, o2);
180: /* First verification */
181: CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2);
182: DMPlexOrientPoint(dm2, 0, o3);
183: DMPlexCheckFaces(dm2, 0);
184: DMPlexEqual(dm1, dm2, &equal);
185: if (!equal) {
186: DMViewFromOptions(dm1, NULL, "-error_dm_view");
187: DMViewFromOptions(dm2, NULL, "-error_dm_view");
188: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Cayley table error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != %" PetscInt_FMT, DMPolytopeTypes[ct], o1, o2, o3);
189: }
190: /* Second verification */
191: DMPlexGetCone(dm1, 0, &cone);
192: DMPolytopeGetOrientation(ct, refcone, cone, &o4);
193: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", o4);
195: DMDestroy(&dm1);
196: DMDestroy(&dm2);
197: }
198: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "\n");
199: }
200: return 0;
201: }
203: static PetscErrorCode VerifyInverse(DM dm, AppCtx *user)
204: {
205: DM dm1, dm2;
206: DMPolytopeType ct;
207: const PetscInt *refcone, *cone;
208: PetscInt No, o, oi, o2;
209: PetscBool equal;
210: const char *name;
213: if (!user->genArr) return 0;
214: PetscObjectGetName((PetscObject)dm, &name);
215: DMPlexGetCellType(dm, 0, &ct);
216: DMPlexGetCone(dm, 0, &refcone);
217: No = DMPolytopeTypeGetNumArrangments(ct) / 2;
218: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "Inverse table for %s\n", DMPolytopeTypes[ct]);
219: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
220: if (ignoreOrnt(user, o)) continue;
221: oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o);
222: CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm1);
223: DMPlexOrientPoint(dm1, 0, o);
224: DMPlexCheckFaces(dm1, 0);
225: DMPlexOrientPoint(dm1, 0, oi);
226: DMPlexCheckFaces(dm1, 0);
227: /* First verification */
228: CreateMesh(PetscObjectComm((PetscObject)dm), user, &dm2);
229: DMPlexEqual(dm1, dm2, &equal);
230: if (!equal) {
231: DMViewFromOptions(dm1, NULL, "-error_dm_view");
232: DMViewFromOptions(dm2, NULL, "-error_dm_view");
233: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Inverse error for %s: %" PetscInt_FMT " * %" PetscInt_FMT " != 0", DMPolytopeTypes[ct], o, oi);
234: }
235: /* Second verification */
236: DMPlexGetCone(dm1, 0, &cone);
237: DMPolytopeGetOrientation(ct, refcone, cone, &o2);
238: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", ", oi);
240: DMDestroy(&dm1);
241: DMDestroy(&dm2);
242: }
243: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "\n");
244: return 0;
245: }
247: /* Suppose that point p has the same arrangement as o from canonical, compare the subcells to canonical subcells */
248: static PetscErrorCode CheckSubcells(DM dm, DM odm, PetscInt p, PetscInt o, AppCtx *user)
249: {
250: DMPlexTransform tr, otr;
251: DMPolytopeType ct;
252: DMPolytopeType *rct;
253: const PetscInt *cone, *ornt, *ocone, *oornt;
254: PetscInt *rsize, *rcone, *rornt;
255: PetscInt Nct, n, oi, debug = 0;
258: DMPlexTransformCreate(PetscObjectComm((PetscObject)dm), &tr);
259: DMPlexTransformSetDM(tr, dm);
260: DMPlexTransformSetFromOptions(tr);
261: DMPlexTransformSetUp(tr);
263: DMPlexTransformCreate(PetscObjectComm((PetscObject)odm), &otr);
264: DMPlexTransformSetDM(otr, odm);
265: DMPlexTransformSetFromOptions(otr);
266: DMPlexTransformSetUp(otr);
268: DMPlexGetCellType(dm, p, &ct);
269: DMPlexGetCone(dm, p, &cone);
270: DMPlexGetConeOrientation(dm, p, &ornt);
271: DMPlexGetCone(odm, p, &ocone);
272: DMPlexGetConeOrientation(odm, p, &oornt);
273: oi = DMPolytopeTypeComposeOrientationInv(ct, 0, o);
274: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "Orientation %" PetscInt_FMT "\n", oi);
276: DMPlexTransformCellTransform(tr, ct, p, NULL, &Nct, &rct, &rsize, &rcone, &rornt);
277: for (n = 0; n < Nct; ++n) {
278: DMPolytopeType ctNew = rct[n];
279: PetscInt r, ro;
281: if (debug) PetscPrintf(PETSC_COMM_SELF, " Checking type %s\n", DMPolytopeTypes[ctNew]);
282: for (r = 0; r < rsize[n]; ++r) {
283: const PetscInt *qcone, *qornt, *oqcone, *oqornt;
284: PetscInt pNew, opNew, oo, pr, fo;
285: PetscBool restore = PETSC_TRUE;
287: DMPlexTransformGetTargetPoint(tr, ct, ctNew, p, r, &pNew);
288: DMPlexTransformGetCone(tr, pNew, &qcone, &qornt);
289: if (debug) {
290: PetscInt c;
292: PetscPrintf(PETSC_COMM_SELF, " Checking replica %" PetscInt_FMT " (%" PetscInt_FMT ")\n Original Cone", r, pNew);
293: for (c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", qcone[c], qornt[c]);
294: PetscPrintf(PETSC_COMM_SELF, "\n");
295: }
296: for (ro = 0; ro < rsize[n]; ++ro) {
297: PetscBool found;
299: DMPlexTransformGetTargetPoint(otr, ct, ctNew, p, ro, &opNew);
300: DMPlexTransformGetConeOriented(otr, opNew, o, &oqcone, &oqornt);
301: DMPolytopeMatchOrientation(ctNew, oqcone, qcone, &oo, &found);
302: if (found) break;
303: DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt);
304: }
305: if (debug) {
306: PetscInt c;
308: PetscPrintf(PETSC_COMM_SELF, " Checking transform replica %" PetscInt_FMT " (%" PetscInt_FMT ") (%" PetscInt_FMT ")\n Transform Cone", ro, opNew, o);
309: for (c = 0; c < DMPolytopeTypeGetConeSize(ctNew); ++c) PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT " (%" PetscInt_FMT ")", oqcone[c], oqornt[c]);
310: PetscPrintf(PETSC_COMM_SELF, "\n");
311: PetscPrintf(PETSC_COMM_SELF, " Matched %" PetscInt_FMT "\n", oo);
312: }
313: if (ro == rsize[n]) {
314: /* The tetrahedron has 3 pairs of opposing edges, and any pair can be connected by the interior segment */
315: if (ct == DM_POLYTOPE_TETRAHEDRON) {
316: /* The segment in a tetrahedron does not map into itself under the group action */
317: if (ctNew == DM_POLYTOPE_SEGMENT) {
318: restore = PETSC_FALSE;
319: ro = r;
320: oo = 0;
321: }
322: /* The last four interior faces do not map into themselves under the group action */
323: if (r > 3 && ctNew == DM_POLYTOPE_TRIANGLE) {
324: restore = PETSC_FALSE;
325: ro = r;
326: oo = 0;
327: }
328: /* The last four interior faces do not map into themselves under the group action */
329: if (r > 3 && ctNew == DM_POLYTOPE_TETRAHEDRON) {
330: restore = PETSC_FALSE;
331: ro = r;
332: oo = 0;
333: }
334: }
336: }
337: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ", %" PetscInt_FMT ", ", ro, oo);
338: DMPlexTransformGetSubcellOrientation(tr, ct, p, oi, ctNew, r, 0, &pr, &fo);
339: if (!user->printTable) {
342: }
343: DMPlexTransformRestoreCone(tr, pNew, &qcone, &qornt);
344: if (restore) DMPlexTransformRestoreCone(otr, pNew, &oqcone, &oqornt);
345: }
346: if (user->printTable) PetscPrintf(PETSC_COMM_SELF, "\n");
347: }
348: DMPlexTransformDestroy(&tr);
349: DMPlexTransformDestroy(&otr);
350: return 0;
351: }
353: static PetscErrorCode RefineArrangments(DM dm, AppCtx *user)
354: {
355: DM odm, rdm;
356: DMPolytopeType ct;
357: PetscInt No, o;
358: const char *name;
361: if (!user->refArr) return 0;
362: PetscObjectGetName((PetscObject)dm, &name);
363: DMPlexGetCellType(dm, 0, &ct);
364: No = DMPolytopeTypeGetNumArrangments(ct) / 2;
365: for (o = PetscMax(-No, user->orntBounds[0]); o < PetscMin(No, user->orntBounds[1]); ++o) {
366: if (ignoreOrnt(user, o)) continue;
367: CreateMesh(PetscObjectComm((PetscObject)dm), user, &odm);
368: if (user->initOrnt) ReorientCell(odm, 0, user->initOrnt, PETSC_FALSE);
369: ReorientCell(odm, 0, o, PETSC_TRUE);
370: DMViewFromOptions(odm, NULL, "-orig_dm_view");
371: DMRefine(odm, MPI_COMM_NULL, &rdm);
372: DMSetFromOptions(rdm);
373: PetscPrintf(PetscObjectComm((PetscObject)dm), "%s orientation %" PetscInt_FMT "\n", name, o);
374: DMViewFromOptions(rdm, NULL, "-ref_dm_view");
375: CheckSubcells(dm, odm, 0, o, user);
376: DMDestroy(&odm);
377: DMDestroy(&rdm);
378: }
379: return 0;
380: }
382: int main(int argc, char **argv)
383: {
384: DM dm;
385: AppCtx user;
388: PetscInitialize(&argc, &argv, NULL, help);
389: ProcessOptions(PETSC_COMM_WORLD, &user);
390: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
391: if (user.initOrnt) {
392: ReorientCell(dm, 0, user.initOrnt, PETSC_FALSE);
393: DMViewFromOptions(dm, NULL, "-ornt_dm_view");
394: }
395: GenerateArrangments(dm, &user);
396: VerifyCayleyTable(dm, &user);
397: VerifyInverse(dm, &user);
398: RefineArrangments(dm, &user);
399: DMDestroy(&dm);
400: PetscFinalize();
401: return 0;
402: }
404: /*TEST
406: testset:
407: args: -dm_coord_space 0 -dm_plex_reference_cell_domain -gen_arrangements \
408: -gen_dm_view ::ascii_latex -dm_plex_view_tikzscale 0.5
410: test:
411: suffix: segment
412: args: -dm_plex_cell segment \
413: -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0
415: test:
416: suffix: triangle
417: args: -dm_plex_cell triangle \
418: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
420: test:
421: suffix: quadrilateral
422: args: -dm_plex_cell quadrilateral \
423: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
425: test:
426: suffix: tensor_segment
427: args: -dm_plex_cell tensor_quad \
428: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
430: test:
431: suffix: tetrahedron
432: args: -dm_plex_cell tetrahedron \
433: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
435: test:
436: suffix: hexahedron
437: args: -dm_plex_cell hexahedron \
438: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0 -dm_plex_view_tikzscale 0.3
440: test:
441: suffix: triangular_prism
442: args: -dm_plex_cell triangular_prism \
443: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
445: test:
446: suffix: tensor_triangular_prism
447: args: -dm_plex_cell tensor_triangular_prism \
448: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
450: test:
451: suffix: tensor_quadrilateral_prism
452: args: -dm_plex_cell tensor_quadrilateral_prism \
453: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
455: test:
456: suffix: pyramid
457: args: -dm_plex_cell pyramid \
458: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
460: testset:
461: args: -dm_coord_space 0 -dm_plex_reference_cell_domain -ref_arrangements -dm_plex_check_all \
462: -ref_dm_view ::ascii_latex -dm_plex_view_tikzscale 1.0
464: test:
465: suffix: ref_segment
466: args: -dm_plex_cell segment \
467: -dm_plex_view_numbers_depth 1,0 -dm_plex_view_colors_depth 1,0
469: test:
470: suffix: ref_triangle
471: args: -dm_plex_cell triangle \
472: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
474: test:
475: suffix: ref_quadrilateral
476: args: -dm_plex_cell quadrilateral \
477: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
479: test:
480: suffix: ref_tensor_segment
481: args: -dm_plex_cell tensor_quad \
482: -dm_plex_view_numbers_depth 1,1,0 -dm_plex_view_colors_depth 1,1,0
484: test:
485: suffix: ref_tetrahedron
486: args: -dm_plex_cell tetrahedron \
487: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
489: test:
490: suffix: ref_hexahedron
491: args: -dm_plex_cell hexahedron \
492: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
494: test:
495: suffix: ref_triangular_prism
496: args: -dm_plex_cell triangular_prism \
497: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
499: test:
500: suffix: ref_tensor_triangular_prism
501: args: -dm_plex_cell tensor_triangular_prism \
502: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
504: test:
505: suffix: ref_tensor_quadrilateral_prism
506: args: -dm_plex_cell tensor_quadrilateral_prism \
507: -dm_plex_view_numbers_depth 1,0,0,0 -dm_plex_view_colors_depth 1,0,0,0
509: # ToBox should recreate the coordinate space since the cell shape changes
510: testset:
511: args: -dm_coord_space 0 -dm_plex_transform_type refine_tobox -ref_arrangements -dm_plex_check_all
513: test:
514: suffix: tobox_triangle
515: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \
516: -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
518: test:
519: suffix: tobox_tensor_segment
520: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \
521: -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
523: test:
524: suffix: tobox_tetrahedron
525: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \
526: -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
528: test:
529: suffix: tobox_triangular_prism
530: args: -dm_plex_reference_cell_domain -dm_plex_cell triangular_prism \
531: -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
533: test:
534: suffix: tobox_tensor_triangular_prism
535: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \
536: -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
538: test:
539: suffix: tobox_tensor_quadrilateral_prism
540: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quadrilateral_prism \
541: -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
543: testset:
544: args: -dm_coord_space 0 -dm_plex_transform_type refine_alfeld -ref_arrangements -dm_plex_check_all
546: test:
547: suffix: alfeld_triangle
548: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle \
549: -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
551: test:
552: suffix: alfeld_tetrahedron
553: args: -dm_plex_reference_cell_domain -dm_plex_cell tetrahedron \
554: -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
556: testset:
557: args: -dm_plex_transform_type refine_sbr -ref_arrangements -dm_plex_check_all
559: # This splits edge 1 of the triangle, and reflects about the added edge
560: test:
561: suffix: sbr_triangle_0
562: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -ornts -2,0 \
563: -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
565: # This splits edge 0 of the triangle, and reflects about the added edge
566: test:
567: suffix: sbr_triangle_1
568: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 1 -ornts -3,0 \
569: -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
571: # This splits edge 2 of the triangle, and reflects about the added edge
572: test:
573: suffix: sbr_triangle_2
574: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5 -init_ornt 2 -ornts -1,0 \
575: -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
577: # This splits edges 1 and 2 of the triangle
578: test:
579: suffix: sbr_triangle_3
580: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -ornts 0 \
581: -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
583: # This splits edges 1 and 0 of the triangle
584: test:
585: suffix: sbr_triangle_4
586: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -ornts 0 \
587: -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
589: # This splits edges 0 and 1 of the triangle
590: test:
591: suffix: sbr_triangle_5
592: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 1 -ornts 0 \
593: -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
595: # This splits edges 0 and 2 of the triangle
596: test:
597: suffix: sbr_triangle_6
598: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 1 -ornts 0 \
599: -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
601: # This splits edges 2 and 0 of the triangle
602: test:
603: suffix: sbr_triangle_7
604: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 5,6 -init_ornt 2 -ornts 0 \
605: -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
607: # This splits edges 2 and 1 of the triangle
608: test:
609: suffix: sbr_triangle_8
610: args: -dm_plex_reference_cell_domain -dm_plex_cell triangle -dm_plex_transform_sbr_ref_cell 4,5 -init_ornt 2 -ornts 0 \
611: -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
613: testset:
614: args: -dm_plex_transform_type refine_boundary_layer -dm_plex_transform_bl_splits 2 -ref_arrangements -dm_plex_check_all
616: test:
617: suffix: bl_tensor_segment
618: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_quad \
619: -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
621: # The subcell check is broken because at orientation 3, the internal triangles do not get properly permuted for the check
622: test:
623: suffix: bl_tensor_triangular_prism
624: requires: TODO
625: args: -dm_plex_reference_cell_domain -dm_plex_cell tensor_triangular_prism \
626: -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
628: TEST*/