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