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