Actual source code: ex1.c
petsc-master 2020-08-25
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*/