Actual source code: ex32.c

  1: static char help[] = "Tests for periodic mesh output\n\n";

  3: #include <petscdmplex.h>

  5: PetscErrorCode CheckMesh(DM dm)
  6: {
  7:   PetscReal detJ, J[9];
  8:   PetscReal vol;
  9:   PetscInt  dim, depth, cStart, cEnd, c;

 11:   DMGetDimension(dm, &dim);
 12:   DMPlexGetDepth(dm, &depth);
 13:   DMPlexGetHeightStratum(dm, 0, &cStart, &cEnd);
 14:   for (c = cStart; c < cEnd; ++c) {
 15:     DMPlexComputeCellGeometryFEM(dm, c, NULL, NULL, J, NULL, &detJ);
 17:     if (depth > 1) {
 18:       DMPlexComputeCellGeometryFVM(dm, c, &vol, NULL, NULL);
 20:     }
 21:   }
 22:   return 0;
 23: }

 25: PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm)
 26: {
 27:   DMCreate(comm, dm);
 28:   DMSetType(*dm, DMPLEX);
 29:   DMSetFromOptions(*dm);
 30:   DMViewFromOptions(*dm, NULL, "-dm_view");
 31:   return 0;
 32: }

 34: int main(int argc, char **argv)
 35: {
 36:   DM dm;

 39:   PetscInitialize(&argc, &argv, NULL, help);
 40:   CreateMesh(PETSC_COMM_WORLD, &dm);
 41:   CheckMesh(dm);
 42:   DMDestroy(&dm);
 43:   PetscFinalize();
 44:   return 0;
 45: }

 47: /*TEST

 49:   test:
 50:     suffix: 0
 51:     args: -dm_plex_simplex 0 -dm_plex_box_faces 3,1,0 -dm_plex_box_bd periodic,none -dm_view ::ascii_info_detail
 52:   test:
 53:     suffix: 1
 54:     nsize: 2
 55:     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
 56:   test:
 57:     suffix: 2
 58:     nsize: 2
 59:     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
 60:   test:
 61:     suffix: 3
 62:     nsize: 4
 63:     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
 64:   test:
 65:     suffix: 4
 66:     nsize: 2
 67:     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
 68:   test:
 69:     suffix: 5
 70:     nsize: 2
 71:     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
 72:   # This checks that the SF with extra root for periodic cut still checks
 73:   test:
 74:     suffix: 5_hdf5
 75:     requires: hdf5
 76:     nsize: 2
 77:     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
 78:   test:
 79:     suffix: 6
 80:     nsize: 4
 81:     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

 83: TEST*/