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 } AppCtx; 10 11 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 12 { 13 PetscErrorCode ierr; 14 15 PetscFunctionBeginUser; 16 options->filename[0] = '\0'; 17 options->report = PETSC_FALSE; 18 options->tol = 0.5; 19 options->condLimit = PETSC_DETERMINE; 20 21 ierr = PetscOptionsBegin(comm, "", "Mesh Quality Evaluation Options", "DMPLEX");CHKERRQ(ierr); 22 ierr = PetscOptionsString("-filename", "The mesh file", "ex9.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr); 23 ierr = PetscOptionsBool("-report", "Output a mesh quality report", "ex9.c", options->report, &options->report, NULL);CHKERRQ(ierr); 24 ierr = PetscOptionsReal("-cond_limit", "Condition number limit for cell output", "ex9.c", options->condLimit, &options->condLimit, NULL);CHKERRQ(ierr); 25 ierr = PetscOptionsReal("-orth_qual_atol", "Absolute tolerance for Orthogonal Quality", "ex9.c", options->tol, &options->tol, NULL);CHKERRQ(ierr); 26 ierr = PetscOptionsEnd();CHKERRQ(ierr); 27 PetscFunctionReturn(0); 28 } 29 30 static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 31 { 32 PetscErrorCode ierr; 33 34 PetscFunctionBeginUser; 35 if (user->filename[0]) {ierr = DMPlexCreateFromFile(comm, user->filename, PETSC_TRUE, dm);CHKERRQ(ierr);} 36 else {ierr = DMPlexCreateBoxMesh(comm, 2, PETSC_TRUE, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);CHKERRQ(ierr);} 37 ierr = PetscObjectSetName((PetscObject) *dm, "Mesh");CHKERRQ(ierr); 38 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 39 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 40 PetscFunctionReturn(0); 41 } 42 43 int main(int argc, char **argv) 44 { 45 DM dm; 46 DMLabel OQLabel; 47 Vec OQ; 48 AppCtx ctx; 49 PetscErrorCode ierr; 50 51 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 52 ierr = ProcessOptions(PETSC_COMM_WORLD, &ctx);CHKERRQ(ierr); 53 ierr = CreateMesh(PETSC_COMM_WORLD, &ctx, &dm);CHKERRQ(ierr); 54 ierr = DMPlexCheckCellShape(dm, ctx.report, ctx.condLimit);CHKERRQ(ierr); 55 ierr = DMPlexComputeOrthogonalQuality(dm, NULL, ctx.tol, &OQ, &OQLabel);CHKERRQ(ierr); 56 ierr = VecDestroy(&OQ);CHKERRQ(ierr); 57 ierr = DMDestroy(&dm);CHKERRQ(ierr); 58 ierr = PetscFinalize(); 59 return ierr; 60 } 61 62 /*TEST 63 64 test: 65 suffix: 0 66 requires: exodusii 67 nsize: {{1 2}} 68 args: -dm_distribute -petscpartitioner_type simple -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -report 69 70 test: 71 suffix: 1 72 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -report 73 74 testset: 75 requires: exodusii 76 args: -dm_plex_orthogonal_quality_label_view -dm_plex_orthogonal_quality_vec_view 77 78 test: 79 suffix: box_1 80 nsize: 1 81 args: -dm_plex_box_simplex 0 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0 82 83 test: 84 suffix: box_2 85 nsize: 2 86 args: -dm_distribute -petscpartitioner_type simple -dm_plex_box_simplex 0 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0 87 88 test: 89 suffix: mesh_1 90 nsize: 1 91 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95 92 93 test: 94 suffix: mesh_2 95 nsize: 2 96 args: -dm_distribute -petscpartitioner_type simple -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95 97 TEST*/ 98