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