1 static char help[] = "Tests for periodic mesh output\n\n"; 2 3 #include <petscdmplex.h> 4 5 PetscErrorCode CheckMesh(DM dm) 6 { 7 PetscReal detJ, J[9]; 8 PetscReal vol; 9 PetscInt dim, depth, cStart, cEnd, c; 10 11 PetscFunctionBegin; 12 PetscCall(DMGetDimension(dm, &dim)); 13 PetscCall(DMPlexGetDepth(dm, &depth)); 14 PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd)); 15 for (c = cStart; c < cEnd; ++c) { 16 PetscCall(DMPlexComputeCellGeometryFEM(dm, c, NULL, NULL, J, NULL, &detJ)); 17 PetscCheck(detJ > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mesh cell %" PetscInt_FMT " is inverted, |J| = %g", c, (double)detJ); 18 if (depth > 1) { 19 PetscCall(DMPlexComputeCellGeometryFVM(dm, c, &vol, NULL, NULL)); 20 PetscCheck(vol > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Mesh cell %" PetscInt_FMT " is inverted, vol = %g", c, (double)vol); 21 } 22 } 23 PetscFunctionReturn(PETSC_SUCCESS); 24 } 25 26 PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm) 27 { 28 PetscFunctionBegin; 29 PetscCall(DMCreate(comm, dm)); 30 PetscCall(DMSetType(*dm, DMPLEX)); 31 PetscCall(DMSetFromOptions(*dm)); 32 PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view")); 33 PetscFunctionReturn(PETSC_SUCCESS); 34 } 35 36 int main(int argc, char **argv) 37 { 38 DM dm; 39 40 PetscFunctionBeginUser; 41 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 42 PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm)); 43 PetscCall(CheckMesh(dm)); 44 PetscCall(DMDestroy(&dm)); 45 PetscCall(PetscFinalize()); 46 return 0; 47 } 48 49 /*TEST 50 51 test: 52 suffix: 0 53 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -dm_view ::ascii_info_detail 54 test: 55 suffix: 1 56 nsize: 2 57 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail 58 test: 59 suffix: 2 60 nsize: 2 61 args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail 62 test: 63 suffix: 3 64 nsize: 4 65 args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -petscpartitioner_type simple -dm_view ::ascii_info_detail 66 test: 67 suffix: 4 68 nsize: 2 69 args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail 70 test: 71 suffix: 5 72 nsize: 2 73 args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail 74 # This checks that the SF with extra root for periodic cut still checks 75 test: 76 suffix: 5_hdf5 77 requires: hdf5 78 nsize: 2 79 args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view hdf5:mesh.h5 80 test: 81 suffix: 6 82 nsize: 4 83 args: -dm_plex_simplex 0 -dm_plex_box_faces 6,2,0 -dm_plex_box_bd periodic,none -dm_plex_periodic_cut -petscpartitioner_type simple -dm_view ::ascii_info_detail 84 85 TEST*/ 86