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); 20*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsBool("-report", "Output a mesh quality report", "ex9.c", options->report, &options->report, NULL)); 21*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsReal("-cond_limit", "Condition number limit for cell output", "ex9.c", options->condLimit, &options->condLimit, NULL)); 22*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsReal("-orth_qual_atol", "Absolute tolerance for Orthogonal Quality", "ex9.c", options->tol, &options->tol, NULL)); 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 PetscFunctionBeginUser; 30*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMCreate(comm, dm)); 31*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetType(*dm, DMPLEX)); 32*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMSetFromOptions(*dm)); 33*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMViewFromOptions(*dm, NULL, "-dm_view")); 34c4762a1bSJed Brown PetscFunctionReturn(0); 35c4762a1bSJed Brown } 36c4762a1bSJed Brown 37c4762a1bSJed Brown int main(int argc, char **argv) 38c4762a1bSJed Brown { 39c4762a1bSJed Brown DM dm; 40f108dbd7SJacob Faibussowitsch DMLabel OQLabel; 41f108dbd7SJacob Faibussowitsch Vec OQ; 42c4762a1bSJed Brown AppCtx ctx; 43c4762a1bSJed Brown PetscErrorCode ierr; 44c4762a1bSJed Brown 45c4762a1bSJed Brown ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 46*5f80ce2aSJacob Faibussowitsch CHKERRQ(ProcessOptions(PETSC_COMM_WORLD, &ctx)); 47*5f80ce2aSJacob Faibussowitsch CHKERRQ(CreateMesh(PETSC_COMM_WORLD, &ctx, &dm)); 48*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMPlexCheckCellShape(dm, ctx.report, ctx.condLimit)); 49*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMPlexComputeOrthogonalQuality(dm, NULL, ctx.tol, &OQ, &OQLabel)); 50*5f80ce2aSJacob Faibussowitsch CHKERRQ(VecDestroy(&OQ)); 51*5f80ce2aSJacob Faibussowitsch CHKERRQ(DMDestroy(&dm)); 52c4762a1bSJed Brown ierr = PetscFinalize(); 53c4762a1bSJed Brown return ierr; 54c4762a1bSJed Brown } 55c4762a1bSJed Brown 56c4762a1bSJed Brown /*TEST 57c4762a1bSJed Brown 58c4762a1bSJed Brown test: 59c4762a1bSJed Brown suffix: 0 60c4762a1bSJed Brown requires: exodusii 61c4762a1bSJed Brown nsize: {{1 2}} 62e600fa54SMatthew G. Knepley args: -petscpartitioner_type simple -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -report 63c4762a1bSJed Brown 64c4762a1bSJed Brown test: 65c4762a1bSJed Brown suffix: 1 6630602db0SMatthew G. Knepley args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square.msh -report 67c4762a1bSJed Brown 68f108dbd7SJacob Faibussowitsch testset: 69f108dbd7SJacob Faibussowitsch args: -dm_plex_orthogonal_quality_label_view -dm_plex_orthogonal_quality_vec_view 70f108dbd7SJacob Faibussowitsch 71f108dbd7SJacob Faibussowitsch test: 72f108dbd7SJacob Faibussowitsch suffix: box_1 73f108dbd7SJacob Faibussowitsch nsize: 1 7430602db0SMatthew G. Knepley args: -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0 75f108dbd7SJacob Faibussowitsch 76f108dbd7SJacob Faibussowitsch test: 77f108dbd7SJacob Faibussowitsch suffix: box_2 78f108dbd7SJacob Faibussowitsch nsize: 2 79e600fa54SMatthew G. Knepley args: -petscpartitioner_type simple -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -orth_qual_atol 1.0 80f108dbd7SJacob Faibussowitsch 81f108dbd7SJacob Faibussowitsch test: 82f108dbd7SJacob Faibussowitsch suffix: mesh_1 83f108dbd7SJacob Faibussowitsch nsize: 1 846ed19f2fSJacob Faibussowitsch requires: exodusii 8530602db0SMatthew G. Knepley args: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95 86f108dbd7SJacob Faibussowitsch 87f108dbd7SJacob Faibussowitsch test: 88f108dbd7SJacob Faibussowitsch suffix: mesh_2 89f108dbd7SJacob Faibussowitsch nsize: 2 906ed19f2fSJacob Faibussowitsch requires: exodusii 91e600fa54SMatthew G. Knepley args: -petscpartitioner_type simple -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sevenside-quad-15.exo -orth_qual_atol 0.95 92c4762a1bSJed Brown TEST*/ 93