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:   PetscFunctionBegin;
 16:   options->fname[0] = '\0';
 17:   PetscOptionsBegin(comm, "", "DMPlex View/Load Test Options", "DMPLEX");
 18:   PetscCall(PetscOptionsString("-fname", "The output mesh file", "ex51.c", options->fname, options->fname, sizeof(options->fname), &flg));
 19:   PetscOptionsEnd();
 20:   PetscFunctionReturn(PETSC_SUCCESS);
 21: }

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

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

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

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

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

109: /*TEST

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

122: TEST*/