1 static char help[] = "Demonstrate HDF5/XDMF load-save-reload cycle\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscviewerhdf5.h> 5 #define EX "ex5.c" 6 7 typedef struct { 8 char infile[PETSC_MAX_PATH_LEN]; /* Input mesh filename */ 9 char outfile[PETSC_MAX_PATH_LEN]; /* Dump/reload mesh filename */ 10 PetscViewerFormat informat; /* Input mesh format */ 11 PetscViewerFormat outformat; /* Dump/reload mesh format */ 12 PetscBool redistribute; /* Redistribute the mesh */ 13 PetscInt ntimes; /* How many times do the cycle */ 14 } AppCtx; 15 16 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 17 { 18 PetscBool flg; 19 PetscErrorCode ierr; 20 21 PetscFunctionBeginUser; 22 options->infile[0] = '\0'; 23 options->outfile[0] = '\0'; 24 options->informat = PETSC_VIEWER_HDF5_XDMF; 25 options->outformat = PETSC_VIEWER_HDF5_XDMF; 26 options->redistribute = PETSC_TRUE; 27 options->ntimes = 2; 28 ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 29 ierr = PetscOptionsString("-infile", "The input mesh file", EX, options->infile, options->infile, sizeof(options->infile), &flg);CHKERRQ(ierr); 30 if (!flg) SETERRQ(comm, PETSC_ERR_USER_INPUT, "-infile needs to be specified"); 31 ierr = PetscOptionsString("-outfile", "The output mesh file (by default it's the same as infile)", EX, options->outfile, options->outfile, sizeof(options->outfile), &flg);CHKERRQ(ierr); 32 if (!flg) SETERRQ(comm, PETSC_ERR_USER_INPUT, "-outfile needs to be specified"); 33 ierr = PetscOptionsEnum("-informat", "Input mesh format", EX, PetscViewerFormats, (PetscEnum)options->informat, (PetscEnum*)&options->informat, NULL);CHKERRQ(ierr); 34 ierr = PetscOptionsEnum("-outformat", "Dump/reload mesh format", EX, PetscViewerFormats, (PetscEnum)options->outformat, (PetscEnum*)&options->outformat, NULL);CHKERRQ(ierr); 35 ierr = PetscOptionsBool("-redistribute", "Redistribute the mesh", EX, options->redistribute, &options->redistribute, NULL);CHKERRQ(ierr); 36 ierr = PetscOptionsInt("-ntimes", "How many times do the cycle", EX, options->ntimes, &options->ntimes, NULL);CHKERRQ(ierr); 37 ierr = PetscOptionsEnd(); 38 PetscFunctionReturn(0); 39 }; 40 41 //TODO test DMLabel I/O (not yet working for PETSC_VIEWER_HDF5_XDMF) 42 int main(int argc, char **argv) 43 { 44 AppCtx user; 45 PetscInt i; 46 PetscBool flg; 47 MPI_Comm comm; 48 PetscMPIInt size; 49 PetscErrorCode ierr; 50 const char *filename; 51 PetscViewerFormat format; 52 53 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 54 comm = PETSC_COMM_WORLD; 55 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 56 ierr = ProcessOptions(comm, &user);CHKERRQ(ierr); 57 58 /* Use infile for the initial load */ 59 filename = user.infile; 60 format = user.informat; 61 62 for (i=0; i<user.ntimes; i++) { 63 DM dm; 64 PetscPartitioner part; 65 PetscViewer v; 66 67 ierr = PetscPrintf(comm, "Begin cycle %D\n",i);CHKERRQ(ierr); 68 69 /* Load data from XDMF into dm in parallel */ 70 /* We could also use 71 ierr = DMPlexCreateFromFile(PETSC_COMM_WORLD, user.filename, PETSC_TRUE, &dm);CHKERRQ(ierr); 72 This currently support a few more formats than DMLoad(). 73 */ 74 ierr = PetscViewerHDF5Open(comm, filename, FILE_MODE_READ, &v);CHKERRQ(ierr); 75 ierr = PetscViewerPushFormat(v, format);CHKERRQ(ierr); 76 ierr = DMCreate(comm, &dm);CHKERRQ(ierr); 77 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 78 ierr = DMSetOptionsPrefix(dm,"loaded_");CHKERRQ(ierr); 79 ierr = DMLoad(dm, v);CHKERRQ(ierr); 80 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 81 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 82 ierr = PetscViewerPopFormat(v);CHKERRQ(ierr); 83 ierr = PetscViewerDestroy(&v);CHKERRQ(ierr); 84 85 /* Use outfile for all I/O except the very initial load */ 86 filename = user.outfile; 87 format = user.outformat; 88 89 /* We just test/demonstrate DM is indeed distributed - unneeded in the application code */ 90 ierr = DMPlexIsDistributed(dm, &flg);CHKERRQ(ierr); 91 ierr = PetscPrintf(comm, "Loaded mesh distributed? %s\n", PetscBools[flg]); 92 93 /* Interpolate */ 94 //TODO we want to be able to do this from options in DMSetFromOptions() probably 95 //TODO we want to be able to do this in-place 96 { 97 DM idm; 98 99 ierr = DMPlexInterpolate(dm, &idm);CHKERRQ(ierr); 100 ierr = DMDestroy(&dm);CHKERRQ(ierr); 101 dm = idm; 102 ierr = DMSetOptionsPrefix(dm,"interpolated_");CHKERRQ(ierr); 103 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 104 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 105 } 106 107 /* Redistribute */ 108 //TODO we want to be able to do this from options in DMSetFromOptions() probably 109 if (user.redistribute) { 110 DM dmdist; 111 112 ierr = DMPlexGetPartitioner(dm, &part);CHKERRQ(ierr); 113 ierr = PetscPartitionerSetFromOptions(part);CHKERRQ(ierr); 114 ierr = DMPlexDistribute(dm, 0, NULL, &dmdist);CHKERRQ(ierr); 115 //TODO we want to be able to do this in-place 116 if (dmdist) { 117 ierr = DMDestroy(&dm);CHKERRQ(ierr); 118 dm = dmdist; 119 ierr = DMSetOptionsPrefix(dm,"redistributed_");CHKERRQ(ierr); 120 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 121 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 122 } 123 } 124 125 /* Save redistributed dm to XDMF in parallel and destroy it */ 126 ierr = PetscViewerHDF5Open(comm, user.outfile, FILE_MODE_WRITE, &v);CHKERRQ(ierr); 127 ierr = PetscViewerPushFormat(v, format);CHKERRQ(ierr); 128 ierr = DMView(dm, v);CHKERRQ(ierr); 129 ierr = PetscViewerPopFormat(v);CHKERRQ(ierr); 130 ierr = PetscViewerDestroy(&v);CHKERRQ(ierr); 131 ierr = DMDestroy(&dm);CHKERRQ(ierr); 132 133 ierr = PetscPrintf(comm, "End cycle %D\n--------\n",i);CHKERRQ(ierr); 134 } 135 136 /* Final clean-up */ 137 ierr = PetscFinalize(); 138 return ierr; 139 } 140 141 /*TEST 142 build: 143 requires: hdf5 144 testset: 145 suffix: 0 146 requires: !complex 147 nsize: {{2 4}} 148 args: -infile ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.h5 -informat hdf5_xdmf 149 args: -outfile ex5_dump.h5 -outformat {{hdf5_xdmf hdf5_petsc}separate output} 150 args: -ntimes 3 151 test: 152 # this partitioner should not shuffle anything, it should yield the same partititioning as the XDMF reader - added just for testing 153 suffix: simple 154 args: -petscpartitioner_type simple 155 test: 156 suffix: parmetis 157 requires: parmetis 158 args: -petscpartitioner_type parmetis 159 test: 160 suffix: ptscotch 161 requires: ptscotch 162 args: -petscpartitioner_type ptscotch 163 TEST*/ 164