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