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