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:   PetscFunctionBegin;
150:   options->testNum      = 0;
151:   options->dim          = 2;
152:   options->simplex      = PETSC_TRUE;
153:   options->useGenerator = PETSC_FALSE;

155:   PetscOptionsBegin(comm, "", "Meshing Interpolation Test Options", "DMPLEX");
156:   PetscCall(PetscOptionsBoundedInt("-testnum", "The mesh to create", "ex7.c", options->testNum, &options->testNum, NULL, 0));
157:   PetscCall(PetscOptionsRangeInt("-dim", "The topological mesh dimension", "ex7.c", options->dim, &options->dim, NULL, 1, 3));
158:   PetscCall(PetscOptionsBool("-simplex", "Use simplices if true, otherwise hexes", "ex7.c", options->simplex, &options->simplex, NULL));
159:   PetscCall(PetscOptionsBool("-use_generator", "Use a mesh generator to build the mesh", "ex7.c", options->useGenerator, &options->useGenerator, NULL));
160:   PetscOptionsEnd();
161:   PetscFunctionReturn(PETSC_SUCCESS);
162: }

164: PetscErrorCode CreateSimplex_2D(MPI_Comm comm, DM dm)
165: {
166:   PetscInt    depth = 1, testNum = 0, p;
167:   PetscMPIInt rank;

169:   PetscFunctionBegin;
170:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
171:   if (rank == 0) {
172:     switch (testNum) {
173:     case 0: {
174:       PetscInt    numPoints[2]        = {4, 2};
175:       PetscInt    coneSize[6]         = {3, 3, 0, 0, 0, 0};
176:       PetscInt    cones[6]            = {2, 3, 4, 5, 4, 3};
177:       PetscInt    coneOrientations[6] = {0, 0, 0, 0, 0, 0};
178:       PetscScalar vertexCoords[8]     = {-0.5, 0.5, 0.0, 0.0, 0.0, 1.0, 0.5, 0.5};
179:       PetscInt    markerPoints[8]     = {2, 1, 3, 1, 4, 1, 5, 1};

181:       PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
182:       for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
183:     } break;
184:     default:
185:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
186:     }
187:   } else {
188:     PetscInt numPoints[2] = {0, 0};

190:     PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL));
191:     PetscCall(DMCreateLabel(dm, "marker"));
192:   }
193:   PetscFunctionReturn(PETSC_SUCCESS);
194: }

196: PetscErrorCode CreateSimplex_3D(MPI_Comm comm, DM dm)
197: {
198:   PetscInt    depth = 1, testNum = 0, p;
199:   PetscMPIInt rank;

201:   PetscFunctionBegin;
202:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
203:   if (rank == 0) {
204:     switch (testNum) {
205:     case 0: {
206:       PetscInt    numPoints[2]        = {5, 2};
207:       PetscInt    coneSize[23]        = {4, 4, 0, 0, 0, 0, 0};
208:       PetscInt    cones[8]            = {2, 4, 3, 5, 3, 4, 6, 5};
209:       PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
210:       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};
211:       PetscInt    markerPoints[8]     = {2, 1, 3, 1, 4, 1, 5, 1};

213:       PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
214:       for (p = 0; p < 4; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
215:     } break;
216:     default:
217:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
218:     }
219:   } else {
220:     PetscInt numPoints[2] = {0, 0};

222:     PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL));
223:     PetscCall(DMCreateLabel(dm, "marker"));
224:   }
225:   PetscFunctionReturn(PETSC_SUCCESS);
226: }

228: PetscErrorCode CreateQuad_2D(MPI_Comm comm, PetscInt testNum, DM dm)
229: {
230:   PetscInt    depth = 1, p;
231:   PetscMPIInt rank;

233:   PetscFunctionBegin;
234:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
235:   if (rank == 0) {
236:     switch (testNum) {
237:     case 0: {
238:       PetscInt    numPoints[2]        = {6, 2};
239:       PetscInt    coneSize[8]         = {4, 4, 0, 0, 0, 0, 0, 0};
240:       PetscInt    cones[8]            = {2, 3, 4, 5, 3, 6, 7, 4};
241:       PetscInt    coneOrientations[8] = {0, 0, 0, 0, 0, 0, 0, 0};
242:       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};
243:       PetscInt    markerPoints[12]    = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1};

245:       PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
246:       for (p = 0; p < 6; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
247:     } break;
248:     case 1: {
249:       PetscInt    numPoints[2]        = {5, 2};
250:       PetscInt    coneSize[7]         = {4, 3, 0, 0, 0, 0, 0};
251:       PetscInt    cones[7]            = {2, 3, 4, 5, 3, 6, 4};
252:       PetscInt    coneOrientations[7] = {0, 0, 0, 0, 0, 0, 0};
253:       PetscScalar vertexCoords[10]    = {-0.5, 0.0, 0.0, 0.0, 0.0, 1.0, -0.5, 1.0, 0.5, 0.0};
254:       PetscInt    markerPoints[10]    = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1};

256:       PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
257:       for (p = 0; p < 5; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
258:     } break;
259:     default:
260:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
261:     }
262:   } else {
263:     PetscInt numPoints[2] = {0, 0};

265:     PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL));
266:     PetscCall(DMCreateLabel(dm, "marker"));
267:   }
268:   PetscFunctionReturn(PETSC_SUCCESS);
269: }

271: PetscErrorCode CreateHex_3D(MPI_Comm comm, DM dm)
272: {
273:   PetscInt    depth = 1, testNum = 0, p;
274:   PetscMPIInt rank;

276:   PetscFunctionBegin;
277:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
278:   if (rank == 0) {
279:     switch (testNum) {
280:     case 0: {
281:       PetscInt    numPoints[2]         = {12, 2};
282:       PetscInt    coneSize[14]         = {8, 8, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
283:       PetscInt    cones[16]            = {2, 5, 4, 3, 6, 7, 8, 9, 3, 4, 11, 10, 7, 12, 13, 8};
284:       PetscInt    coneOrientations[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
285:       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};
286:       PetscInt    markerPoints[16]     = {2, 1, 3, 1, 4, 1, 5, 1, 6, 1, 7, 1, 8, 1, 9, 1};

288:       PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, coneSize, cones, coneOrientations, vertexCoords));
289:       for (p = 0; p < 8; ++p) PetscCall(DMSetLabelValue(dm, "marker", markerPoints[p * 2], markerPoints[p * 2 + 1]));
290:     } break;
291:     default:
292:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "No test mesh %" PetscInt_FMT, testNum);
293:     }
294:   } else {
295:     PetscInt numPoints[2] = {0, 0};

297:     PetscCall(DMPlexCreateFromDAG(dm, depth, numPoints, NULL, NULL, NULL, NULL));
298:     PetscCall(DMCreateLabel(dm, "marker"));
299:   }
300:   PetscFunctionReturn(PETSC_SUCCESS);
301: }

303: PetscErrorCode CreateMesh(MPI_Comm comm, PetscInt testNum, AppCtx *user, DM *dm)
304: {
305:   PetscBool useGenerator = user->useGenerator;

307:   PetscFunctionBegin;
308:   PetscCall(DMCreate(comm, dm));
309:   PetscCall(DMSetType(*dm, DMPLEX));
310:   if (!useGenerator) {
311:     PetscInt  dim     = user->dim;
312:     PetscBool simplex = user->simplex;

314:     PetscCall(DMSetDimension(*dm, dim));
315:     switch (dim) {
316:     case 2:
317:       if (simplex) {
318:         PetscCall(CreateSimplex_2D(comm, *dm));
319:       } else {
320:         PetscCall(CreateQuad_2D(comm, testNum, *dm));
321:       }
322:       break;
323:     case 3:
324:       if (simplex) {
325:         PetscCall(CreateSimplex_3D(comm, *dm));
326:       } else {
327:         PetscCall(CreateHex_3D(comm, *dm));
328:       }
329:       break;
330:     default:
331:       SETERRQ(comm, PETSC_ERR_ARG_OUTOFRANGE, "Cannot make meshes for dimension %" PetscInt_FMT, dim);
332:     }
333:   }
334:   PetscCall(DMSetFromOptions(*dm));
335:   PetscCall(PetscObjectSetName((PetscObject)*dm, "Interpolated Mesh"));
336:   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
337:   PetscFunctionReturn(PETSC_SUCCESS);
338: }

340: int main(int argc, char **argv)
341: {
342:   DM     dm;
343:   AppCtx user; /* user-defined work context */

345:   PetscFunctionBeginUser;
346:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
347:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
348:   PetscCall(CreateMesh(PETSC_COMM_WORLD, user.testNum, &user, &dm));
349:   PetscCall(DMDestroy(&dm));
350:   PetscCall(PetscFinalize());
351:   return 0;
352: }

354: /*TEST
355:   testset:
356:     args: -dm_plex_interpolate_pre -dm_plex_check_all

358:     # Two cell test meshes 0-7
359:     test:
360:       suffix: 0
361:       args: -dim 2 -dm_view ascii::ascii_info_detail
362:     test:
363:       suffix: 1
364:       nsize: 2
365:       args: -petscpartitioner_type simple -dim 2 -dm_view ascii::ascii_info_detail
366:     test:
367:       suffix: 2
368:       args: -dim 2 -simplex 0 -dm_view ascii::ascii_info_detail
369:     test:
370:       suffix: 3
371:       nsize: 2
372:       args: -petscpartitioner_type simple -dim 2 -simplex 0 -dm_view ascii::ascii_info_detail
373:     test:
374:       suffix: 4
375:       args: -dim 3 -dm_view ascii::ascii_info_detail
376:     test:
377:       suffix: 5
378:       nsize: 2
379:       args: -petscpartitioner_type simple -dim 3 -dm_view ascii::ascii_info_detail
380:     test:
381:       suffix: 6
382:       args: -dim 3 -simplex 0 -dm_view ascii::ascii_info_detail
383:     test:
384:       suffix: 7
385:       nsize: 2
386:       args: -petscpartitioner_type simple -dim 3 -simplex 0 -dm_view ascii::ascii_info_detail
387:     # 2D Hybrid Mesh 8
388:     test:
389:       suffix: 8
390:       args: -dim 2 -simplex 0 -testnum 1 -dm_view ascii::ascii_info_detail

392:   testset:
393:     args: -dm_plex_check_all

395:     # Reference cells
396:     test:
397:       suffix: 12
398:       args: -use_generator -dm_plex_reference_cell_domain -dm_plex_cell pyramid -dm_plex_check_all
399:     # TetGen meshes 9-10
400:     test:
401:       suffix: 9
402:       requires: triangle
403:       args: -use_generator -dm_view ascii::ascii_info_detail -dm_coord_space 0
404:     test:
405:       suffix: 10
406:       requires: ctetgen
407:       args: -use_generator -dm_plex_dim 3 -dm_view ascii::ascii_info_detail -dm_coord_space 0
408:     # Cubit meshes 11
409:     test:
410:       suffix: 11
411:       requires: exodusii
412:       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

414: TEST*/