1 static char help[] = "Evaluate the shape quality of a mesh\n\n"; 2 3 #include <petscdmplex.h> 4 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; 12 13 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 14 { 15 PetscErrorCode ierr; 16 17 PetscFunctionBeginUser; 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; 27 28 ierr = PetscOptionsBegin(comm, "", "Mesh Quality Evaluation Options", "DMPLEX");CHKERRQ(ierr); 29 ierr = PetscOptionsString("-filename", "The mesh file", "ex9.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr); 30 ierr = PetscOptionsBool("-report", "Output a mesh quality report", "ex9.c", options->report, &options->report, NULL);CHKERRQ(ierr); 31 ierr = PetscOptionsReal("-cond_limit", "Condition number limit for cell output", "ex9.c", options->condLimit, &options->condLimit, NULL);CHKERRQ(ierr); 32 ierr = PetscOptionsReal("-orth_qual_atol", "Absolute tolerance for Orthogonal Quality", "ex9.c", options->tol, &options->tol, NULL);CHKERRQ(ierr); 33 ierr = PetscOptionsBool("-interpolate", "Interpolate mesh", "ex9.c", options->interpolate, &options->interpolate, NULL);CHKERRQ(ierr); 34 ierr = PetscOptionsBool("-distribute", "Distribute mesh", "ex9.c", options->distribute, &options->distribute, NULL);CHKERRQ(ierr); 35 ierr = PetscOptionsBool("-simplex", "Create simplex mesh elements", "ex9.c", options->simplex, &options->simplex, NULL);CHKERRQ(ierr); 36 ierr = PetscOptionsInt("-overlap", "Number of overlap levels for dm", "ex9.c", options->overlap, &options->overlap, NULL);CHKERRQ(ierr); 37 ierr = PetscOptionsInt("-dim", "Dimension of mesh if generated", "ex9.c", options->dim, &options->dim, NULL);CHKERRQ(ierr); 38 ierr = PetscOptionsEnd();CHKERRQ(ierr); 39 PetscFunctionReturn(0); 40 } 41 42 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 43 { 44 PetscErrorCode ierr; 45 46 PetscFunctionBeginUser; 47 if (user->filename[0]) { 48 ierr = DMPlexCreateFromFile(comm, user->filename, PETSC_TRUE, dm);CHKERRQ(ierr); 49 } else { 50 ierr = DMPlexCreateBoxMesh(comm, user->dim, user->simplex, NULL, NULL, NULL, NULL, user->interpolate, dm);CHKERRQ(ierr); 51 ierr = PetscObjectSetName((PetscObject) *dm, "Mesh");CHKERRQ(ierr); 52 } 53 ierr = DMSetUp(*dm);CHKERRQ(ierr); 54 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 55 if (user->distribute) { 56 DM dmDist = NULL; 57 ierr = DMPlexDistribute(*dm, user->overlap, NULL, &dmDist);CHKERRQ(ierr); 58 if (dmDist) { 59 ierr = DMDestroy(dm);CHKERRQ(ierr); 60 *dm = dmDist; 61 } 62 } 63 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 64 PetscFunctionReturn(0); 65 } 66 67 int main(int argc, char **argv) 68 { 69 DM dm; 70 DMLabel OQLabel; 71 Vec OQ; 72 AppCtx ctx; 73 PetscErrorCode ierr; 74 75 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 76 ierr = ProcessOptions(PETSC_COMM_WORLD, &ctx);CHKERRQ(ierr); 77 ierr = CreateMesh(PETSC_COMM_WORLD, &ctx, &dm);CHKERRQ(ierr); 78 ierr = DMPlexCheckCellShape(dm, ctx.report, ctx.condLimit);CHKERRQ(ierr); 79 ierr = DMPlexComputeOrthogonalQuality(dm, NULL, ctx.tol, &OQ, &OQLabel);CHKERRQ(ierr); 80 ierr = VecDestroy(&OQ);CHKERRQ(ierr); 81 ierr = DMDestroy(&dm);CHKERRQ(ierr); 82 ierr = PetscFinalize(); 83 return ierr; 84 } 85 86 /*TEST 87 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 93 94 test: 95 suffix: 1 96 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -report 97 98 testset: 99 requires: exodusii 100 args: -dm_plex_orthogonal_quality_label_view -dm_plex_orthogonal_quality_vec_view 101 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 106 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 111 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 116 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*/ 122