Actual source code: ex51.c

  1: static char help[] = "Tests save/load plex with distribution in HDF5.\n\n";

  3: #include <petscdmplex.h>
  4: #include <petscsf.h>
  5: #include <petsclayouthdf5.h>

  7: typedef struct {
  8:   char fname[PETSC_MAX_PATH_LEN]; /* Output mesh filename */
  9: } AppCtx;

 11: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
 12: {
 13:   PetscBool flg;

 15:   options->fname[0] = '\0';
 16:   PetscOptionsBegin(comm, "", "DMPlex View/Load Test Options", "DMPLEX");
 17:   PetscOptionsString("-fname", "The output mesh file", "ex51.c", options->fname, options->fname, sizeof(options->fname), &flg);
 18:   PetscOptionsEnd();
 19:   return 0;
 20: }

 22: int main(int argc, char **argv)
 23: {
 24:   const char        exampleDMPlexName[]       = "exampleDMPlex";
 25:   const char        exampleDistributionName[] = "exampleDistribution";
 26:   PetscViewerFormat format                    = PETSC_VIEWER_HDF5_PETSC;
 27:   AppCtx            user;

 29:   PetscInitialize(&argc, &argv, NULL, help);
 30:   ProcessOptions(PETSC_COMM_WORLD, &user);
 31:   /* Save */
 32:   {
 33:     DM          dm;
 34:     PetscViewer viewer;

 36:     PetscViewerHDF5Open(PETSC_COMM_WORLD, user.fname, FILE_MODE_WRITE, &viewer);
 37:     /* Save exampleDMPlex */
 38:     {
 39:       DM               pdm;
 40:       PetscPartitioner part;
 41:       const PetscInt   faces[2] = {6, 1};
 42:       PetscSF          sf;
 43:       PetscInt         overlap = 1;

 45:       DMPlexCreateBoxMesh(PETSC_COMM_WORLD, 2, PETSC_FALSE, faces, NULL, NULL, NULL, PETSC_TRUE, &dm);
 46:       DMPlexGetPartitioner(dm, &part);
 47:       PetscPartitionerSetFromOptions(part);
 48:       DMPlexDistribute(dm, overlap, &sf, &pdm);
 49:       if (pdm) {
 50:         DMDestroy(&dm);
 51:         dm = pdm;
 52:       }
 53:       PetscSFDestroy(&sf);
 54:       PetscObjectSetName((PetscObject)dm, exampleDMPlexName);
 55:       DMPlexDistributionSetName(dm, exampleDistributionName);
 56:       PetscViewerPushFormat(viewer, format);
 57:       DMPlexTopologyView(dm, viewer);
 58:       DMPlexLabelsView(dm, viewer);
 59:       PetscViewerPopFormat(viewer);
 60:     }
 61:     /* Save coordinates */
 62:     PetscViewerPushFormat(viewer, format);
 63:     DMPlexCoordinatesView(dm, viewer);
 64:     PetscViewerPopFormat(viewer);
 65:     PetscObjectSetName((PetscObject)dm, "Save: DM");
 66:     DMViewFromOptions(dm, NULL, "-dm_view");
 67:     DMDestroy(&dm);
 68:     PetscViewerDestroy(&viewer);
 69:   }
 70:   /* Load */
 71:   {
 72:     DM          dm;
 73:     PetscSF     sfXC;
 74:     PetscViewer viewer;

 76:     PetscViewerHDF5Open(PETSC_COMM_WORLD, user.fname, FILE_MODE_READ, &viewer);
 77:     /* Load exampleDMPlex */
 78:     DMCreate(PETSC_COMM_WORLD, &dm);
 79:     DMSetType(dm, DMPLEX);
 80:     PetscObjectSetName((PetscObject)dm, exampleDMPlexName);
 81:     DMPlexDistributionSetName(dm, exampleDistributionName);
 82:     /* sfXC: X -> C                         */
 83:     /* X: set of globalPointNumbers, [0, N) */
 84:     /* C: loaded in-memory plex             */
 85:     PetscViewerPushFormat(viewer, format);
 86:     DMPlexTopologyLoad(dm, viewer, &sfXC);
 87:     PetscViewerPopFormat(viewer);
 88:     /* Do not distribute (Already distributed just like the saved plex) */
 89:     /* Load labels */
 90:     DMPlexLabelsLoad(dm, viewer, sfXC);
 91:     /* Load coordinates */
 92:     PetscViewerPushFormat(viewer, format);
 93:     DMPlexCoordinatesLoad(dm, viewer, sfXC);
 94:     PetscViewerPopFormat(viewer);
 95:     DMSetFromOptions(dm);
 96:     /* Print the exact same plex as the saved one */
 97:     PetscObjectSetName((PetscObject)dm, "Load: DM");
 98:     DMViewFromOptions(dm, NULL, "-dm_view");
 99:     PetscSFDestroy(&sfXC);
100:     DMDestroy(&dm);
101:     PetscViewerDestroy(&viewer);
102:   }
103:   /* Finalize */
104:   PetscFinalize();
105:   return 0;
106: }

108: /*TEST

110:   build:
111:     requires: hdf5
112:     requires: !complex
113:   test:
114:     suffix: 0
115:     requires: parmetis
116:     nsize: {{1 2 4}separate output}
117:     args: -fname ex51_dump.h5 -dm_view ascii::ascii_info_detail
118:     args: -petscpartitioner_type parmetis
119:     args: -dm_plex_view_hdf5_storage_version 2.1.0

121: TEST*/