xref: /petsc/src/dm/impls/plex/tests/ex7.c (revision b698fc57f0bea7237255b29c1b77df0acc362ffd)
1 static char help[] = "Tests for mesh interpolation\n\n";
2 
3 /* TODO
4 */
5 
6 #include <petscdmplex.h>
7 
8 /* List of test meshes
9 
10 Triangle
11 --------
12 Test 0:
13 Two triangles sharing a face
14 
15         4
16       / | \
17      /  |  \
18     /   |   \
19    2  0 | 1  5
20     \   |   /
21      \  |  /
22       \ | /
23         3
24 
25 should become
26 
27         4
28       / | \
29      8  |  9
30     /   |   \
31    2  0 7 1  5
32     \   |   /
33      6  |  10
34       \ | /
35         3
36 
37 Tetrahedron
38 -----------
39 Test 0:
40 Two tets sharing a face
41 
42  cell   5 _______    cell
43  0    / | \      \       1
44      /  |  \      \
45     /   |   \      \
46    2----|----4-----6
47     \   |   /      /
48      \  |  /     /
49       \ | /      /
50         3-------
51 
52 should become
53 
54  cell   5 _______    cell
55  0    / | \      \       1
56      /  |  \      \
57    17   |   18 13  22
58    / 8 19 10 \      \
59   /     |     \      \
60  2---14-|------4--20--6
61   \     |     /      /
62    \ 9  | 7  /      /
63    16   |   15 11  21
64      \  |  /      /
65       \ | /      /
66         3-------
67 
68 
69 Quadrilateral
70 -------------
71 Test 0:
72 Two quads sharing a face
73 
74    5-------4-------7
75    |       |       |
76    |   0   |   1   |
77    |       |       |
78    2-------3-------6
79 
80 should become
81 
82    5--10---4--14---7
83    |       |       |
84   11   0   9   1  13
85    |       |       |
86    2---8---3--12---6
87 
88 Test 1:
89 A quad and a triangle sharing a face
90 
91    5-------4
92    |       | \
93    |   0   |  \
94    |       | 1 \
95    2-------3----6
96 
97 should become
98 
99    5---9---4
100    |       | \
101   10   0   8  12
102    |       | 1 \
103    2---7---3-11-6
104 
105 Hexahedron
106 ----------
107 Test 0:
108 Two hexes sharing a face
109 
110 cell   9-------------8-------------13 cell
111 0     /|            /|            /|     1
112      / |   15      / |   21      / |
113     /  |          /  |          /  |
114    6-------------7-------------12  |
115    |   |     18  |   |     24  |   |
116    |   |         |   |         |   |
117    |19 |         |17 |         |23 |
118    |   |  16     |   |   22    |   |
119    |   5---------|---4---------|---11
120    |  /          |  /          |  /
121    | /   14      | /    20     | /
122    |/            |/            |/
123    2-------------3-------------10
124 
125 should become
126 
127 cell   9-----31------8-----42------13 cell
128 0     /|            /|            /|     1
129     32 |   15      30|   21      41|
130     /  |          /  |          /  |
131    6-----29------7-----40------12  |
132    |   |     17  |   |     23  |   |
133    |  35         |  36         |   44
134    |19 |         |18 |         |24 |
135   34   |  16    33   |   22   43   |
136    |   5-----26--|---4-----37--|---11
137    |  /          |  /          |  /
138    | 25   14     | 27    20    | 38
139    |/            |/            |/
140    2-----28------3-----39------10
141 
142 */
143 
144 typedef struct {
145   DM        dm;
146   PetscInt  testNum;                      /* Indicates the mesh to create */
147   PetscInt  dim;                          /* The topological mesh dimension */
148   PetscBool cellSimplex;                  /* Use simplices or hexes */
149   PetscBool useGenerator;                 /* Construct mesh with a mesh generator */
150   PetscBool refCell;                      /* Construct reference cell by type */
151   PetscBool femGeometry;                  /* Flag to compute FEM geometry */
152   char      filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
153 } AppCtx;
154 
155 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
156 {
157   PetscErrorCode ierr;
158 
159   PetscFunctionBegin;
160   options->testNum      = 0;
161   options->dim          = 2;
162   options->cellSimplex  = PETSC_TRUE;
163   options->useGenerator = PETSC_FALSE;
164   options->refCell      = PETSC_FALSE;
165   options->femGeometry  = PETSC_TRUE;
166   options->filename[0]  = '\0';
167 
168   ierr = PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");CHKERRQ(ierr);
169   ierr = PetscOptionsBoundedInt("-testnum", "The mesh to create", "ex7.c", options->testNum, &options->testNum, NULL,0);CHKERRQ(ierr);
170   ierr = PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex7.c", options->dim, &options->dim, NULL,1,3);CHKERRQ(ierr);
171   ierr = PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex7.c", options->cellSimplex, &options->cellSimplex, NULL);CHKERRQ(ierr);
172   ierr = PetscOptionsBool("-use_generator", "Use a mesh generator to build the mesh", "ex7.c", options->useGenerator, &options->useGenerator, NULL);CHKERRQ(ierr);
173   ierr = PetscOptionsBool("-ref_cell", "Create a reference cell", "ex7.c", options->refCell, &options->refCell, NULL);CHKERRQ(ierr);
174   ierr = PetscOptionsBool("-fem_geometry", "Flag to compute FEM geometry", "ex7.c", options->femGeometry, &options->femGeometry, NULL);CHKERRQ(ierr);
175   ierr = PetscOptionsString("-filename", "The mesh file", "ex7.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr);
176   ierr = PetscOptionsEnd();
177   PetscFunctionReturn(0);
178 }
179 
180 PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM dm)
181 {
182   PetscInt       depth = 1, testNum  = 0, p;
183   PetscMPIInt    rank;
184   PetscErrorCode ierr;
185 
186   PetscFunctionBegin;
187   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
188   if (!rank) {
189     switch (testNum) {
190     case 0:
191     {
192       PetscInt    numPoints[2]        = {4, 2};
193       PetscInt    coneSize[6]         = {3, 3, 0, 0, 0, 0};
194       PetscInt    cones[6]            = {2, 3, 4,  5, 4, 3};
195       PetscInt    coneOrientations[6] = {0, 0, 0,  0, 0, 0};
196       PetscScalar vertexCoords[8]     = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5};
197       PetscInt    markerPoints[8]     = {2, 1, 3, 1, 4, 1, 5, 1};
198 
199       ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
200       for (p = 0; p < 4; ++p) {
201         ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);
202       }
203     }
204     break;
205     default:
206       SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum);
207     }
208   } else {
209     PetscInt numPoints[2] = {0, 0};
210 
211     ierr = DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr);
212   }
213   PetscFunctionReturn(0);
214 }
215 
216 PetscErrorCode CreateSimplex_3D(MPI_Comm comm, DM dm)
217 {
218   PetscInt       depth = 1, testNum  = 0, p;
219   PetscMPIInt    rank;
220   PetscErrorCode ierr;
221 
222   PetscFunctionBegin;
223   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
224   if (!rank) {
225     switch (testNum) {
226     case 0:
227     {
228       PetscInt    numPoints[2]        = {5, 2};
229       PetscInt    coneSize[23]        = {4, 4, 0, 0, 0, 0, 0};
230       PetscInt    cones[8]            = {2, 4, 3, 5,  3, 4, 6, 5};
231       PetscInt    coneOrientations[8] = {0, 0, 0, 0,  0, 0, 0, 0};
232       PetscScalar vertexCoords[15]    = {0.0, 0.0, -0.5,  0.0, -0.5, 0.0,  1.0, 0.0, 0.0,  0.0, 0.5, 0.0,  0.0, 0.0, 0.5};
233       PetscInt    markerPoints[8]     = {2, 1, 3, 1, 4, 1, 5, 1};
234 
235       ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
236       for (p = 0; p < 4; ++p) {
237         ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);
238       }
239     }
240     break;
241     default:
242       SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum);
243     }
244   } else {
245     PetscInt numPoints[2] = {0, 0};
246 
247     ierr = DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr);
248   }
249   PetscFunctionReturn(0);
250 }
251 
252 PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscInt testNum, DM dm)
253 {
254   PetscInt       depth = 1, p;
255   PetscMPIInt    rank;
256   PetscErrorCode ierr;
257 
258   PetscFunctionBegin;
259   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
260   if (!rank) {
261     switch (testNum) {
262     case 0:
263     {
264       PetscInt    numPoints[2]        = {6, 2};
265       PetscInt    coneSize[8]         = {4, 4, 0, 0, 0, 0, 0, 0};
266       PetscInt    cones[8]            = {2, 3, 4, 5,  3, 6, 7, 4};
267       PetscInt    coneOrientations[8] = {0, 0, 0, 0,  0, 0, 0, 0};
268       PetscScalar vertexCoords[12]    = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0, 0.5, 1.0};
269       PetscInt    markerPoints[12]    = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1};
270 
271       ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
272       for (p = 0; p < 6; ++p) {
273         ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);
274       }
275     }
276     break;
277     case 1:
278     {
279       PetscInt    numPoints[2]        = {5, 2};
280       PetscInt    coneSize[7]         = {4, 3, 0, 0, 0, 0, 0};
281       PetscInt    cones[7]            = {2, 3, 4, 5,  3, 6, 4};
282       PetscInt    coneOrientations[7] = {0, 0, 0, 0,  0, 0, 0};
283       PetscScalar vertexCoords[10]    = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0};
284       PetscInt    markerPoints[10]    = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};
285 
286       ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
287       for (p = 0; p < 5; ++p) {
288         ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);
289       }
290     }
291     break;
292     default:
293       SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum);
294     }
295   } else {
296     PetscInt numPoints[2] = {0, 0};
297 
298     ierr = DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr);
299   }
300   PetscFunctionReturn(0);
301 }
302 
303 PetscErrorCode CreateHex_3D(MPI_Comm comm, DM dm)
304 {
305   PetscInt       depth = 1, testNum  = 0, p;
306   PetscMPIInt    rank;
307   PetscErrorCode ierr;
308 
309   PetscFunctionBegin;
310   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
311   if (!rank) {
312     switch (testNum) {
313     case 0:
314     {
315       PetscInt    numPoints[2]         = {12, 2};
316       PetscInt    coneSize[14]         = {8, 8, 0,0,0,0,0,0,0,0,0,0,0,0};
317       PetscInt    cones[16]            = {2,5,4,3,6,7,8,9,  3,4,11,10,7,12,13,8};
318       PetscInt    coneOrientations[16] = {0,0,0,0,0,0,0,0,  0,0, 0, 0,0, 0, 0,0};
319       PetscScalar vertexCoords[36]     = {-0.5,0.0,0.0, 0.0,0.0,0.0, 0.0,1.0,0.0, -0.5,1.0,0.0,
320                                           -0.5,0.0,1.0, 0.0,0.0,1.0, 0.0,1.0,1.0, -0.5,1.0,1.0,
321                                            0.5,0.0,0.0, 0.5,1.0,0.0, 0.5,0.0,1.0,  0.5,1.0,1.0};
322       PetscInt    markerPoints[16]     = {2,1,3,1,4,1,5,1,6,1,7,1,8,1,9,1};
323 
324       ierr = DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);CHKERRQ(ierr);
325       for (p = 0; p < 8; ++p) {
326         ierr = DMSetLabelValue(dm, "marker", markerPoints[p*2], markerPoints[p*2+1]);CHKERRQ(ierr);
327       }
328     }
329     break;
330     default:
331       SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %d", testNum);
332     }
333   } else {
334     PetscInt numPoints[2] = {0, 0};
335 
336     ierr = DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);CHKERRQ(ierr);
337   }
338   PetscFunctionReturn(0);
339 }
340 
341 PetscErrorCode CheckMesh(DM dm, AppCtx *user)
342 {
343   PetscReal      detJ, J[9], refVol = 1.0;
344   PetscReal      vol;
345   PetscInt       dim, depth, d, cStart, cEnd, c;
346   PetscErrorCode ierr;
347 
348   PetscFunctionBegin;
349   ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr);
350   ierr = DMPlexGetDepth(dm, &depth);CHKERRQ(ierr);
351   for (d = 0; d < dim; ++d) {
352     refVol *= 2.0;
353   }
354   ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr);
355   for (c = cStart; c < cEnd; ++c) {
356     if (user->femGeometry) {
357       ierr = DMPlexComputeCellGeometryFEM(dm, c, NULL, NULL, J, NULL, &detJ);CHKERRQ(ierr);
358       if (detJ <= 0.0) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mesh cell %D is inverted, |J| = %g", c, (double) detJ);
359       ierr = PetscInfo1(dm, "FEM Volume: %g\n", (double) detJ*refVol);CHKERRQ(ierr);
360     }
361     if (depth > 1) {
362       ierr = DMPlexComputeCellGeometryFVM(dm, c, &vol, NULL, NULL);CHKERRQ(ierr);
363       if (vol <= 0.0) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mesh cell %D is inverted, vol = %g", c, (double) vol);
364       ierr = PetscInfo1(dm, "FVM Volume: %g\n", (double) vol);CHKERRQ(ierr);
365     }
366   }
367   PetscFunctionReturn(0);
368 }
369 
370 PetscErrorCode CompareCones(DM dm, DM idm)
371 {
372   PetscInt       cStart, cEnd, c, vStart, vEnd, v;
373   PetscErrorCode ierr;
374 
375   PetscFunctionBegin;
376   ierr = DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);CHKERRQ(ierr);
377   ierr = DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);CHKERRQ(ierr);
378   for (c = cStart; c < cEnd; ++c) {
379     const PetscInt *cone;
380     PetscInt       *points = NULL, numPoints, p, numVertices = 0, coneSize;
381 
382     ierr = DMPlexGetCone(dm, c, &cone);CHKERRQ(ierr);
383     ierr = DMPlexGetConeSize(dm, c, &coneSize);CHKERRQ(ierr);
384     ierr = DMPlexGetTransitiveClosure(idm, c, PETSC_TRUE, &numPoints, &points);CHKERRQ(ierr);
385     for (p = 0; p < numPoints*2; p += 2) {
386       const PetscInt point = points[p];
387       if ((point >= vStart) && (point < vEnd)) points[numVertices++] = point;
388     }
389     if (numVertices != coneSize) SETERRQ3(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "In cell %d, cone size %d != %d vertices in closure", c, coneSize, numVertices);
390     for (v = 0; v < numVertices; ++v) {
391       if (cone[v] != points[v]) SETERRQ4(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "In cell %d, cone point %d is %d != %d vertex in closure", c, v, cone[v], points[v]);
392     }
393     ierr = DMPlexRestoreTransitiveClosure(idm, c, PETSC_TRUE, &numPoints, &points);CHKERRQ(ierr);
394   }
395   PetscFunctionReturn(0);
396 }
397 
398 PetscErrorCode CreateMesh(MPI_Comm comm, PetscInt testNum, AppCtx *user, DM *dm)
399 {
400   PetscInt       dim          = user->dim;
401   PetscBool      cellSimplex  = user->cellSimplex;
402   PetscBool      useGenerator = user->useGenerator;
403   const char    *filename     = user->filename;
404   size_t         len;
405   PetscMPIInt    rank;
406   PetscErrorCode ierr;
407 
408   PetscFunctionBegin;
409   ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr);
410   ierr = PetscStrlen(filename, &len);CHKERRQ(ierr);
411   if (len) {
412     ierr = DMPlexCreateFromFile(comm, filename, PETSC_FALSE, dm);CHKERRQ(ierr);
413     ierr = DMGetDimension(*dm, &dim);CHKERRQ(ierr);
414   } else if (useGenerator) {
415     ierr = DMPlexCreateBoxMesh(comm, dim, cellSimplex, NULL, NULL, NULL, NULL, PETSC_FALSE, dm);CHKERRQ(ierr);
416   } else if (user->refCell) {
417     ierr = DMPlexCreateReferenceCellByType(comm, DM_POLYTOPE_UNKNOWN, dm);CHKERRQ(ierr);
418   } else {
419     ierr = DMCreate(comm, dm);CHKERRQ(ierr);
420     ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr);
421     ierr = DMSetDimension(*dm, dim);CHKERRQ(ierr);
422     switch (dim) {
423     case 2:
424       if (cellSimplex) {
425         ierr = CreateSimplex_2D(comm, *dm);CHKERRQ(ierr);
426       } else {
427         ierr = CreateQuad_2D(comm, testNum, *dm);CHKERRQ(ierr);
428       }
429       break;
430     case 3:
431       if (cellSimplex) {
432         ierr = CreateSimplex_3D(comm, *dm);CHKERRQ(ierr);
433       } else {
434         ierr = CreateHex_3D(comm, *dm);CHKERRQ(ierr);
435       }
436       break;
437     default:
438       SETERRQ1(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %d", dim);
439     }
440   }
441   {
442     DMPlexInterpolatedFlag interpolated;
443 
444     ierr = CheckMesh(*dm, user);CHKERRQ(ierr);
445     ierr = DMPlexIsInterpolated(*dm, &interpolated);CHKERRQ(ierr);
446     if (interpolated == DMPLEX_INTERPOLATED_FULL) {
447       DM udm;
448 
449       ierr = DMPlexUninterpolate(*dm, &udm);CHKERRQ(ierr);
450       ierr = CompareCones(udm, *dm);CHKERRQ(ierr);
451       ierr = DMDestroy(&udm);CHKERRQ(ierr);
452     } else {
453       DM idm;
454 
455       ierr = DMPlexInterpolate(*dm, &idm);CHKERRQ(ierr);
456       ierr = CompareCones(*dm, idm);CHKERRQ(ierr);
457       ierr = DMDestroy(dm);CHKERRQ(ierr);
458       *dm  = idm;
459     }
460   }
461   {
462     DM               distributedMesh = NULL;
463     PetscPartitioner part;
464 
465     ierr = DMPlexGetPartitioner(*dm, &part);CHKERRQ(ierr);
466     ierr = PetscPartitionerSetFromOptions(part);CHKERRQ(ierr);
467     /* Distribute mesh over processes */
468     ierr = DMPlexDistribute(*dm, 0, NULL, &distributedMesh);CHKERRQ(ierr);
469     if (distributedMesh) {
470       ierr = DMViewFromOptions(distributedMesh, NULL, "-dm_view");CHKERRQ(ierr);
471       ierr = DMDestroy(dm);CHKERRQ(ierr);
472       *dm  = distributedMesh;
473     }
474   }
475   ierr = PetscObjectSetName((PetscObject) *dm, "Interpolated Mesh");CHKERRQ(ierr);
476   ierr = DMSetFromOptions(*dm);CHKERRQ(ierr);
477   ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr);
478   user->dm = *dm;
479   PetscFunctionReturn(0);
480 }
481 
482 int main(int argc, char **argv)
483 {
484   AppCtx         user;                 /* user-defined work context */
485   PetscErrorCode ierr;
486 
487   ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
488   ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr);
489   ierr = CreateMesh(PETSC_COMM_WORLD, user.testNum, &user, &user.dm);CHKERRQ(ierr);
490   ierr = CheckMesh(user.dm, &user);CHKERRQ(ierr);
491   ierr = DMDestroy(&user.dm);CHKERRQ(ierr);
492   ierr = PetscFinalize();
493   return ierr;
494 }
495 
496 /*TEST
497 
498   # Two cell test meshes 0-7
499   test:
500     suffix: 0
501     args: -dim 2 -dm_view ascii::ascii_info_detail
502   test:
503     suffix: 1
504     nsize: 2
505     args: -petscpartitioner_type simple -dim 2 -dm_view ascii::ascii_info_detail
506   test:
507     suffix: 2
508     args: -dim 2 -cell_simplex 0 -dm_view ascii::ascii_info_detail
509   test:
510     suffix: 3
511     nsize: 2
512     args: -petscpartitioner_type simple -dim 2 -cell_simplex 0 -dm_view ascii::ascii_info_detail
513   test:
514     suffix: 4
515     args: -dim 3 -dm_view ascii::ascii_info_detail
516   test:
517     suffix: 5
518     nsize: 2
519     args: -petscpartitioner_type simple -dim 3 -dm_view ascii::ascii_info_detail
520   test:
521     suffix: 6
522     args: -dim 3 -cell_simplex 0 -dm_view ascii::ascii_info_detail
523   test:
524     suffix: 7
525     nsize: 2
526     args: -petscpartitioner_type simple -dim 3 -cell_simplex 0 -dm_view ascii::ascii_info_detail
527   # 2D Hybrid Mesh 8
528   test:
529     suffix: 8
530     args: -dim 2 -cell_simplex 0 -testnum 1 -dm_view ascii::ascii_info_detail
531   # Reference cells
532   test:
533     suffix: 12
534     args: -ref_cell -dm_plex_ref_type pyramid -fem_geometry 0 -dm_plex_check_symmetry -dm_plex_check_skeleton -dm_plex_check_faces
535   # TetGen meshes 9-10
536   test:
537     suffix: 9
538     requires: triangle
539     args: -dim 2 -use_generator -dm_view ascii::ascii_info_detail
540   test:
541     suffix: 10
542     requires: ctetgen
543     args: -dim 3 -use_generator -dm_view ascii::ascii_info_detail
544   # Cubit meshes 11
545   test:
546     suffix: 11
547     requires: exodusii
548     args: -dim 3 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -dm_view ascii::ascii_info_detail
549 
550 TEST*/
551