Actual source code: ex1.c

petsc-master 2020-08-25
Report Typos and Errors
  1: static char help[] = "Tests various DMPlex routines to construct, refine and distribute a mesh.\n\n";

  3:  #include <petscdmplex.h>
  4:  #include <petscsf.h>

  6: typedef enum {BOX, CYLINDER, SPHERE, BALL} DomainShape;
  7: enum {STAGE_LOAD, STAGE_DISTRIBUTE, STAGE_REFINE, STAGE_OVERLAP};

  9: typedef struct {
 10:   DM            dm;                /* REQUIRED in order to use SNES evaluation functions */
 11:   PetscInt      debug;             /* The debugging level */
 12:   PetscLogEvent createMeshEvent;
 13:   PetscLogStage stages[4];
 14:   /* Domain and mesh definition */
 15:   PetscInt      dim;                             /* The topological mesh dimension */
 16:   PetscBool     interpolate;                     /* Generate intermediate mesh elements */
 17:   PetscReal     refinementLimit;                 /* The largest allowable cell volume */
 18:   PetscBool     cellSimplex;                     /* Use simplices or hexes */
 19:   PetscBool     cellWedge;                       /* Use wedges */
 20:   DomainShape   domainShape;                     /* Shape of the region to be meshed */
 21:   PetscInt      *domainBoxSizes;                 /* Sizes of the box mesh */
 22:   PetscReal     *domainBoxL,*domainBoxU;         /* Lower left, upper right corner of the box mesh */
 23:   DMBoundaryType periodicity[3];                 /* The domain periodicity */
 24:   char          filename[PETSC_MAX_PATH_LEN];    /* Import mesh from file */
 25:   char          bdfilename[PETSC_MAX_PATH_LEN];  /* Import mesh boundary from file */
 26:   char          extfilename[PETSC_MAX_PATH_LEN]; /* Import 2D mesh to be extruded from file */
 27:   PetscBool     testPartition;                   /* Use a fixed partitioning for testing */
 28:   PetscInt      overlap;                         /* The cell overlap to use during partitioning */
 29:   PetscReal     extrude_thickness;               /* Thickness of extrusion */
 30:   PetscInt      extrude_layers;                  /* Layers to be extruded */
 31:   PetscBool     extrude_hfirst;                  /* New numbering: height first? */
 32:   PetscBool     testp4est[2];
 33:   PetscBool     redistribute;
 34:   PetscBool     final_ref;                       /* Run refinement at the end */
 35:   PetscBool     final_diagnostics;               /* Run diagnostics on the final mesh */
 36: } AppCtx;

 38: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 39: {
 40:   const char       *dShapes[4] = {"box", "cylinder", "sphere", "ball"};
 41:   PetscInt         shape, bd, n;
 42:   static PetscInt  domainBoxSizes[3] = {1,1,1};
 43:   static PetscReal domainBoxL[3] = {0.,0.,0.};
 44:   static PetscReal domainBoxU[3] = {1.,1.,1.};
 45:   PetscBool        flg;
 46:   PetscErrorCode   ierr;

 49:   options->debug             = 0;
 50:   options->dim               = 2;
 51:   options->interpolate       = PETSC_FALSE;
 52:   options->refinementLimit   = 0.0;
 53:   options->cellSimplex       = PETSC_TRUE;
 54:   options->cellWedge         = PETSC_FALSE;
 55:   options->domainShape       = BOX;
 56:   options->domainBoxSizes    = NULL;
 57:   options->domainBoxL        = NULL;
 58:   options->domainBoxU        = NULL;
 59:   options->periodicity[0]    = DM_BOUNDARY_NONE;
 60:   options->periodicity[1]    = DM_BOUNDARY_NONE;
 61:   options->periodicity[2]    = DM_BOUNDARY_NONE;
 62:   options->filename[0]       = '\0';
 63:   options->bdfilename[0]     = '\0';
 64:   options->extfilename[0]    = '\0';
 65:   options->testPartition     = PETSC_FALSE;
 66:   options->overlap           = 0;
 67:   options->extrude_layers    = 0;
 68:   options->extrude_thickness = 0.1;
 69:   options->extrude_hfirst    = PETSC_TRUE;
 70:   options->testp4est[0]      = PETSC_FALSE;
 71:   options->testp4est[1]      = PETSC_FALSE;
 72:   options->redistribute      = PETSC_FALSE;
 73:   options->final_ref         = PETSC_FALSE;
 74:   options->final_diagnostics = PETSC_TRUE;

 76:   PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");
 77:   PetscOptionsBoundedInt("-debug", "The debugging level", "ex1.c", options->debug, &options->debug, NULL,0);
 78:   PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex1.c", options->dim, &options->dim, NULL,1,3);
 79:   PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex1.c", options->interpolate, &options->interpolate, NULL);
 80:   PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex1.c", options->refinementLimit, &options->refinementLimit, NULL);
 81:   PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex1.c", options->cellSimplex, &options->cellSimplex, NULL);
 82:   PetscOptionsBool("-cell_wedge", "Use wedges if true", "ex1.c", options->cellWedge, &options->cellWedge, NULL);
 83:   shape = options->domainShape;
 84:   PetscOptionsEList("-domain_shape","The shape of the domain","ex1.c", dShapes, 4, dShapes[options->domainShape], &shape, NULL);
 85:   options->domainShape = (DomainShape) shape;
 86:   PetscOptionsIntArray("-domain_box_sizes","The sizes of the box domain","ex1.c", domainBoxSizes, (n=3,&n), &flg);
 87:   if (flg) { options->domainShape = BOX; options->domainBoxSizes = domainBoxSizes;}
 88:   PetscOptionsRealArray("-domain_box_ll","Coordinates of the lower left corner of the box domain","ex1.c", domainBoxL, (n=3,&n), &flg);
 89:   if (flg) { options->domainBoxL = domainBoxL;}
 90:   PetscOptionsRealArray("-domain_box_ur","Coordinates of the upper right corner of the box domain","ex1.c", domainBoxU, (n=3,&n), &flg);
 91:   if (flg) { options->domainBoxU = domainBoxU;}
 92:   bd = options->periodicity[0];
 93:   PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex1.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
 94:   options->periodicity[0] = (DMBoundaryType) bd;
 95:   bd = options->periodicity[1];
 96:   PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex1.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
 97:   options->periodicity[1] = (DMBoundaryType) bd;
 98:   bd = options->periodicity[2];
 99:   PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex1.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
100:   options->periodicity[2] = (DMBoundaryType) bd;
101:   PetscOptionsString("-filename", "The mesh file", "ex1.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);
102:   PetscOptionsString("-bd_filename", "The mesh boundary file", "ex1.c", options->bdfilename, options->bdfilename, sizeof(options->bdfilename), NULL);
103:   PetscOptionsString("-ext_filename", "The 2D mesh file to be extruded", "ex1.c", options->extfilename, options->extfilename, sizeof(options->extfilename), &flg);
104:   if (flg) options->extrude_layers = 2;
105:   PetscOptionsBoundedInt("-ext_layers", "The number of layers to extrude", "ex1.c", options->extrude_layers, &options->extrude_layers, NULL,0);
106:   PetscOptionsReal("-ext_thickness", "The thickness of the layer to be extruded", "ex1.c", options->extrude_thickness, &options->extrude_thickness, NULL);
107:   PetscOptionsBool("-ext_hfirst", "Order the cells in the height first", "ex1.c", options->extrude_hfirst, &options->extrude_hfirst, NULL);
108:   PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex1.c", options->testPartition, &options->testPartition, NULL);
109:   PetscOptionsBoundedInt("-overlap", "The cell overlap for partitioning", "ex1.c", options->overlap, &options->overlap, NULL,0);
110:   PetscOptionsBool("-test_p4est_seq", "Test p4est with sequential base DM", "ex1.c", options->testp4est[0], &options->testp4est[0], NULL);
111:   PetscOptionsBool("-test_p4est_par", "Test p4est with parallel base DM", "ex1.c", options->testp4est[1], &options->testp4est[1], NULL);
112:   PetscOptionsBool("-test_redistribute", "Test redistribution", "ex1.c", options->redistribute, &options->redistribute, NULL);
113:   PetscOptionsBool("-final_ref", "Run uniform refinement on the final mesh", "ex1.c", options->final_ref, &options->final_ref, NULL);
114:   PetscOptionsBool("-final_diagnostics", "Run diagnostics on the final mesh", "ex1.c", options->final_diagnostics, &options->final_diagnostics, NULL);
115:   PetscOptionsEnd();

117:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
118:   PetscLogStageRegister("MeshLoad",       &options->stages[STAGE_LOAD]);
119:   PetscLogStageRegister("MeshDistribute", &options->stages[STAGE_DISTRIBUTE]);
120:   PetscLogStageRegister("MeshRefine",     &options->stages[STAGE_REFINE]);
121:   PetscLogStageRegister("MeshOverlap",    &options->stages[STAGE_OVERLAP]);
122:   return(0);
123: }

125: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
126: {
127:   PetscInt       dim                  = user->dim;
128:   PetscBool      interpolate          = user->interpolate;
129:   PetscReal      refinementLimit      = user->refinementLimit;
130:   PetscBool      cellSimplex          = user->cellSimplex;
131:   PetscBool      cellWedge            = user->cellWedge;
132:   const char    *filename             = user->filename;
133:   const char    *bdfilename           = user->bdfilename;
134:   const char    *extfilename          = user->extfilename;
135:   PetscBool      testp4est_seq        = user->testp4est[0];
136:   PetscBool      testp4est_par        = user->testp4est[1];
137:   PetscInt       triSizes_n2[2]       = {4, 4};
138:   PetscInt       triPoints_n2[8]      = {3, 5, 6, 7, 0, 1, 2, 4};
139:   PetscInt       triSizes_n8[8]       = {1, 1, 1, 1, 1, 1, 1, 1};
140:   PetscInt       triPoints_n8[8]      = {0, 1, 2, 3, 4, 5, 6, 7};
141:   PetscInt       quadSizes[2]         = {2, 2};
142:   PetscInt       quadPoints[4]        = {2, 3, 0, 1};
143:   PetscInt       gmshSizes_n3[3]      = {14, 14, 14};
144:   PetscInt       gmshPoints_n3[42]    = {1, 2,  4,  5,  9, 10, 11, 15, 16, 20, 21, 27, 28, 29,
145:                                          3, 8, 12, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40,
146:                                          0, 6,  7, 13, 14, 17, 18, 19, 22, 23, 24, 25, 26, 41};
147:   PetscInt       fluentSizes_n3[3]    = {50, 50, 50};
148:   PetscInt       fluentPoints_n3[150] = { 5,  6,  7,  8, 12, 14, 16,  34,  36,  37,  38,  39,  40,  41,  42,  43,  44,  45,  46,  48,  50,  51,  80,  81,  89,
149:                                          91, 93, 94, 95, 96, 97, 98,  99, 100, 101, 104, 121, 122, 124, 125, 126, 127, 128, 129, 131, 133, 143, 144, 145, 147,
150:                                           1,  3,  4,  9, 10, 17, 18,  19,  24,  25,  26,  27,  28,  29,  30,  31,  32,  33,  35,  47,  61,  71,  72,  73,  74,
151:                                          75, 76, 77, 78, 79, 86, 87,  88,  90,  92, 113, 115, 116, 117, 118, 119, 120, 123, 138, 140, 141, 142, 146, 148, 149,
152:                                           0,  2, 11, 13, 15, 20, 21,  22,  23,  49,  52,  53,  54,  55,  56,  57,  58,  59,  60,  62,  63,  64,  65,  66,  67,
153:                                          68, 69, 70, 82, 83, 84, 85, 102, 103, 105, 106, 107, 108, 109, 110, 111, 112, 114, 130, 132, 134, 135, 136, 137, 139};
154:   size_t         len, bdlen, extlen;
155:   PetscMPIInt    rank, size;

159:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
160:   MPI_Comm_rank(comm, &rank);
161:   MPI_Comm_size(comm, &size);
162:   PetscStrlen(filename, &len);
163:   PetscStrlen(bdfilename, &bdlen);
164:   PetscStrlen(extfilename, &extlen);
165:   PetscLogStagePush(user->stages[STAGE_LOAD]);
166:   if (len) {
167:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
168:   } else if (bdlen) {
169:     DM boundary;

171:     DMPlexCreateFromFile(comm, bdfilename, interpolate, &boundary);
172:     DMPlexGenerate(boundary, NULL, interpolate, dm);
173:     DMDestroy(&boundary);
174:   } else if (extlen) {
175:     DM edm;

177:     DMPlexCreateFromFile(comm, extfilename, interpolate, &edm);
178:     DMPlexExtrude(edm, user->extrude_layers, user->extrude_thickness, user->extrude_hfirst, NULL, interpolate, dm);
179:     DMDestroy(&edm);
180:   } else {
181:     switch (user->domainShape) {
182:     case BOX:
183:       if (cellWedge) {
184:         if (dim != 3) SETERRQ1(comm, PETSC_ERR_ARG_WRONG, "Dimension must be 3 for a wedge mesh, not %D", dim);
185:         DMPlexCreateWedgeBoxMesh(comm, user->domainBoxSizes, user->domainBoxL, user->domainBoxU, user->periodicity, PETSC_FALSE, interpolate, dm);
186:       } else {
187:         DMPlexCreateBoxMesh(comm, dim, cellSimplex, user->domainBoxSizes, user->domainBoxL, user->domainBoxU, user->periodicity, interpolate, dm);
188:       }
189:       break;
190:     case CYLINDER:
191:       if (cellSimplex) SETERRQ(comm, PETSC_ERR_ARG_WRONG, "Cannot mesh a cylinder with simplices");
192:       if (dim != 3)    SETERRQ1(comm, PETSC_ERR_ARG_WRONG, "Dimension must be 3 for a cylinder mesh, not %D", dim);
193:       if (cellWedge) {
194:         DMPlexCreateWedgeCylinderMesh(comm, 6, interpolate, dm);
195:       } else {
196:         DMPlexCreateHexCylinderMesh(comm, 1, user->periodicity[2], dm);
197:       }
198:       break;
199:     case SPHERE:
200:       DMPlexCreateSphereMesh(comm, dim, cellSimplex, 1.0, dm);
201:       break;
202:     case BALL:
203:       DMPlexCreateBallMesh(comm, dim, 1.0, dm);
204:       break;
205:     default: SETERRQ1(comm, PETSC_ERR_ARG_WRONG, "Unknown domain shape %D", user->domainShape);
206:     }
207:   }
208:   if (!extlen && user->extrude_layers > 0) {
209:     DM edm;

211:     DMPlexExtrude(*dm, user->extrude_layers, user->extrude_thickness, user->extrude_hfirst, NULL, interpolate, &edm);
212:     DMDestroy(dm);
213:     *dm  = edm;
214:   }
215:   DMLocalizeCoordinates(*dm); /* needed for periodic */
216:   DMViewFromOptions(*dm,NULL,"-init_dm_view");
217:   DMGetDimension(*dm,&dim);

219:   if (testp4est_seq) {
220: #if defined(PETSC_HAVE_P4EST)
221:     DM dmConv = NULL;

223:     DMPlexCheckSymmetry(*dm);
224:     DMPlexCheckSkeleton(*dm, 0);
225:     DMPlexCheckFaces(*dm, 0);
226:     DMPlexCheckGeometry(*dm);
227:     DMPlexCheckPointSF(*dm);
228:     DMPlexCheckInterfaceCones(*dm);
229:     DMPlexSetRefinementUniform(*dm, PETSC_TRUE);
230:     DMPlexSetCellRefinerType(*dm, DM_REFINER_TO_BOX);
231:     DMRefine(*dm, PETSC_COMM_WORLD, &dmConv);
232:     if (dmConv) {
233:       DMDestroy(dm);
234:       *dm  = dmConv;
235:     }
236:     DMViewFromOptions(*dm,NULL,"-initref_dm_view");
237:     DMPlexCheckSymmetry(*dm);
238:     DMPlexCheckSkeleton(*dm, 0);
239:     DMPlexCheckFaces(*dm, 0);
240:     DMPlexCheckGeometry(*dm);
241:     DMPlexCheckPointSF(*dm);
242:     DMPlexCheckInterfaceCones(*dm);
243:     user->cellSimplex = PETSC_FALSE;

245:     DMConvert(*dm,dim == 2 ? DMP4EST : DMP8EST,&dmConv);
246:     if (dmConv) {
247:       PetscObjectSetOptionsPrefix((PetscObject) dmConv, "conv_seq_1_");
248:       DMSetFromOptions(dmConv);
249:       DMDestroy(dm);
250:       *dm  = dmConv;
251:     }
252:     PetscObjectSetOptionsPrefix((PetscObject) *dm, "conv_seq_1_");
253:     DMSetUp(*dm);
254:     DMViewFromOptions(*dm, NULL, "-dm_view");
255:     DMConvert(*dm,DMPLEX,&dmConv);
256:     if (dmConv) {
257:       PetscObjectSetOptionsPrefix((PetscObject) dmConv, "conv_seq_2_");
258:       DMSetFromOptions(dmConv);
259:       DMDestroy(dm);
260:       *dm  = dmConv;
261:     }
262:     PetscObjectSetOptionsPrefix((PetscObject) *dm, "conv_seq_2_");
263:     DMViewFromOptions(*dm, NULL, "-dm_view");
264:     PetscObjectSetOptionsPrefix((PetscObject) *dm, NULL);
265: #else
266:     SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Recompile with --download-p4est");
267: #endif
268:   }

270:   PetscLogStagePop();
271:   if (!testp4est_seq) {
272:     DM refinedMesh     = NULL;
273:     DM distributedMesh = NULL;

275:     if (user->testPartition) {
276:       const PetscInt  *sizes = NULL;
277:       const PetscInt  *points = NULL;
278:       PetscPartitioner part;

280:       if (!rank) {
281:         if (dim == 2 && cellSimplex && size == 2) {
282:            sizes = triSizes_n2; points = triPoints_n2;
283:         } else if (dim == 2 && cellSimplex && size == 8) {
284:           sizes = triSizes_n8; points = triPoints_n8;
285:         } else if (dim == 2 && !cellSimplex && size == 2) {
286:           sizes = quadSizes; points = quadPoints;
287:         } else if (dim == 2 && size == 3) {
288:           PetscInt Nc;

290:           DMPlexGetHeightStratum(*dm, 0, NULL, &Nc);
291:           if (Nc == 42) { /* Gmsh 3 & 4 */
292:             sizes = gmshSizes_n3; points = gmshPoints_n3;
293:           } else if (Nc == 150) { /* Fluent 1 */
294:             sizes = fluentSizes_n3; points = fluentPoints_n3;
295:           } else if (Nc == 42) { /* Med 1 */
296:           } else if (Nc == 161) { /* Med 3 */
297:           }
298:         }
299:       }
300:       DMPlexGetPartitioner(*dm, &part);
301:       PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);
302:       PetscPartitionerShellSetPartition(part, size, sizes, points);
303:     } else {
304:       PetscPartitioner part;

306:       DMPlexGetPartitioner(*dm,&part);
307:       PetscPartitionerSetFromOptions(part);
308:     }
309:     /* Distribute mesh over processes */
310:     PetscLogStagePush(user->stages[STAGE_DISTRIBUTE]);
311:     DMViewFromOptions(*dm, NULL, "-dm_pre_dist_view");
312:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
313:     if (distributedMesh) {
314:       DMDestroy(dm);
315:       *dm  = distributedMesh;
316:     }
317:     PetscLogStagePop();
318:     DMViewFromOptions(*dm, NULL, "-distributed_dm_view");
319:     /* Refine mesh using a volume constraint */
320:     PetscLogStagePush(user->stages[STAGE_REFINE]);
321:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
322:     DMPlexSetRefinementLimit(*dm, refinementLimit);
323:     DMRefine(*dm, comm, &refinedMesh);
324:     if (refinedMesh) {
325:       DMDestroy(dm);
326:       *dm  = refinedMesh;
327:     }
328:     PetscLogStagePop();
329:   }
330:   PetscLogStagePush(user->stages[STAGE_REFINE]);
331:   DMSetFromOptions(*dm);
332:   PetscLogStagePop();

334:   if (testp4est_par) {
335: #if defined(PETSC_HAVE_P4EST)
336:     DM dmConv = NULL;

338:     DMViewFromOptions(*dm, NULL, "-dm_tobox_view");
339:     DMPlexSetRefinementUniform(*dm, PETSC_TRUE);
340:     DMPlexSetCellRefinerType(*dm, DM_REFINER_TO_BOX);
341:     DMRefine(*dm, PETSC_COMM_WORLD, &dmConv);
342:     if (dmConv) {
343:       DMDestroy(dm);
344:       *dm  = dmConv;
345:     }
346:     user->cellSimplex = PETSC_FALSE;
347:     DMViewFromOptions(*dm, NULL, "-dm_tobox_view");
348:     DMPlexCheckSymmetry(*dm);
349:     DMPlexCheckSkeleton(*dm, 0);
350:     DMPlexCheckFaces(*dm, 0);
351:     DMPlexCheckGeometry(*dm);
352:     DMPlexCheckPointSF(*dm);
353:     DMPlexCheckInterfaceCones(*dm);

355:     DMConvert(*dm,dim == 2 ? DMP4EST : DMP8EST,&dmConv);
356:     if (dmConv) {
357:       PetscObjectSetOptionsPrefix((PetscObject) dmConv, "conv_par_1_");
358:       DMSetFromOptions(dmConv);
359:       DMDestroy(dm);
360:       *dm  = dmConv;
361:     }
362:     PetscObjectSetOptionsPrefix((PetscObject) *dm, "conv_par_1_");
363:     DMSetUp(*dm);
364:     DMViewFromOptions(*dm, NULL, "-dm_view");
365:     DMConvert(*dm, DMPLEX, &dmConv);
366:     if (dmConv) {
367:       PetscObjectSetOptionsPrefix((PetscObject) dmConv, "conv_par_2_");
368:       DMSetFromOptions(dmConv);
369:       DMDestroy(dm);
370:       *dm  = dmConv;
371:     }
372:     PetscObjectSetOptionsPrefix((PetscObject) *dm, "conv_par_2_");
373:     DMViewFromOptions(*dm, NULL, "-dm_view");
374:     PetscObjectSetOptionsPrefix((PetscObject) *dm, NULL);
375: #else
376:     SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Recompile with --download-p4est");
377: #endif
378:   }

380:   /* test redistribution of an already distributed mesh */
381:   if (user->redistribute) {
382:     DM       distributedMesh;
383:     PetscSF  sf;
384:     PetscInt nranks;

386:     DMViewFromOptions(*dm, NULL, "-dm_pre_redist_view");
387:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
388:     if (distributedMesh) {
389:       DMGetPointSF(distributedMesh, &sf);
390:       PetscSFSetUp(sf);
391:       DMGetNeighbors(distributedMesh, &nranks, NULL);
392:       MPI_Allreduce(MPI_IN_PLACE, &nranks, 1, MPIU_INT, MPI_MIN, PetscObjectComm((PetscObject)*dm));
393:       PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)*dm)), "Minimum number of neighbors: %D\n", nranks);
394:       DMDestroy(dm);
395:       *dm  = distributedMesh;
396:     }
397:     DMViewFromOptions(*dm, NULL, "-dm_post_redist_view");
398:   }

400:   if (user->overlap) {
401:     DM overlapMesh = NULL;

403:     /* Add the overlap to refined mesh */
404:     PetscLogStagePush(user->stages[STAGE_OVERLAP]);
405:     DMViewFromOptions(*dm, NULL, "-dm_pre_overlap_view");
406:     DMPlexDistributeOverlap(*dm, user->overlap, NULL, &overlapMesh);
407:     if (overlapMesh) {
408:       PetscInt overlap;
409:       DMPlexGetOverlap(overlapMesh, &overlap);
410:       PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD, "Overlap: %D\n", overlap);
411:       DMDestroy(dm);
412:       *dm = overlapMesh;
413:     }
414:     DMViewFromOptions(*dm, NULL, "-dm_post_overlap_view");
415:     PetscLogStagePop();
416:   }
417:   if (user->final_ref) {
418:     DM refinedMesh = NULL;

420:     DMPlexSetRefinementUniform(*dm, PETSC_TRUE);
421:     DMRefine(*dm, comm, &refinedMesh);
422:     if (refinedMesh) {
423:       DMDestroy(dm);
424:       *dm  = refinedMesh;
425:     }
426:   }

428:   PetscObjectSetName((PetscObject) *dm, "Simplicial Mesh");
429:   DMViewFromOptions(*dm, NULL, "-dm_view");
430:   if (user->final_diagnostics) {
431:     DMPlexInterpolatedFlag interpolated;
432:     PetscInt  dim, depth;

434:     DMGetDimension(*dm, &dim);
435:     DMPlexGetDepth(*dm, &depth);
436:     DMPlexIsInterpolatedCollective(*dm, &interpolated);

438:     DMPlexCheckSymmetry(*dm);
439:     if (interpolated == DMPLEX_INTERPOLATED_FULL) {
440:       DMPlexCheckFaces(*dm, 0);
441:     }
442:     DMPlexCheckSkeleton(*dm, 0);
443:     DMPlexCheckGeometry(*dm);
444:   }
445:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
446:   user->dm = *dm;
447:   return(0);
448: }

450: int main(int argc, char **argv)
451: {
452:   AppCtx         user;                 /* user-defined work context */

455:   PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
456:   ProcessOptions(PETSC_COMM_WORLD, &user);
457:   CreateMesh(PETSC_COMM_WORLD, &user, &user.dm);
458:   DMDestroy(&user.dm);
459:   PetscFinalize();
460:   return ierr;
461: }

463: /*TEST

465:   # CTetGen 0-1
466:   test:
467:     suffix: 0
468:     requires: ctetgen
469:     args: -dim 3 -ctetgen_verbose 4 -dm_view ascii::ascii_info_detail -info :~sys
470:   test:
471:     suffix: 1
472:     requires: ctetgen
473:     args: -dim 3 -ctetgen_verbose 4 -refinement_limit 0.0625 -dm_view ascii::ascii_info_detail -info :~sys


476:   # 2D LaTex and ASCII output 2-9
477:   test:
478:     suffix: 2
479:     requires: triangle
480:     args: -dim 2 -dm_view ascii::ascii_latex
481:   test:
482:     suffix: 3
483:     requires: triangle
484:     args: -dim 2 -dm_refine 1 -interpolate 1 -dm_view ascii::ascii_info_detail
485:   test:
486:     suffix: 4
487:     requires: triangle
488:     nsize: 2
489:     args: -dim 2 -dm_refine 1 -interpolate 1 -test_partition -dm_view ascii::ascii_info_detail
490:   test:
491:     suffix: 5
492:     requires: triangle
493:     nsize: 2
494:     args: -dim 2 -dm_refine 1 -interpolate 1 -test_partition -dm_view ascii::ascii_latex
495:   test:
496:     suffix: 6
497:     args: -dim 2 -cell_simplex 0 -interpolate -dm_view ascii::ascii_info_detail
498:   test:
499:     suffix: 7
500:     args: -dim 2 -cell_simplex 0 -interpolate -dm_refine 1 -dm_view ascii::ascii_info_detail
501:   test:
502:     suffix: 8
503:     nsize: 2
504:     args: -dim 2 -cell_simplex 0 -interpolate -dm_refine 1 -interpolate 1 -test_partition -dm_view ascii::ascii_latex

506:   # 1D ASCII output
507:   test:
508:     suffix: 1d_0
509:     args: -dim 1 -domain_shape box -dm_view ascii::ascii_info_detail
510:   test:
511:     suffix: 1d_1
512:     args: -dim 1 -domain_shape box -dm_refine 2 -dm_view ascii::ascii_info_detail
513:   test:
514:     suffix: 1d_2
515:     args: -dim 1 -domain_box_sizes 5 -x_periodicity periodic -dm_view ascii::ascii_info_detail -dm_plex_check_all

517:   # Parallel refinement tests with overlap
518:   test:
519:     suffix: refine_overlap_1d
520:     nsize: 2
521:     args: -dim 1 -domain_box_sizes 4 -dm_refine 1 -overlap {{0 1 2}separate output} -petscpartitioner_type simple -dm_view ascii::ascii_info
522:   test:
523:     suffix: refine_overlap_2d
524:     requires: triangle
525:     nsize: {{2 8}separate output}
526:     args: -dim 2 -cell_simplex 1 -dm_refine 1 -interpolate 1 -test_partition -overlap {{0 1 2}separate output} -dm_view ascii::ascii_info

528:   # Parallel simple partitioner tests
529:   test:
530:     suffix: part_simple_0
531:     requires: triangle
532:     nsize: 2
533:     args: -dim 2 -cell_simplex 1 -dm_refine 0 -interpolate 0 -petscpartitioner_type simple -partition_view -dm_view ascii::ascii_info_detail
534:   test:
535:     suffix: part_simple_1
536:     requires: triangle
537:     nsize: 8
538:     args: -dim 2 -cell_simplex 1 -dm_refine 1 -interpolate 1 -petscpartitioner_type simple -partition_view -dm_view ascii::ascii_info_detail

540:   # Parallel partitioner tests
541:   test:
542:     suffix: part_parmetis_0
543:     requires: parmetis
544:     nsize: 2
545:     args: -dim 2 -cell_simplex 0 -dm_refine 1 -interpolate 1 -petscpartitioner_type parmetis -dm_view -petscpartitioner_view -test_redistribute -dm_plex_csr_via_mat {{0 1}} -dm_pre_redist_view ::load_balance -dm_post_redist_view ::load_balance -petscpartitioner_view_graph
546:   test:
547:     suffix: part_ptscotch_0
548:     requires: ptscotch
549:     nsize: 2
550:     args: -dim 2 -cell_simplex 0 -dm_refine 0 -interpolate 1 -petscpartitioner_type ptscotch -petscpartitioner_view -petscpartitioner_ptscotch_strategy quality -test_redistribute -dm_plex_csr_via_mat {{0 1}} -dm_pre_redist_view ::load_balance -dm_post_redist_view ::load_balance -petscpartitioner_view_graph
551:   test:
552:     suffix: part_ptscotch_1
553:     requires: ptscotch
554:     nsize: 8
555:     args: -dim 2 -cell_simplex 0 -dm_refine 1 -interpolate 1 -petscpartitioner_type ptscotch -petscpartitioner_view -petscpartitioner_ptscotch_imbalance 0.1

557:   # CGNS reader tests 10-11 (need to find smaller test meshes)
558:   test:
559:     suffix: cgns_0
560:     requires: cgns
561:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/tut21.cgns -interpolate 1 -dm_view

563:   # Gmsh mesh reader tests
564:   test:
565:     suffix: gmsh_0
566:     requires: !single
567:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh -interpolate 1 -dm_view
568:   test:
569:     suffix: gmsh_1
570:     requires: !single
571:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -interpolate 1 -dm_view
572:   test:
573:     suffix: gmsh_2
574:     requires: !single
575:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin.msh -interpolate 1 -dm_view
576:   test:
577:     suffix: gmsh_3
578:     nsize: 3
579:     requires: !single
580:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -test_partition -interpolate 1 -dm_view
581:   test:
582:     suffix: gmsh_4
583:     nsize: 3
584:     requires: !single
585:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin.msh -test_partition -interpolate 1 -dm_view
586:   test:
587:     suffix: gmsh_5
588:     requires: !single
589:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_quad.msh -interpolate 1 -dm_view
590:   # TODO: it seems the mesh is not a valid gmsh (inverted cell)
591:   test:
592:     suffix: gmsh_6
593:     requires: !single
594:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin_physnames.msh -interpolate 1 -dm_view -final_diagnostics 0
595:   test:
596:     suffix: gmsh_7
597:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere_bin.msh -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
598:   test:
599:     suffix: gmsh_8
600:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere.msh -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
601:   testset:
602:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic_bin.msh -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
603:     test:
604:       suffix: gmsh_9
605:     test:
606:       suffix: gmsh_9_periodic_0
607:       args: -dm_plex_gmsh_periodic 0
608:   testset:
609:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
610:     test:
611:       suffix: gmsh_10
612:     test:
613:       suffix: gmsh_10_periodic_0
614:       args: -dm_plex_gmsh_periodic 0
615:   testset:
616:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all -dm_refine 1
617:     test:
618:       suffix: gmsh_11
619:     test:
620:       suffix: gmsh_11_periodic_0
621:       args: -dm_plex_gmsh_periodic 0
622:   # TODO: it seems the mesh is not a valid gmsh (inverted cell)
623:   test:
624:     suffix: gmsh_12
625:     nsize: 4
626:     requires: !single mpiio
627:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin_physnames.msh -viewer_binary_mpiio -petscpartitioner_type simple -interpolate 1 -dm_view -final_diagnostics 0
628:   test:
629:     suffix: gmsh_13_hybs2t
630:     nsize: 4
631:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_triquad.msh -petscpartitioner_type simple -interpolate 1 -dm_view -dm_refine 1 -dm_plex_cell_refiner tobox -dm_plex_check_all
632:   test:
633:     suffix: gmsh_14_ext
634:     requires: !single
635:     args: -ext_layers 2 -ext_thickness 1.5 -ext_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin.msh -dm_view -interpolate -dm_plex_check_all
636:   test:
637:     suffix: gmsh_14_ext_s2t
638:     requires: !single
639:     args: -ext_layers 2 -ext_thickness 1.5 -ext_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_bin.msh -dm_view -interpolate -dm_plex_check_all -dm_refine 1 -dm_plex_cell_refiner tobox
640:   test:
641:     suffix: gmsh_15_hyb3d
642:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh -dm_view -interpolate -dm_plex_check_all
643:   test:
644:     suffix: gmsh_15_hyb3d_vtk
645:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh -dm_view vtk: -interpolate -dm_plex_gmsh_hybrid -dm_plex_check_all
646:   test:
647:     suffix: gmsh_15_hyb3d_s2t
648:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh -dm_view -interpolate -dm_plex_check_all -dm_refine 1 -dm_plex_cell_refiner tobox
649:   test:
650:     suffix: gmsh_16_spheresurface
651:     nsize : 4
652:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple
653:   test:
654:     suffix: gmsh_16_spheresurface_s2t
655:     nsize : 4
656:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -dm_refine 1 -dm_plex_cell_refiner tobox -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple
657:   test:
658:     suffix: gmsh_16_spheresurface_extruded
659:     nsize : 4
660:     args: -ext_layers 3 -ext_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple
661:   test:
662:     suffix: gmsh_16_spheresurface_extruded_s2t
663:     nsize : 4
664:     args: -ext_layers 3 -ext_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -dm_refine 1 -dm_plex_cell_refiner tobox -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple
665:   test:
666:     suffix: gmsh_17_hyb3d_interp_ascii
667:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_hexwedge.msh -dm_view -interpolate -dm_plex_check_all
668:   test:
669:     suffix: exodus_17_hyb3d_interp_ascii
670:     requires: exodusii
671:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_hexwedge.exo -dm_view -interpolate -dm_plex_check_all

673:   # Legacy Gmsh v22/v40 ascii/binary reader tests
674:   testset:
675:     output_file: output/ex1_gmsh_3d_legacy.out
676:     args: -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
677:     test:
678:       suffix: gmsh_3d_ascii_v22
679:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii.msh2
680:     test:
681:       suffix: gmsh_3d_ascii_v40
682:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii.msh4
683:     test:
684:       suffix: gmsh_3d_binary_v22
685:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary.msh2
686:     test:
687:       suffix: gmsh_3d_binary_v40
688:       requires: long64
689:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary.msh4

691:   # Gmsh v41 ascii/binary reader tests
692:   testset: # 32bit mesh, sequential
693:     args: -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
694:     output_file: output/ex1_gmsh_3d_32.out
695:     test:
696:       suffix: gmsh_3d_ascii_v41_32
697:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii-32.msh
698:     test:
699:       suffix: gmsh_3d_binary_v41_32
700:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-32.msh
701:     test:
702:       suffix: gmsh_3d_binary_v41_32_mpiio
703:       requires: define(PETSC_HAVE_MPIIO)
704:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-32.msh -viewer_binary_mpiio
705:   testset:  # 32bit mesh, parallel
706:     args:  -petscpartitioner_type simple -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
707:     nsize: 2
708:     output_file: output/ex1_gmsh_3d_32_np2.out
709:     test:
710:       suffix: gmsh_3d_ascii_v41_32_np2
711:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii-32.msh
712:     test:
713:       suffix: gmsh_3d_binary_v41_32_np2
714:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-32.msh
715:     test:
716:       suffix: gmsh_3d_binary_v41_32_np2_mpiio
717:       requires: define(PETSC_HAVE_MPIIO)
718:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-32.msh -viewer_binary_mpiio
719:   testset: # 64bit mesh, sequential
720:     args: -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
721:     output_file: output/ex1_gmsh_3d_64.out
722:     test:
723:       suffix: gmsh_3d_ascii_v41_64
724:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii-64.msh
725:     test:
726:       suffix: gmsh_3d_binary_v41_64
727:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-64.msh
728:     test:
729:       suffix: gmsh_3d_binary_v41_64_mpiio
730:       requires: define(PETSC_HAVE_MPIIO)
731:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-64.msh -viewer_binary_mpiio
732:   testset:  # 64bit mesh, parallel
733:     args:  -petscpartitioner_type simple -dm_view ::ascii_info_detail -interpolate -dm_plex_check_all
734:     nsize: 2
735:     output_file: output/ex1_gmsh_3d_64_np2.out
736:     test:
737:       suffix: gmsh_3d_ascii_v41_64_np2
738:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-ascii-64.msh
739:     test:
740:       suffix: gmsh_3d_binary_v41_64_np2
741:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-64.msh
742:     test:
743:       suffix: gmsh_3d_binary_v41_64_np2_mpiio
744:       requires: define(PETSC_HAVE_MPIIO)
745:       args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/gmsh-3d-binary-64.msh -viewer_binary_mpiio

747:   # Fluent mesh reader tests
748:   # TODO: Geometry checks fail
749:   test:
750:     suffix: fluent_0
751:     requires: !complex
752:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.cas -interpolate 1 -dm_view -final_diagnostics 0
753:   test:
754:     suffix: fluent_1
755:     nsize: 3
756:     requires: !complex
757:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.cas -interpolate 1 -test_partition -dm_view -final_diagnostics 0
758:   test:
759:     suffix: fluent_2
760:     requires: !complex
761:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_5tets_ascii.cas -interpolate 1 -dm_view -final_diagnostics 0
762:   test:
763:     suffix: fluent_3
764:     requires: !complex
765:     TODO: Fails on non-linux: fseek(), fileno() ? https://gitlab.com/petsc/petsc/merge_requests/2206#note_238166382
766:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cube_5tets.cas -interpolate 1 -dm_view -final_diagnostics 0

768:   # Med mesh reader tests, including parallel file reads
769:   test:
770:     suffix: med_0
771:     requires: med
772:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.med -interpolate 1 -dm_view
773:   test:
774:     suffix: med_1
775:     requires: med
776:     nsize: 3
777:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.med -interpolate 1 -petscpartitioner_type simple -dm_view
778:   test:
779:     suffix: med_2
780:     requires: med
781:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cylinder.med -interpolate 1 -dm_view
782:   test:
783:     suffix: med_3
784:     requires: med
785:     TODO: MED
786:     nsize: 3
787:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/cylinder.med -interpolate 1 -petscpartitioner_type simple -dm_view

789:   # Test shape quality
790:   test:
791:     suffix: test_shape
792:     requires: ctetgen
793:     args: -dim 3 -interpolate -dm_refine_hierarchy 3 -dm_plex_check_all -dm_plex_check_cell_shape

795:   # Test simplex to tensor conversion
796:   test:
797:     suffix: s2t2
798:     requires: triangle
799:     args: -dim 2 -dm_refine 1 -interpolate -dm_plex_cell_refiner tobox -refinement_limit 0.0625 -dm_view ascii::ascii_info_detail

801:   test:
802:     suffix: s2t3
803:     requires: ctetgen
804:     args: -dim 3 -dm_refine 1 -interpolate -dm_plex_cell_refiner tobox -refinement_limit 0.0625 -dm_view ascii::ascii_info_detail

806:   # Test domain shapes
807:   test:
808:     suffix: cylinder
809:     args: -dim 3 -cell_simplex 0 -interpolate -domain_shape cylinder -dm_plex_check_all -dm_view

811:   test:
812:     suffix: cylinder_per
813:     args: -dim 3 -cell_simplex 0 -interpolate -domain_shape cylinder -z_periodicity periodic -dm_plex_check_all -dm_view

815:   test:
816:     suffix: cylinder_wedge
817:     args: -dim 3 -cell_simplex 0 -interpolate -cell_wedge -domain_shape cylinder -dm_view vtk: -dm_plex_check_all

819:   test:
820:     suffix: cylinder_wedge_int
821:     output_file: output/ex1_cylinder_wedge.out
822:     args: -dim 3 -cell_simplex 0 -interpolate -cell_wedge -domain_shape cylinder -dm_view vtk: -dm_plex_check_all

824:   test:
825:     suffix: box_2d
826:     args: -dim 2 -cell_simplex 0 -interpolate -domain_shape box -dm_refine 2 -dm_plex_check_all -dm_view

828:   test:
829:     suffix: box_2d_per
830:     args: -dim 2 -cell_simplex 0 -interpolate -domain_shape box -dm_refine 2 -dm_plex_check_all -dm_view

832:   test:
833:     suffix: box_2d_per_unint
834:     args: -dim 2 -cell_simplex 0 -interpolate 0 -domain_shape box -domain_box_sizes 3,3 -dm_plex_check_all -dm_view ::ascii_info_detail

836:   test:
837:     suffix: box_3d
838:     args: -dim 3 -cell_simplex 0 -interpolate -domain_shape box -dm_refine 3 -dm_plex_check_all -dm_view

840:   test:
841:     requires: triangle
842:     suffix: box_wedge
843:     args: -dim 3 -cell_simplex 0 -interpolate -cell_wedge -domain_shape box -dm_view vtk: -dm_plex_check_all

845:   testset:
846:     requires: triangle
847:     args: -dim 3 -cell_simplex 0 -interpolate -cell_wedge -domain_shape box -domain_box_sizes 2,3,1 -dm_view -dm_plex_check_all -dm_refine 1 -dm_plex_cell_refiner tobox
848:     test:
849:       suffix: box_wedge_s2t
850:     test:
851:       nsize: 3
852:       args: -petscpartitioner_type simple
853:       suffix: box_wedge_s2t_parallel

855:   # Test GLVis output
856:   test:
857:     suffix: glvis_2d_tet
858:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_gmsh_periodic 0 -dm_view glvis:

860:   test:
861:     suffix: glvis_2d_tet_per
862:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary 0

864:   test:
865:     suffix: glvis_2d_tet_per_mfem
866:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -viewer_glvis_dm_plex_enable_boundary -viewer_glvis_dm_plex_enable_mfem -dm_view glvis: -interpolate

868:   test:
869:     suffix: glvis_2d_quad
870:     args: -dim 2 -cell_simplex 0 -interpolate -domain_shape box -domain_box_sizes 3,3 -dm_view glvis:

872:   test:
873:     suffix: glvis_2d_quad_per
874:     args: -dim 2 -cell_simplex 0 -interpolate -domain_shape box -domain_box_sizes 3,3 -x_periodicity periodic -y_periodicity periodic -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary

876:   test:
877:     suffix: glvis_2d_quad_per_mfem
878:     args: -dim 2 -cell_simplex 0 -interpolate -domain_shape box -domain_box_sizes 3,3 -x_periodicity periodic -y_periodicity periodic -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary -viewer_glvis_dm_plex_enable_mfem

880:   test:
881:     suffix: glvis_3d_tet
882:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere_bin.msh -dm_plex_gmsh_periodic 0 -dm_view glvis:

884:   test:
885:     suffix: glvis_3d_tet_per
886:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere_bin.msh -dm_view glvis: -interpolate -viewer_glvis_dm_plex_enable_boundary

888:   test:
889:     suffix: glvis_3d_tet_per_mfem
890:     TODO: broken
891:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere_bin.msh -viewer_glvis_dm_plex_enable_mfem -dm_view glvis: -interpolate

893:   test:
894:     suffix: glvis_3d_hex
895:     args: -dim 3 -cell_simplex 0 -interpolate -domain_shape box -domain_box_sizes 3,3,3 -dm_view glvis:

897:   test:
898:     suffix: glvis_3d_hex_per
899:     args: -dim 3 -cell_simplex 0 -interpolate -domain_shape box -domain_box_sizes 3,3,3 -x_periodicity periodic -y_periodicity periodic -z_periodicity periodic -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary 0

901:   test:
902:     suffix: glvis_3d_hex_per_mfem
903:     args: -dim 3 -cell_simplex 0 -domain_shape box -domain_box_sizes 3,3,3 -x_periodicity periodic -y_periodicity periodic -z_periodicity periodic -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary -viewer_glvis_dm_plex_enable_mfem -interpolate

905:   # Test P4EST
906:   testset:
907:     requires: p4est
908:     args: -interpolate -dm_view -test_p4est_seq -conv_seq_2_dm_plex_check_all -conv_seq_1_dm_forest_minimum_refinement 1
909:     test:
910:       suffix: p4est_periodic
911:       args: -dim 2 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -domain_box_sizes 3,5 -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 2 -conv_seq_1_dm_p4est_refine_pattern hash
912:     test:
913:       suffix: p4est_periodic_3d
914:       args: -dim 3 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -z_periodicity -domain_box_sizes 3,5,4 -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 2 -conv_seq_1_dm_p4est_refine_pattern hash
915:     test:
916:       suffix: p4est_gmsh_periodic
917:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
918:     test:
919:       suffix: p4est_gmsh_surface
920:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3
921:     test:
922:       suffix: p4est_gmsh_surface_parallel
923:       nsize: 2
924:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -petscpartitioner_type simple -dm_view ::load_balance
925:     test:
926:       suffix: p4est_hyb_2d
927:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_triquad.msh
928:     test:
929:       suffix: p4est_hyb_3d
930:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh
931:     test:
932:       requires: ctetgen
933:       suffix: p4est_s2t_bugfaces_3d
934:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 0 -dim 3 -domain_box_sizes 1,1 -cell_simplex
935:     test:
936:       suffix: p4est_bug_overlapsf
937:       nsize: 3
938:       args: -dim 3 -cell_simplex 0 -domain_box_sizes 2,2,1 -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash  -petscpartitioner_type simple
939:     test:
940:       suffix: p4est_redistribute
941:       nsize: 3
942:       args: -dim 3 -cell_simplex 0 -domain_box_sizes 2,2,1 -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash  -petscpartitioner_type simple -test_redistribute -dm_plex_csr_via_mat {{0 1}} -dm_view ::load_balance
943:     test:
944:       suffix: p4est_gmsh_s2t_3d
945:       args: -conv_seq_1_dm_forest_initial_refinement 1 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
946:     test:
947:       suffix: p4est_gmsh_s2t_3d_hash
948:       args: -conv_seq_1_dm_forest_initial_refinement 1 -conv_seq_1_dm_forest_maximum_refinement 2 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
949:     test:
950:       requires: long_runtime
951:       suffix: p4est_gmsh_periodic_3d
952:       args: -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 1 -conv_seq_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere.msh

954:   testset:
955:     requires: p4est
956:     nsize: 6
957:     args: -interpolate -test_p4est_par -conv_par_2_dm_plex_check_all -conv_par_1_dm_forest_minimum_refinement 1 -conv_par_1_dm_forest_partition_overlap 0
958:     test:
959:       TODO: interface cones do not conform
960:       suffix: p4est_par_periodic
961:       args: -dim 2 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -domain_box_sizes 3,5 -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash
962:     test:
963:       TODO: interface cones do not conform
964:       suffix: p4est_par_periodic_3d
965:       args: -dim 3 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -z_periodicity periodic -domain_box_sizes 3,5,4 -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash
966:     test:
967:       TODO: interface cones do not conform
968:       suffix: p4est_par_gmsh_periodic
969:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
970:     test:
971:       suffix: p4est_par_gmsh_surface
972:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3
973:     test:
974:       suffix: p4est_par_gmsh_s2t_3d
975:       args: -conv_par_1_dm_forest_initial_refinement 1 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
976:     test:
977:       TODO: interface cones do not conform
978:       suffix: p4est_par_gmsh_s2t_3d_hash
979:       args: -conv_par_1_dm_forest_initial_refinement 1 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
980:     test:
981:       requires: long_runtime
982:       suffix: p4est_par_gmsh_periodic_3d
983:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere.msh

985:   testset:
986:     requires: p4est
987:     nsize: 6
988:     args: -interpolate -test_p4est_par -conv_par_2_dm_plex_check_all -conv_par_1_dm_forest_minimum_refinement 1 -conv_par_1_dm_forest_partition_overlap 1 -petscpartitioner_type simple
989:     test:
990:       suffix: p4est_par_ovl_periodic
991:       args: -dim 2 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -domain_box_sizes 3,5 -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash
992:     #TODO Mesh cell 201 is inverted, vol = 0. (FVM Volume. Is it correct? -> Diagnostics disabled)
993:     test:
994:       suffix: p4est_par_ovl_periodic_3d
995:       args: -dim 3 -domain_shape box -cell_simplex 0 -x_periodicity periodic -y_periodicity periodic -z_periodicity -domain_box_sizes 3,5,4 -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash -final_diagnostics 0
996:     test:
997:       suffix: p4est_par_ovl_gmsh_periodic
998:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
999:     test:
1000:       suffix: p4est_par_ovl_gmsh_surface
1001:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3
1002:     test:
1003:       suffix: p4est_par_ovl_gmsh_s2t_3d
1004:       args: -conv_par_1_dm_forest_initial_refinement 1 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
1005:     test:
1006:       suffix: p4est_par_ovl_gmsh_s2t_3d_hash
1007:       args: -conv_par_1_dm_forest_initial_refinement 1 -conv_par_1_dm_forest_maximum_refinement 2 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/doublet-tet.msh
1008:     test:
1009:       requires: long_runtime
1010:       suffix: p4est_par_ovl_gmsh_periodic_3d
1011:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/mesh-3d-box-innersphere.msh
1012:     test:
1013:       suffix: p4est_par_ovl_hyb_2d
1014:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_triquad.msh
1015:     test:
1016:       suffix: p4est_par_ovl_hyb_3d
1017:       args: -conv_par_1_dm_forest_initial_refinement 0 -conv_par_1_dm_forest_maximum_refinement 1 -conv_par_1_dm_p4est_refine_pattern hash -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh

1019:   test:
1020:     TODO: broken
1021:     requires: p4est
1022:     nsize: 2
1023:     suffix: p4est_bug_labels_noovl
1024:     args: -interpolate -test_p4est_seq -dm_plex_check_all -dm_forest_minimum_refinement 0 -dm_forest_partition_overlap 1 -dim 2 -domain_shape box -cell_simplex 0 -domain_box_sizes 3,3 -dm_forest_initial_refinement 0 -dm_forest_maximum_refinement 2 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -dm_forest_print_label_error

1026:   test:
1027:     requires: p4est
1028:     nsize: 2
1029:     suffix: p4est_bug_distribute_overlap
1030:     args: -interpolate -test_p4est_seq -conv_seq_2_dm_plex_check_all -conv_seq_1_dm_forest_minimum_refinement 0 -conv_seq_1_dm_forest_partition_overlap 0 -dim 2 -domain_shape box -cell_simplex 0 -domain_box_sizes 3,3 -conv_seq_1_dm_forest_initial_refinement 0 -conv_seq_1_dm_forest_maximum_refinement 2 -conv_seq_1_dm_p4est_refine_pattern hash -petscpartitioner_type simple -overlap 1 -dm_view ::load_balance
1031:     args: -dm_post_overlap_view

1033:   test:
1034:     suffix: glvis_2d_hyb
1035:     args: -dim 2 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_triquad.msh -interpolate -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary -petscpartitioner_type simple

1037:   test:
1038:     suffix: glvis_3d_hyb
1039:     args: -dim 3 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_tetwedge.msh -interpolate -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary -petscpartitioner_type simple

1041:   test:
1042:     suffix: glvis_3d_hyb_s2t
1043:     args: -dim 3 -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_3d_cube.msh -interpolate -dm_view glvis: -viewer_glvis_dm_plex_enable_boundary -petscpartitioner_type simple -dm_refine 1 -dm_plex_cell_refiner tobox -dm_plex_check_all

1045:   test:
1046:     suffix: ref_alfeld2d_0
1047:     requires: triangle
1048:     args: -dim 2 -domain_shape box -cell_simplex 1 -domain_box_sizes 5 -dm_view -interpolate -dm_plex_check_all -dm_refine 1 -dm_plex_cell_refiner alfeld2d -final_diagnostics
1049:   test:
1050:     suffix: ref_alfeld3d_0
1051:     requires: ctetgen
1052:     args: -dim 3 -domain_shape box -cell_simplex 1 -domain_box_sizes 5 -dm_view -interpolate -dm_plex_check_all -dm_refine 1 -dm_plex_cell_refiner alfeld3d -final_diagnostics

1054:   # Boundary layer refiners
1055:   test:
1056:     suffix: ref_bl_1
1057:     args: -dim 1 -domain_shape box -cell_simplex 0 -domain_box_sizes 5 -dm_view -interpolate -dm_plex_check_all 0 -dm_refine 1 -dm_plex_cell_refiner boundarylayer -ext_layers 2 -final_diagnostics -dm_plex_refine_boundarylayer_splits 3 -ext_hfirst {{0 1}}
1058:   test:
1059:     suffix: ref_bl_2_tri
1060:     requires: triangle
1061:     args: -dim 2 -domain_shape box -cell_simplex 1 -domain_box_sizes 5 -dm_view -interpolate -dm_plex_check_all 0 -dm_refine 1 -dm_plex_cell_refiner boundarylayer -ext_layers 3 -final_diagnostics -dm_plex_refine_boundarylayer_splits 4 -ext_hfirst {{0 1}}
1062:   test:
1063:     suffix: ref_bl_3_quad
1064:     args: -dim 2 -domain_shape box -cell_simplex 0 -domain_box_sizes 5 -dm_view -interpolate -dm_plex_check_all 0 -dm_refine 1 -dm_plex_cell_refiner boundarylayer -ext_layers 3 -final_diagnostics -dm_plex_refine_boundarylayer_splits 4 -ext_hfirst {{0 1}}
1065:   test:
1066:     suffix: ref_bl_spheresurface_extruded
1067:     nsize : 4
1068:     args: -ext_layers 3 -ext_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/surfacesphere_bin.msh -dm_plex_gmsh_spacedim 3 -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple -ext_hfirst {{0 1}separate output} -final_diagnostics -dm_refine 1 -dm_plex_cell_refiner boundarylayer -dm_plex_refine_boundarylayer_splits 2
1069:   test:
1070:     suffix: ref_bl_3d_hyb
1071:     nsize : 4
1072:     args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/hybrid_3d_cube.msh -dm_plex_check_all -dm_view -interpolate -petscpartitioner_type simple -final_diagnostics -dm_refine 1 -dm_plex_cell_refiner boundarylayer -dm_plex_refine_boundarylayer_splits 4 -dm_plex_refine_boundarylayer_progression 3.1

1074:   test:
1075:     suffix: sphere_0
1076:     args: -dim 2 -domain_shape sphere -dm_plex_check_all -dm_view ::ascii_info_detail

1078:   test:
1079:     suffix: sphere_1
1080:     args: -dim 2 -domain_shape sphere -dm_plex_check_all -dm_refine 2 -dm_view

1082:   test:
1083:     suffix: sphere_2
1084:     args: -dim 2 -cell_simplex 0 -domain_shape sphere -dm_plex_check_all -dm_view ::ascii_info_detail

1086:   test:
1087:     suffix: sphere_3
1088:     args: -dim 2 -cell_simplex 0 -domain_shape sphere -dm_plex_check_all -dm_refine 2 -dm_view

1090:   test:
1091:     suffix: ball_0
1092:     requires: ctetgen
1093:     args: -dim 3 -domain_shape ball -dm_plex_check_all -dm_view

1095:   test:
1096:     suffix: ball_1
1097:     requires: ctetgen
1098:     args: -dim 3 -domain_shape ball -bd_dm_refine 2 -dm_plex_check_all -dm_view

1100: TEST*/