Actual source code: ex7.c
1: static char help[] = "Tests for mesh interpolation\n\n";
3: #include <petscdmplex.h>
5: /* List of test meshes
7: Triangle
8: --------
9: Test 0:
10: Two triangles sharing a face
12: 4
13: / | \
14: / | \
15: / | \
16: 2 0 | 1 5
17: \ | /
18: \ | /
19: \ | /
20: 3
22: should become
24: 4
25: / | \
26: 8 | 9
27: / | \
28: 2 0 7 1 5
29: \ | /
30: 6 | 10
31: \ | /
32: 3
34: Tetrahedron
35: -----------
36: Test 0:
37: Two tets sharing a face
39: cell 5 _______ cell
40: 0 / | \ \ 1
41: / | \ \
42: / | \ \
43: 2----|----4-----6
44: \ | / /
45: \ | / /
46: \ | / /
47: 3-------
49: should become
51: cell 5 _______ cell
52: 0 / | \ \ 1
53: / | \ \
54: 17 | 18 13 22
55: / 8 19 10 \ \
56: / | \ \
57: 2---14-|------4--20--6
58: \ | / /
59: \ 9 | 7 / /
60: 16 | 15 11 21
61: \ | / /
62: \ | / /
63: 3-------
65: Quadrilateral
66: -------------
67: Test 0:
68: Two quads sharing a face
70: 5-------4-------7
71: | | |
72: | 0 | 1 |
73: | | |
74: 2-------3-------6
76: should become
78: 5--10---4--14---7
79: | | |
80: 11 0 9 1 13
81: | | |
82: 2---8---3--12---6
84: Test 1:
85: A quad and a triangle sharing a face
87: 5-------4
88: | | \
89: | 0 | \
90: | | 1 \
91: 2-------3----6
93: should become
95: 5---9---4
96: | | \
97: 10 0 8 12
98: | | 1 \
99: 2---7---3-11-6
101: Hexahedron
102: ----------
103: Test 0:
104: Two hexes sharing a face
106: cell 9-------------8-------------13 cell
107: 0 /| /| /| 1
108: / | 15 / | 21 / |
109: / | / | / |
110: 6-------------7-------------12 |
111: | | 18 | | 24 | |
112: | | | | | |
113: |19 | |17 | |23 |
114: | | 16 | | 22 | |
115: | 5---------|---4---------|---11
116: | / | / | /
117: | / 14 | / 20 | /
118: |/ |/ |/
119: 2-------------3-------------10
121: should become
123: cell 9-----31------8-----42------13 cell
124: 0 /| /| /| 1
125: 32 | 15 30| 21 41|
126: / | / | / |
127: 6-----29------7-----40------12 |
128: | | 17 | | 23 | |
129: | 35 | 36 | 44
130: |19 | |18 | |24 |
131: 34 | 16 33 | 22 43 |
132: | 5-----26--|---4-----37--|---11
133: | / | / | /
134: | 25 14 | 27 20 | 38
135: |/ |/ |/
136: 2-----28------3-----39------10
138: */
140: typedef struct {
141: PetscInt testNum; /* Indicates the mesh to create */
142: PetscInt dim; /* The topological mesh dimension */
143: PetscBool simplex; /* Use simplices or hexes */
144: PetscBool useGenerator; /* Construct mesh with a mesh generator */
145: } AppCtx;
147: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
148: {
149: options->testNum = 0;
150: options->dim = 2;
151: options->simplex = PETSC_TRUE;
152: options->useGenerator = PETSC_FALSE;
154: PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");
155: PetscOptionsBoundedInt("-testnum", "The mesh to create", "ex7.c", options->testNum, &options->testNum, NULL, 0);
156: PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex7.c", options->dim, &options->dim, NULL, 1, 3);
157: PetscOptionsBool("-simplex", "Use simplices if true, otherwise hexes", "ex7.c", options->simplex, &options->simplex, NULL);
158: PetscOptionsBool("-use_generator", "Use a mesh generator to build the mesh", "ex7.c", options->useGenerator, &options->useGenerator, NULL);
159: PetscOptionsEnd();
160: return 0;
161: }
163: PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM dm)
164: {
165: PetscInt depth = 1, testNum = 0, p;
166: PetscMPIInt rank;
168: MPI_Comm_rank(comm, &rank);
169: if (rank == 0) {
170: switch (testNum) {
171: case 0: {
172: PetscInt numPoints[2] = {4, 2};
173: PetscInt coneSize[6] = {3, 3, 0, 0, 0, 0};
174: PetscInt cones[6] = {2, 3, 4, 5, 4, 3};
175: PetscInt coneOrientations[6] = {0, 0, 0, 0, 0, 0};
176: PetscScalar vertexCoords[8] = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5};
177: PetscInt markerPoints[8] = {2, 1, 3, 1, 4, 1, 5, 1};
179: DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
180: for (p = 0; p < 4; ++p) DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]);
181: } break;
182: default:
183: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
184: }
185: } else {
186: PetscInt numPoints[2] = {0, 0};
188: DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);
189: DMCreateLabel(dm, "marker");
190: }
191: return 0;
192: }
194: PetscErrorCode CreateSimplex_3D(MPI_Comm comm, DM dm)
195: {
196: PetscInt depth = 1, testNum = 0, p;
197: PetscMPIInt rank;
199: MPI_Comm_rank(comm, &rank);
200: if (rank == 0) {
201: switch (testNum) {
202: case 0: {
203: PetscInt numPoints[2] = {5, 2};
204: PetscInt coneSize[23] = {4, 4, 0, 0, 0, 0, 0};
205: PetscInt cones[8] = {2, 4, 3, 5, 3, 4, 6, 5};
206: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
207: 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};
208: PetscInt markerPoints[8] = {2, 1, 3, 1, 4, 1, 5, 1};
210: DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
211: for (p = 0; p < 4; ++p) DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]);
212: } break;
213: default:
214: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
215: }
216: } else {
217: PetscInt numPoints[2] = {0, 0};
219: DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);
220: DMCreateLabel(dm, "marker");
221: }
222: return 0;
223: }
225: PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscInt testNum, DM dm)
226: {
227: PetscInt depth = 1, p;
228: PetscMPIInt rank;
230: MPI_Comm_rank(comm, &rank);
231: if (rank == 0) {
232: switch (testNum) {
233: case 0: {
234: PetscInt numPoints[2] = {6, 2};
235: PetscInt coneSize[8] = {4, 4, 0, 0, 0, 0, 0, 0};
236: PetscInt cones[8] = {2, 3, 4, 5, 3, 6, 7, 4};
237: PetscInt coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
238: PetscScalar vertexCoords[12] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0, 0.5, 1.0};
239: PetscInt markerPoints[12] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1};
241: DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
242: for (p = 0; p < 6; ++p) DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]);
243: } break;
244: case 1: {
245: PetscInt numPoints[2] = {5, 2};
246: PetscInt coneSize[7] = {4, 3, 0, 0, 0, 0, 0};
247: PetscInt cones[7] = {2, 3, 4, 5, 3, 6, 4};
248: PetscInt coneOrientations[7] = {0, 0, 0, 0, 0, 0, 0};
249: PetscScalar vertexCoords[10] = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0};
250: PetscInt markerPoints[10] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};
252: DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
253: for (p = 0; p < 5; ++p) DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]);
254: } break;
255: default:
256: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
257: }
258: } else {
259: PetscInt numPoints[2] = {0, 0};
261: DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);
262: DMCreateLabel(dm, "marker");
263: }
264: return 0;
265: }
267: PetscErrorCode CreateHex_3D(MPI_Comm comm, DM dm)
268: {
269: PetscInt depth = 1, testNum = 0, p;
270: PetscMPIInt rank;
272: MPI_Comm_rank(comm, &rank);
273: if (rank == 0) {
274: switch (testNum) {
275: case 0: {
276: PetscInt numPoints[2] = {12, 2};
277: PetscInt coneSize[14] = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
278: PetscInt cones[16] = {2, 5, 4, 3, 6, 7, 8, 9, 3, 4, 11, 10, 7, 12, 13, 8};
279: PetscInt coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
280: PetscScalar vertexCoords[36] = {-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, 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};
281: PetscInt markerPoints[16] = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1, 8, 1, 9, 1};
283: DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords);
284: for (p = 0; p < 8; ++p) DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]);
285: } break;
286: default:
287: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
288: }
289: } else {
290: PetscInt numPoints[2] = {0, 0};
292: DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL);
293: DMCreateLabel(dm, "marker");
294: }
295: return 0;
296: }
298: PetscErrorCode CreateMesh(MPI_Comm comm, PetscInt testNum, AppCtx *user, DM *dm)
299: {
300: PetscBool useGenerator = user->useGenerator;
302: DMCreate(comm, dm);
303: DMSetType(*dm, DMPLEX);
304: if (!useGenerator) {
305: PetscInt dim = user->dim;
306: PetscBool simplex = user->simplex;
308: DMSetDimension(*dm, dim);
309: switch (dim) {
310: case 2:
311: if (simplex) {
312: CreateSimplex_2D(comm, *dm);
313: } else {
314: CreateQuad_2D(comm, testNum, *dm);
315: }
316: break;
317: case 3:
318: if (simplex) {
319: CreateSimplex_3D(comm, *dm);
320: } else {
321: CreateHex_3D(comm, *dm);
322: }
323: break;
324: default:
325: SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, dim);
326: }
327: }
328: DMSetFromOptions(*dm);
329: PetscObjectSetName((PetscObject)*dm, "Interpolated Mesh");
330: DMViewFromOptions(*dm, NULL, "-dm_view");
331: return 0;
332: }
334: int main(int argc, char **argv)
335: {
336: DM dm;
337: AppCtx user; /* user-defined work context */
340: PetscInitialize(&argc, &argv, NULL, help);
341: ProcessOptions(PETSC_COMM_WORLD, &user);
342: CreateMesh(PETSC_COMM_WORLD, user.testNum, &user, &dm);
343: DMDestroy(&dm);
344: PetscFinalize();
345: return 0;
346: }
348: /*TEST
349: testset:
350: args: -dm_plex_interpolate_pre -dm_plex_check_all
352: # Two cell test meshes 0-7
353: test:
354: suffix: 0
355: args: -dim 2 -dm_view ascii::ascii_info_detail
356: test:
357: suffix: 1
358: nsize: 2
359: args: -petscpartitioner_type simple -dim 2 -dm_view ascii::ascii_info_detail
360: test:
361: suffix: 2
362: args: -dim 2 -simplex 0 -dm_view ascii::ascii_info_detail
363: test:
364: suffix: 3
365: nsize: 2
366: args: -petscpartitioner_type simple -dim 2 -simplex 0 -dm_view ascii::ascii_info_detail
367: test:
368: suffix: 4
369: args: -dim 3 -dm_view ascii::ascii_info_detail
370: test:
371: suffix: 5
372: nsize: 2
373: args: -petscpartitioner_type simple -dim 3 -dm_view ascii::ascii_info_detail
374: test:
375: suffix: 6
376: args: -dim 3 -simplex 0 -dm_view ascii::ascii_info_detail
377: test:
378: suffix: 7
379: nsize: 2
380: args: -petscpartitioner_type simple -dim 3 -simplex 0 -dm_view ascii::ascii_info_detail
381: # 2D Hybrid Mesh 8
382: test:
383: suffix: 8
384: args: -dim 2 -simplex 0 -testnum 1 -dm_view ascii::ascii_info_detail
386: testset:
387: args: -dm_plex_check_all
389: # Reference cells
390: test:
391: suffix: 12
392: args: -use_generator -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_plex_check_all
393: # TetGen meshes 9-10
394: test:
395: suffix: 9
396: requires: triangle
397: args: -use_generator -dm_view ascii::ascii_info_detail -dm_coord_space 0
398: test:
399: suffix: 10
400: requires: ctetgen
401: args: -use_generator -dm_plex_dim 3 -dm_view ascii::ascii_info_detail -dm_coord_space 0
402: # Cubit meshes 11
403: test:
404: suffix: 11
405: requires: exodusii
406: args: -use_generator -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -dm_view ascii::ascii_info_detail -dm_coord_space 0
408: TEST*/