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