1 static char help[] = "Tests save/load plex with distribution in HDF5.\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscsf.h> 5 #include <petsclayouthdf5.h> 6 7 typedef struct { 8 char fname[PETSC_MAX_PATH_LEN]; /* Output mesh filename */ 9 } AppCtx; 10 11 PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) { 12 PetscBool flg; 13 14 PetscFunctionBegin; 15 options->fname[0] = '\0'; 16 PetscOptionsBegin(comm, "", "DMPlex View/Load Test Options", "DMPLEX"); 17 PetscCall(PetscOptionsString("-fname", "The output mesh file", "ex51.c", options->fname, options->fname, sizeof(options->fname), &flg)); 18 PetscOptionsEnd(); 19 PetscFunctionReturn(0); 20 } 21 22 int main(int argc, char **argv) { 23 const char exampleDMPlexName[] = "exampleDMPlex"; 24 const char exampleDistributionName[] = "exampleDistribution"; 25 PetscViewerFormat format = PETSC_VIEWER_HDF5_PETSC; 26 AppCtx user; 27 28 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 29 PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user)); 30 /* Save */ 31 { 32 DM dm; 33 PetscViewer viewer; 34 35 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, user.fname, FILE_MODE_WRITE, &viewer)); 36 /* Save exampleDMPlex */ 37 { 38 DM pdm; 39 PetscPartitioner part; 40 const PetscInt faces[2] = {6, 1}; 41 PetscSF sf; 42 PetscInt overlap = 1; 43 44 PetscCall(DMPlexCreateBoxMesh(PETSC_COMM_WORLD, 2, PETSC_FALSE, faces, NULL, NULL, NULL, PETSC_TRUE, &dm)); 45 PetscCall(DMPlexGetPartitioner(dm, &part)); 46 PetscCall(PetscPartitionerSetFromOptions(part)); 47 PetscCall(DMPlexDistribute(dm, overlap, &sf, &pdm)); 48 if (pdm) { 49 PetscCall(DMDestroy(&dm)); 50 dm = pdm; 51 } 52 PetscCall(PetscSFDestroy(&sf)); 53 PetscCall(PetscObjectSetName((PetscObject)dm, exampleDMPlexName)); 54 PetscCall(DMPlexDistributionSetName(dm, exampleDistributionName)); 55 PetscCall(PetscViewerPushFormat(viewer, format)); 56 PetscCall(DMPlexTopologyView(dm, viewer)); 57 PetscCall(DMPlexLabelsView(dm, viewer)); 58 PetscCall(PetscViewerPopFormat(viewer)); 59 } 60 /* Save coordinates */ 61 PetscCall(PetscViewerPushFormat(viewer, format)); 62 PetscCall(DMPlexCoordinatesView(dm, viewer)); 63 PetscCall(PetscViewerPopFormat(viewer)); 64 PetscCall(PetscObjectSetName((PetscObject)dm, "Save: DM")); 65 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 66 PetscCall(DMDestroy(&dm)); 67 PetscCall(PetscViewerDestroy(&viewer)); 68 } 69 /* Load */ 70 { 71 DM dm; 72 PetscSF sfXC; 73 PetscViewer viewer; 74 75 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, user.fname, FILE_MODE_READ, &viewer)); 76 /* Load exampleDMPlex */ 77 PetscCall(DMCreate(PETSC_COMM_WORLD, &dm)); 78 PetscCall(DMSetType(dm, DMPLEX)); 79 PetscCall(PetscObjectSetName((PetscObject)dm, exampleDMPlexName)); 80 PetscCall(DMPlexDistributionSetName(dm, exampleDistributionName)); 81 /* sfXC: X -> C */ 82 /* X: set of globalPointNumbers, [0, N) */ 83 /* C: loaded in-memory plex */ 84 PetscCall(PetscViewerPushFormat(viewer, format)); 85 PetscCall(DMPlexTopologyLoad(dm, viewer, &sfXC)); 86 PetscCall(PetscViewerPopFormat(viewer)); 87 /* Do not distribute (Already distributed just like the saved plex) */ 88 /* Load labels */ 89 PetscCall(DMPlexLabelsLoad(dm, viewer, sfXC)); 90 /* Load coordinates */ 91 PetscCall(PetscViewerPushFormat(viewer, format)); 92 PetscCall(DMPlexCoordinatesLoad(dm, viewer, sfXC)); 93 PetscCall(PetscViewerPopFormat(viewer)); 94 PetscCall(DMSetFromOptions(dm)); 95 /* Print the exact same plex as the saved one */ 96 PetscCall(PetscObjectSetName((PetscObject)dm, "Load: DM")); 97 PetscCall(DMViewFromOptions(dm, NULL, "-dm_view")); 98 PetscCall(PetscSFDestroy(&sfXC)); 99 PetscCall(DMDestroy(&dm)); 100 PetscCall(PetscViewerDestroy(&viewer)); 101 } 102 /* Finalize */ 103 PetscCall(PetscFinalize()); 104 return 0; 105 } 106 107 /*TEST 108 109 build: 110 requires: hdf5 111 requires: !complex 112 test: 113 suffix: 0 114 requires: parmetis 115 nsize: {{1 2 4}separate output} 116 args: -fname ex51_dump.h5 -dm_view ascii::ascii_info_detail 117 args: -petscpartitioner_type parmetis 118 args: -dm_plex_view_hdf5_storage_version 2.0.0 119 120 TEST*/ 121