Actual source code: ex4.c
1: static char help[] = "Tests for refinement of meshes created by hand\n\n";
3: #include <petscdmplex.h>
5: typedef struct {
6: PetscInt debug; /* The debugging level */
7: PetscInt dim; /* The topological mesh dimension */
8: PetscBool cellHybrid; /* Use a hybrid mesh */
9: PetscBool cellSimplex; /* Use simplices or hexes */
10: PetscBool testPartition; /* Use a fixed partitioning for testing */
11: PetscInt testNum; /* The particular mesh to test */
12: PetscBool uninterpolate; /* Uninterpolate the mesh at the end */
13: PetscBool reinterpolate; /* Reinterpolate the mesh at the end */
14: } AppCtx;
16: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
17: {
18: options->debug = 0;
19: options->dim = 2;
20: options->cellHybrid = PETSC_TRUE;
21: options->cellSimplex = PETSC_TRUE;
22: options->testPartition = PETSC_TRUE;
23: options->testNum = 0;
24: options->uninterpolate = PETSC_FALSE;
25: options->reinterpolate = PETSC_FALSE;
27: PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");
28: PetscOptionsBoundedInt("-debug", "The debugging level", "ex4.c", options->debug, &options->debug, NULL, 0);
29: PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex4.c", options->dim, &options->dim, NULL, 1, 3);
30: PetscOptionsBool("-cell_hybrid", "Use a hybrid mesh", "ex4.c", options->cellHybrid, &options->cellHybrid, NULL);
31: PetscOptionsBool("-cell_simplex", "Use simplices if true, otherwise hexes", "ex4.c", options->cellSimplex, &options->cellSimplex, NULL);
32: PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex4.c", options->testPartition, &options->testPartition, NULL);
33: PetscOptionsBoundedInt("-test_num", "The particular mesh to test", "ex4.c", options->testNum, &options->testNum, NULL, 0);
34: PetscOptionsBool("-uninterpolate", "Uninterpolate the mesh at the end", "ex4.c", options->uninterpolate, &options->uninterpolate, NULL);
35: PetscOptionsBool("-reinterpolate", "Reinterpolate the mesh at the end", "ex4.c", options->reinterpolate, &options->reinterpolate, NULL);
36: PetscOptionsEnd();
37: return 0;
38: }
40: /* Two segments
42: 2-------0-------3-------1-------4
44: become
46: 4---0---7---1---5---2---8---3---6
48: */
49: PetscErrorCode CreateSimplex_1D(MPI_Comm comm, DM *dm)
50: {
51: PetscInt depth = 1;
52: PetscMPIInt rank;
54: MPI_Comm_rank(comm, &rank);
55: if (rank == 0) {
56: PetscInt numPoints[2] = {3, 2};
57: PetscInt coneSize[5] = {2, 2, 0, 0, 0};
58: PetscInt cones[4] = {2, 3, 3, 4};
59: PetscInt coneOrientations[16] = {0, 0, 0, 0};
60: PetscScalar vertexCoords[3] = {-1.0, 0.0, 1.0};
62: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
63: } else {
64: PetscInt numPoints[2] = {0, 0};
66: DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
67: }
68: return 0;
69: }
71: /* Two triangles
72: 4
73: / | \
74: 8 | 10
75: / | \
76: 2 0 7 1 5
77: \ | /
78: 6 | 9
79: \ | /
80: 3
82: Becomes
83: 10
84: / | \
85: 21 | 26
86: / | \
87: 14 2 20 4 16
88: /|\ | /|\
89: 22 | 28 | 32 | 25
90: / | \ | / | 6\
91: 8 29 3 13 7 31 11
92: \0 | / | \ | /
93: 17 | 27 | 30 | 24
94: \|/ | \|/
95: 12 1 19 5 15
96: \ | /
97: 18 | 23
98: \ | /
99: 9
100: */
101: PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM *dm)
102: {
103: PetscInt depth = 2;
104: PetscMPIInt rank;
106: MPI_Comm_rank(comm, &rank);
107: if (rank == 0) {
108: PetscInt numPoints[3] = {4, 5, 2};
109: PetscInt coneSize[11] = {3, 3, 0, 0, 0, 0, 2, 2, 2, 2, 2};
110: PetscInt cones[16] = {6, 7, 8, 7, 9, 10, 2, 3, 3, 4, 4, 2, 3, 5, 5, 4};
111: PetscInt coneOrientations[16] = {0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
112: PetscScalar vertexCoords[8] = {-0.5, 0.0, 0.0, -0.5, 0.0, 0.5, 0.5, 0.0};
114: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
115: } else {
116: PetscInt numPoints[3] = {0, 0, 0};
118: DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
119: }
120: return 0;
121: }
123: /* Two triangles separated by a zero-volume cell with 4 vertices/2 edges
124: 5--16--8
125: / | | \
126: 11 | | 12
127: / | | \
128: 3 0 10 2 14 1 6
129: \ | | /
130: 9 | | 13
131: \ | | /
132: 4--15--7
133: */
134: PetscErrorCode CreateSimplexHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
135: {
136: DM idm, hdm = NULL;
137: DMLabel faultLabel, hybridLabel;
138: PetscInt p;
139: PetscMPIInt rank;
141: MPI_Comm_rank(comm, &rank);
142: if (rank == 0) {
143: switch (testNum) {
144: case 0: {
145: PetscInt numPoints[2] = {4, 2};
146: PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0};
147: PetscInt cones[6] = {2, 3, 4, 5, 4, 3};
148: PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0};
149: PetscScalar vertexCoords[8] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, 1.0, 0.5};
150: PetscInt faultPoints[2] = {3, 4};
152: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
153: for (p = 0; p < 2; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
154: } break;
155: case 1: {
156: PetscInt numPoints[2] = {5, 4};
157: PetscInt coneSize[9] = {3, 3, 3, 3, 0, 0, 0, 0, 0};
158: PetscInt cones[12] = {4, 5, 6, 6, 7, 4, 6, 5, 8, 6, 8, 7};
159: PetscInt coneOrientations[12] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
160: PetscScalar vertexCoords[10] = {-1.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 1.0, 0.0};
161: PetscInt faultPoints[3] = {5, 6, 7};
163: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
164: for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
165: } break;
166: default:
167: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
168: }
169: DMPlexInterpolate(*dm, &idm);
170: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
171: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
172: DMSetFromOptions(idm);
173: DMViewFromOptions(idm, NULL, "-dm_view");
174: DMGetLabel(*dm, "fault", &faultLabel);
175: DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
176: DMLabelDestroy(&hybridLabel);
177: DMDestroy(&idm);
178: DMDestroy(dm);
179: *dm = hdm;
180: } else {
181: PetscInt numPoints[2] = {0, 0};
183: DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
184: DMPlexInterpolate(*dm, &idm);
185: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
186: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
187: DMSetFromOptions(idm);
188: DMViewFromOptions(idm, NULL, "-dm_view");
189: DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
190: DMDestroy(&idm);
191: DMDestroy(dm);
192: *dm = hdm;
193: }
194: return 0;
195: }
197: /* Two quadrilaterals
199: 5----10-----4----14-----7
200: | | |
201: | | |
202: | | |
203: 11 0 9 1 13
204: | | |
205: | | |
206: | | |
207: 2-----8-----3----12-----6
208: */
209: PetscErrorCode CreateTensorProduct_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
210: {
211: PetscInt depth = 2;
212: PetscMPIInt rank;
214: MPI_Comm_rank(comm, &rank);
215: if (rank == 0) {
216: PetscInt numPoints[3] = {6, 7, 2};
217: PetscInt coneSize[15] = {4, 4, 0, 0, 0, 0, 0, 0, 2, 2, 2, 2, 2, 2, 2};
218: PetscInt cones[22] = {8, 9, 10, 11, 12, 13, 14, 9, 2, 3, 3, 4, 4, 5, 5, 2, 3, 6, 6, 7, 7, 4};
219: PetscInt coneOrientations[22] = {0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
220: PetscScalar vertexCoords[12] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5};
222: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
223: } else {
224: PetscInt numPoints[3] = {0, 0, 0};
226: DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
227: }
228: return 0;
229: }
231: PetscErrorCode CreateTensorProductHybrid_2D(MPI_Comm comm, PetscInt testNum, DM *dm)
232: {
233: DM idm, hdm = NULL;
234: DMLabel faultLabel, hybridLabel;
235: PetscInt p;
236: PetscMPIInt rank;
238: MPI_Comm_rank(comm, &rank);
239: if (rank == 0) {
240: PetscInt numPoints[2] = {6, 2};
241: PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0};
242: PetscInt cones[8] = {
243: 2, 3, 4, 5, 3, 6, 7, 4,
244: };
245: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
246: PetscScalar vertexCoords[12] = {-1.0, -0.5, 0.0, -0.5, 0.0, 0.5, -1.0, 0.5, 1.0, -0.5, 1.0, 0.5};
247: PetscInt faultPoints[2] = {3, 4};
249: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
250: for (p = 0; p < 2; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
251: DMPlexInterpolate(*dm, &idm);
252: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
253: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
254: DMSetFromOptions(idm);
255: DMViewFromOptions(idm, NULL, "-dm_view");
256: DMGetLabel(*dm, "fault", &faultLabel);
257: DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
258: DMLabelDestroy(&hybridLabel);
259: } else {
260: PetscInt numPoints[3] = {0, 0, 0};
262: DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
263: DMPlexInterpolate(*dm, &idm);
264: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
265: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
266: DMSetFromOptions(idm);
267: DMViewFromOptions(idm, NULL, "-dm_view");
268: DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
269: }
270: DMDestroy(&idm);
271: DMDestroy(dm);
272: *dm = hdm;
273: return 0;
274: }
276: /* Two tetrahedrons
278: cell 5 5______ cell
279: 0 / | \ |\ \ 1
280: 17 | 18 | 18 13 21
281: /8 19 10\ 19 \ \
282: 2-14-|----4 | 4--22--6
283: \ 9 | 7 / |10 / /
284: 16 | 15 | 15 12 20
285: \ | / |/ /
286: 3 3------
287: */
288: PetscErrorCode CreateSimplex_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
289: {
290: DM idm;
291: PetscInt depth = 3;
292: PetscMPIInt rank;
294: MPI_Comm_rank(comm, &rank);
295: if (rank == 0) {
296: switch (testNum) {
297: case 0: {
298: PetscInt numPoints[4] = {5, 9, 7, 2};
299: PetscInt coneSize[23] = {4, 4, 0, 0, 0, 0, 0, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2};
300: PetscInt cones[47] = {7, 8, 9, 10, 10, 11, 12, 13, 14, 15, 16, 17, 18, 14, 16, 19, 17, 15, 18, 19, 20, 21, 19, 15, 22, 20, 18, 21, 22, 2, 4, 4, 3, 3, 2, 2, 5, 5, 4, 3, 5, 3, 6, 6, 5, 4, 6};
301: PetscInt coneOrientations[47] = {0, 0, 0, 0, -2, 0, 0, 0, 0, 0, 0, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, -1, -1, 0, -1, -1, -1, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
302: PetscScalar vertexCoords[15] = {0.0, 0.0, -0.5, 0.0, -0.5, 0.0, 1.0, 0.0, 0.0, 0.0, 0.5, 0.0, 0.0, 0.0, 0.5};
304: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
305: } break;
306: case 1: {
307: PetscInt numPoints[2] = {5, 2};
308: PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0};
309: PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6};
310: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
311: PetscScalar vertexCoords[15] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0};
313: depth = 1;
314: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
315: DMPlexInterpolate(*dm, &idm);
316: DMViewFromOptions(idm, NULL, "-in_dm_view");
317: DMDestroy(dm);
318: *dm = idm;
319: } break;
320: case 2: {
321: PetscInt numPoints[2] = {4, 1};
322: PetscInt coneSize[5] = {4, 0, 0, 0, 0};
323: PetscInt cones[4] = {2, 3, 4, 1};
324: PetscInt coneOrientations[4] = {0, 0, 0, 0};
325: PetscScalar vertexCoords[12] = {0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0};
327: depth = 1;
328: DMPlexCreateFromDAG(*dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
329: DMPlexInterpolate(*dm, &idm);
330: DMViewFromOptions(idm, NULL, "-in_dm_view");
331: DMDestroy(dm);
332: *dm = idm;
333: } break;
334: default:
335: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
336: }
337: } else {
338: PetscInt numPoints[4] = {0, 0, 0, 0};
340: DMPlexCreateFromDAG(*dm, depth, numPoints, NULL, NULL, NULL, NULL);
341: switch (testNum) {
342: case 1:
343: DMPlexInterpolate(*dm, &idm);
344: DMViewFromOptions(idm, NULL, "-in_dm_view");
345: DMDestroy(dm);
346: *dm = idm;
347: break;
348: }
349: }
350: return 0;
351: }
353: /* Two tetrahedrons separated by a zero-volume cell with 6 vertices
355: cell 6 ___33___10______ cell
356: 0 / | \ |\ \ 1
357: 21 | 23 | 29 27
358: /12 24 14\ 30 \ \
359: 3-20-|----5--32-|---9--26--7
360: \ 13| 11/ |18 / /
361: 19 | 22 | 28 25
362: \ | / |/ /
363: 4----31----8------
364: cell 2
365: */
366: PetscErrorCode CreateSimplexHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
367: {
368: DM idm, hdm = NULL;
369: DMLabel faultLabel, hybridLabel;
370: PetscInt p;
371: PetscMPIInt rank;
373: MPI_Comm_rank(comm, &rank);
374: if (rank == 0) {
375: switch (testNum) {
376: case 0: {
377: PetscInt numPoints[2] = {5, 2};
378: PetscInt coneSize[7] = {4, 4, 0, 0, 0, 0, 0};
379: PetscInt cones[8] = {4, 3, 5, 2, 5, 3, 4, 6};
380: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
381: PetscScalar vertexCoords[15] = {-1.0, 0.0, 0.0, 0.0, -1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 0.0, 0.0};
382: PetscInt faultPoints[3] = {3, 4, 5};
384: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
385: for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
386: } break;
387: case 1: {
388: /* Tets 0,3,5 and 1,2,4 */
389: PetscInt numPoints[2] = {9, 6};
390: PetscInt coneSize[15] = {4, 4, 4, 4, 4, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0};
391: PetscInt cones[24] = {7, 8, 9, 6, 11, 13, 9, 14, 10, 11, 13, 9, 9, 10, 11, 7, 9, 14, 13, 12, 7, 8, 11, 9};
392: PetscInt coneOrientations[24] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
393: PetscScalar vertexCoords[27] = {-2.0, -1.0, 0.0, -2.0, 0.0, 0.0, -2.0, 0.0, 1.0, 0.0, -1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1.0, 2.0, -1.0, 0.0, 2.0, 0.0, 0.0, 2.0, 0.0, 1.0};
394: PetscInt faultPoints[3] = {9, 10, 11};
396: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
397: for (p = 0; p < 3; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
398: } break;
399: default:
400: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
401: }
402: DMPlexInterpolate(*dm, &idm);
403: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
404: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
405: DMSetFromOptions(idm);
406: DMViewFromOptions(idm, NULL, "-dm_view");
407: DMGetLabel(*dm, "fault", &faultLabel);
408: DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, &hybridLabel, NULL, NULL, &hdm);
409: DMLabelDestroy(&hybridLabel);
410: DMDestroy(&idm);
411: DMDestroy(dm);
412: *dm = hdm;
413: } else {
414: PetscInt numPoints[4] = {0, 0, 0, 0};
416: DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
417: DMPlexInterpolate(*dm, &idm);
418: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
419: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
420: DMSetFromOptions(idm);
421: DMViewFromOptions(idm, NULL, "-dm_view");
422: DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
423: DMDestroy(&idm);
424: DMDestroy(dm);
425: *dm = hdm;
426: }
427: return 0;
428: }
430: PetscErrorCode CreateTensorProduct_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
431: {
432: DM idm;
433: PetscMPIInt rank;
435: MPI_Comm_rank(comm, &rank);
436: if (rank == 0) {
437: switch (testNum) {
438: case 0: {
439: PetscInt numPoints[2] = {12, 2};
440: PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
441: PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8};
442: PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
443: PetscScalar vertexCoords[36] = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5};
445: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
446: } break;
447: case 1: {
448: PetscInt numPoints[2] = {8, 1};
449: PetscInt coneSize[9] = {8, 0, 0, 0, 0, 0, 0, 0, 0};
450: PetscInt cones[8] = {1, 2, 3, 4, 5, 6, 7, 8};
451: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
452: PetscScalar vertexCoords[24] = {-1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, -1.0, 1.0, -1.0, -1.0, -1.0, -1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0};
454: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
455: } break;
456: default:
457: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
458: }
459: } else {
460: PetscInt numPoints[4] = {0, 0, 0, 0};
462: DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
463: }
464: DMPlexInterpolate(*dm, &idm);
465: DMViewFromOptions(idm, NULL, "-in_dm_view");
466: DMDestroy(dm);
467: *dm = idm;
468: return 0;
469: }
471: PetscErrorCode CreateTensorProductHybrid_3D(MPI_Comm comm, PetscInt testNum, DM *dm)
472: {
473: DM idm, hdm = NULL;
474: DMLabel faultLabel;
475: PetscInt p;
476: PetscMPIInt rank;
478: MPI_Comm_rank(comm, &rank);
479: if (rank == 0) {
480: switch (testNum) {
481: case 0: {
482: PetscInt numPoints[2] = {12, 2};
483: PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
484: PetscInt cones[16] = {2, 3, 4, 5, 6, 7, 8, 9, 5, 4, 10, 11, 7, 12, 13, 8};
485: PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
486: PetscScalar vertexCoords[36] = {-1.0, -0.5, -0.5, -1.0, 0.5, -0.5, 0.0, 0.5, -0.5, 0.0, -0.5, -0.5, -1.0, -0.5, 0.5, 0.0, -0.5, 0.5, 0.0, 0.5, 0.5, -1.0, 0.5, 0.5, 1.0, 0.5, -0.5, 1.0, -0.5, -0.5, 1.0, -0.5, 0.5, 1.0, 0.5, 0.5};
487: PetscInt faultPoints[4] = {2, 3, 5, 6};
489: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
490: for (p = 0; p < 4; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
491: } break;
492: case 1: {
493: PetscInt numPoints[2] = {30, 7};
494: PetscInt coneSize[37] = {8, 8, 8, 8, 8, 8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
495: PetscInt cones[56] = {8, 21, 20, 7, 13, 12, 23, 24, 14, 15, 10, 9, 13, 8, 21, 24, 15, 16, 11, 10, 24, 21, 22, 25, 30, 29, 28, 21, 35, 24, 33, 34, 24, 21, 30, 35, 25, 36, 31, 22, 27, 20, 21, 28, 32, 33, 24, 23, 15, 24, 13, 14, 19, 18, 17, 26};
496: PetscInt coneOrientations[56] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
497: PetscScalar vertexCoords[90] = {-2.0, -2.0, -2.0, -2.0, -1.0, -2.0, -3.0, 0.0, -2.0, -2.0, 1.0, -2.0, -2.0, 2.0, -2.0, -2.0, -2.0, 0.0, -2.0, -1.0, 0.0, -3.0, 0.0, 0.0, -2.0, 1.0, 0.0, -2.0, 2.0, 0.0,
498: -2.0, -1.0, 2.0, -3.0, 0.0, 2.0, -2.0, 1.0, 2.0, 0.0, -2.0, -2.0, 0.0, 0.0, -2.0, 0.0, 2.0, -2.0, 0.0, -2.0, 0.0, 0.0, 0.0, 0.0, 0.0, 2.0, 0.0, 0.0, 0.0, 2.0,
499: 2.0, -2.0, -2.0, 2.0, -1.0, -2.0, 3.0, 0.0, -2.0, 2.0, 1.0, -2.0, 2.0, 2.0, -2.0, 2.0, -2.0, 0.0, 2.0, -1.0, 0.0, 3.0, 0.0, 0.0, 2.0, 1.0, 0.0, 2.0, 2.0, 0.0};
500: PetscInt faultPoints[6] = {20, 21, 22, 23, 24, 25};
502: DMPlexCreateFromDAG(*dm, 1, numPoints, coneSize, cones, coneOrientations, vertexCoords);
503: for (p = 0; p < 6; ++p) DMSetLabelValue(*dm, "fault", faultPoints[p], 1);
504: } break;
505: default:
506: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
507: }
508: DMPlexInterpolate(*dm, &idm);
509: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
510: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
511: DMSetFromOptions(idm);
512: DMViewFromOptions(idm, NULL, "-dm_view");
513: DMGetLabel(*dm, "fault", &faultLabel);
514: DMPlexCreateHybridMesh(idm, faultLabel, NULL, 0, NULL, NULL, NULL, &hdm);
515: DMDestroy(&idm);
516: DMDestroy(dm);
517: *dm = hdm;
518: } else {
519: PetscInt numPoints[4] = {0, 0, 0, 0};
521: DMPlexCreateFromDAG(*dm, 1, numPoints, NULL, NULL, NULL, NULL);
522: DMPlexInterpolate(*dm, &idm);
523: PetscObjectSetOptionsPrefix((PetscObject)idm, "in_");
524: DMPlexDistributeSetDefault(idm, PETSC_FALSE);
525: DMSetFromOptions(idm);
526: DMViewFromOptions(idm, NULL, "-dm_view");
527: DMPlexCreateHybridMesh(idm, NULL, NULL, 0, NULL, NULL, NULL, &hdm);
528: DMDestroy(&idm);
529: DMDestroy(dm);
530: *dm = hdm;
531: }
532: return 0;
533: }
535: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
536: {
537: PetscInt dim = user->dim;
538: PetscBool cellHybrid = user->cellHybrid;
539: PetscBool cellSimplex = user->cellSimplex;
540: PetscMPIInt rank, size;
542: MPI_Comm_rank(comm, &rank);
543: MPI_Comm_size(comm, &size);
544: DMCreate(comm, dm);
545: DMSetType(*dm, DMPLEX);
546: DMSetDimension(*dm, dim);
547: switch (dim) {
548: case 1:
550: CreateSimplex_1D(comm, dm);
551: break;
552: case 2:
553: if (cellSimplex) {
554: if (cellHybrid) {
555: CreateSimplexHybrid_2D(comm, user->testNum, dm);
556: } else {
557: CreateSimplex_2D(comm, dm);
558: }
559: } else {
560: if (cellHybrid) {
561: CreateTensorProductHybrid_2D(comm, user->testNum, dm);
562: } else {
563: CreateTensorProduct_2D(comm, user->testNum, dm);
564: }
565: }
566: break;
567: case 3:
568: if (cellSimplex) {
569: if (cellHybrid) {
570: CreateSimplexHybrid_3D(comm, user->testNum, dm);
571: } else {
572: CreateSimplex_3D(comm, user->testNum, dm);
573: }
574: } else {
575: if (cellHybrid) {
576: CreateTensorProductHybrid_3D(comm, user->testNum, dm);
577: } else {
578: CreateTensorProduct_3D(comm, user->testNum, dm);
579: }
580: }
581: break;
582: default:
583: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, dim);
584: }
585: if (user->testPartition && size > 1) {
586: PetscPartitioner part;
587: PetscInt *sizes = NULL;
588: PetscInt *points = NULL;
590: if (rank == 0) {
591: if (dim == 2 && cellSimplex && !cellHybrid && size == 2) {
592: switch (user->testNum) {
593: case 0: {
594: PetscInt triSizes_p2[2] = {1, 1};
595: PetscInt triPoints_p2[2] = {0, 1};
597: PetscMalloc2(2, &sizes, 2, &points);
598: PetscArraycpy(sizes, triSizes_p2, 2);
599: PetscArraycpy(points, triPoints_p2, 2);
600: break;
601: }
602: default:
603: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular mesh on 2 procs", user->testNum);
604: }
605: } else if (dim == 2 && cellSimplex && cellHybrid && size == 2) {
606: switch (user->testNum) {
607: case 0: {
608: PetscInt triSizes_p2[2] = {1, 2};
609: PetscInt triPoints_p2[3] = {0, 1, 2};
611: PetscMalloc2(2, &sizes, 3, &points);
612: PetscArraycpy(sizes, triSizes_p2, 2);
613: PetscArraycpy(points, triPoints_p2, 3);
614: break;
615: }
616: default:
617: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for triangular hybrid mesh on 2 procs", user->testNum);
618: }
619: } else if (dim == 2 && !cellSimplex && !cellHybrid && size == 2) {
620: switch (user->testNum) {
621: case 0: {
622: PetscInt quadSizes_p2[2] = {1, 1};
623: PetscInt quadPoints_p2[2] = {0, 1};
625: PetscMalloc2(2, &sizes, 2, &points);
626: PetscArraycpy(sizes, quadSizes_p2, 2);
627: PetscArraycpy(points, quadPoints_p2, 2);
628: break;
629: }
630: default:
631: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral mesh on 2 procs", user->testNum);
632: }
633: } else if (dim == 2 && !cellSimplex && cellHybrid && size == 2) {
634: switch (user->testNum) {
635: case 0: {
636: PetscInt quadSizes_p2[2] = {1, 2};
637: PetscInt quadPoints_p2[3] = {0, 1, 2};
639: PetscMalloc2(2, &sizes, 3, &points);
640: PetscArraycpy(sizes, quadSizes_p2, 2);
641: PetscArraycpy(points, quadPoints_p2, 3);
642: break;
643: }
644: default:
645: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for quadrilateral hybrid mesh on 2 procs", user->testNum);
646: }
647: } else if (dim == 3 && cellSimplex && !cellHybrid && size == 2) {
648: switch (user->testNum) {
649: case 0: {
650: PetscInt tetSizes_p2[2] = {1, 1};
651: PetscInt tetPoints_p2[2] = {0, 1};
653: PetscMalloc2(2, &sizes, 2, &points);
654: PetscArraycpy(sizes, tetSizes_p2, 2);
655: PetscArraycpy(points, tetPoints_p2, 2);
656: break;
657: }
658: case 1: {
659: PetscInt tetSizes_p2[2] = {1, 1};
660: PetscInt tetPoints_p2[2] = {0, 1};
662: PetscMalloc2(2, &sizes, 2, &points);
663: PetscArraycpy(sizes, tetSizes_p2, 2);
664: PetscArraycpy(points, tetPoints_p2, 2);
665: break;
666: }
667: default:
668: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for tetrahedral mesh on 2 procs", user->testNum);
669: }
670: } else if (dim == 3 && cellSimplex && cellHybrid && size == 2) {
671: switch (user->testNum) {
672: case 0: {
673: PetscInt tetSizes_p2[2] = {1, 2};
674: PetscInt tetPoints_p2[3] = {0, 1, 2};
676: PetscMalloc2(2, &sizes, 3, &points);
677: PetscArraycpy(sizes, tetSizes_p2, 2);
678: PetscArraycpy(points, tetPoints_p2, 3);
679: break;
680: }
681: case 1: {
682: PetscInt tetSizes_p2[2] = {3, 4};
683: PetscInt tetPoints_p2[7] = {0, 3, 5, 1, 2, 4, 6};
685: PetscMalloc2(2, &sizes, 7, &points);
686: PetscArraycpy(sizes, tetSizes_p2, 2);
687: PetscArraycpy(points, tetPoints_p2, 7);
688: break;
689: }
690: default:
691: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for tetrahedral hybrid mesh on 2 procs", user->testNum);
692: }
693: } else if (dim == 3 && !cellSimplex && !cellHybrid && size == 2) {
694: switch (user->testNum) {
695: case 0: {
696: PetscInt hexSizes_p2[2] = {1, 1};
697: PetscInt hexPoints_p2[2] = {0, 1};
699: PetscMalloc2(2, &sizes, 2, &points);
700: PetscArraycpy(sizes, hexSizes_p2, 2);
701: PetscArraycpy(points, hexPoints_p2, 2);
702: break;
703: }
704: default:
705: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral mesh on 2 procs", user->testNum);
706: }
707: } else if (dim == 3 && !cellSimplex && cellHybrid && size == 2) {
708: switch (user->testNum) {
709: case 0: {
710: PetscInt hexSizes_p2[2] = {1, 1};
711: PetscInt hexPoints_p2[2] = {0, 1};
713: PetscMalloc2(2, &sizes, 2, &points);
714: PetscArraycpy(sizes, hexSizes_p2, 2);
715: PetscArraycpy(points, hexPoints_p2, 2);
716: break;
717: }
718: case 1: {
719: PetscInt hexSizes_p2[2] = {5, 4};
720: PetscInt hexPoints_p2[9] = {3, 4, 5, 7, 8, 0, 1, 2, 6};
722: PetscMalloc2(2, &sizes, 9, &points);
723: PetscArraycpy(sizes, hexSizes_p2, 2);
724: PetscArraycpy(points, hexPoints_p2, 9);
725: break;
726: }
727: default:
728: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test number %" PetscInt_FMT " for hexahedral hybrid mesh on 2 procs", user->testNum);
729: }
730: } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Could not find matching test partition");
731: }
732: DMPlexGetPartitioner(*dm, &part);
733: PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);
734: PetscPartitionerShellSetPartition(part, size, sizes, points);
735: PetscFree2(sizes, points);
736: } else {
737: PetscPartitioner part;
739: DMPlexGetPartitioner(*dm, &part);
740: PetscPartitionerSetFromOptions(part);
741: }
742: {
743: DM pdm = NULL;
745: DMPlexDistribute(*dm, 0, NULL, &pdm);
746: if (pdm) {
747: DMViewFromOptions(pdm, NULL, "-dm_view");
748: DMDestroy(dm);
749: *dm = pdm;
750: }
751: }
752: DMPlexDistributeSetDefault(*dm, PETSC_FALSE);
753: DMViewFromOptions(*dm, NULL, "-dm_view_pre");
754: DMSetFromOptions(*dm);
755: if (user->uninterpolate || user->reinterpolate) {
756: DM udm = NULL;
758: DMPlexUninterpolate(*dm, &udm);
759: DMPlexCopyCoordinates(*dm, udm);
760: DMDestroy(dm);
761: *dm = udm;
762: }
763: if (user->reinterpolate) {
764: DM idm = NULL;
766: DMPlexInterpolate(*dm, &idm);
767: DMPlexCopyCoordinates(*dm, idm);
768: DMDestroy(dm);
769: *dm = idm;
770: }
771: DMPlexDistributeSetDefault(*dm, PETSC_FALSE);
772: PetscObjectSetName((PetscObject)*dm, "Hybrid Mesh");
773: DMViewFromOptions(*dm, NULL, "-dm_view");
774: PetscObjectSetOptionsPrefix((PetscObject)*dm, "hyb_");
775: DMSetFromOptions(*dm);
776: return 0;
777: }
779: int main(int argc, char **argv)
780: {
781: DM dm;
782: AppCtx user; /* user-defined work context */
785: PetscInitialize(&argc, &argv, NULL, help);
786: ProcessOptions(PETSC_COMM_WORLD, &user);
787: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
788: DMDestroy(&dm);
789: PetscFinalize();
790: return 0;
791: }
793: /*TEST
795: # 1D Simplex 29-31
796: testset:
797: args: -dim 1 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
798: test:
799: suffix: 29
800: test:
801: suffix: 30
802: args: -dm_refine 1
803: test:
804: suffix: 31
805: args: -dm_refine 5
807: # 2D Simplex 0-3
808: testset:
809: args: -dim 2 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
810: test:
811: suffix: 0
812: test:
813: suffix: 1
814: args: -dm_refine 1
815: test:
816: suffix: 2
817: nsize: 2
818: test:
819: suffix: 3
820: nsize: 2
821: args: -dm_refine 1
822: test:
823: suffix: 32
824: args: -dm_refine 1 -uninterpolate
825: test:
826: suffix: 33
827: nsize: 2
828: args: -dm_refine 1 -uninterpolate
829: test:
830: suffix: 34
831: nsize: 2
832: args: -dm_refine 3 -uninterpolate
834: # 2D Hybrid Simplex 4-7
835: testset:
836: args: -dim 2 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
837: test:
838: suffix: 4
839: test:
840: suffix: 5
841: args: -dm_refine 1
842: test:
843: suffix: 6
844: nsize: 2
845: test:
846: suffix: 7
847: nsize: 2
848: args: -dm_refine 1
849: test:
850: suffix: 24
851: args: -test_num 1 -dm_refine 1
853: # 2D Quad 12-13
854: testset:
855: args: -dim 2 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
856: test:
857: suffix: 12
858: args: -dm_refine 1
859: test:
860: suffix: 13
861: nsize: 2
862: args: -dm_refine 1
864: # 2D Hybrid Quad 27-28
865: testset:
866: args: -dim 2 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
867: test:
868: suffix: 27
869: args: -dm_refine 1
870: test:
871: suffix: 28
872: nsize: 2
873: args: -dm_refine 1
875: # 3D Simplex 8-11
876: testset:
877: args: -dim 3 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
878: test:
879: suffix: 8
880: args: -dm_refine 1
881: test:
882: suffix: 9
883: nsize: 2
884: args: -dm_refine 1
885: test:
886: suffix: 10
887: args: -test_num 1 -dm_refine 1
888: test:
889: suffix: 11
890: nsize: 2
891: args: -test_num 1 -dm_refine 1
892: test:
893: suffix: 25
894: args: -test_num 2 -dm_refine 2
896: # 3D Hybrid Simplex 16-19
897: testset:
898: args: -dim 3 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
899: test:
900: suffix: 16
901: args: -dm_refine 1
902: test:
903: suffix: 17
904: nsize: 2
905: args: -dm_refine 1
906: test:
907: suffix: 18
908: args: -test_num 1 -dm_refine 1
909: test:
910: suffix: 19
911: nsize: 2
912: args: -test_num 1 -dm_refine 1
914: # 3D Hex 14-15
915: testset:
916: args: -dim 3 -cell_simplex 0 -cell_hybrid 0 -hyb_dm_plex_check_all -dm_plex_check_all
917: test:
918: suffix: 14
919: args: -dm_refine 1
920: test:
921: suffix: 15
922: nsize: 2
923: args: -dm_refine 1
924: test:
925: suffix: 26
926: args: -test_num 1 -dm_refine 2
928: # 3D Hybrid Hex 20-23
929: testset:
930: args: -dim 3 -cell_simplex 0 -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
931: test:
932: suffix: 20
933: args: -dm_refine 1
934: test:
935: suffix: 21
936: nsize: 2
937: args: -dm_refine 1
938: test:
939: suffix: 22
940: args: -test_num 1 -dm_refine 1
941: test:
942: suffix: 23
943: nsize: 2
944: args: -test_num 1 -dm_refine 1
946: # Hybrid interpolation
947: # TODO Setup new tests (like -reinterpolate) that interpolate hybrid cells
948: testset:
949: nsize: 2
950: args: -test_partition 0 -petscpartitioner_type simple -dm_view -hyb_dm_plex_check_all -in_dm_plex_check_all -dm_plex_check_all
951: test:
952: suffix: hybint_2d_0
953: args: -dim 2 -dm_refine 2
954: test:
955: suffix: hybint_2d_1
956: args: -dim 2 -dm_refine 2 -test_num 1
957: test:
958: suffix: hybint_3d_0
959: args: -dim 3 -dm_refine 1
960: test:
961: suffix: hybint_3d_1
962: args: -dim 3 -dm_refine 1 -test_num 1
964: TEST*/