1 static char help[] = "Load and save the mesh and fields to HDF5 and ExodusII\n\n"; 2 3 #include <petscdmplex.h> 4 #include <petscviewerhdf5.h> 5 #include <petscsf.h> 6 7 typedef struct { 8 PetscBool compare; /* Compare the meshes using DMPlexEqual() */ 9 PetscBool distribute; /* Distribute the mesh */ 10 PetscBool interpolate; /* Generate intermediate mesh elements */ 11 char filename[PETSC_MAX_PATH_LEN]; /* Mesh filename */ 12 PetscViewerFormat format; /* Format to write and read */ 13 PetscBool second_write_read; /* Write and read for the 2nd time */ 14 PetscBool use_low_level_functions; /* Use low level functions for viewing and loading */ 15 } AppCtx; 16 17 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 18 { 19 PetscErrorCode ierr; 20 21 PetscFunctionBeginUser; 22 options->compare = PETSC_FALSE; 23 options->distribute = PETSC_TRUE; 24 options->interpolate = PETSC_FALSE; 25 options->filename[0] = '\0'; 26 options->format = PETSC_VIEWER_DEFAULT; 27 options->second_write_read = PETSC_FALSE; 28 options->use_low_level_functions = PETSC_FALSE; 29 30 ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 31 ierr = PetscOptionsBool("-compare", "Compare the meshes using DMPlexEqual()", "ex55.c", options->compare, &options->compare, NULL);CHKERRQ(ierr); 32 ierr = PetscOptionsBool("-distribute", "Distribute the mesh", "ex55.c", options->distribute, &options->distribute, NULL);CHKERRQ(ierr); 33 ierr = PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex55.c", options->interpolate, &options->interpolate, NULL);CHKERRQ(ierr); 34 ierr = PetscOptionsString("-filename", "The mesh file", "ex55.c", options->filename, options->filename, sizeof(options->filename), NULL);CHKERRQ(ierr); 35 ierr = PetscOptionsEnum("-format", "Format to write and read", "ex55.c", PetscViewerFormats, (PetscEnum)options->format, (PetscEnum*)&options->format, NULL);CHKERRQ(ierr); 36 ierr = PetscOptionsBool("-second_write_read", "Write and read for the 2nd time", "ex55.c", options->second_write_read, &options->second_write_read, NULL);CHKERRQ(ierr); 37 ierr = PetscOptionsBool("-use_low_level_functions", "Use low level functions for viewing and loading", "ex55.c", options->use_low_level_functions, &options->use_low_level_functions, NULL);CHKERRQ(ierr); 38 ierr = PetscOptionsEnd();CHKERRQ(ierr); 39 PetscFunctionReturn(0); 40 }; 41 42 static PetscErrorCode DMPlexWriteAndReadHDF5(DM dm, const char filename[], const char prefix[], AppCtx user, DM *dm_new) 43 { 44 DM dmnew; 45 PetscViewer v; 46 PetscErrorCode ierr; 47 48 PetscFunctionBeginUser; 49 ierr = PetscViewerHDF5Open(PetscObjectComm((PetscObject) dm), filename, FILE_MODE_WRITE, &v);CHKERRQ(ierr); 50 ierr = PetscViewerPushFormat(v, user.format);CHKERRQ(ierr); 51 if (user.use_low_level_functions) { 52 ierr = DMPlexTopologyView(dm, v);CHKERRQ(ierr); 53 ierr = DMPlexCoordinatesView(dm, v);CHKERRQ(ierr); 54 ierr = DMPlexLabelsView(dm, v);CHKERRQ(ierr); 55 } else { 56 ierr = DMView(dm, v);CHKERRQ(ierr); 57 } 58 59 ierr = PetscViewerFileSetMode(v, FILE_MODE_READ);CHKERRQ(ierr); 60 ierr = DMCreate(PETSC_COMM_WORLD, &dmnew);CHKERRQ(ierr); 61 ierr = DMSetType(dmnew, DMPLEX);CHKERRQ(ierr); 62 ierr = DMSetOptionsPrefix(dmnew, prefix);CHKERRQ(ierr); 63 if (user.use_low_level_functions) { 64 ierr = DMPlexTopologyLoad(dmnew, v, NULL);CHKERRQ(ierr); 65 ierr = DMPlexCoordinatesLoad(dmnew, v);CHKERRQ(ierr); 66 ierr = DMPlexLabelsLoad(dmnew, v);CHKERRQ(ierr); 67 } else { 68 ierr = DMLoad(dmnew, v);CHKERRQ(ierr); 69 } 70 ierr = PetscObjectSetName((PetscObject)dmnew,"Mesh_new");CHKERRQ(ierr); 71 72 ierr = PetscViewerPopFormat(v);CHKERRQ(ierr); 73 ierr = PetscViewerDestroy(&v);CHKERRQ(ierr); 74 *dm_new = dmnew; 75 PetscFunctionReturn(0); 76 } 77 78 int main(int argc, char **argv) 79 { 80 DM dm, dmnew; 81 PetscPartitioner part; 82 AppCtx user; 83 PetscBool flg; 84 PetscErrorCode ierr; 85 86 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 87 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 88 ierr = DMPlexCreateFromFile(PETSC_COMM_WORLD, user.filename, user.interpolate, &dm);CHKERRQ(ierr); 89 ierr = DMSetOptionsPrefix(dm,"orig_");CHKERRQ(ierr); 90 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 91 92 if (user.distribute) { 93 DM dmdist; 94 95 ierr = DMPlexGetPartitioner(dm, &part);CHKERRQ(ierr); 96 ierr = PetscPartitionerSetFromOptions(part);CHKERRQ(ierr); 97 ierr = DMPlexDistribute(dm, 0, NULL, &dmdist);CHKERRQ(ierr); 98 if (dmdist) { 99 ierr = DMDestroy(&dm);CHKERRQ(ierr); 100 dm = dmdist; 101 } 102 } 103 104 ierr = DMSetOptionsPrefix(dm,NULL);CHKERRQ(ierr); 105 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 106 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 107 108 ierr = DMPlexWriteAndReadHDF5(dm, "dmdist.h5", "new_", user, &dmnew);CHKERRQ(ierr); 109 110 if (user.second_write_read) { 111 ierr = DMDestroy(&dm);CHKERRQ(ierr); 112 dm = dmnew; 113 ierr = DMPlexWriteAndReadHDF5(dm, "dmdist.h5", "new_", user, &dmnew);CHKERRQ(ierr); 114 } 115 116 ierr = DMViewFromOptions(dmnew, NULL, "-dm_view");CHKERRQ(ierr); 117 /* TODO: Is it still true? */ 118 /* The NATIVE format for coordiante viewing is killing parallel output, since we have a local vector. Map it to global, and it will work. */ 119 120 /* This currently makes sense only for sequential meshes. */ 121 if (user.compare) { 122 ierr = DMPlexEqual(dmnew, dm, &flg);CHKERRQ(ierr); 123 if (flg) {ierr = PetscPrintf(PETSC_COMM_WORLD,"DMs equal\n");CHKERRQ(ierr);} 124 else {ierr = PetscPrintf(PETSC_COMM_WORLD,"DMs are not equal\n");CHKERRQ(ierr);} 125 } 126 127 ierr = DMDestroy(&dm);CHKERRQ(ierr); 128 ierr = DMDestroy(&dmnew);CHKERRQ(ierr); 129 ierr = PetscFinalize(); 130 return ierr; 131 } 132 133 /*TEST 134 build: 135 requires: hdf5 136 # Idempotence of saving/loading 137 # Have to replace Exodus file, which is creating uninterpolated edges 138 test: 139 suffix: 0 140 requires: exodusii broken 141 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/Rect-tri3.exo -dm_view ascii::ascii_info_detail 142 args: -format hdf5_petsc -compare 143 test: 144 suffix: 1 145 requires: exodusii parmetis !defined(PETSC_USE_64BIT_INDICES) broken 146 nsize: 2 147 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/Rect-tri3.exo -dm_view ascii::ascii_info_detail 148 args: -petscpartitioner_type parmetis 149 args: -format hdf5_petsc -new_dm_view ascii::ascii_info_detail 150 testset: 151 requires: exodusii 152 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo 153 args: -petscpartitioner_type simple 154 args: -dm_view ascii::ascii_info_detail 155 args: -new_dm_view ascii::ascii_info_detail 156 test: 157 suffix: 2 158 nsize: {{1 2 4 8}separate output} 159 args: -format {{default hdf5_petsc}separate output} 160 args: -interpolate {{0 1}separate output} 161 test: 162 suffix: 2a 163 nsize: {{1 2 4 8}separate output} 164 args: -format {{hdf5_xdmf hdf5_viz}separate output} 165 test: 166 suffix: 3 167 requires: exodusii 168 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo -compare 169 170 # Load HDF5 file in XDMF format in parallel, write, read dm1, write, read dm2, and compare dm1 and dm2 171 test: 172 suffix: 4 173 requires: !complex 174 nsize: {{1 2 3 4 8}} 175 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.h5 176 args: -dm_plex_create_from_hdf5_xdmf -distribute 0 -format hdf5_xdmf -second_write_read -compare 177 178 # Use low level functions, DMPlexTopologyView()/Load(), DMPlexCoordinatesView()/Load(), and DMPlexLabelsView()/Load() 179 # Output must be the same as ex55_2_nsize-2_format-hdf5_petsc_interpolate-0.out 180 test: 181 suffix: 5 182 requires: exodusii 183 nsize: 2 184 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/blockcylinder-50.exo 185 args: -petscpartitioner_type simple 186 args: -dm_view ascii::ascii_info_detail 187 args: -new_dm_view ascii::ascii_info_detail 188 args: -format hdf5_petsc -use_low_level_functions 189 190 testset: 191 # the same data and settings as dm_impls_plex_tests-ex18_9% 192 requires: hdf5 !complex datafilespath 193 nsize: {{1 2 4}} 194 args: -dm_plex_check_symmetry -dm_plex_check_skeleton -dm_plex_check_geometry 195 args: -filename ${DATAFILESPATH}/meshes/cube-hexahedra-refined.h5 -dm_plex_create_from_hdf5_xdmf -dm_plex_hdf5_topology_path /cells -dm_plex_hdf5_geometry_path /coordinates 196 args: -format hdf5_xdmf -second_write_read -compare 197 test: 198 suffix: 9_hdf5_seqload 199 args: -distribute -petscpartitioner_type simple 200 args: -interpolate {{0 1}} 201 args: -dm_plex_hdf5_force_sequential 202 test: 203 suffix: 9_hdf5_seqload_metis 204 requires: parmetis 205 args: -distribute -petscpartitioner_type parmetis 206 args: -interpolate 1 207 args: -dm_plex_hdf5_force_sequential 208 test: 209 suffix: 9_hdf5 210 args: -interpolate 1 -petscpartitioner_type simple 211 test: 212 suffix: 9_hdf5_repart 213 requires: parmetis 214 args: -distribute -petscpartitioner_type parmetis 215 args: -interpolate 1 216 test: 217 TODO: Parallel partitioning of uninterpolated meshes not supported 218 suffix: 9_hdf5_repart_ppu 219 requires: parmetis 220 args: -distribute -petscpartitioner_type parmetis 221 args: -interpolate 0 222 223 # reproduce PetscSFView() crash - fixed, left as regression test 224 test: 225 suffix: new_dm_view 226 requires: exodusii 227 nsize: 2 228 args: -filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/TwoQuads.exo -new_dm_view ascii:ex5_new_dm_view.log:ascii_info_detail 229 TEST*/ 230