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