Actual source code: ex4.c

  1: static char help[] = "Tests for refinement of meshes created by hand\n\n";

  3: #include <petscdmplex.h>

  5: typedef struct {
  6:   PetscInt  debug;         /* The debugging level */
  7:   PetscInt  dim;           /* The topological mesh dimension */
  8:   PetscBool cellHybrid;    /* Use a hybrid mesh */
  9:   PetscBool cellSimplex;   /* Use simplices or hexes */
 10:   PetscBool testPartition; /* Use a fixed partitioning for testing */
 11:   PetscInt  testNum;       /* The particular mesh to test */
 12:   PetscBool uninterpolate; /* Uninterpolate the mesh at the end */
 13:   PetscBool reinterpolate; /* Reinterpolate the mesh at the end */
 14: } AppCtx;

 16: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 17: {
 18:   options->debug         = 0;
 19:   options->dim           = 2;
 20:   options->cellHybrid    = PETSC_TRUE;
 21:   options->cellSimplex   = PETSC_TRUE;
 22:   options->testPartition = PETSC_TRUE;
 23:   options->testNum       = 0;
 24:   options->uninterpolate = PETSC_FALSE;
 25:   options->reinterpolate = PETSC_FALSE;

 27:   PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");
 28:   PetscOptionsBoundedInt("-debug", "The debugging level", "ex4.c", options->debug, &options->debug, NULL, 0);
 29:   PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex4.c", options->dim, &options->dim, NULL, 1, 3);
 30:   PetscOptionsBool("-cell_hybrid", "Use a hybrid mesh", "ex4.c", options->cellHybrid, &options->cellHybrid, NULL);
 31:   PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex4.c", options->cellSimplex, &options->cellSimplex, NULL);
 32:   PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex4.c", options->testPartition, &options->testPartition, NULL);
 33:   PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex4.c", options->testNum, &options->testNum, NULL, 0);
 34:   PetscOptionsBool("-uninterpolate", "Uninterpolate the mesh at the end", "ex4.c", options->uninterpolate, &options->uninterpolate, NULL);
 35:   PetscOptionsBool("-reinterpolate", "Reinterpolate the mesh at the end", "ex4.c", options->reinterpolate, &options->reinterpolate, NULL);
 36:   PetscOptionsEnd();
 37:   return 0;
 38: }

 40: /* Two segments

 42:   2-------0-------3-------1-------4

 44: become

 46:   4---0---7---1---5---2---8---3---6

 48: */
 49: PetscErrorCode CreateSimplex_1D(MPI_Comm comm, DM *dm)
 50: {
 51:   PetscInt    depth = 1;
 52:   PetscMPIInt rank;

 54:   MPI_Comm_rank(comm, &rank);
 55:   if (rank == 0) {
 56:     PetscInt    numPoints[2]         = {3, 2};
 57:     PetscInt    coneSize[5]          = {2, 2, 0, 0, 0};
 58:     PetscInt    cones[4]             = {2, 3, 3, 4};
 59:     PetscInt    coneOrientations[16] = {0, 0, 0, 0};
 60:     PetscScalar vertexCoords[3]      = {-1.0, 0.0, 1.0};

 62:     DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
 63:   } else {
 64:     PetscInt numPoints[2] = {0, 0};

 66:     DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
 67:   }
 68:   return 0;
 69: }

 71: /* Two triangles
 72:         4
 73:       / | \
 74:      8  |  10
 75:     /   |   \
 76:    2  0 7  1 5
 77:     \   |   /
 78:      6  |  9
 79:       \ | /
 80:         3

 82: Becomes
 83:            10
 84:           / | \
 85:         21  |  26
 86:         /   |   \
 87:       14 2 20 4  16
 88:       /|\   |   /|\
 89:     22 | 28 | 32 | 25
 90:     /  |  \ | /  | 6\
 91:    8  29 3 13  7 31  11
 92:     \0 |  / | \  |  /
 93:     17 | 27 | 30 | 24
 94:       \|/   |   \|/
 95:       12 1 19 5  15
 96:         \   |   /
 97:         18  |  23
 98:           \ | /
 99:             9
100: */
101: PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM *dm)
102: {
103:   PetscInt    depth = 2;
104:   PetscMPIInt rank;

106:   MPI_Comm_rank(comm, &rank);
107:   if (rank == 0) {
108:     PetscInt    numPoints[3]         = {4, 5, 2};
109:     PetscInt    coneSize[11]         = {3, 3, 0, 0, 0, 0, 2, 2, 2, 2, 2};
110:     PetscInt    cones[16]            = {6, 7, 8, 7, 9, 10, 2, 3, 3, 4, 4, 2, 3, 5, 5, 4};
111:     PetscInt    coneOrientations[16] = {0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
112:     PetscScalar vertexCoords[8]      = {-0.5, 0.0, 0.0, -0.5, 0.0, 0.5, 0.5, 0.0};

114:     DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
115:   } else {
116:     PetscInt numPoints[3] = {0, 0, 0};

118:     DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
119:   }
120:   return 0;
121: }

123: /* Two triangles separated by a zero-volume cell with 4 vertices/2 edges
124:         5--16--8
125:       / |      | \
126:     11  |      |  12
127:     /   |      |   \
128:    3  0 10  2 14 1  6
129:     \   |      |   /
130:      9  |      |  13
131:       \ |      | /
132:         4--15--7
133: */
134: PetscErrorCode CreateSimplexHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
135: {
136:   DM          idm, hdm = NULL;
137:   DMLabel     faultLabel, hybridLabel;
138:   PetscInt    p;
139:   PetscMPIInt rank;

141:   MPI_Comm_rank(comm, &rank);
142:   if (rank == 0) {
143:     switch (testNum) {
144:     case 0: {
145:       PetscInt    numPoints[2]        = {4, 2};
146:       PetscInt    coneSize[6]         = {3, 3, 0, 0, 0, 0};
147:       PetscInt    cones[6]            = {2, 3, 4, 5, 4, 3};
148:       PetscInt    coneOrientations[6] = {0, 0, 0, 0, 0, 0};
149:       PetscScalar vertexCoords[8]     = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, 1.0, 0.5};
150:       PetscInt    faultPoints[2]      = {3, 4};

152:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
153:       for (p = 0; p < 2; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
154:     } break;
155:     case 1: {
156:       PetscInt    numPoints[2]         = {5, 4};
157:       PetscInt    coneSize[9]          = {3, 3, 3, 3, 0, 0, 0, 0, 0};
158:       PetscInt    cones[12]            = {4, 5, 6, 6, 7, 4, 6, 5, 8, 6, 8, 7};
159:       PetscInt    coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
160:       PetscScalar vertexCoords[10]     = {-1.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0};
161:       PetscInt    faultPoints[3]       = {5, 6, 7};

163:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
164:       for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
165:     } break;
166:     default:
167:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
168:     }
169:     DMPlexInterpolate(*dm, &idm);
170:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
171:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
172:     DMSetFromOptions(idm);
173:     DMViewFromOptions(idm, NULL, "-dm_view");
174:     DMGetLabel(*dm, "fault", &faultLabel);
175:     DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
176:     DMLabelDestroy(&hybridLabel);
177:     DMDestroy(&idm);
178:     DMDestroy(dm);
179:     *dm = hdm;
180:   } else {
181:     PetscInt numPoints[2] = {0, 0};

183:     DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
184:     DMPlexInterpolate(*dm, &idm);
185:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
186:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
187:     DMSetFromOptions(idm);
188:     DMViewFromOptions(idm, NULL, "-dm_view");
189:     DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
190:     DMDestroy(&idm);
191:     DMDestroy(dm);
192:     *dm = hdm;
193:   }
194:   return 0;
195: }

197: /* Two quadrilaterals

199:   5----10-----4----14-----7
200:   |           |           |
201:   |           |           |
202:   |           |           |
203:  11     0     9     1     13
204:   |           |           |
205:   |           |           |
206:   |           |           |
207:   2-----8-----3----12-----6
208: */
209: PetscErrorCode CreateTensorProduct_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
210: {
211:   PetscInt    depth = 2;
212:   PetscMPIInt rank;

214:   MPI_Comm_rank(comm, &rank);
215:   if (rank == 0) {
216:     PetscInt    numPoints[3]         = {6, 7, 2};
217:     PetscInt    coneSize[15]         = {4, 4, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2};
218:     PetscInt    cones[22]            = {8, 9, 10, 11, 12, 13, 14, 9, 2, 3, 3, 4, 4, 5, 5, 2, 3, 6, 6, 7, 7, 4};
219:     PetscInt    coneOrientations[22] = {0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
220:     PetscScalar vertexCoords[12]     = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5};

222:     DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
223:   } else {
224:     PetscInt numPoints[3] = {0, 0, 0};

226:     DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
227:   }
228:   return 0;
229: }

231: PetscErrorCode CreateTensorProductHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
232: {
233:   DM          idm, hdm = NULL;
234:   DMLabel     faultLabel, hybridLabel;
235:   PetscInt    p;
236:   PetscMPIInt rank;

238:   MPI_Comm_rank(comm, &rank);
239:   if (rank == 0) {
240:     PetscInt numPoints[2] = {6, 2};
241:     PetscInt coneSize[8]  = {4, 4, 0, 0, 0, 0, 0, 0};
242:     PetscInt cones[8]     = {
243:       2, 3, 4, 5, 3, 6, 7, 4,
244:     };
245:     PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
246:     PetscScalar vertexCoords[12]    = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5};
247:     PetscInt    faultPoints[2]      = {3, 4};

249:     DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
250:     for (p = 0; p < 2; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
251:     DMPlexInterpolate(*dm, &idm);
252:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
253:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
254:     DMSetFromOptions(idm);
255:     DMViewFromOptions(idm, NULL, "-dm_view");
256:     DMGetLabel(*dm, "fault", &faultLabel);
257:     DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
258:     DMLabelDestroy(&hybridLabel);
259:   } else {
260:     PetscInt numPoints[3] = {0, 0, 0};

262:     DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
263:     DMPlexInterpolate(*dm, &idm);
264:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
265:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
266:     DMSetFromOptions(idm);
267:     DMViewFromOptions(idm, NULL, "-dm_view");
268:     DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
269:   }
270:   DMDestroy(&idm);
271:   DMDestroy(dm);
272:   *dm = hdm;
273:   return 0;
274: }

276: /* Two tetrahedrons

278:  cell   5          5______    cell
279:  0    / | \        |\      \     1
280:     17  |  18      | 18 13  21
281:     /8 19 10\     19  \      \
282:    2-14-|----4     |   4--22--6
283:     \ 9 | 7 /      |10 /      /
284:     16  |  15      | 15  12 20
285:       \ | /        |/      /
286:         3          3------
287: */
288: PetscErrorCode CreateSimplex_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
289: {
290:   DM          idm;
291:   PetscInt    depth = 3;
292:   PetscMPIInt rank;

294:   MPI_Comm_rank(comm, &rank);
295:   if (rank == 0) {
296:     switch (testNum) {
297:     case 0: {
298:       PetscInt    numPoints[4]         = {5, 9, 7, 2};
299:       PetscInt    coneSize[23]         = {4, 4, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2};
300:       PetscInt    cones[47]            = {7, 8, 9, 10, 10, 11, 12, 13, 14, 15, 16, 17, 18, 14, 16, 19, 17, 15, 18, 19, 20, 21, 19, 15, 22, 20, 18, 21, 22, 2, 4, 4, 3, 3, 2, 2, 5, 5, 4, 3, 5, 3, 6, 6, 5, 4, 6};
301:       PetscInt    coneOrientations[47] = {0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
302:       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};

304:       DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
305:     } break;
306:     case 1: {
307:       PetscInt    numPoints[2]        = {5, 2};
308:       PetscInt    coneSize[7]         = {4, 4, 0, 0, 0, 0, 0};
309:       PetscInt    cones[8]            = {4, 3, 5, 2, 5, 3, 4, 6};
310:       PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
311:       PetscScalar vertexCoords[15]    = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0};

313:       depth = 1;
314:       DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
315:       DMPlexInterpolate(*dm, &idm);
316:       DMViewFromOptions(idm, NULL, "-in_dm_view");
317:       DMDestroy(dm);
318:       *dm = idm;
319:     } break;
320:     case 2: {
321:       PetscInt    numPoints[2]        = {4, 1};
322:       PetscInt    coneSize[5]         = {4, 0, 0, 0, 0};
323:       PetscInt    cones[4]            = {2, 3, 4, 1};
324:       PetscInt    coneOrientations[4] = {0, 0, 0, 0};
325:       PetscScalar vertexCoords[12]    = {0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0};

327:       depth = 1;
328:       DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
329:       DMPlexInterpolate(*dm, &idm);
330:       DMViewFromOptions(idm, NULL, "-in_dm_view");
331:       DMDestroy(dm);
332:       *dm = idm;
333:     } break;
334:     default:
335:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
336:     }
337:   } else {
338:     PetscInt numPoints[4] = {0, 0, 0, 0};

340:     DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
341:     switch (testNum) {
342:     case 1:
343:       DMPlexInterpolate(*dm, &idm);
344:       DMViewFromOptions(idm, NULL, "-in_dm_view");
345:       DMDestroy(dm);
346:       *dm = idm;
347:       break;
348:     }
349:   }
350:   return 0;
351: }

353: /* Two tetrahedrons separated by a zero-volume cell with 6 vertices

355:  cell   6 ___33___10______    cell
356:  0    / | \        |\      \     1
357:     21  |  23      | 29     27
358:     /12 24 14\    30  \      \
359:    3-20-|----5--32-|---9--26--7
360:     \ 13| 11/      |18 /      /
361:     19  |  22      | 28     25
362:       \ | /        |/      /
363:         4----31----8------
364:          cell 2
365: */
366: PetscErrorCode CreateSimplexHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
367: {
368:   DM          idm, hdm = NULL;
369:   DMLabel     faultLabel, hybridLabel;
370:   PetscInt    p;
371:   PetscMPIInt rank;

373:   MPI_Comm_rank(comm, &rank);
374:   if (rank == 0) {
375:     switch (testNum) {
376:     case 0: {
377:       PetscInt    numPoints[2]        = {5, 2};
378:       PetscInt    coneSize[7]         = {4, 4, 0, 0, 0, 0, 0};
379:       PetscInt    cones[8]            = {4, 3, 5, 2, 5, 3, 4, 6};
380:       PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
381:       PetscScalar vertexCoords[15]    = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0};
382:       PetscInt    faultPoints[3]      = {3, 4, 5};

384:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
385:       for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
386:     } break;
387:     case 1: {
388:       /* Tets 0,3,5 and 1,2,4 */
389:       PetscInt    numPoints[2]         = {9, 6};
390:       PetscInt    coneSize[15]         = {4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0};
391:       PetscInt    cones[24]            = {7, 8, 9, 6, 11, 13, 9, 14, 10, 11, 13, 9, 9, 10, 11, 7, 9, 14, 13, 12, 7, 8, 11, 9};
392:       PetscInt    coneOrientations[24] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
393:       PetscScalar vertexCoords[27]     = {-2.0, -1.0, 0.0, -2.0, 0.0, 0.0, -2.0, 0.0, 1.0, 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 2.0, -1.0, 0.0, 2.0, 0.0, 0.0, 2.0, 0.0, 1.0};
394:       PetscInt    faultPoints[3]       = {9, 10, 11};

396:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
397:       for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
398:     } break;
399:     default:
400:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
401:     }
402:     DMPlexInterpolate(*dm, &idm);
403:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
404:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
405:     DMSetFromOptions(idm);
406:     DMViewFromOptions(idm, NULL, "-dm_view");
407:     DMGetLabel(*dm, "fault", &faultLabel);
408:     DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
409:     DMLabelDestroy(&hybridLabel);
410:     DMDestroy(&idm);
411:     DMDestroy(dm);
412:     *dm = hdm;
413:   } else {
414:     PetscInt numPoints[4] = {0, 0, 0, 0};

416:     DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
417:     DMPlexInterpolate(*dm, &idm);
418:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
419:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
420:     DMSetFromOptions(idm);
421:     DMViewFromOptions(idm, NULL, "-dm_view");
422:     DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
423:     DMDestroy(&idm);
424:     DMDestroy(dm);
425:     *dm = hdm;
426:   }
427:   return 0;
428: }

430: PetscErrorCode CreateTensorProduct_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
431: {
432:   DM          idm;
433:   PetscMPIInt rank;

435:   MPI_Comm_rank(comm, &rank);
436:   if (rank == 0) {
437:     switch (testNum) {
438:     case 0: {
439:       PetscInt    numPoints[2]         = {12, 2};
440:       PetscInt    coneSize[14]         = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
441:       PetscInt    cones[16]            = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8};
442:       PetscInt    coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
443:       PetscScalar vertexCoords[36]     = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5};

445:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
446:     } break;
447:     case 1: {
448:       PetscInt    numPoints[2]        = {8, 1};
449:       PetscInt    coneSize[9]         = {8, 0, 0, 0, 0, 0, 0, 0, 0};
450:       PetscInt    cones[8]            = {1, 2, 3, 4, 5, 6, 7, 8};
451:       PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
452:       PetscScalar vertexCoords[24]    = {-1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0};

454:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
455:     } break;
456:     default:
457:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
458:     }
459:   } else {
460:     PetscInt numPoints[4] = {0, 0, 0, 0};

462:     DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
463:   }
464:   DMPlexInterpolate(*dm, &idm);
465:   DMViewFromOptions(idm, NULL, "-in_dm_view");
466:   DMDestroy(dm);
467:   *dm = idm;
468:   return 0;
469: }

471: PetscErrorCode CreateTensorProductHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
472: {
473:   DM          idm, hdm = NULL;
474:   DMLabel     faultLabel;
475:   PetscInt    p;
476:   PetscMPIInt rank;

478:   MPI_Comm_rank(comm, &rank);
479:   if (rank == 0) {
480:     switch (testNum) {
481:     case 0: {
482:       PetscInt    numPoints[2]         = {12, 2};
483:       PetscInt    coneSize[14]         = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
484:       PetscInt    cones[16]            = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8};
485:       PetscInt    coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
486:       PetscScalar vertexCoords[36]     = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5};
487:       PetscInt    faultPoints[4]       = {2, 3, 5, 6};

489:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
490:       for (p = 0; p < 4; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
491:     } break;
492:     case 1: {
493:       PetscInt numPoints[2] = {30, 7};
494:       PetscInt coneSize[37] = {8, 8, 8, 8, 8, 8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
495:       PetscInt cones[56] = {8, 21, 20, 7, 13, 12, 23, 24, 14, 15, 10, 9, 13, 8, 21, 24, 15, 16, 11, 10, 24, 21, 22, 25, 30, 29, 28, 21, 35, 24, 33, 34, 24, 21, 30, 35, 25, 36, 31, 22, 27, 20, 21, 28, 32, 33, 24, 23, 15, 24, 13, 14, 19, 18, 17, 26};
496:       PetscInt coneOrientations[56] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
497:       PetscScalar vertexCoords[90]  = {-2.0, -2.0, -2.0, -2.0, -1.0, -2.0, -3.0, 0.0, -2.0, -2.0, 1.0,  -2.0, -2.0, 2.0, -2.0, -2.0, -2.0, 0.0,  -2.0, -1.0, 0.0, -3.0, 0.0, 0.0, -2.0, 1.0, 0.0, -2.0, 2.0, 0.0,
498:                                        -2.0, -1.0, 2.0,  -3.0, 0.0,  2.0,  -2.0, 1.0, 2.0,  0.0,  -2.0, -2.0, 0.0,  0.0, -2.0, 0.0,  2.0,  -2.0, 0.0,  -2.0, 0.0, 0.0,  0.0, 0.0, 0.0,  2.0, 0.0, 0.0,  0.0, 2.0,
499:                                        2.0,  -2.0, -2.0, 2.0,  -1.0, -2.0, 3.0,  0.0, -2.0, 2.0,  1.0,  -2.0, 2.0,  2.0, -2.0, 2.0,  -2.0, 0.0,  2.0,  -1.0, 0.0, 3.0,  0.0, 0.0, 2.0,  1.0, 0.0, 2.0,  2.0, 0.0};
500:       PetscInt    faultPoints[6]    = {20, 21, 22, 23, 24, 25};

502:       DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
503:       for (p = 0; p < 6; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
504:     } break;
505:     default:
506:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
507:     }
508:     DMPlexInterpolate(*dm, &idm);
509:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
510:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
511:     DMSetFromOptions(idm);
512:     DMViewFromOptions(idm, NULL, "-dm_view");
513:     DMGetLabel(*dm, "fault", &faultLabel);
514:     DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, NULL, NULL, NULL, &hdm);
515:     DMDestroy(&idm);
516:     DMDestroy(dm);
517:     *dm = hdm;
518:   } else {
519:     PetscInt numPoints[4] = {0, 0, 0, 0};

521:     DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
522:     DMPlexInterpolate(*dm, &idm);
523:     PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
524:     DMPlexDistributeSetDefault(idm, PETSC_FALSE);
525:     DMSetFromOptions(idm);
526:     DMViewFromOptions(idm, NULL, "-dm_view");
527:     DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
528:     DMDestroy(&idm);
529:     DMDestroy(dm);
530:     *dm = hdm;
531:   }
532:   return 0;
533: }

535: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
536: {
537:   PetscInt    dim         = user->dim;
538:   PetscBool   cellHybrid  = user->cellHybrid;
539:   PetscBool   cellSimplex = user->cellSimplex;
540:   PetscMPIInt rank, size;

542:   MPI_Comm_rank(comm, &rank);
543:   MPI_Comm_size(comm, &size);
544:   DMCreate(comm, dm);
545:   DMSetType(*dm, DMPLEX);
546:   DMSetDimension(*dm, dim);
547:   switch (dim) {
548:   case 1:
550:     CreateSimplex_1D(comm, dm);
551:     break;
552:   case 2:
553:     if (cellSimplex) {
554:       if (cellHybrid) {
555:         CreateSimplexHybrid_2D(comm, user->testNum, dm);
556:       } else {
557:         CreateSimplex_2D(comm, dm);
558:       }
559:     } else {
560:       if (cellHybrid) {
561:         CreateTensorProductHybrid_2D(comm, user->testNum, dm);
562:       } else {
563:         CreateTensorProduct_2D(comm, user->testNum, dm);
564:       }
565:     }
566:     break;
567:   case 3:
568:     if (cellSimplex) {
569:       if (cellHybrid) {
570:         CreateSimplexHybrid_3D(comm, user->testNum, dm);
571:       } else {
572:         CreateSimplex_3D(comm, user->testNum, dm);
573:       }
574:     } else {
575:       if (cellHybrid) {
576:         CreateTensorProductHybrid_3D(comm, user->testNum, dm);
577:       } else {
578:         CreateTensorProduct_3D(comm, user->testNum, dm);
579:       }
580:     }
581:     break;
582:   default:
583:     SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, dim);
584:   }
585:   if (user->testPartition && size > 1) {
586:     PetscPartitioner part;
587:     PetscInt        *sizes  = NULL;
588:     PetscInt        *points = NULL;

590:     if (rank == 0) {
591:       if (dim == 2 && cellSimplex && !cellHybrid && size == 2) {
592:         switch (user->testNum) {
593:         case 0: {
594:           PetscInt triSizes_p2[2]  = {1, 1};
595:           PetscInt triPoints_p2[2] = {0, 1};

597:           PetscMalloc2(2, &sizes, 2, &points);
598:           PetscArraycpy(sizes, triSizes_p2, 2);
599:           PetscArraycpy(points, triPoints_p2, 2);
600:           break;
601:         }
602:         default:
603:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 2 procs", user->testNum);
604:         }
605:       } else if (dim == 2 && cellSimplex && cellHybrid && size == 2) {
606:         switch (user->testNum) {
607:         case 0: {
608:           PetscInt triSizes_p2[2]  = {1, 2};
609:           PetscInt triPoints_p2[3] = {0, 1, 2};

611:           PetscMalloc2(2, &sizes, 3, &points);
612:           PetscArraycpy(sizes, triSizes_p2, 2);
613:           PetscArraycpy(points, triPoints_p2, 3);
614:           break;
615:         }
616:         default:
617:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular hybrid mesh on 2 procs", user->testNum);
618:         }
619:       } else if (dim == 2 && !cellSimplex && !cellHybrid && size == 2) {
620:         switch (user->testNum) {
621:         case 0: {
622:           PetscInt quadSizes_p2[2]  = {1, 1};
623:           PetscInt quadPoints_p2[2] = {0, 1};

625:           PetscMalloc2(2, &sizes, 2, &points);
626:           PetscArraycpy(sizes, quadSizes_p2, 2);
627:           PetscArraycpy(points, quadPoints_p2, 2);
628:           break;
629:         }
630:         default:
631:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral mesh on 2 procs", user->testNum);
632:         }
633:       } else if (dim == 2 && !cellSimplex && cellHybrid && size == 2) {
634:         switch (user->testNum) {
635:         case 0: {
636:           PetscInt quadSizes_p2[2]  = {1, 2};
637:           PetscInt quadPoints_p2[3] = {0, 1, 2};

639:           PetscMalloc2(2, &sizes, 3, &points);
640:           PetscArraycpy(sizes, quadSizes_p2, 2);
641:           PetscArraycpy(points, quadPoints_p2, 3);
642:           break;
643:         }
644:         default:
645:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral hybrid mesh on 2 procs", user->testNum);
646:         }
647:       } else if (dim == 3 && cellSimplex && !cellHybrid && size == 2) {
648:         switch (user->testNum) {
649:         case 0: {
650:           PetscInt tetSizes_p2[2]  = {1, 1};
651:           PetscInt tetPoints_p2[2] = {0, 1};

653:           PetscMalloc2(2, &sizes, 2, &points);
654:           PetscArraycpy(sizes, tetSizes_p2, 2);
655:           PetscArraycpy(points, tetPoints_p2, 2);
656:           break;
657:         }
658:         case 1: {
659:           PetscInt tetSizes_p2[2]  = {1, 1};
660:           PetscInt tetPoints_p2[2] = {0, 1};

662:           PetscMalloc2(2, &sizes, 2, &points);
663:           PetscArraycpy(sizes, tetSizes_p2, 2);
664:           PetscArraycpy(points, tetPoints_p2, 2);
665:           break;
666:         }
667:         default:
668:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for tetrahedral mesh on 2 procs", user->testNum);
669:         }
670:       } else if (dim == 3 && cellSimplex && cellHybrid && size == 2) {
671:         switch (user->testNum) {
672:         case 0: {
673:           PetscInt tetSizes_p2[2]  = {1, 2};
674:           PetscInt tetPoints_p2[3] = {0, 1, 2};

676:           PetscMalloc2(2, &sizes, 3, &points);
677:           PetscArraycpy(sizes, tetSizes_p2, 2);
678:           PetscArraycpy(points, tetPoints_p2, 3);
679:           break;
680:         }
681:         case 1: {
682:           PetscInt tetSizes_p2[2]  = {3, 4};
683:           PetscInt tetPoints_p2[7] = {0, 3, 5, 1, 2, 4, 6};

685:           PetscMalloc2(2, &sizes, 7, &points);
686:           PetscArraycpy(sizes, tetSizes_p2, 2);
687:           PetscArraycpy(points, tetPoints_p2, 7);
688:           break;
689:         }
690:         default:
691:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for tetrahedral hybrid mesh on 2 procs", user->testNum);
692:         }
693:       } else if (dim == 3 && !cellSimplex && !cellHybrid && size == 2) {
694:         switch (user->testNum) {
695:         case 0: {
696:           PetscInt hexSizes_p2[2]  = {1, 1};
697:           PetscInt hexPoints_p2[2] = {0, 1};

699:           PetscMalloc2(2, &sizes, 2, &points);
700:           PetscArraycpy(sizes, hexSizes_p2, 2);
701:           PetscArraycpy(points, hexPoints_p2, 2);
702:           break;
703:         }
704:         default:
705:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral mesh on 2 procs", user->testNum);
706:         }
707:       } else if (dim == 3 && !cellSimplex && cellHybrid && size == 2) {
708:         switch (user->testNum) {
709:         case 0: {
710:           PetscInt hexSizes_p2[2]  = {1, 1};
711:           PetscInt hexPoints_p2[2] = {0, 1};

713:           PetscMalloc2(2, &sizes, 2, &points);
714:           PetscArraycpy(sizes, hexSizes_p2, 2);
715:           PetscArraycpy(points, hexPoints_p2, 2);
716:           break;
717:         }
718:         case 1: {
719:           PetscInt hexSizes_p2[2]  = {5, 4};
720:           PetscInt hexPoints_p2[9] = {3, 4, 5, 7, 8, 0, 1, 2, 6};

722:           PetscMalloc2(2, &sizes, 9, &points);
723:           PetscArraycpy(sizes, hexSizes_p2, 2);
724:           PetscArraycpy(points, hexPoints_p2, 9);
725:           break;
726:         }
727:         default:
728:           SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral hybrid mesh on 2 procs", user->testNum);
729:         }
730:       } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition");
731:     }
732:     DMPlexGetPartitioner(*dm, &part);
733:     PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);
734:     PetscPartitionerShellSetPartition(part, size, sizes, points);
735:     PetscFree2(sizes, points);
736:   } else {
737:     PetscPartitioner part;

739:     DMPlexGetPartitioner(*dm, &part);
740:     PetscPartitionerSetFromOptions(part);
741:   }
742:   {
743:     DM pdm = NULL;

745:     DMPlexDistribute(*dm, 0, NULL, &pdm);
746:     if (pdm) {
747:       DMViewFromOptions(pdm, NULL, "-dm_view");
748:       DMDestroy(dm);
749:       *dm = pdm;
750:     }
751:   }
752:   DMPlexDistributeSetDefault(*dm, PETSC_FALSE);
753:   DMViewFromOptions(*dm, NULL, "-dm_view_pre");
754:   DMSetFromOptions(*dm);
755:   if (user->uninterpolate || user->reinterpolate) {
756:     DM udm = NULL;

758:     DMPlexUninterpolate(*dm, &udm);
759:     DMPlexCopyCoordinates(*dm, udm);
760:     DMDestroy(dm);
761:     *dm = udm;
762:   }
763:   if (user->reinterpolate) {
764:     DM idm = NULL;

766:     DMPlexInterpolate(*dm, &idm);
767:     DMPlexCopyCoordinates(*dm, idm);
768:     DMDestroy(dm);
769:     *dm = idm;
770:   }
771:   DMPlexDistributeSetDefault(*dm, PETSC_FALSE);
772:   PetscObjectSetName((PetscObject)*dm, "Hybrid Mesh");
773:   DMViewFromOptions(*dm, NULL, "-dm_view");
774:   PetscObjectSetOptionsPrefix((PetscObject)*dm, "hyb_");
775:   DMSetFromOptions(*dm);
776:   return 0;
777: }

779: int main(int argc, char **argv)
780: {
781:   DM     dm;
782:   AppCtx user; /* user-defined work context */

785:   PetscInitialize(&argc, &argv, NULL, help);
786:   ProcessOptions(PETSC_COMM_WORLD, &user);
787:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
788:   DMDestroy(&dm);
789:   PetscFinalize();
790:   return 0;
791: }

793: /*TEST

795:   # 1D Simplex 29-31
796:   testset:
797:     args: -dim 1 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
798:     test:
799:       suffix: 29
800:     test:
801:       suffix: 30
802:       args: -dm_refine 1
803:     test:
804:       suffix: 31
805:       args: -dm_refine 5

807:   # 2D Simplex 0-3
808:   testset:
809:     args: -dim 2 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
810:     test:
811:       suffix: 0
812:     test:
813:       suffix: 1
814:       args: -dm_refine 1
815:     test:
816:       suffix: 2
817:       nsize: 2
818:     test:
819:       suffix: 3
820:       nsize: 2
821:       args: -dm_refine 1
822:     test:
823:       suffix: 32
824:       args: -dm_refine 1 -uninterpolate
825:     test:
826:       suffix: 33
827:       nsize: 2
828:       args: -dm_refine 1 -uninterpolate
829:     test:
830:       suffix: 34
831:       nsize: 2
832:       args: -dm_refine 3 -uninterpolate

834:   # 2D Hybrid Simplex 4-7
835:   testset:
836:     args: -dim 2 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
837:     test:
838:       suffix: 4
839:     test:
840:       suffix: 5
841:       args: -dm_refine 1
842:     test:
843:       suffix: 6
844:       nsize: 2
845:     test:
846:       suffix: 7
847:       nsize: 2
848:       args: -dm_refine 1
849:     test:
850:       suffix: 24
851:       args: -test_num 1 -dm_refine 1

853:   # 2D Quad 12-13
854:   testset:
855:     args: -dim 2 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
856:     test:
857:       suffix: 12
858:       args: -dm_refine 1
859:     test:
860:       suffix: 13
861:       nsize: 2
862:       args: -dm_refine 1

864:   # 2D Hybrid Quad 27-28
865:   testset:
866:     args: -dim 2 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
867:     test:
868:       suffix: 27
869:       args: -dm_refine 1
870:     test:
871:       suffix: 28
872:       nsize: 2
873:       args: -dm_refine 1

875:   # 3D Simplex 8-11
876:   testset:
877:     args: -dim 3 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
878:     test:
879:       suffix: 8
880:       args: -dm_refine 1
881:     test:
882:       suffix: 9
883:       nsize: 2
884:       args: -dm_refine 1
885:     test:
886:       suffix: 10
887:       args: -test_num 1 -dm_refine 1
888:     test:
889:       suffix: 11
890:       nsize: 2
891:       args: -test_num 1 -dm_refine 1
892:     test:
893:       suffix: 25
894:       args: -test_num 2 -dm_refine 2

896:   # 3D Hybrid Simplex 16-19
897:   testset:
898:     args: -dim 3 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
899:     test:
900:       suffix: 16
901:       args: -dm_refine 1
902:     test:
903:       suffix: 17
904:       nsize: 2
905:       args: -dm_refine 1
906:     test:
907:       suffix: 18
908:       args: -test_num 1 -dm_refine 1
909:     test:
910:       suffix: 19
911:       nsize: 2
912:       args: -test_num 1 -dm_refine 1

914:   # 3D Hex 14-15
915:   testset:
916:     args: -dim 3 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
917:     test:
918:       suffix: 14
919:       args: -dm_refine 1
920:     test:
921:       suffix: 15
922:       nsize: 2
923:      args: -dm_refine 1
924:     test:
925:       suffix: 26
926:       args: -test_num 1 -dm_refine 2

928:   # 3D Hybrid Hex 20-23
929:   testset:
930:     args: -dim 3 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
931:     test:
932:       suffix: 20
933:       args: -dm_refine 1
934:     test:
935:       suffix: 21
936:       nsize: 2
937:       args: -dm_refine 1
938:     test:
939:       suffix: 22
940:       args: -test_num 1 -dm_refine 1
941:     test:
942:       suffix: 23
943:       nsize: 2
944:       args: -test_num 1 -dm_refine 1

946:   # Hybrid interpolation
947:   #   TODO Setup new tests (like -reinterpolate) that interpolate hybrid cells
948:   testset:
949:     nsize: 2
950:     args: -test_partition 0 -petscpartitioner_type simple -dm_view -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
951:     test:
952:       suffix: hybint_2d_0
953:       args: -dim 2 -dm_refine 2
954:     test:
955:       suffix: hybint_2d_1
956:       args: -dim 2 -dm_refine 2 -test_num 1
957:     test:
958:       suffix: hybint_3d_0
959:       args: -dim 3 -dm_refine 1
960:     test:
961:       suffix: hybint_3d_1
962:       args: -dim 3 -dm_refine 1 -test_num 1

964: TEST*/