Actual source code: ex18.c
1: static char help[] = "Tests for parallel mesh loading and parallel topological interpolation\n\n";
3: #include <petsc/private/dmpleximpl.h>
4: /* List of test meshes
6: Network
7: -------
8: Test 0 (2 ranks):
10: network=0:
11: ---------
12: cell 0 cell 1 cell 2 nCells-1 (edge)
13: 0 ------ 1 ------ 2 ------ 3 -- -- v -- -- nCells (vertex)
15: vertex distribution:
16: rank 0: 0 1
17: rank 1: 2 3 ... nCells
18: cell(edge) distribution:
19: rank 0: 0 1
20: rank 1: 2 ... nCells-1
22: network=1:
23: ---------
24: v2
25: ^
26: |
27: cell 2
28: |
29: v0 --cell 0--> v3--cell 1--> v1
31: vertex distribution:
32: rank 0: 0 1 3
33: rank 1: 2
34: cell(edge) distribution:
35: rank 0: 0 1
36: rank 1: 2
38: example:
39: mpiexec -n 2 ./ex18 -distribute 1 -dim 1 -orig_dm_view -dist_dm_view -dist_dm_view -petscpartitioner_type parmetis -ncells 50
41: Triangle
42: --------
43: Test 0 (2 ranks):
44: Two triangles sharing a face
46: 2
47: / | \
48: / | \
49: / | \
50: 0 0 | 1 3
51: \ | /
52: \ | /
53: \ | /
54: 1
56: vertex distribution:
57: rank 0: 0 1
58: rank 1: 2 3
59: cell distribution:
60: rank 0: 0
61: rank 1: 1
63: Test 1 (3 ranks):
64: Four triangles partitioned across 3 ranks
66: 0 _______ 3
67: | \ / |
68: | \ 1 / |
69: | \ / |
70: | 0 2 2 |
71: | / \ |
72: | / 3 \ |
73: | / \ |
74: 1 ------- 4
76: vertex distribution:
77: rank 0: 0 1
78: rank 1: 2 3
79: rank 2: 4
80: cell distribution:
81: rank 0: 0
82: rank 1: 1
83: rank 2: 2 3
85: Test 2 (3 ranks):
86: Four triangles partitioned across 3 ranks
88: 1 _______ 3
89: | \ / |
90: | \ 1 / |
91: | \ / |
92: | 0 0 2 |
93: | / \ |
94: | / 3 \ |
95: | / \ |
96: 2 ------- 4
98: vertex distribution:
99: rank 0: 0 1
100: rank 1: 2 3
101: rank 2: 4
102: cell distribution:
103: rank 0: 0
104: rank 1: 1
105: rank 2: 2 3
107: Tetrahedron
108: -----------
109: Test 0:
110: Two tets sharing a face
112: cell 3 _______ cell
113: 0 / | \ \ 1
114: / | \ \
115: / | \ \
116: 0----|----4-----2
117: \ | / /
118: \ | / /
119: \ | / /
120: 1-------
121: y
122: | x
123: |/
124: *----z
126: vertex distribution:
127: rank 0: 0 1
128: rank 1: 2 3 4
129: cell distribution:
130: rank 0: 0
131: rank 1: 1
133: Quadrilateral
134: -------------
135: Test 0 (2 ranks):
136: Two quads sharing a face
138: 3-------2-------5
139: | | |
140: | 0 | 1 |
141: | | |
142: 0-------1-------4
144: vertex distribution:
145: rank 0: 0 1 2
146: rank 1: 3 4 5
147: cell distribution:
148: rank 0: 0
149: rank 1: 1
151: TODO Test 1:
152: A quad and a triangle sharing a face
154: 5-------4
155: | | \
156: | 0 | \
157: | | 1 \
158: 2-------3----6
160: Hexahedron
161: ----------
162: Test 0 (2 ranks):
163: Two hexes sharing a face
165: cell 7-------------6-------------11 cell
166: 0 /| /| /| 1
167: / | F1 / | F7 / |
168: / | / | / |
169: 4-------------5-------------10 |
170: | | F4 | | F10 | |
171: | | | | | |
172: |F5 | |F3 | |F9 |
173: | | F2 | | F8 | |
174: | 3---------|---2---------|---9
175: | / | / | /
176: | / F0 | / F6 | /
177: |/ |/ |/
178: 0-------------1-------------8
180: vertex distribution:
181: rank 0: 0 1 2 3 4 5
182: rank 1: 6 7 8 9 10 11
183: cell distribution:
184: rank 0: 0
185: rank 1: 1
187: */
189: typedef enum {
190: NONE,
191: CREATE,
192: AFTER_CREATE,
193: AFTER_DISTRIBUTE
194: } InterpType;
196: typedef struct {
197: PetscInt debug; /* The debugging level */
198: PetscInt testNum; /* Indicates the mesh to create */
199: PetscInt dim; /* The topological mesh dimension */
200: PetscBool cellSimplex; /* Use simplices or hexes */
201: PetscBool distribute; /* Distribute the mesh */
202: InterpType interpolate; /* Interpolate the mesh before or after DMPlexDistribute() */
203: PetscBool useGenerator; /* Construct mesh with a mesh generator */
204: PetscBool testOrientIF; /* Test for different original interface orientations */
205: PetscBool testHeavy; /* Run the heavy PointSF test */
206: PetscBool customView; /* Show results of DMPlexIsInterpolated() etc. */
207: PetscInt ornt[2]; /* Orientation of interface on rank 0 and rank 1 */
208: PetscInt faces[3]; /* Number of faces per dimension for generator */
209: PetscScalar coords[128];
210: PetscReal coordsTol;
211: PetscInt ncoords;
212: PetscInt pointsToExpand[128];
213: PetscInt nPointsToExpand;
214: PetscBool testExpandPointsEmpty;
215: char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
216: } AppCtx;
218: struct _n_PortableBoundary {
219: Vec coordinates;
220: PetscInt depth;
221: PetscSection *sections;
222: };
223: typedef struct _n_PortableBoundary *PortableBoundary;
225: static PetscLogStage stage[3];
227: static PetscErrorCode DMPlexCheckPointSFHeavy(DM, PortableBoundary);
228: static PetscErrorCode DMPlexSetOrientInterface_Private(DM, PetscBool);
229: static PetscErrorCode DMPlexGetExpandedBoundary_Private(DM, PortableBoundary *);
230: static PetscErrorCode DMPlexExpandedConesToFaces_Private(DM, IS, PetscSection, IS *);
232: static PetscErrorCode PortableBoundaryDestroy(PortableBoundary *bnd)
233: {
234: PetscInt d;
236: PetscFunctionBegin;
237: if (!*bnd) PetscFunctionReturn(PETSC_SUCCESS);
238: PetscCall(VecDestroy(&(*bnd)->coordinates));
239: for (d = 0; d < (*bnd)->depth; d++) PetscCall(PetscSectionDestroy(&(*bnd)->sections[d]));
240: PetscCall(PetscFree((*bnd)->sections));
241: PetscCall(PetscFree(*bnd));
242: PetscFunctionReturn(PETSC_SUCCESS);
243: }
245: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
246: {
247: const char *interpTypes[4] = {"none", "create", "after_create", "after_distribute"};
248: PetscInt interp = NONE, dim;
249: PetscBool flg1, flg2;
251: PetscFunctionBegin;
252: options->debug = 0;
253: options->testNum = 0;
254: options->dim = 2;
255: options->cellSimplex = PETSC_TRUE;
256: options->distribute = PETSC_FALSE;
257: options->interpolate = NONE;
258: options->useGenerator = PETSC_FALSE;
259: options->testOrientIF = PETSC_FALSE;
260: options->testHeavy = PETSC_TRUE;
261: options->customView = PETSC_FALSE;
262: options->testExpandPointsEmpty = PETSC_FALSE;
263: options->ornt[0] = 0;
264: options->ornt[1] = 0;
265: options->faces[0] = 2;
266: options->faces[1] = 2;
267: options->faces[2] = 2;
268: options->filename[0] = '\0';
269: options->coordsTol = PETSC_DEFAULT;
271: PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");
272: PetscCall(PetscOptionsBoundedInt("-debug", "The debugging level", "ex18.c", options->debug, &options->debug, NULL, 0));
273: PetscCall(PetscOptionsBoundedInt("-testnum", "The mesh to create", "ex18.c", options->testNum, &options->testNum, NULL, 0));
274: PetscCall(PetscOptionsBool("-cell_simplex", "Generate simplices if true, otherwise hexes", "ex18.c", options->cellSimplex, &options->cellSimplex, NULL));
275: PetscCall(PetscOptionsBool("-distribute", "Distribute the mesh", "ex18.c", options->distribute, &options->distribute, NULL));
276: PetscCall(PetscOptionsEList("-interpolate", "Type of mesh interpolation (none, create, after_create, after_distribute)", "ex18.c", interpTypes, 4, interpTypes[options->interpolate], &interp, NULL));
277: options->interpolate = (InterpType)interp;
278: PetscCheck(options->distribute || options->interpolate != AFTER_DISTRIBUTE, comm, PETSC_ERR_SUP, "-interpolate after_distribute needs -distribute 1");
279: PetscCall(PetscOptionsBool("-use_generator", "Use a mesh generator to build the mesh", "ex18.c", options->useGenerator, &options->useGenerator, NULL));
280: options->ncoords = 128;
281: PetscCall(PetscOptionsScalarArray("-view_vertices_from_coords", "Print DAG points corresponding to vertices with given coordinates", "ex18.c", options->coords, &options->ncoords, NULL));
282: PetscCall(PetscOptionsReal("-view_vertices_from_coords_tol", "Tolerance for -view_vertices_from_coords", "ex18.c", options->coordsTol, &options->coordsTol, NULL));
283: options->nPointsToExpand = 128;
284: PetscCall(PetscOptionsIntArray("-test_expand_points", "Expand given array of DAG point using DMPlexGetConeRecursive() and print results", "ex18.c", options->pointsToExpand, &options->nPointsToExpand, NULL));
285: if (options->nPointsToExpand) PetscCall(PetscOptionsBool("-test_expand_points_empty", "For -test_expand_points, rank 0 will have empty input array", "ex18.c", options->testExpandPointsEmpty, &options->testExpandPointsEmpty, NULL));
286: PetscCall(PetscOptionsBool("-test_heavy", "Run the heavy PointSF test", "ex18.c", options->testHeavy, &options->testHeavy, NULL));
287: PetscCall(PetscOptionsBool("-custom_view", "Custom DMPlex view", "ex18.c", options->customView, &options->customView, NULL));
288: PetscCall(PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex18.c", options->dim, &options->dim, &flg1, 1, 3));
289: dim = 3;
290: PetscCall(PetscOptionsIntArray("-faces", "Number of faces per dimension", "ex18.c", options->faces, &dim, &flg2));
291: if (flg2) {
292: PetscCheck(!flg1 || dim == options->dim, comm, PETSC_ERR_ARG_OUTOFRANGE, "specified -dim %" PetscInt_FMT " is not equal to length %" PetscInt_FMT " of -faces (note that -dim can be omitted)", options->dim, dim);
293: options->dim = dim;
294: }
295: PetscCall(PetscOptionsString("-filename", "The mesh file", "ex18.c", options->filename, options->filename, sizeof(options->filename), NULL));
296: PetscCall(PetscOptionsBoundedInt("-rotate_interface_0", "Rotation (relative orientation) of interface on rank 0; implies -interpolate create -distribute 0", "ex18.c", options->ornt[0], &options->ornt[0], &options->testOrientIF, 0));
297: PetscCall(PetscOptionsBoundedInt("-rotate_interface_1", "Rotation (relative orientation) of interface on rank 1; implies -interpolate create -distribute 0", "ex18.c", options->ornt[1], &options->ornt[1], &flg2, 0));
298: PetscCheck(flg2 == options->testOrientIF, comm, PETSC_ERR_ARG_OUTOFRANGE, "neither or both -rotate_interface_0 -rotate_interface_1 must be set");
299: if (options->testOrientIF) {
300: PetscInt i;
301: for (i = 0; i < 2; i++) {
302: if (options->ornt[i] >= 10) options->ornt[i] = -(options->ornt[i] - 10); /* 11 12 13 become -1 -2 -3 */
303: }
304: options->filename[0] = 0;
305: options->useGenerator = PETSC_FALSE;
306: options->dim = 3;
307: options->cellSimplex = PETSC_TRUE;
308: options->interpolate = CREATE;
309: options->distribute = PETSC_FALSE;
310: }
311: PetscOptionsEnd();
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
315: static PetscErrorCode CreateMesh_1D(MPI_Comm comm, PetscBool interpolate, AppCtx *user, DM *dm)
316: {
317: PetscInt testNum = user->testNum;
318: PetscMPIInt rank, size;
319: PetscInt numCorners = 2, i;
320: PetscInt numCells, numVertices, network;
321: PetscInt *cells;
322: PetscReal *coords;
324: PetscFunctionBegin;
325: PetscCallMPI(MPI_Comm_rank(comm, &rank));
326: PetscCallMPI(MPI_Comm_size(comm, &size));
327: PetscCheck(size <= 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for <=2 processes", testNum);
329: numCells = 3;
330: PetscCall(PetscOptionsGetInt(NULL, NULL, "-ncells", &numCells, NULL));
331: PetscCheck(numCells >= 3, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test ncells %" PetscInt_FMT " must >=3", numCells);
333: if (size == 1) {
334: numVertices = numCells + 1;
335: PetscCall(PetscMalloc2(2 * numCells, &cells, 2 * numVertices, &coords));
336: for (i = 0; i < numCells; i++) {
337: cells[2 * i] = i;
338: cells[2 * i + 1] = i + 1;
339: coords[2 * i] = i;
340: coords[2 * i + 1] = i + 1;
341: }
343: PetscCall(DMPlexCreateFromCellListPetsc(comm, user->dim, numCells, numVertices, numCorners, PETSC_FALSE, cells, user->dim, coords, dm));
344: PetscCall(PetscFree2(cells, coords));
345: PetscFunctionReturn(PETSC_SUCCESS);
346: }
348: network = 0;
349: PetscCall(PetscOptionsGetInt(NULL, NULL, "-network_case", &network, NULL));
350: if (network == 0) {
351: switch (rank) {
352: case 0: {
353: numCells = 2;
354: numVertices = numCells;
355: PetscCall(PetscMalloc2(2 * numCells, &cells, 2 * numCells, &coords));
356: cells[0] = 0;
357: cells[1] = 1;
358: cells[2] = 1;
359: cells[3] = 2;
360: coords[0] = 0.;
361: coords[1] = 1.;
362: coords[2] = 1.;
363: coords[3] = 2.;
364: } break;
365: case 1: {
366: numCells -= 2;
367: numVertices = numCells + 1;
368: PetscCall(PetscMalloc2(2 * numCells, &cells, 2 * numCells, &coords));
369: for (i = 0; i < numCells; i++) {
370: cells[2 * i] = 2 + i;
371: cells[2 * i + 1] = 2 + i + 1;
372: coords[2 * i] = 2 + i;
373: coords[2 * i + 1] = 2 + i + 1;
374: }
375: } break;
376: default:
377: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
378: }
379: } else { /* network_case = 1 */
380: /* ----------------------- */
381: switch (rank) {
382: case 0: {
383: numCells = 2;
384: numVertices = 3;
385: PetscCall(PetscMalloc2(2 * numCells, &cells, 2 * numCells, &coords));
386: cells[0] = 0;
387: cells[1] = 3;
388: cells[2] = 3;
389: cells[3] = 1;
390: } break;
391: case 1: {
392: numCells = 1;
393: numVertices = 1;
394: PetscCall(PetscMalloc2(2 * numCells, &cells, 2 * numCells, &coords));
395: cells[0] = 3;
396: cells[1] = 2;
397: } break;
398: default:
399: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
400: }
401: }
402: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, PETSC_FALSE, cells, user->dim, coords, NULL, NULL, dm));
403: PetscCall(PetscFree2(cells, coords));
404: PetscFunctionReturn(PETSC_SUCCESS);
405: }
407: static PetscErrorCode CreateSimplex_2D(MPI_Comm comm, PetscBool interpolate, AppCtx *user, DM *dm)
408: {
409: PetscInt testNum = user->testNum, p;
410: PetscMPIInt rank, size;
412: PetscFunctionBegin;
413: PetscCallMPI(MPI_Comm_rank(comm, &rank));
414: PetscCallMPI(MPI_Comm_size(comm, &size));
415: switch (testNum) {
416: case 0:
417: PetscCheck(size == 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 2 processes", testNum);
418: switch (rank) {
419: case 0: {
420: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
421: const PetscInt cells[3] = {0, 1, 2};
422: PetscReal coords[4] = {-0.5, 0.5, 0.0, 0.0};
423: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
425: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
426: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
427: } break;
428: case 1: {
429: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
430: const PetscInt cells[3] = {1, 3, 2};
431: PetscReal coords[4] = {0.0, 1.0, 0.5, 0.5};
432: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
434: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
435: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
436: } break;
437: default:
438: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
439: }
440: break;
441: case 1:
442: PetscCheck(size == 3, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 3 processes", testNum);
443: switch (rank) {
444: case 0: {
445: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
446: const PetscInt cells[3] = {0, 1, 2};
447: PetscReal coords[4] = {0.0, 1.0, 0.0, 0.0};
448: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
450: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
451: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
452: } break;
453: case 1: {
454: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
455: const PetscInt cells[3] = {0, 2, 3};
456: PetscReal coords[4] = {0.5, 0.5, 1.0, 1.0};
457: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
459: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
460: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
461: } break;
462: case 2: {
463: const PetscInt numCells = 2, numVertices = 1, numCorners = 3;
464: const PetscInt cells[6] = {2, 4, 3, 2, 1, 4};
465: PetscReal coords[2] = {1.0, 0.0};
466: PetscInt markerPoints[10] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};
468: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
469: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
470: } break;
471: default:
472: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
473: }
474: break;
475: case 2:
476: PetscCheck(size == 3, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 3 processes", testNum);
477: switch (rank) {
478: case 0: {
479: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
480: const PetscInt cells[3] = {1, 2, 0};
481: PetscReal coords[4] = {0.5, 0.5, 0.0, 1.0};
482: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
484: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
485: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
486: } break;
487: case 1: {
488: const PetscInt numCells = 1, numVertices = 2, numCorners = 3;
489: const PetscInt cells[3] = {1, 0, 3};
490: PetscReal coords[4] = {0.0, 0.0, 1.0, 1.0};
491: PetscInt markerPoints[6] = {1, 1, 2, 1, 3, 1};
493: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
494: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
495: } break;
496: case 2: {
497: const PetscInt numCells = 2, numVertices = 1, numCorners = 3;
498: const PetscInt cells[6] = {0, 4, 3, 0, 2, 4};
499: PetscReal coords[2] = {1.0, 0.0};
500: PetscInt markerPoints[10] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};
502: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
503: for (p = 0; p < 3; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
504: } break;
505: default:
506: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
507: }
508: break;
509: default:
510: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
511: }
512: PetscFunctionReturn(PETSC_SUCCESS);
513: }
515: static PetscErrorCode CreateSimplex_3D(MPI_Comm comm, PetscBool interpolate, AppCtx *user, DM *dm)
516: {
517: PetscInt testNum = user->testNum, p;
518: PetscMPIInt rank, size;
520: PetscFunctionBegin;
521: PetscCallMPI(MPI_Comm_rank(comm, &rank));
522: PetscCallMPI(MPI_Comm_size(comm, &size));
523: switch (testNum) {
524: case 0:
525: PetscCheck(size == 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 2 processes", testNum);
526: switch (rank) {
527: case 0: {
528: const PetscInt numCells = 1, numVertices = 2, numCorners = 4;
529: const PetscInt cells[4] = {0, 2, 1, 3};
530: PetscReal coords[6] = {0.0, 0.0, -0.5, 0.0, -0.5, 0.0};
531: PetscInt markerPoints[8] = {1, 1, 2, 1, 3, 1, 4, 1};
533: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
534: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
535: } break;
536: case 1: {
537: const PetscInt numCells = 1, numVertices = 3, numCorners = 4;
538: const PetscInt cells[4] = {1, 2, 4, 3};
539: PetscReal coords[9] = {1.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.5};
540: PetscInt markerPoints[8] = {1, 1, 2, 1, 3, 1, 4, 1};
542: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
543: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
544: } break;
545: default:
546: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
547: }
548: break;
549: default:
550: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
551: }
552: if (user->testOrientIF) {
553: PetscInt ifp[] = {8, 6};
555: PetscCall(PetscObjectSetName((PetscObject)*dm, "Mesh before orientation"));
556: PetscCall(DMViewFromOptions(*dm, NULL, "-before_orientation_dm_view"));
557: /* rotate interface face ifp[rank] by given orientation ornt[rank] */
558: PetscCall(DMPlexOrientPoint(*dm, ifp[rank], user->ornt[rank]));
559: PetscCall(DMViewFromOptions(*dm, NULL, "-before_orientation_dm_view"));
560: PetscCall(DMPlexCheckFaces(*dm, 0));
561: PetscCall(DMPlexOrientInterface_Internal(*dm));
562: PetscCall(PetscPrintf(comm, "Orientation test PASSED\n"));
563: }
564: PetscFunctionReturn(PETSC_SUCCESS);
565: }
567: static PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscBool interpolate, AppCtx *user, DM *dm)
568: {
569: PetscInt testNum = user->testNum, p;
570: PetscMPIInt rank, size;
572: PetscFunctionBegin;
573: PetscCallMPI(MPI_Comm_rank(comm, &rank));
574: PetscCallMPI(MPI_Comm_size(comm, &size));
575: switch (testNum) {
576: case 0:
577: PetscCheck(size == 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 2 processes", testNum);
578: switch (rank) {
579: case 0: {
580: const PetscInt numCells = 1, numVertices = 3, numCorners = 4;
581: const PetscInt cells[4] = {0, 1, 2, 3};
582: PetscReal coords[6] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0};
583: PetscInt markerPoints[4 * 2] = {1, 1, 2, 1, 3, 1, 4, 1};
585: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
586: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
587: } break;
588: case 1: {
589: const PetscInt numCells = 1, numVertices = 3, numCorners = 4;
590: const PetscInt cells[4] = {1, 4, 5, 2};
591: PetscReal coords[6] = {-0.5, 1.0, 0.5, 0.0, 0.5, 1.0};
592: PetscInt markerPoints[4 * 2] = {1, 1, 2, 1, 3, 1, 4, 1};
594: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
595: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
596: } break;
597: default:
598: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
599: }
600: break;
601: default:
602: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
603: }
604: PetscFunctionReturn(PETSC_SUCCESS);
605: }
607: static PetscErrorCode CreateHex_3D(MPI_Comm comm, PetscBool interpolate, AppCtx *user, DM *dm)
608: {
609: PetscInt testNum = user->testNum, p;
610: PetscMPIInt rank, size;
612: PetscFunctionBegin;
613: PetscCallMPI(MPI_Comm_rank(comm, &rank));
614: PetscCallMPI(MPI_Comm_size(comm, &size));
615: switch (testNum) {
616: case 0:
617: PetscCheck(size == 2, comm, PETSC_ERR_ARG_OUTOFRANGE, "Test mesh %" PetscInt_FMT " only for 2 processes", testNum);
618: switch (rank) {
619: case 0: {
620: const PetscInt numCells = 1, numVertices = 6, numCorners = 8;
621: const PetscInt cells[8] = {0, 3, 2, 1, 4, 5, 6, 7};
622: PetscReal coords[6 * 3] = {-0.5, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 0.0, -0.5, 1.0, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 1.0};
623: PetscInt markerPoints[8 * 2] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1, 8, 1, 9, 1};
625: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
626: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
627: } break;
628: case 1: {
629: const PetscInt numCells = 1, numVertices = 6, numCorners = 8;
630: const PetscInt cells[8] = {1, 2, 9, 8, 5, 10, 11, 6};
631: PetscReal coords[6 * 3] = {0.0, 1.0, 1.0, -0.5, 1.0, 1.0, 0.5, 0.0, 0.0, 0.5, 1.0, 0.0, 0.5, 0.0, 1.0, 0.5, 1.0, 1.0};
632: PetscInt markerPoints[8 * 2] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1, 8, 1, 9, 1};
634: PetscCall(DMPlexCreateFromCellListParallelPetsc(comm, user->dim, numCells, numVertices, PETSC_DECIDE, numCorners, interpolate, cells, user->dim, coords, NULL, NULL, dm));
635: for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(*dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
636: } break;
637: default:
638: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh for rank %d", rank);
639: }
640: break;
641: default:
642: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
643: }
644: PetscFunctionReturn(PETSC_SUCCESS);
645: }
647: static PetscErrorCode CustomView(DM dm, PetscViewer v)
648: {
649: DMPlexInterpolatedFlag interpolated;
650: PetscBool distributed;
652: PetscFunctionBegin;
653: PetscCall(DMPlexIsDistributed(dm, &distributed));
654: PetscCall(DMPlexIsInterpolatedCollective(dm, &interpolated));
655: PetscCall(PetscViewerASCIIPrintf(v, "DMPlexIsDistributed: %s\n", PetscBools[distributed]));
656: PetscCall(PetscViewerASCIIPrintf(v, "DMPlexIsInterpolatedCollective: %s\n", DMPlexInterpolatedFlags[interpolated]));
657: PetscFunctionReturn(PETSC_SUCCESS);
658: }
660: static PetscErrorCode CreateMeshFromFile(MPI_Comm comm, AppCtx *user, DM *dm, DM *serialDM)
661: {
662: const char *filename = user->filename;
663: PetscBool testHeavy = user->testHeavy;
664: PetscBool interpCreate = user->interpolate == CREATE ? PETSC_TRUE : PETSC_FALSE;
665: PetscBool distributed = PETSC_FALSE;
667: PetscFunctionBegin;
668: *serialDM = NULL;
669: if (testHeavy && interpCreate) PetscCall(DMPlexSetOrientInterface_Private(NULL, PETSC_FALSE));
670: PetscCall(PetscLogStagePush(stage[0]));
671: PetscCall(DMPlexCreateFromFile(comm, filename, "ex18_plex", interpCreate, dm)); /* with DMPlexOrientInterface_Internal() call skipped so that PointSF issues are left to DMPlexCheckPointSFHeavy() */
672: PetscCall(PetscLogStagePop());
673: if (testHeavy && interpCreate) PetscCall(DMPlexSetOrientInterface_Private(NULL, PETSC_TRUE));
674: PetscCall(DMPlexIsDistributed(*dm, &distributed));
675: PetscCall(PetscPrintf(comm, "DMPlexCreateFromFile produced %s mesh.\n", distributed ? "distributed" : "serial"));
676: if (testHeavy && distributed) {
677: PetscCall(PetscOptionsSetValue(NULL, "-dm_plex_hdf5_force_sequential", NULL));
678: PetscCall(DMPlexCreateFromFile(comm, filename, "ex18_plex", interpCreate, serialDM));
679: PetscCall(DMPlexIsDistributed(*serialDM, &distributed));
680: PetscCheck(!distributed, comm, PETSC_ERR_PLIB, "unable to create a serial DM from file");
681: }
682: PetscCall(DMGetDimension(*dm, &user->dim));
683: PetscFunctionReturn(PETSC_SUCCESS);
684: }
686: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
687: {
688: PetscPartitioner part;
689: PortableBoundary boundary = NULL;
690: DM serialDM = NULL;
691: PetscBool cellSimplex = user->cellSimplex;
692: PetscBool useGenerator = user->useGenerator;
693: PetscBool interpCreate = user->interpolate == CREATE ? PETSC_TRUE : PETSC_FALSE;
694: PetscBool interpSerial = user->interpolate == AFTER_CREATE ? PETSC_TRUE : PETSC_FALSE;
695: PetscBool interpParallel = user->interpolate == AFTER_DISTRIBUTE ? PETSC_TRUE : PETSC_FALSE;
696: PetscBool testHeavy = user->testHeavy;
697: PetscMPIInt rank;
699: PetscFunctionBegin;
700: PetscCallMPI(MPI_Comm_rank(comm, &rank));
701: if (user->filename[0]) {
702: PetscCall(CreateMeshFromFile(comm, user, dm, &serialDM));
703: } else if (useGenerator) {
704: PetscCall(PetscLogStagePush(stage[0]));
705: PetscCall(DMPlexCreateBoxMesh(comm, user->dim, cellSimplex, user->faces, NULL, NULL, NULL, interpCreate, 0, PETSC_TRUE, dm));
706: PetscCall(PetscLogStagePop());
707: } else {
708: PetscCall(PetscLogStagePush(stage[0]));
709: switch (user->dim) {
710: case 1:
711: PetscCall(CreateMesh_1D(comm, interpCreate, user, dm));
712: break;
713: case 2:
714: if (cellSimplex) {
715: PetscCall(CreateSimplex_2D(comm, interpCreate, user, dm));
716: } else {
717: PetscCall(CreateQuad_2D(comm, interpCreate, user, dm));
718: }
719: break;
720: case 3:
721: if (cellSimplex) {
722: PetscCall(CreateSimplex_3D(comm, interpCreate, user, dm));
723: } else {
724: PetscCall(CreateHex_3D(comm, interpCreate, user, dm));
725: }
726: break;
727: default:
728: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, user->dim);
729: }
730: PetscCall(PetscLogStagePop());
731: }
732: PetscCheck(user->ncoords % user->dim == 0, comm, PETSC_ERR_ARG_OUTOFRANGE, "length of coordinates array %" PetscInt_FMT " must be divisible by spatial dimension %" PetscInt_FMT, user->ncoords, user->dim);
733: PetscCall(PetscObjectSetName((PetscObject)*dm, "Original Mesh"));
734: PetscCall(DMViewFromOptions(*dm, NULL, "-orig_dm_view"));
736: if (interpSerial) {
737: DM idm;
739: if (testHeavy) PetscCall(DMPlexSetOrientInterface_Private(*dm, PETSC_FALSE));
740: PetscCall(PetscLogStagePush(stage[2]));
741: PetscCall(DMPlexInterpolate(*dm, &idm)); /* with DMPlexOrientInterface_Internal() call skipped so that PointSF issues are left to DMPlexCheckPointSFHeavy() */
742: PetscCall(PetscLogStagePop());
743: if (testHeavy) PetscCall(DMPlexSetOrientInterface_Private(*dm, PETSC_TRUE));
744: PetscCall(DMDestroy(dm));
745: *dm = idm;
746: PetscCall(PetscObjectSetName((PetscObject)*dm, "Interpolated Mesh"));
747: PetscCall(DMViewFromOptions(*dm, NULL, "-intp_dm_view"));
748: }
750: /* Set partitioner options */
751: PetscCall(DMPlexGetPartitioner(*dm, &part));
752: if (part) {
753: PetscCall(PetscPartitionerSetType(part, PETSCPARTITIONERSIMPLE));
754: PetscCall(PetscPartitionerSetFromOptions(part));
755: }
757: if (user->customView) PetscCall(CustomView(*dm, PETSC_VIEWER_STDOUT_(comm)));
758: if (testHeavy) {
759: PetscBool distributed;
761: PetscCall(DMPlexIsDistributed(*dm, &distributed));
762: if (!serialDM && !distributed) {
763: serialDM = *dm;
764: PetscCall(PetscObjectReference((PetscObject)*dm));
765: }
766: if (serialDM) PetscCall(DMPlexGetExpandedBoundary_Private(serialDM, &boundary));
767: if (boundary) {
768: /* check DM which has been created in parallel and already interpolated */
769: PetscCall(DMPlexCheckPointSFHeavy(*dm, boundary));
770: }
771: /* Orient interface because it could be deliberately skipped above. It is idempotent. */
772: PetscCall(DMPlexOrientInterface_Internal(*dm));
773: }
774: if (user->distribute) {
775: DM pdm = NULL;
777: /* Redistribute mesh over processes using that partitioner */
778: PetscCall(PetscLogStagePush(stage[1]));
779: PetscCall(DMPlexDistribute(*dm, 0, NULL, &pdm));
780: PetscCall(PetscLogStagePop());
781: if (pdm) {
782: PetscCall(DMDestroy(dm));
783: *dm = pdm;
784: PetscCall(PetscObjectSetName((PetscObject)*dm, "Redistributed Mesh"));
785: PetscCall(DMViewFromOptions(*dm, NULL, "-dist_dm_view"));
786: }
788: if (interpParallel) {
789: DM idm;
791: if (testHeavy) PetscCall(DMPlexSetOrientInterface_Private(*dm, PETSC_FALSE));
792: PetscCall(PetscLogStagePush(stage[2]));
793: PetscCall(DMPlexInterpolate(*dm, &idm)); /* with DMPlexOrientInterface_Internal() call skipped so that PointSF issues are left to DMPlexCheckPointSFHeavy() */
794: PetscCall(PetscLogStagePop());
795: if (testHeavy) PetscCall(DMPlexSetOrientInterface_Private(*dm, PETSC_TRUE));
796: PetscCall(DMDestroy(dm));
797: *dm = idm;
798: PetscCall(PetscObjectSetName((PetscObject)*dm, "Interpolated Redistributed Mesh"));
799: PetscCall(DMViewFromOptions(*dm, NULL, "-intp_dm_view"));
800: }
801: }
802: if (testHeavy) {
803: if (boundary) PetscCall(DMPlexCheckPointSFHeavy(*dm, boundary));
804: /* Orient interface because it could be deliberately skipped above. It is idempotent. */
805: PetscCall(DMPlexOrientInterface_Internal(*dm));
806: }
808: PetscCall(PetscObjectSetName((PetscObject)*dm, "Parallel Mesh"));
809: PetscCall(DMPlexDistributeSetDefault(*dm, PETSC_FALSE));
810: PetscCall(DMSetFromOptions(*dm));
811: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
813: if (user->customView) PetscCall(CustomView(*dm, PETSC_VIEWER_STDOUT_(comm)));
814: PetscCall(DMDestroy(&serialDM));
815: PetscCall(PortableBoundaryDestroy(&boundary));
816: PetscFunctionReturn(PETSC_SUCCESS);
817: }
819: #define ps2d(number) ((double)PetscRealPart(number))
820: static inline PetscErrorCode coord2str(char buf[], size_t len, PetscInt dim, const PetscScalar coords[], PetscReal tol)
821: {
822: PetscFunctionBegin;
823: PetscCheck(dim <= 3, PETSC_COMM_SELF, PETSC_ERR_SUP, "dim must be less than or equal 3");
824: if (tol >= 1e-3) {
825: switch (dim) {
826: case 1:
827: PetscCall(PetscSNPrintf(buf, len, "(%12.3f)", ps2d(coords[0])));
828: case 2:
829: PetscCall(PetscSNPrintf(buf, len, "(%12.3f, %12.3f)", ps2d(coords[0]), ps2d(coords[1])));
830: case 3:
831: PetscCall(PetscSNPrintf(buf, len, "(%12.3f, %12.3f, %12.3f)", ps2d(coords[0]), ps2d(coords[1]), ps2d(coords[2])));
832: }
833: } else {
834: switch (dim) {
835: case 1:
836: PetscCall(PetscSNPrintf(buf, len, "(%12.6f)", ps2d(coords[0])));
837: case 2:
838: PetscCall(PetscSNPrintf(buf, len, "(%12.6f, %12.6f)", ps2d(coords[0]), ps2d(coords[1])));
839: case 3:
840: PetscCall(PetscSNPrintf(buf, len, "(%12.6f, %12.6f, %12.6f)", ps2d(coords[0]), ps2d(coords[1]), ps2d(coords[2])));
841: }
842: }
843: PetscFunctionReturn(PETSC_SUCCESS);
844: }
846: static PetscErrorCode ViewVerticesFromCoords(DM dm, Vec coordsVec, PetscReal tol, PetscViewer viewer)
847: {
848: PetscInt dim, i, npoints;
849: IS pointsIS;
850: const PetscInt *points;
851: const PetscScalar *coords;
852: char coordstr[128];
853: MPI_Comm comm;
854: PetscMPIInt rank;
856: PetscFunctionBegin;
857: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
858: PetscCallMPI(MPI_Comm_rank(comm, &rank));
859: PetscCall(DMGetDimension(dm, &dim));
860: PetscCall(PetscViewerASCIIPushSynchronized(viewer));
861: PetscCall(DMPlexFindVertices(dm, coordsVec, tol, &pointsIS));
862: PetscCall(ISGetIndices(pointsIS, &points));
863: PetscCall(ISGetLocalSize(pointsIS, &npoints));
864: PetscCall(VecGetArrayRead(coordsVec, &coords));
865: for (i = 0; i < npoints; i++) {
866: PetscCall(coord2str(coordstr, sizeof(coordstr), dim, &coords[i * dim], tol));
867: if (rank == 0 && i) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "-----\n"));
868: PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "[%d] %s --> points[%" PetscInt_FMT "] = %" PetscInt_FMT "\n", rank, coordstr, i, points[i]));
869: PetscCall(PetscViewerFlush(viewer));
870: }
871: PetscCall(PetscViewerASCIIPopSynchronized(viewer));
872: PetscCall(VecRestoreArrayRead(coordsVec, &coords));
873: PetscCall(ISRestoreIndices(pointsIS, &points));
874: PetscCall(ISDestroy(&pointsIS));
875: PetscFunctionReturn(PETSC_SUCCESS);
876: }
878: static PetscErrorCode TestExpandPoints(DM dm, AppCtx *user)
879: {
880: IS is;
881: PetscSection *sects;
882: IS *iss;
883: PetscInt d, depth;
884: PetscMPIInt rank;
885: PetscViewer viewer = PETSC_VIEWER_STDOUT_WORLD, sviewer;
887: PetscFunctionBegin;
888: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
889: if (user->testExpandPointsEmpty && rank == 0) {
890: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, PETSC_USE_POINTER, &is));
891: } else {
892: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, user->nPointsToExpand, user->pointsToExpand, PETSC_USE_POINTER, &is));
893: }
894: PetscCall(DMPlexGetConeRecursive(dm, is, &depth, &iss, §s));
895: PetscCall(PetscViewerGetSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
896: PetscCall(PetscViewerASCIIPrintf(sviewer, "[%d] ==========================\n", rank));
897: for (d = depth - 1; d >= 0; d--) {
898: IS checkIS;
899: PetscBool flg;
901: PetscCall(PetscViewerASCIIPrintf(sviewer, "depth %" PetscInt_FMT " ---------------\n", d));
902: PetscCall(PetscSectionView(sects[d], sviewer));
903: PetscCall(ISView(iss[d], sviewer));
904: /* check reverse operation */
905: if (d < depth - 1) {
906: PetscCall(DMPlexExpandedConesToFaces_Private(dm, iss[d], sects[d], &checkIS));
907: PetscCall(ISEqualUnsorted(checkIS, iss[d + 1], &flg));
908: PetscCheck(flg, PetscObjectComm((PetscObject)checkIS), PETSC_ERR_PLIB, "DMPlexExpandedConesToFaces_Private produced wrong IS");
909: PetscCall(ISDestroy(&checkIS));
910: }
911: }
912: PetscCall(PetscViewerRestoreSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
913: PetscCall(DMPlexRestoreConeRecursive(dm, is, &depth, &iss, §s));
914: PetscCall(ISDestroy(&is));
915: PetscFunctionReturn(PETSC_SUCCESS);
916: }
918: static PetscErrorCode DMPlexExpandedConesToFaces_Private(DM dm, IS is, PetscSection section, IS *newis)
919: {
920: PetscInt n, n1, ncone, numCoveredPoints, o, p, q, start, end;
921: const PetscInt *coveredPoints;
922: const PetscInt *arr, *cone;
923: PetscInt *newarr;
925: PetscFunctionBegin;
926: PetscCall(ISGetLocalSize(is, &n));
927: PetscCall(PetscSectionGetStorageSize(section, &n1));
928: PetscCall(PetscSectionGetChart(section, &start, &end));
929: PetscCheck(n == n1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "IS size = %" PetscInt_FMT " != %" PetscInt_FMT " = section storage size", n, n1);
930: PetscCall(ISGetIndices(is, &arr));
931: PetscCall(PetscMalloc1(end - start, &newarr));
932: for (q = start; q < end; q++) {
933: PetscCall(PetscSectionGetDof(section, q, &ncone));
934: PetscCall(PetscSectionGetOffset(section, q, &o));
935: cone = &arr[o];
936: if (ncone == 1) {
937: numCoveredPoints = 1;
938: p = cone[0];
939: } else {
940: PetscInt i;
941: p = PETSC_INT_MAX;
942: for (i = 0; i < ncone; i++)
943: if (cone[i] < 0) {
944: p = -1;
945: break;
946: }
947: if (p >= 0) {
948: PetscCall(DMPlexGetJoin(dm, ncone, cone, &numCoveredPoints, &coveredPoints));
949: PetscCheck(numCoveredPoints <= 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "more than one covered points for section point %" PetscInt_FMT, q);
950: if (numCoveredPoints) p = coveredPoints[0];
951: else p = -2;
952: PetscCall(DMPlexRestoreJoin(dm, ncone, cone, &numCoveredPoints, &coveredPoints));
953: }
954: }
955: newarr[q - start] = p;
956: }
957: PetscCall(ISRestoreIndices(is, &arr));
958: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, end - start, newarr, PETSC_OWN_POINTER, newis));
959: PetscFunctionReturn(PETSC_SUCCESS);
960: }
962: static PetscErrorCode DMPlexExpandedVerticesToFaces_Private(DM dm, IS boundary_expanded_is, PetscInt depth, PetscSection sections[], IS *boundary_is)
963: {
964: PetscInt d;
965: IS is, newis;
967: PetscFunctionBegin;
968: is = boundary_expanded_is;
969: PetscCall(PetscObjectReference((PetscObject)is));
970: for (d = 0; d < depth - 1; ++d) {
971: PetscCall(DMPlexExpandedConesToFaces_Private(dm, is, sections[d], &newis));
972: PetscCall(ISDestroy(&is));
973: is = newis;
974: }
975: *boundary_is = is;
976: PetscFunctionReturn(PETSC_SUCCESS);
977: }
979: #define CHKERRQI(incall, ierr) \
980: do { \
981: if (ierr) incall = PETSC_FALSE; \
982: } while (0)
984: static PetscErrorCode DMLabelViewFromOptionsOnComm_Private(DMLabel label, const char optionname[], MPI_Comm comm)
985: {
986: PetscViewer viewer;
987: PetscBool flg;
988: static PetscBool incall = PETSC_FALSE;
989: PetscViewerFormat format;
991: PetscFunctionBegin;
992: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
993: incall = PETSC_TRUE;
994: CHKERRQI(incall, PetscOptionsCreateViewer(comm, ((PetscObject)label)->options, ((PetscObject)label)->prefix, optionname, &viewer, &format, &flg));
995: if (flg) {
996: CHKERRQI(incall, PetscViewerPushFormat(viewer, format));
997: CHKERRQI(incall, DMLabelView(label, viewer));
998: CHKERRQI(incall, PetscViewerPopFormat(viewer));
999: CHKERRQI(incall, PetscViewerDestroy(&viewer));
1000: }
1001: incall = PETSC_FALSE;
1002: PetscFunctionReturn(PETSC_SUCCESS);
1003: }
1005: /* TODO: this is hotfixing DMLabelGetStratumIS() - it should be fixed systematically instead */
1006: static inline PetscErrorCode DMLabelGetStratumISOnComm_Private(DMLabel label, PetscInt value, MPI_Comm comm, IS *is)
1007: {
1008: IS tmpis;
1010: PetscFunctionBegin;
1011: PetscCall(DMLabelGetStratumIS(label, value, &tmpis));
1012: if (!tmpis) PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 0, NULL, PETSC_USE_POINTER, &tmpis));
1013: PetscCall(ISOnComm(tmpis, comm, PETSC_COPY_VALUES, is));
1014: PetscCall(ISDestroy(&tmpis));
1015: PetscFunctionReturn(PETSC_SUCCESS);
1016: }
1018: /* currently only for simple PetscSection without fields or constraints */
1019: static PetscErrorCode PetscSectionReplicate_Private(MPI_Comm comm, PetscMPIInt rootrank, PetscSection sec0, PetscSection *secout)
1020: {
1021: PetscSection sec;
1022: PetscInt chart[2], p;
1023: PetscInt *dofarr;
1024: PetscMPIInt rank;
1026: PetscFunctionBegin;
1027: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1028: if (rank == rootrank) PetscCall(PetscSectionGetChart(sec0, &chart[0], &chart[1]));
1029: PetscCallMPI(MPI_Bcast(chart, 2, MPIU_INT, rootrank, comm));
1030: PetscCall(PetscMalloc1(chart[1] - chart[0], &dofarr));
1031: if (rank == rootrank) {
1032: for (p = chart[0]; p < chart[1]; p++) PetscCall(PetscSectionGetDof(sec0, p, &dofarr[p - chart[0]]));
1033: }
1034: PetscCallMPI(MPI_Bcast(dofarr, (PetscMPIInt)(chart[1] - chart[0]), MPIU_INT, rootrank, comm));
1035: PetscCall(PetscSectionCreate(comm, &sec));
1036: PetscCall(PetscSectionSetChart(sec, chart[0], chart[1]));
1037: for (p = chart[0]; p < chart[1]; p++) PetscCall(PetscSectionSetDof(sec, p, dofarr[p - chart[0]]));
1038: PetscCall(PetscSectionSetUp(sec));
1039: PetscCall(PetscFree(dofarr));
1040: *secout = sec;
1041: PetscFunctionReturn(PETSC_SUCCESS);
1042: }
1044: static PetscErrorCode DMPlexExpandedVerticesCoordinatesToFaces_Private(DM ipdm, PortableBoundary bnd, IS *face_is)
1045: {
1046: IS faces_expanded_is;
1048: PetscFunctionBegin;
1049: PetscCall(DMPlexFindVertices(ipdm, bnd->coordinates, 0.0, &faces_expanded_is));
1050: PetscCall(DMPlexExpandedVerticesToFaces_Private(ipdm, faces_expanded_is, bnd->depth, bnd->sections, face_is));
1051: PetscCall(ISDestroy(&faces_expanded_is));
1052: PetscFunctionReturn(PETSC_SUCCESS);
1053: }
1055: /* hack disabling DMPlexOrientInterface() call in DMPlexInterpolate() via -dm_plex_interpolate_orient_interfaces option */
1056: static PetscErrorCode DMPlexSetOrientInterface_Private(DM dm, PetscBool enable)
1057: {
1058: PetscOptions options = NULL;
1059: const char *prefix = NULL;
1060: const char opt[] = "-dm_plex_interpolate_orient_interfaces";
1061: char prefix_opt[512];
1062: PetscBool flg, set;
1063: static PetscBool wasSetTrue = PETSC_FALSE;
1065: PetscFunctionBegin;
1066: if (dm) {
1067: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
1068: options = ((PetscObject)dm)->options;
1069: }
1070: PetscCall(PetscStrncpy(prefix_opt, "-", sizeof(prefix_opt)));
1071: PetscCall(PetscStrlcat(prefix_opt, prefix, sizeof(prefix_opt)));
1072: PetscCall(PetscStrlcat(prefix_opt, &opt[1], sizeof(prefix_opt)));
1073: PetscCall(PetscOptionsGetBool(options, prefix, opt, &flg, &set));
1074: if (!enable) {
1075: if (set && flg) wasSetTrue = PETSC_TRUE;
1076: PetscCall(PetscOptionsSetValue(options, prefix_opt, "0"));
1077: } else if (set && !flg) {
1078: if (wasSetTrue) {
1079: PetscCall(PetscOptionsSetValue(options, prefix_opt, "1"));
1080: } else {
1081: /* default is PETSC_TRUE */
1082: PetscCall(PetscOptionsClearValue(options, prefix_opt));
1083: }
1084: wasSetTrue = PETSC_FALSE;
1085: }
1086: if (PetscDefined(USE_DEBUG)) {
1087: PetscCall(PetscOptionsGetBool(options, prefix, opt, &flg, &set));
1088: PetscCheck(!set || flg == enable, PetscObjectComm((PetscObject)dm), PETSC_ERR_PLIB, "PetscOptionsSetValue did not have the desired effect");
1089: }
1090: PetscFunctionReturn(PETSC_SUCCESS);
1091: }
1093: /* get coordinate description of the whole-domain boundary */
1094: static PetscErrorCode DMPlexGetExpandedBoundary_Private(DM dm, PortableBoundary *boundary)
1095: {
1096: PortableBoundary bnd0, bnd;
1097: MPI_Comm comm;
1098: DM idm;
1099: DMLabel label;
1100: PetscInt d;
1101: const char boundaryName[] = "DMPlexDistributeInterpolateMarkInterface_boundary";
1102: IS boundary_is;
1103: IS *boundary_expanded_iss;
1104: PetscMPIInt rootrank = 0;
1105: PetscMPIInt rank, size;
1106: PetscInt value = 1;
1107: DMPlexInterpolatedFlag intp;
1108: PetscBool flg;
1110: PetscFunctionBegin;
1111: PetscCall(PetscNew(&bnd));
1112: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
1113: PetscCallMPI(MPI_Comm_rank(comm, &rank));
1114: PetscCallMPI(MPI_Comm_size(comm, &size));
1115: PetscCall(DMPlexIsDistributed(dm, &flg));
1116: PetscCheck(!flg, comm, PETSC_ERR_ARG_WRONG, "serial DM (all points on one rank) needed");
1118: /* interpolate serial DM if not yet interpolated */
1119: PetscCall(DMPlexIsInterpolatedCollective(dm, &intp));
1120: if (intp == DMPLEX_INTERPOLATED_FULL) {
1121: idm = dm;
1122: PetscCall(PetscObjectReference((PetscObject)dm));
1123: } else {
1124: PetscCall(DMPlexInterpolate(dm, &idm));
1125: PetscCall(DMViewFromOptions(idm, NULL, "-idm_view"));
1126: }
1128: /* mark whole-domain boundary of the serial DM */
1129: PetscCall(DMLabelCreate(PETSC_COMM_SELF, boundaryName, &label));
1130: PetscCall(DMAddLabel(idm, label));
1131: PetscCall(DMPlexMarkBoundaryFaces(idm, value, label));
1132: PetscCall(DMLabelViewFromOptionsOnComm_Private(label, "-idm_boundary_view", comm));
1133: PetscCall(DMLabelGetStratumIS(label, value, &boundary_is));
1135: /* translate to coordinates */
1136: PetscCall(PetscNew(&bnd0));
1137: PetscCall(DMGetCoordinatesLocalSetUp(idm));
1138: if (rank == rootrank) {
1139: PetscCall(DMPlexGetConeRecursive(idm, boundary_is, &bnd0->depth, &boundary_expanded_iss, &bnd0->sections));
1140: PetscCall(DMGetCoordinatesLocalTuple(dm, boundary_expanded_iss[0], NULL, &bnd0->coordinates));
1141: /* self-check */
1142: {
1143: IS is0;
1144: PetscCall(DMPlexExpandedVerticesCoordinatesToFaces_Private(idm, bnd0, &is0));
1145: PetscCall(ISEqual(is0, boundary_is, &flg));
1146: PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "DMPlexExpandedVerticesCoordinatesToFaces_Private produced a wrong IS");
1147: PetscCall(ISDestroy(&is0));
1148: }
1149: } else {
1150: PetscCall(VecCreateFromOptions(PETSC_COMM_SELF, NULL, 1, 0, 0, &bnd0->coordinates));
1151: }
1153: {
1154: Vec tmp;
1155: VecScatter sc;
1156: IS xis;
1157: PetscInt n;
1159: /* just convert seq vectors to mpi vector */
1160: PetscCall(VecGetLocalSize(bnd0->coordinates, &n));
1161: PetscCallMPI(MPI_Bcast(&n, 1, MPIU_INT, rootrank, comm));
1162: if (rank == rootrank) {
1163: PetscCall(VecCreateFromOptions(comm, NULL, 1, n, n, &tmp));
1164: } else {
1165: PetscCall(VecCreateFromOptions(comm, NULL, 1, 0, n, &tmp));
1166: }
1167: PetscCall(VecCopy(bnd0->coordinates, tmp));
1168: PetscCall(VecDestroy(&bnd0->coordinates));
1169: bnd0->coordinates = tmp;
1171: /* replicate coordinates from root rank to all ranks */
1172: PetscCall(VecCreateFromOptions(comm, NULL, 1, n, n * size, &bnd->coordinates));
1173: PetscCall(ISCreateStride(comm, n, 0, 1, &xis));
1174: PetscCall(VecScatterCreate(bnd0->coordinates, xis, bnd->coordinates, NULL, &sc));
1175: PetscCall(VecScatterBegin(sc, bnd0->coordinates, bnd->coordinates, INSERT_VALUES, SCATTER_FORWARD));
1176: PetscCall(VecScatterEnd(sc, bnd0->coordinates, bnd->coordinates, INSERT_VALUES, SCATTER_FORWARD));
1177: PetscCall(VecScatterDestroy(&sc));
1178: PetscCall(ISDestroy(&xis));
1179: }
1180: bnd->depth = bnd0->depth;
1181: PetscCallMPI(MPI_Bcast(&bnd->depth, 1, MPIU_INT, rootrank, comm));
1182: PetscCall(PetscMalloc1(bnd->depth, &bnd->sections));
1183: for (d = 0; d < bnd->depth; d++) PetscCall(PetscSectionReplicate_Private(comm, rootrank, (rank == rootrank) ? bnd0->sections[d] : NULL, &bnd->sections[d]));
1185: if (rank == rootrank) PetscCall(DMPlexRestoreConeRecursive(idm, boundary_is, &bnd0->depth, &boundary_expanded_iss, &bnd0->sections));
1186: PetscCall(PortableBoundaryDestroy(&bnd0));
1187: PetscCall(DMRemoveLabelBySelf(idm, &label, PETSC_TRUE));
1188: PetscCall(DMLabelDestroy(&label));
1189: PetscCall(ISDestroy(&boundary_is));
1190: PetscCall(DMDestroy(&idm));
1191: *boundary = bnd;
1192: PetscFunctionReturn(PETSC_SUCCESS);
1193: }
1195: /* get faces of inter-partition interface */
1196: static PetscErrorCode DMPlexGetInterfaceFaces_Private(DM ipdm, IS boundary_faces_is, IS *interface_faces_is)
1197: {
1198: MPI_Comm comm;
1199: DMLabel label;
1200: IS part_boundary_faces_is;
1201: const char partBoundaryName[] = "DMPlexDistributeInterpolateMarkInterface_partBoundary";
1202: PetscInt value = 1;
1203: DMPlexInterpolatedFlag intp;
1205: PetscFunctionBegin;
1206: PetscCall(PetscObjectGetComm((PetscObject)ipdm, &comm));
1207: PetscCall(DMPlexIsInterpolatedCollective(ipdm, &intp));
1208: PetscCheck(intp == DMPLEX_INTERPOLATED_FULL, comm, PETSC_ERR_ARG_WRONG, "only for fully interpolated DMPlex");
1210: /* get ipdm partition boundary (partBoundary) */
1211: {
1212: PetscSF sf;
1214: PetscCall(DMLabelCreate(PETSC_COMM_SELF, partBoundaryName, &label));
1215: PetscCall(DMAddLabel(ipdm, label));
1216: PetscCall(DMGetPointSF(ipdm, &sf));
1217: PetscCall(DMSetPointSF(ipdm, NULL));
1218: PetscCall(DMPlexMarkBoundaryFaces(ipdm, value, label));
1219: PetscCall(DMSetPointSF(ipdm, sf));
1220: PetscCall(DMLabelViewFromOptionsOnComm_Private(label, "-ipdm_part_boundary_view", comm));
1221: PetscCall(DMLabelGetStratumISOnComm_Private(label, value, comm, &part_boundary_faces_is));
1222: PetscCall(DMRemoveLabelBySelf(ipdm, &label, PETSC_TRUE));
1223: PetscCall(DMLabelDestroy(&label));
1224: }
1226: /* remove ipdm whole-domain boundary (boundary_faces_is) from ipdm partition boundary (part_boundary_faces_is), resulting just in inter-partition interface */
1227: PetscCall(ISDifference(part_boundary_faces_is, boundary_faces_is, interface_faces_is));
1228: PetscCall(ISDestroy(&part_boundary_faces_is));
1229: PetscFunctionReturn(PETSC_SUCCESS);
1230: }
1232: /* compute inter-partition interface including edges and vertices */
1233: static PetscErrorCode DMPlexComputeCompleteInterface_Private(DM ipdm, IS interface_faces_is, IS *interface_is)
1234: {
1235: DMLabel label;
1236: PetscInt value = 1;
1237: const char interfaceName[] = "DMPlexDistributeInterpolateMarkInterface_interface";
1238: DMPlexInterpolatedFlag intp;
1239: MPI_Comm comm;
1241: PetscFunctionBegin;
1242: PetscCall(PetscObjectGetComm((PetscObject)ipdm, &comm));
1243: PetscCall(DMPlexIsInterpolatedCollective(ipdm, &intp));
1244: PetscCheck(intp == DMPLEX_INTERPOLATED_FULL, comm, PETSC_ERR_ARG_WRONG, "only for fully interpolated DMPlex");
1246: PetscCall(DMLabelCreate(PETSC_COMM_SELF, interfaceName, &label));
1247: PetscCall(DMAddLabel(ipdm, label));
1248: PetscCall(DMLabelSetStratumIS(label, value, interface_faces_is));
1249: PetscCall(DMLabelViewFromOptionsOnComm_Private(label, "-interface_faces_view", comm));
1250: PetscCall(DMPlexLabelComplete(ipdm, label));
1251: PetscCall(DMLabelViewFromOptionsOnComm_Private(label, "-interface_view", comm));
1252: PetscCall(DMLabelGetStratumISOnComm_Private(label, value, comm, interface_is));
1253: PetscCall(PetscObjectSetName((PetscObject)*interface_is, "interface_is"));
1254: PetscCall(ISViewFromOptions(*interface_is, NULL, "-interface_is_view"));
1255: PetscCall(DMRemoveLabelBySelf(ipdm, &label, PETSC_TRUE));
1256: PetscCall(DMLabelDestroy(&label));
1257: PetscFunctionReturn(PETSC_SUCCESS);
1258: }
1260: static PetscErrorCode PointSFGetOutwardInterfacePoints(PetscSF sf, IS *is)
1261: {
1262: PetscInt n;
1263: const PetscInt *arr;
1265: PetscFunctionBegin;
1266: PetscCall(PetscSFGetGraph(sf, NULL, &n, &arr, NULL));
1267: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)sf), n, arr, PETSC_USE_POINTER, is));
1268: PetscFunctionReturn(PETSC_SUCCESS);
1269: }
1271: static PetscErrorCode PointSFGetInwardInterfacePoints(PetscSF sf, IS *is)
1272: {
1273: PetscInt n;
1274: const PetscInt *rootdegree;
1275: PetscInt *arr;
1277: PetscFunctionBegin;
1278: PetscCall(PetscSFSetUp(sf));
1279: PetscCall(PetscSFComputeDegreeBegin(sf, &rootdegree));
1280: PetscCall(PetscSFComputeDegreeEnd(sf, &rootdegree));
1281: PetscCall(PetscSFComputeMultiRootOriginalNumbering(sf, rootdegree, &n, &arr));
1282: PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)sf), n, arr, PETSC_OWN_POINTER, is));
1283: PetscFunctionReturn(PETSC_SUCCESS);
1284: }
1286: static PetscErrorCode PointSFGetInterfacePoints_Private(PetscSF pointSF, IS *is)
1287: {
1288: IS pointSF_out_is, pointSF_in_is;
1290: PetscFunctionBegin;
1291: PetscCall(PointSFGetOutwardInterfacePoints(pointSF, &pointSF_out_is));
1292: PetscCall(PointSFGetInwardInterfacePoints(pointSF, &pointSF_in_is));
1293: PetscCall(ISExpand(pointSF_out_is, pointSF_in_is, is));
1294: PetscCall(ISDestroy(&pointSF_out_is));
1295: PetscCall(ISDestroy(&pointSF_in_is));
1296: PetscFunctionReturn(PETSC_SUCCESS);
1297: }
1299: #define CHKERRMY(ierr) PetscCheck(!ierr, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PointSF is wrong. Unable to show details!")
1301: static PetscErrorCode ViewPointsWithType_Internal(DM dm, IS pointsIS, PetscViewer v)
1302: {
1303: DMLabel label;
1304: PetscSection coordsSection;
1305: Vec coordsVec;
1306: PetscScalar *coordsScalar;
1307: PetscInt coneSize, depth, dim, i, p, npoints;
1308: const PetscInt *points;
1310: PetscFunctionBegin;
1311: PetscCall(DMGetDimension(dm, &dim));
1312: PetscCall(DMGetCoordinateSection(dm, &coordsSection));
1313: PetscCall(DMGetCoordinatesLocal(dm, &coordsVec));
1314: PetscCall(VecGetArray(coordsVec, &coordsScalar));
1315: PetscCall(ISGetLocalSize(pointsIS, &npoints));
1316: PetscCall(ISGetIndices(pointsIS, &points));
1317: PetscCall(DMPlexGetDepthLabel(dm, &label));
1318: PetscCall(PetscViewerASCIIPushTab(v));
1319: for (i = 0; i < npoints; i++) {
1320: p = points[i];
1321: PetscCall(DMLabelGetValue(label, p, &depth));
1322: if (!depth) {
1323: PetscInt n, o;
1324: char coordstr[128];
1326: PetscCall(PetscSectionGetDof(coordsSection, p, &n));
1327: PetscCall(PetscSectionGetOffset(coordsSection, p, &o));
1328: PetscCall(coord2str(coordstr, sizeof(coordstr), n, &coordsScalar[o], 1.0));
1329: PetscCall(PetscViewerASCIISynchronizedPrintf(v, "vertex %" PetscInt_FMT " w/ coordinates %s\n", p, coordstr));
1330: } else {
1331: char entityType[16];
1333: switch (depth) {
1334: case 1:
1335: PetscCall(PetscStrncpy(entityType, "edge", sizeof(entityType)));
1336: break;
1337: case 2:
1338: PetscCall(PetscStrncpy(entityType, "face", sizeof(entityType)));
1339: break;
1340: case 3:
1341: PetscCall(PetscStrncpy(entityType, "cell", sizeof(entityType)));
1342: break;
1343: default:
1344: SETERRQ(PetscObjectComm((PetscObject)v), PETSC_ERR_SUP, "Only for depth <= 3");
1345: }
1346: if (depth == dim && dim < 3) PetscCall(PetscStrlcat(entityType, " (cell)", sizeof(entityType)));
1347: PetscCall(PetscViewerASCIISynchronizedPrintf(v, "%s %" PetscInt_FMT "\n", entityType, p));
1348: }
1349: PetscCall(DMPlexGetConeSize(dm, p, &coneSize));
1350: if (coneSize) {
1351: const PetscInt *cone;
1352: IS coneIS;
1354: PetscCall(DMPlexGetCone(dm, p, &cone));
1355: PetscCall(ISCreateGeneral(PETSC_COMM_SELF, coneSize, cone, PETSC_USE_POINTER, &coneIS));
1356: PetscCall(ViewPointsWithType_Internal(dm, coneIS, v));
1357: PetscCall(ISDestroy(&coneIS));
1358: }
1359: }
1360: PetscCall(PetscViewerASCIIPopTab(v));
1361: PetscCall(VecRestoreArray(coordsVec, &coordsScalar));
1362: PetscCall(ISRestoreIndices(pointsIS, &points));
1363: PetscFunctionReturn(PETSC_SUCCESS);
1364: }
1366: static PetscErrorCode ViewPointsWithType(DM dm, IS points, PetscViewer v)
1367: {
1368: PetscBool flg;
1369: PetscInt npoints;
1370: PetscMPIInt rank;
1372: PetscFunctionBegin;
1373: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &flg));
1374: PetscCheck(flg, PetscObjectComm((PetscObject)v), PETSC_ERR_SUP, "Only for ASCII viewer");
1375: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)v), &rank));
1376: PetscCall(PetscViewerASCIIPushSynchronized(v));
1377: PetscCall(ISGetLocalSize(points, &npoints));
1378: if (npoints) {
1379: PetscCall(PetscViewerASCIISynchronizedPrintf(v, "[%d] --------\n", rank));
1380: PetscCall(ViewPointsWithType_Internal(dm, points, v));
1381: }
1382: PetscCall(PetscViewerFlush(v));
1383: PetscCall(PetscViewerASCIIPopSynchronized(v));
1384: PetscFunctionReturn(PETSC_SUCCESS);
1385: }
1387: static PetscErrorCode DMPlexComparePointSFWithInterface_Private(DM ipdm, IS interface_is)
1388: {
1389: PetscSF pointsf;
1390: IS pointsf_is;
1391: PetscBool flg;
1392: MPI_Comm comm;
1393: PetscMPIInt size;
1395: PetscFunctionBegin;
1396: PetscCall(PetscObjectGetComm((PetscObject)ipdm, &comm));
1397: PetscCallMPI(MPI_Comm_size(comm, &size));
1398: PetscCall(DMGetPointSF(ipdm, &pointsf));
1399: if (pointsf) {
1400: PetscInt nroots;
1401: PetscCall(PetscSFGetGraph(pointsf, &nroots, NULL, NULL, NULL));
1402: if (nroots < 0) pointsf = NULL; /* uninitialized SF */
1403: }
1404: if (!pointsf) {
1405: PetscInt N = 0;
1406: if (interface_is) PetscCall(ISGetSize(interface_is, &N));
1407: PetscCheck(!N, comm, PETSC_ERR_PLIB, "interface_is should be NULL or empty for PointSF being NULL");
1408: PetscFunctionReturn(PETSC_SUCCESS);
1409: }
1411: /* get PointSF points as IS pointsf_is */
1412: PetscCall(PointSFGetInterfacePoints_Private(pointsf, &pointsf_is));
1414: /* compare pointsf_is with interface_is */
1415: PetscCall(ISEqual(interface_is, pointsf_is, &flg));
1416: PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &flg, 1, MPI_C_BOOL, MPI_LAND, comm));
1417: if (!flg) {
1418: IS pointsf_extra_is, pointsf_missing_is;
1419: PetscViewer errv = PETSC_VIEWER_STDERR_(comm);
1420: CHKERRMY(ISDifference(interface_is, pointsf_is, &pointsf_missing_is));
1421: CHKERRMY(ISDifference(pointsf_is, interface_is, &pointsf_extra_is));
1422: CHKERRMY(PetscViewerASCIIPrintf(errv, "Points missing in PointSF:\n"));
1423: CHKERRMY(ViewPointsWithType(ipdm, pointsf_missing_is, errv));
1424: CHKERRMY(PetscViewerASCIIPrintf(errv, "Extra points in PointSF:\n"));
1425: CHKERRMY(ViewPointsWithType(ipdm, pointsf_extra_is, errv));
1426: CHKERRMY(ISDestroy(&pointsf_extra_is));
1427: CHKERRMY(ISDestroy(&pointsf_missing_is));
1428: SETERRQ(comm, PETSC_ERR_PLIB, "PointSF is wrong! See details above.");
1429: }
1430: PetscCall(ISDestroy(&pointsf_is));
1431: PetscFunctionReturn(PETSC_SUCCESS);
1432: }
1434: /* remove faces & edges from label, leave just vertices */
1435: static PetscErrorCode DMPlexISFilterVertices_Private(DM dm, IS points)
1436: {
1437: PetscInt vStart, vEnd;
1438: MPI_Comm comm;
1440: PetscFunctionBegin;
1441: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
1442: PetscCall(DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd));
1443: PetscCall(ISGeneralFilter(points, vStart, vEnd));
1444: PetscFunctionReturn(PETSC_SUCCESS);
1445: }
1447: /*
1448: DMPlexCheckPointSFHeavy - Thoroughly test that the PointSF after parallel DMPlexInterpolate() includes exactly all interface points.
1450: Collective
1452: Input Parameter:
1453: . dm - The DMPlex object
1455: Notes:
1456: The input DMPlex must be serial (one partition has all points, the other partitions have no points).
1457: This is a heavy test which involves DMPlexInterpolate() if the input DM is not interpolated yet, and depends on having a representation of the whole-domain boundary (PortableBoundary), which can be obtained only by DMPlexGetExpandedBoundary_Private() (which involves DMPlexInterpolate() of a sequential DM).
1458: This is mainly intended for debugging/testing purposes.
1460: Algorithm:
1461: 1. boundary faces of the serial version of the whole mesh are found using DMPlexMarkBoundaryFaces()
1462: 2. boundary faces are translated into vertices using DMPlexGetConeRecursive() and these are translated into coordinates - this description (aka PortableBoundary) is completely independent of partitioning and point numbering
1463: 3. the mesh is distributed or loaded in parallel
1464: 4. boundary faces of the distributed mesh are reconstructed from PortableBoundary using DMPlexFindVertices()
1465: 5. partition boundary faces of the parallel mesh are found using DMPlexMarkBoundaryFaces()
1466: 6. partition interfaces are computed as set difference of partition boundary faces minus the reconstructed boundary
1467: 7. check that interface covered by PointSF (union of inward and outward points) is equal to the partition interface for each rank, otherwise print the difference and throw an error
1469: Level: developer
1471: .seealso: DMGetPointSF(), DMPlexCheckSymmetry(), DMPlexCheckSkeleton(), DMPlexCheckFaces()
1472: */
1473: static PetscErrorCode DMPlexCheckPointSFHeavy(DM dm, PortableBoundary bnd)
1474: {
1475: DM ipdm = NULL;
1476: IS boundary_faces_is, interface_faces_is, interface_is;
1477: DMPlexInterpolatedFlag intp;
1478: MPI_Comm comm;
1480: PetscFunctionBegin;
1481: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
1483: PetscCall(DMPlexIsInterpolatedCollective(dm, &intp));
1484: if (intp == DMPLEX_INTERPOLATED_FULL) {
1485: ipdm = dm;
1486: } else {
1487: /* create temporary interpolated DM if input DM is not interpolated */
1488: PetscCall(DMPlexSetOrientInterface_Private(dm, PETSC_FALSE));
1489: PetscCall(DMPlexInterpolate(dm, &ipdm)); /* with DMPlexOrientInterface_Internal() call skipped so that PointSF issues are left to DMPlexComparePointSFWithInterface_Private() below */
1490: PetscCall(DMPlexSetOrientInterface_Private(dm, PETSC_TRUE));
1491: }
1492: PetscCall(DMViewFromOptions(ipdm, NULL, "-ipdm_view"));
1494: /* recover ipdm whole-domain boundary faces from the expanded vertices coordinates */
1495: PetscCall(DMPlexExpandedVerticesCoordinatesToFaces_Private(ipdm, bnd, &boundary_faces_is));
1496: /* get inter-partition interface faces (interface_faces_is)*/
1497: PetscCall(DMPlexGetInterfaceFaces_Private(ipdm, boundary_faces_is, &interface_faces_is));
1498: /* compute inter-partition interface including edges and vertices (interface_is) */
1499: PetscCall(DMPlexComputeCompleteInterface_Private(ipdm, interface_faces_is, &interface_is));
1500: /* destroy immediate ISs */
1501: PetscCall(ISDestroy(&boundary_faces_is));
1502: PetscCall(ISDestroy(&interface_faces_is));
1504: /* for uninterpolated case, keep just vertices in interface */
1505: if (!intp) {
1506: PetscCall(DMPlexISFilterVertices_Private(ipdm, interface_is));
1507: PetscCall(DMDestroy(&ipdm));
1508: }
1510: /* compare PointSF with the boundary reconstructed from coordinates */
1511: PetscCall(DMPlexComparePointSFWithInterface_Private(dm, interface_is));
1512: PetscCall(PetscPrintf(comm, "DMPlexCheckPointSFHeavy PASSED\n"));
1513: PetscCall(ISDestroy(&interface_is));
1514: PetscFunctionReturn(PETSC_SUCCESS);
1515: }
1517: int main(int argc, char **argv)
1518: {
1519: DM dm;
1520: AppCtx user;
1522: PetscFunctionBeginUser;
1523: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
1524: PetscCall(PetscLogStageRegister("create", &stage[0]));
1525: PetscCall(PetscLogStageRegister("distribute", &stage[1]));
1526: PetscCall(PetscLogStageRegister("interpolate", &stage[2]));
1527: PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
1528: PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
1529: if (user.nPointsToExpand) PetscCall(TestExpandPoints(dm, &user));
1530: if (user.ncoords) {
1531: Vec coords;
1533: PetscCall(VecCreateSeqWithArray(PETSC_COMM_SELF, user.ncoords, user.ncoords, user.coords, &coords));
1534: PetscCall(ViewVerticesFromCoords(dm, coords, user.coordsTol, PETSC_VIEWER_STDOUT_WORLD));
1535: PetscCall(VecDestroy(&coords));
1536: }
1537: PetscCall(DMDestroy(&dm));
1538: PetscCall(PetscFinalize());
1539: return 0;
1540: }
1542: /*TEST
1544: testset:
1545: nsize: 2
1546: args: -dm_view ascii::ascii_info_detail
1547: args: -dm_plex_check_all
1548: test:
1549: suffix: 1_tri_dist0
1550: args: -distribute 0 -interpolate {{none create}separate output}
1551: test:
1552: suffix: 1_tri_dist1
1553: args: -distribute 1 -interpolate {{none create after_distribute}separate output}
1554: test:
1555: suffix: 1_quad_dist0
1556: args: -cell_simplex 0 -distribute 0 -interpolate {{none create}separate output}
1557: test:
1558: suffix: 1_quad_dist1
1559: args: -cell_simplex 0 -distribute 1 -interpolate {{none create after_distribute}separate output}
1560: test:
1561: suffix: 1_1d_dist1
1562: args: -dim 1 -distribute 1
1564: testset:
1565: nsize: 3
1566: args: -testnum 1 -interpolate create
1567: args: -dm_plex_check_all
1568: test:
1569: suffix: 2
1570: args: -dm_view ascii::ascii_info_detail
1571: test:
1572: suffix: 2a
1573: args: -dm_plex_check_cones_conform_on_interfaces_verbose
1574: test:
1575: suffix: 2b
1576: args: -test_expand_points 0,1,2,5,6
1577: test:
1578: suffix: 2c
1579: args: -test_expand_points 0,1,2,5,6 -test_expand_points_empty
1581: testset:
1582: # the same as 1% for 3D
1583: nsize: 2
1584: args: -dim 3 -dm_view ascii::ascii_info_detail
1585: args: -dm_plex_check_all
1586: test:
1587: suffix: 4_tet_dist0
1588: args: -distribute 0 -interpolate {{none create}separate output}
1589: test:
1590: suffix: 4_tet_dist1
1591: args: -distribute 1 -interpolate {{none create after_distribute}separate output}
1592: test:
1593: suffix: 4_hex_dist0
1594: args: -cell_simplex 0 -distribute 0 -interpolate {{none create}separate output}
1595: test:
1596: suffix: 4_hex_dist1
1597: args: -cell_simplex 0 -distribute 1 -interpolate {{none create after_distribute}separate output}
1599: test:
1600: # the same as 4_tet_dist0 but test different initial orientations
1601: suffix: 4_tet_test_orient
1602: nsize: 2
1603: args: -dim 3 -distribute 0
1604: args: -dm_plex_check_all
1605: args: -rotate_interface_0 {{0 1 2 11 12 13}}
1606: args: -rotate_interface_1 {{0 1 2 11 12 13}}
1608: testset:
1609: requires: exodusii
1610: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/TwoQuads.exo
1611: args: -dm_view ascii::ascii_info_detail
1612: args: -dm_plex_check_all
1613: args: -custom_view
1614: test:
1615: suffix: 5_seq
1616: nsize: 1
1617: args: -distribute 0 -interpolate {{none create}separate output}
1618: test:
1619: # Detail viewing in a non-distributed mesh is broken because the DMLabelView() is collective, but the label is not shared
1620: suffix: 5_dist0
1621: nsize: 2
1622: args: -distribute 0 -interpolate {{none create}separate output} -dm_view
1623: test:
1624: suffix: 5_dist1
1625: nsize: 2
1626: args: -distribute 1 -interpolate {{none create after_distribute}separate output}
1628: testset:
1629: nsize: {{1 2 4}}
1630: args: -use_generator
1631: args: -dm_plex_check_all
1632: args: -distribute -interpolate none
1633: test:
1634: suffix: 6_tri
1635: requires: triangle
1636: args: -faces {{2,2 1,3 7,4}} -cell_simplex 1 -dm_generator triangle
1637: test:
1638: suffix: 6_quad
1639: args: -faces {{2,2 1,3 7,4}} -cell_simplex 0
1640: test:
1641: suffix: 6_tet
1642: requires: ctetgen
1643: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 1 -dm_generator ctetgen
1644: test:
1645: suffix: 6_hex
1646: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 0
1647: testset:
1648: nsize: {{1 2 4}}
1649: args: -use_generator
1650: args: -dm_plex_check_all
1651: args: -distribute -interpolate create
1652: test:
1653: suffix: 6_int_tri
1654: requires: triangle
1655: args: -faces {{2,2 1,3 7,4}} -cell_simplex 1 -dm_generator triangle
1656: test:
1657: suffix: 6_int_quad
1658: args: -faces {{2,2 1,3 7,4}} -cell_simplex 0
1659: test:
1660: suffix: 6_int_tet
1661: requires: ctetgen
1662: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 1 -dm_generator ctetgen
1663: test:
1664: suffix: 6_int_hex
1665: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 0
1666: testset:
1667: nsize: {{2 4}}
1668: args: -use_generator
1669: args: -dm_plex_check_all
1670: args: -distribute -interpolate after_distribute
1671: test:
1672: suffix: 6_parint_tri
1673: requires: triangle
1674: args: -faces {{2,2 1,3 7,4}} -cell_simplex 1 -dm_generator triangle
1675: test:
1676: suffix: 6_parint_quad
1677: args: -faces {{2,2 1,3 7,4}} -cell_simplex 0
1678: test:
1679: suffix: 6_parint_tet
1680: requires: ctetgen
1681: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 1 -dm_generator ctetgen
1682: test:
1683: suffix: 6_parint_hex
1684: args: -faces {{2,2,2 1,3,5 3,4,7}} -cell_simplex 0
1686: testset: # 7 EXODUS
1687: requires: exodusii
1688: args: -dm_plex_check_all
1689: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo
1690: args: -distribute
1691: test: # seq load, simple partitioner
1692: suffix: 7_exo
1693: nsize: {{1 2 4 5}}
1694: args: -interpolate none
1695: test: # seq load, seq interpolation, simple partitioner
1696: suffix: 7_exo_int_simple
1697: nsize: {{1 2 4 5}}
1698: args: -interpolate create
1699: test: # seq load, seq interpolation, metis partitioner
1700: suffix: 7_exo_int_metis
1701: requires: parmetis
1702: nsize: {{2 4 5}}
1703: args: -interpolate create
1704: args: -petscpartitioner_type parmetis
1705: test: # seq load, simple partitioner, par interpolation
1706: suffix: 7_exo_simple_int
1707: nsize: {{2 4 5}}
1708: args: -interpolate after_distribute
1709: test: # seq load, metis partitioner, par interpolation
1710: suffix: 7_exo_metis_int
1711: requires: parmetis
1712: nsize: {{2 4 5}}
1713: args: -interpolate after_distribute
1714: args: -petscpartitioner_type parmetis
1716: testset: # 7 HDF5 SEQUANTIAL LOAD
1717: requires: hdf5 !complex
1718: args: -dm_plex_check_all
1719: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.h5 -dm_plex_create_from_hdf5_xdmf
1720: args: -dm_plex_hdf5_force_sequential
1721: args: -distribute
1722: test: # seq load, simple partitioner
1723: suffix: 7_seq_hdf5_simple
1724: nsize: {{1 2 4 5}}
1725: args: -interpolate none
1726: test: # seq load, seq interpolation, simple partitioner
1727: suffix: 7_seq_hdf5_int_simple
1728: nsize: {{1 2 4 5}}
1729: args: -interpolate after_create
1730: test: # seq load, seq interpolation, metis partitioner
1731: nsize: {{2 4 5}}
1732: suffix: 7_seq_hdf5_int_metis
1733: requires: parmetis
1734: args: -interpolate after_create
1735: args: -petscpartitioner_type parmetis
1736: test: # seq load, simple partitioner, par interpolation
1737: suffix: 7_seq_hdf5_simple_int
1738: nsize: {{2 4 5}}
1739: args: -interpolate after_distribute
1740: test: # seq load, metis partitioner, par interpolation
1741: nsize: {{2 4 5}}
1742: suffix: 7_seq_hdf5_metis_int
1743: requires: parmetis
1744: args: -interpolate after_distribute
1745: args: -petscpartitioner_type parmetis
1747: testset: # 7 HDF5 PARALLEL LOAD
1748: requires: hdf5 !complex
1749: nsize: {{2 4 5}}
1750: args: -dm_plex_check_all
1751: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.h5 -dm_plex_create_from_hdf5_xdmf
1752: test: # par load
1753: suffix: 7_par_hdf5
1754: args: -interpolate none
1755: test: # par load, par interpolation
1756: suffix: 7_par_hdf5_int
1757: args: -interpolate after_create
1758: test: # par load, parmetis repartitioner
1759: TODO: Parallel partitioning of uninterpolated meshes not supported
1760: suffix: 7_par_hdf5_parmetis
1761: requires: parmetis
1762: args: -distribute -petscpartitioner_type parmetis
1763: args: -interpolate none
1764: test: # par load, par interpolation, parmetis repartitioner
1765: suffix: 7_par_hdf5_int_parmetis
1766: requires: parmetis
1767: args: -distribute -petscpartitioner_type parmetis
1768: args: -interpolate after_create
1769: test: # par load, parmetis partitioner, par interpolation
1770: TODO: Parallel partitioning of uninterpolated meshes not supported
1771: suffix: 7_par_hdf5_parmetis_int
1772: requires: parmetis
1773: args: -distribute -petscpartitioner_type parmetis
1774: args: -interpolate after_distribute
1776: test:
1777: suffix: 7_hdf5_hierarch
1778: requires: hdf5 ptscotch !complex
1779: nsize: {{2 3 4}separate output}
1780: args: -distribute
1781: args: -interpolate after_create
1782: args: -petscpartitioner_type matpartitioning -petscpartitioner_view ::ascii_info
1783: args: -mat_partitioning_type hierarch -mat_partitioning_hierarchical_nfineparts 2
1784: args: -mat_partitioning_hierarchical_coarseparttype ptscotch -mat_partitioning_hierarchical_fineparttype ptscotch
1786: test:
1787: suffix: 8
1788: requires: hdf5 !complex
1789: nsize: 4
1790: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.h5 -dm_plex_create_from_hdf5_xdmf
1791: args: -distribute 0 -interpolate after_create
1792: args: -view_vertices_from_coords 0.,1.,0.,-0.5,1.,0.,0.583,-0.644,0.,-2.,-2.,-2. -view_vertices_from_coords_tol 1e-3
1793: args: -dm_plex_check_all
1794: args: -custom_view
1796: testset: # 9 HDF5 SEQUANTIAL LOAD
1797: requires: hdf5 !complex datafilespath
1798: args: -dm_plex_check_all
1799: args: -filename ${DATAFILESPATH}/meshes/cube-hexahedra-refined.h5 -dm_plex_create_from_hdf5_xdmf -dm_plex_hdf5_topology_path /cells -dm_plex_hdf5_geometry_path /coordinates
1800: args: -dm_plex_hdf5_force_sequential
1801: args: -distribute
1802: test: # seq load, simple partitioner
1803: suffix: 9_seq_hdf5_simple
1804: nsize: {{1 2 4 5}}
1805: args: -interpolate none
1806: test: # seq load, seq interpolation, simple partitioner
1807: suffix: 9_seq_hdf5_int_simple
1808: nsize: {{1 2 4 5}}
1809: args: -interpolate after_create
1810: test: # seq load, seq interpolation, metis partitioner
1811: nsize: {{2 4 5}}
1812: suffix: 9_seq_hdf5_int_metis
1813: requires: parmetis
1814: args: -interpolate after_create
1815: args: -petscpartitioner_type parmetis
1816: test: # seq load, simple partitioner, par interpolation
1817: suffix: 9_seq_hdf5_simple_int
1818: nsize: {{2 4 5}}
1819: args: -interpolate after_distribute
1820: test: # seq load, simple partitioner, par interpolation
1821: # This is like 9_seq_hdf5_simple_int but testing error output of DMPlexCheckPointSFHeavy().
1822: # Once 9_seq_hdf5_simple_int gets fixed, this one gets broken.
1823: # We can then provide an intentionally broken mesh instead.
1824: TODO: This test is broken because PointSF is fixed.
1825: suffix: 9_seq_hdf5_simple_int_err
1826: nsize: 4
1827: args: -interpolate after_distribute
1828: filter: sed -e "/PETSC ERROR/,$$d"
1829: test: # seq load, metis partitioner, par interpolation
1830: nsize: {{2 4 5}}
1831: suffix: 9_seq_hdf5_metis_int
1832: requires: parmetis
1833: args: -interpolate after_distribute
1834: args: -petscpartitioner_type parmetis
1836: testset: # 9 HDF5 PARALLEL LOAD
1837: requires: hdf5 !complex datafilespath
1838: nsize: {{2 4 5}}
1839: args: -dm_plex_check_all
1840: args: -filename ${DATAFILESPATH}/meshes/cube-hexahedra-refined.h5 -dm_plex_create_from_hdf5_xdmf -dm_plex_hdf5_topology_path /cells -dm_plex_hdf5_geometry_path /coordinates
1841: test: # par load
1842: suffix: 9_par_hdf5
1843: args: -interpolate none
1844: test: # par load, par interpolation
1845: suffix: 9_par_hdf5_int
1846: args: -interpolate after_create
1847: test: # par load, parmetis repartitioner
1848: TODO: Parallel partitioning of uninterpolated meshes not supported
1849: suffix: 9_par_hdf5_parmetis
1850: requires: parmetis
1851: args: -distribute -petscpartitioner_type parmetis
1852: args: -interpolate none
1853: test: # par load, par interpolation, parmetis repartitioner
1854: suffix: 9_par_hdf5_int_parmetis
1855: requires: parmetis
1856: args: -distribute -petscpartitioner_type parmetis
1857: args: -interpolate after_create
1858: test: # par load, parmetis partitioner, par interpolation
1859: TODO: Parallel partitioning of uninterpolated meshes not supported
1860: suffix: 9_par_hdf5_parmetis_int
1861: requires: parmetis
1862: args: -distribute -petscpartitioner_type parmetis
1863: args: -interpolate after_distribute
1865: testset: # 10 HDF5 PARALLEL LOAD
1866: requires: hdf5 !complex datafilespath
1867: nsize: {{2 4 7}}
1868: args: -dm_plex_check_all
1869: args: -filename ${DATAFILESPATH}/meshes/cube-hexahedra-refined2.h5 -dm_plex_create_from_hdf5_xdmf -dm_plex_hdf5_topology_path /topo -dm_plex_hdf5_geometry_path /geom
1870: test: # par load, par interpolation
1871: suffix: 10_par_hdf5_int
1872: args: -interpolate after_create
1873: TEST*/