Actual source code: ex9.c
petsc-master 2020-08-25
1: static char help[] = "Evaluate the shape quality of a mesh\n\n";
3: #include <petscdmplex.h>
5: typedef struct {
6: char filename[PETSC_MAX_PATH_LEN]; /* Import mesh from file */
7: PetscBool report; /* Print a quality report */
8: PetscReal condLimit, tol; /* Condition number limit for cell output */
9: PetscBool interpolate, distribute, simplex;
10: PetscInt dim, overlap;
11: } AppCtx;
13: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
14: {
18: options->filename[0] = '\0';
19: options->report = PETSC_FALSE;
20: options->interpolate = PETSC_FALSE;
21: options->distribute = PETSC_FALSE;
22: options->simplex = PETSC_FALSE;
23: options->overlap = 0;
24: options->dim = 2;
25: options->tol = 0.5;
26: options->condLimit = PETSC_DETERMINE;
28: PetscOptionsBegin(comm, "", "Mesh Quality Evaluation Options", "DMPLEX");
29: PetscOptionsString("-filename", "The mesh file", "ex9.c", options->filename, options->filename, sizeof(options->filename), NULL);
30: PetscOptionsBool("-report", "Output a mesh quality report", "ex9.c", options->report, &options->report, NULL);
31: PetscOptionsReal("-cond_limit", "Condition number limit for cell output", "ex9.c", options->condLimit, &options->condLimit, NULL);
32: PetscOptionsReal("-orth_qual_atol", "Absolute tolerance for Orthogonal Quality", "ex9.c", options->tol, &options->tol, NULL);
33: PetscOptionsBool("-interpolate", "Interpolate mesh", "ex9.c", options->interpolate, &options->interpolate, NULL);
34: PetscOptionsBool("-distribute", "Distribute mesh", "ex9.c", options->distribute, &options->distribute, NULL);
35: PetscOptionsBool("-simplex", "Create simplex mesh elements", "ex9.c", options->simplex, &options->simplex, NULL);
36: PetscOptionsInt("-overlap", "Number of overlap levels for dm", "ex9.c", options->overlap, &options->overlap, NULL);
37: PetscOptionsInt("-dim", "Dimension of mesh if generated", "ex9.c", options->dim, &options->dim, NULL);
38: PetscOptionsEnd();
39: return(0);
40: }
42: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
43: {
44: PetscErrorCode ierr;
47: if (user->filename[0]) {
48: DMPlexCreateFromFile(comm, user->filename, PETSC_TRUE, dm);
49: } else {
50: DMPlexCreateBoxMesh(comm, user->dim, user->simplex, NULL, NULL, NULL, NULL, user->interpolate, dm);
51: PetscObjectSetName((PetscObject) *dm, "Mesh");
52: }
53: DMSetUp(*dm);
54: DMSetFromOptions(*dm);
55: if (user->distribute) {
56: DM dmDist = NULL;
57: DMPlexDistribute(*dm, user->overlap, NULL, &dmDist);
58: if (dmDist) {
59: DMDestroy(dm);
60: *dm = dmDist;
61: }
62: }
63: DMViewFromOptions(*dm, NULL, "-dm_view");
64: return(0);
65: }
67: int main(int argc, char **argv)
68: {
69: DM dm;
70: DMLabel OQLabel;
71: Vec OQ;
72: AppCtx ctx;
75: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
76: ProcessOptions(PETSC_COMM_WORLD, &ctx);
77: CreateMesh(PETSC_COMM_WORLD, &ctx, &dm);
78: DMPlexCheckCellShape(dm, ctx.report, ctx.condLimit);
79: DMPlexComputeOrthogonalQuality(dm, NULL, ctx.tol, &OQ, &OQLabel);
80: VecDestroy(&OQ);
81: DMDestroy(&dm);
82: PetscFinalize();
83: return ierr;
84: }
86: /*TEST
88: test:
89: suffix: 0
90: requires: exodusii
91: nsize: {{1 2}}
92: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -report
94: test:
95: suffix: 1
96: args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -report
98: testset:
99: requires: exodusii
100: args: -dm_plex_orthogonal_quality_label_view -dm_plex_orthogonal_quality_vec_view
102: test:
103: suffix: box_1
104: nsize: 1
105: args: -interpolate -distribute -dim 2 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0
107: test:
108: suffix: box_2
109: nsize: 2
110: args: -interpolate -distribute -dim 2 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0
112: test:
113: suffix: mesh_1
114: nsize: 1
115: args: -interpolate -distribute -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95
117: test:
118: suffix: mesh_2
119: nsize: 2
120: args: -interpolate -distribute -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95
121: TEST*/