1 #include <petscdmplex.h> 2 #include <petscviewerhdf5.h> 3 #include <petscsf.h> 4 5 static const char help[] = "Load and save the mesh to the native HDF5 format\n\n"; 6 static const char EX[] = "ex56.c"; 7 static const char LABEL_NAME[] = "BoundaryVertices"; 8 static const PetscInt LABEL_VALUE = 12345; 9 typedef struct { 10 MPI_Comm comm; 11 const char *meshname; /* Mesh name */ 12 PetscBool compare; /* Compare the meshes using DMPlexEqual() and DMCompareLabels() */ 13 PetscBool compare_labels; /* Compare labels in the meshes using DMCompareLabels() */ 14 PetscBool compare_boundary; /* Check label I/O via boundary vertex coordinates */ 15 PetscBool compare_pre_post; /* Compare labels loaded before distribution with those loaded after distribution */ 16 char outfile[PETSC_MAX_PATH_LEN]; /* Output file */ 17 PetscBool use_low_level_functions; /* Use low level functions for viewing and loading */ 18 //TODO This is meant as temporary option; can be removed once we have full parallel loading in place 19 PetscBool distribute_after_topo_load; /* Distribute topology right after DMPlexTopologyLoad(), if use_low_level_functions=true */ 20 PetscBool verbose; 21 } AppCtx; 22 23 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options) 24 { 25 PetscErrorCode ierr; 26 27 PetscFunctionBeginUser; 28 options->comm = comm; 29 options->compare = PETSC_FALSE; 30 options->compare_labels = PETSC_FALSE; 31 options->compare_boundary = PETSC_FALSE; 32 options->compare_pre_post = PETSC_FALSE; 33 options->outfile[0] = '\0'; 34 options->use_low_level_functions = PETSC_FALSE; 35 options->distribute_after_topo_load = PETSC_FALSE; 36 options->verbose = PETSC_FALSE; 37 38 ierr = PetscOptionsBegin(comm, "", "Meshing Problem Options", "DMPLEX");CHKERRQ(ierr); 39 ierr = PetscOptionsBool("-compare", "Compare the meshes using DMPlexEqual() and DMCompareLabels()", EX, options->compare, &options->compare, NULL);CHKERRQ(ierr); 40 ierr = PetscOptionsBool("-compare_labels", "Compare labels in the meshes using DMCompareLabels()", "ex55.c", options->compare_labels, &options->compare_labels, NULL);CHKERRQ(ierr); 41 ierr = PetscOptionsBool("-compare_boundary", "Check label I/O via boundary vertex coordinates", "ex55.c", options->compare_boundary, &options->compare_boundary, NULL);CHKERRQ(ierr); 42 ierr = PetscOptionsBool("-compare_pre_post", "Compare labels loaded before distribution with those loaded after distribution", "ex55.c", options->compare_pre_post, &options->compare_pre_post, NULL);CHKERRQ(ierr); 43 ierr = PetscOptionsString("-outfile", "Output mesh file", EX, options->outfile, options->outfile, sizeof(options->outfile), NULL);CHKERRQ(ierr); 44 ierr = PetscOptionsBool("-use_low_level_functions", "Use low level functions for viewing and loading", EX, options->use_low_level_functions, &options->use_low_level_functions, NULL);CHKERRQ(ierr); 45 ierr = PetscOptionsBool("-distribute_after_topo_load", "Distribute topology right after DMPlexTopologyLoad(), if use_low_level_functions=true", EX, options->distribute_after_topo_load, &options->distribute_after_topo_load, NULL);CHKERRQ(ierr); 46 ierr = PetscOptionsBool("-verbose", "Verbose printing", EX, options->verbose, &options->verbose, NULL);CHKERRQ(ierr); 47 ierr = PetscOptionsEnd();CHKERRQ(ierr); 48 PetscFunctionReturn(0); 49 }; 50 51 static PetscErrorCode CreateMesh(AppCtx *options, DM *newdm) 52 { 53 DM dm; 54 PetscErrorCode ierr; 55 56 PetscFunctionBeginUser; 57 ierr = DMCreate(options->comm, &dm);CHKERRQ(ierr); 58 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 59 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 60 ierr = PetscObjectGetName((PetscObject)dm, &options->meshname);CHKERRQ(ierr); 61 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 62 *newdm = dm; 63 PetscFunctionReturn(0); 64 } 65 66 static PetscErrorCode SaveMesh(AppCtx *options, DM dm) 67 { 68 PetscViewer v; 69 PetscErrorCode ierr; 70 71 PetscFunctionBeginUser; 72 ierr = PetscViewerHDF5Open(PetscObjectComm((PetscObject) dm), options->outfile, FILE_MODE_WRITE, &v);CHKERRQ(ierr); 73 if (options->use_low_level_functions) { 74 ierr = DMPlexTopologyView(dm, v);CHKERRQ(ierr); 75 ierr = DMPlexCoordinatesView(dm, v);CHKERRQ(ierr); 76 ierr = DMPlexLabelsView(dm, v);CHKERRQ(ierr); 77 } else { 78 ierr = DMView(dm, v);CHKERRQ(ierr); 79 } 80 ierr = PetscViewerDestroy(&v);CHKERRQ(ierr); 81 PetscFunctionReturn(0); 82 } 83 84 typedef enum {NONE=0, PRE_DIST=1, POST_DIST=2} AuxObjLoadMode; 85 86 static PetscErrorCode LoadMeshLowLevel(AppCtx *options, PetscViewer v, PetscBool explicitDistribute, AuxObjLoadMode mode, DM *newdm) 87 { 88 DM dm; 89 PetscSF sfXC; 90 PetscErrorCode ierr; 91 92 PetscFunctionBeginUser; 93 ierr = DMCreate(options->comm, &dm);CHKERRQ(ierr); 94 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 95 ierr = PetscObjectSetName((PetscObject) dm, options->meshname);CHKERRQ(ierr); 96 ierr = DMPlexTopologyLoad(dm, v, &sfXC);CHKERRQ(ierr); 97 if (mode == PRE_DIST) { 98 ierr = DMPlexCoordinatesLoad(dm, v, sfXC);CHKERRQ(ierr); 99 ierr = DMPlexLabelsLoad(dm, v, sfXC);CHKERRQ(ierr); 100 } 101 if (explicitDistribute) { 102 DM dmdist; 103 PetscSF sfXB = sfXC, sfBC; 104 105 ierr = DMPlexDistribute(dm, 0, &sfBC, &dmdist);CHKERRQ(ierr); 106 if (dmdist) { 107 const char *name; 108 109 ierr = PetscObjectGetName((PetscObject) dm, &name);CHKERRQ(ierr); 110 ierr = PetscObjectSetName((PetscObject) dmdist, name);CHKERRQ(ierr); 111 ierr = PetscSFCompose(sfXB, sfBC, &sfXC);CHKERRQ(ierr); 112 ierr = PetscSFDestroy(&sfXB);CHKERRQ(ierr); 113 ierr = PetscSFDestroy(&sfBC);CHKERRQ(ierr); 114 ierr = DMDestroy(&dm);CHKERRQ(ierr); 115 dm = dmdist; 116 } 117 } 118 if (mode == POST_DIST) { 119 ierr = DMPlexCoordinatesLoad(dm, v, sfXC);CHKERRQ(ierr); 120 ierr = DMPlexLabelsLoad(dm, v, sfXC);CHKERRQ(ierr); 121 } 122 ierr = PetscSFDestroy(&sfXC);CHKERRQ(ierr); 123 *newdm = dm; 124 PetscFunctionReturn(0); 125 } 126 127 static PetscErrorCode LoadMesh(AppCtx *options, DM *dmnew) 128 { 129 DM dm; 130 PetscViewer v; 131 PetscErrorCode ierr; 132 133 PetscFunctionBeginUser; 134 ierr = PetscViewerHDF5Open(options->comm, options->outfile, FILE_MODE_READ, &v);CHKERRQ(ierr); 135 if (options->use_low_level_functions) { 136 if (options->compare_pre_post) { 137 DM dm0; 138 139 ierr = LoadMeshLowLevel(options, v, PETSC_TRUE, PRE_DIST, &dm0);CHKERRQ(ierr); 140 ierr = LoadMeshLowLevel(options, v, PETSC_TRUE, POST_DIST, &dm);CHKERRQ(ierr); 141 ierr = DMCompareLabels(dm0, dm, NULL, NULL);CHKERRQ(ierr); 142 ierr = DMDestroy(&dm0);CHKERRQ(ierr); 143 } else { 144 ierr = LoadMeshLowLevel(options, v, options->distribute_after_topo_load, POST_DIST, &dm);CHKERRQ(ierr); 145 } 146 } else { 147 ierr = DMCreate(options->comm, &dm);CHKERRQ(ierr); 148 ierr = DMSetType(dm, DMPLEX);CHKERRQ(ierr); 149 ierr = PetscObjectSetName((PetscObject) dm, options->meshname);CHKERRQ(ierr); 150 ierr = DMLoad(dm, v);CHKERRQ(ierr); 151 } 152 ierr = PetscViewerDestroy(&v);CHKERRQ(ierr); 153 154 ierr = DMSetOptionsPrefix(dm, "load_");CHKERRQ(ierr); 155 ierr = DMSetFromOptions(dm);CHKERRQ(ierr); 156 ierr = DMViewFromOptions(dm, NULL, "-dm_view");CHKERRQ(ierr); 157 *dmnew = dm; 158 PetscFunctionReturn(0); 159 } 160 161 static PetscErrorCode CompareMeshes(AppCtx *options, DM dm0, DM dm1) 162 { 163 PetscBool flg; 164 PetscErrorCode ierr; 165 166 PetscFunctionBeginUser; 167 if (options->compare) { 168 ierr = DMPlexEqual(dm0, dm1, &flg);CHKERRQ(ierr); 169 PetscCheckFalse(!flg,options->comm, PETSC_ERR_ARG_INCOMP, "DMs are not equal"); 170 ierr = PetscPrintf(options->comm,"DMs equal\n");CHKERRQ(ierr); 171 } 172 if (options->compare_labels) { 173 ierr = DMCompareLabels(dm0, dm1, NULL, NULL);CHKERRQ(ierr); 174 ierr = PetscPrintf(options->comm,"DMLabels equal\n");CHKERRQ(ierr); 175 } 176 PetscFunctionReturn(0); 177 } 178 179 static PetscErrorCode MarkBoundaryVertices(DM dm, PetscInt value, DMLabel *label) 180 { 181 DMLabel l; 182 IS points; 183 const PetscInt *idx; 184 PetscInt i, n; 185 PetscErrorCode ierr; 186 187 PetscFunctionBeginUser; 188 ierr = DMLabelCreate(PetscObjectComm((PetscObject)dm), LABEL_NAME, &l);CHKERRQ(ierr); 189 ierr = DMPlexMarkBoundaryFaces(dm, value, l);CHKERRQ(ierr); 190 ierr = DMPlexLabelComplete(dm, l);CHKERRQ(ierr); 191 ierr = DMLabelGetStratumIS(l, value, &points);CHKERRQ(ierr); 192 193 ierr = ISGetLocalSize(points, &n);CHKERRQ(ierr); 194 ierr = ISGetIndices(points, &idx);CHKERRQ(ierr); 195 for (i=0; i<n; i++) { 196 const PetscInt p = idx[i]; 197 PetscInt d; 198 199 ierr = DMPlexGetPointDepth(dm, p, &d);CHKERRQ(ierr); 200 if (d != 0) { 201 ierr = DMLabelClearValue(l, p, value);CHKERRQ(ierr); 202 } 203 } 204 ierr = ISRestoreIndices(points, &idx);CHKERRQ(ierr); 205 ierr = ISDestroy(&points);CHKERRQ(ierr); 206 *label = l; 207 PetscFunctionReturn(0); 208 } 209 210 static PetscErrorCode VertexCoordinatesToAll(DM dm, IS vertices, Vec *allCoords) 211 { 212 Vec coords, allCoords_; 213 VecScatter sc; 214 MPI_Comm comm; 215 PetscErrorCode ierr; 216 217 PetscFunctionBeginUser; 218 ierr = PetscObjectGetComm((PetscObject)dm, &comm);CHKERRQ(ierr); 219 ierr = DMGetCoordinatesLocalSetUp(dm);CHKERRQ(ierr); 220 if (vertices) { 221 ierr = DMGetCoordinatesLocalTuple(dm, vertices, NULL, &coords);CHKERRQ(ierr); 222 } else { 223 ierr = VecCreateSeq(PETSC_COMM_SELF, 0, &coords);CHKERRQ(ierr); 224 } 225 { 226 PetscInt n; 227 Vec mpivec; 228 229 ierr = VecGetLocalSize(coords, &n);CHKERRQ(ierr); 230 ierr = VecCreateMPI(comm, n, PETSC_DECIDE, &mpivec);CHKERRQ(ierr); 231 ierr = VecCopy(coords, mpivec);CHKERRQ(ierr); 232 ierr = VecDestroy(&coords);CHKERRQ(ierr); 233 coords = mpivec; 234 } 235 236 ierr = VecScatterCreateToAll(coords, &sc, &allCoords_);CHKERRQ(ierr); 237 ierr = VecScatterBegin(sc,coords,allCoords_,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 238 ierr = VecScatterEnd(sc,coords,allCoords_,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 239 ierr = VecScatterDestroy(&sc);CHKERRQ(ierr); 240 ierr = VecDestroy(&coords);CHKERRQ(ierr); 241 *allCoords = allCoords_; 242 PetscFunctionReturn(0); 243 } 244 245 /* Compute boundary label, remember boundary vertices using coordinates, save and load label, check it is defined on the original boundary vertices */ 246 static PetscErrorCode DMAddBoundaryLabel_GetCoordinateRepresentation(DM dm, Vec *allCoords) 247 { 248 DMLabel label; 249 IS vertices; 250 PetscErrorCode ierr; 251 252 PetscFunctionBeginUser; 253 ierr = MarkBoundaryVertices(dm, LABEL_VALUE, &label);CHKERRQ(ierr); 254 ierr = DMLabelGetStratumIS(label, LABEL_VALUE, &vertices);CHKERRQ(ierr); 255 ierr = VertexCoordinatesToAll(dm, vertices, allCoords);CHKERRQ(ierr); 256 ierr = DMAddLabel(dm, label);CHKERRQ(ierr); 257 ierr = DMLabelDestroy(&label);CHKERRQ(ierr); 258 ierr = ISDestroy(&vertices);CHKERRQ(ierr); 259 PetscFunctionReturn(0); 260 } 261 262 static PetscErrorCode DMGetBoundaryLabel_CompareWithCoordinateRepresentation(AppCtx *user, DM dm, Vec allCoords) 263 { 264 DMLabel label; 265 IS pointsIS; 266 const PetscInt *points; 267 PetscInt i, n; 268 PetscBool fail = PETSC_FALSE; 269 MPI_Comm comm; 270 PetscMPIInt rank; 271 PetscErrorCode ierr; 272 273 PetscFunctionBeginUser; 274 ierr = PetscObjectGetComm((PetscObject)dm, &comm);CHKERRQ(ierr); 275 ierr = MPI_Comm_rank(comm, &rank);CHKERRMPI(ierr); 276 ierr = DMGetLabel(dm, LABEL_NAME, &label);CHKERRQ(ierr); 277 PetscCheckFalse(!label,PETSC_COMM_SELF, PETSC_ERR_PLIB, "Label \"%s\" was not loaded", LABEL_NAME); 278 { 279 PetscInt pStart, pEnd; 280 281 ierr = DMPlexGetChart(dm, &pStart, &pEnd);CHKERRQ(ierr); 282 ierr = DMLabelCreateIndex(label, pStart, pEnd);CHKERRQ(ierr); 283 } 284 ierr = DMPlexFindVertices(dm, allCoords, 0.0, &pointsIS);CHKERRQ(ierr); 285 ierr = ISGetIndices(pointsIS, &points);CHKERRQ(ierr); 286 ierr = ISGetLocalSize(pointsIS, &n);CHKERRQ(ierr); 287 if (user->verbose) {ierr = DMLabelView(label, PETSC_VIEWER_STDOUT_(comm));CHKERRQ(ierr);} 288 for (i=0; i<n; i++) { 289 const PetscInt p = points[i]; 290 PetscBool has; 291 PetscInt v; 292 293 if (p < 0) continue; 294 ierr = DMLabelHasPoint(label, p, &has);CHKERRQ(ierr); 295 if (!has) { 296 ierr = PetscSynchronizedFPrintf(comm, PETSC_STDERR, "[%d] Label does not have point %D\n", rank, p);CHKERRQ(ierr); 297 fail = PETSC_TRUE; 298 continue; 299 } 300 ierr = DMLabelGetValue(label, p, &v);CHKERRQ(ierr); 301 if (v != LABEL_VALUE) { 302 ierr = PetscSynchronizedFPrintf(comm, PETSC_STDERR, "Point %D has bad value %D", p, v);CHKERRQ(ierr); 303 fail = PETSC_TRUE; 304 continue; 305 } 306 if (user->verbose) {ierr = PetscSynchronizedPrintf(comm, "[%d] OK point %D\n", rank, p);CHKERRQ(ierr);} 307 } 308 ierr = PetscSynchronizedFlush(comm, PETSC_STDOUT);CHKERRQ(ierr); 309 ierr = PetscSynchronizedFlush(comm, PETSC_STDERR);CHKERRQ(ierr); 310 ierr = ISRestoreIndices(pointsIS, &points);CHKERRQ(ierr); 311 ierr = ISDestroy(&pointsIS);CHKERRQ(ierr); 312 ierr = MPI_Allreduce(MPI_IN_PLACE, &fail, 1, MPIU_BOOL, MPI_LOR, comm);CHKERRMPI(ierr); 313 PetscCheckFalse(fail,comm, PETSC_ERR_PLIB, "Label \"%s\" was not loaded correctly - see details above", LABEL_NAME); 314 PetscFunctionReturn(0); 315 } 316 317 int main(int argc, char **argv) 318 { 319 DM dm, dmnew; 320 AppCtx user; 321 Vec allCoords = NULL; 322 PetscErrorCode ierr; 323 324 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 325 ierr = ProcessOptions(PETSC_COMM_WORLD, &user);CHKERRQ(ierr); 326 ierr = CreateMesh(&user, &dm);CHKERRQ(ierr); 327 if (user.compare_boundary) { 328 ierr = DMAddBoundaryLabel_GetCoordinateRepresentation(dm, &allCoords);CHKERRQ(ierr); 329 } 330 ierr = SaveMesh(&user, dm);CHKERRQ(ierr); 331 ierr = LoadMesh(&user, &dmnew);CHKERRQ(ierr); 332 ierr = CompareMeshes(&user, dm, dmnew);CHKERRQ(ierr); 333 if (user.compare_boundary) { 334 ierr = DMGetBoundaryLabel_CompareWithCoordinateRepresentation(&user, dmnew, allCoords);CHKERRQ(ierr); 335 } 336 ierr = DMDestroy(&dm);CHKERRQ(ierr); 337 ierr = DMDestroy(&dmnew);CHKERRQ(ierr); 338 ierr = VecDestroy(&allCoords);CHKERRQ(ierr); 339 ierr = PetscFinalize(); 340 return ierr; 341 } 342 343 //TODO we can -compare once the new parallel topology format is in place 344 /*TEST 345 build: 346 requires: hdf5 347 348 # load old format, save in new format, reload, distribute 349 testset: 350 suffix: 1 351 requires: !complex datafilespath 352 args: -dm_plex_name plex 353 args: -dm_plex_check_all -dm_plex_view_hdf5_storage_version 2.0.0 354 args: -dm_plex_interpolate 355 args: -load_dm_plex_check_all 356 args: -use_low_level_functions {{0 1}} -compare_boundary 357 args: -outfile ex56_1.h5 358 nsize: {{1 3}} 359 test: 360 suffix: a 361 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/annulus-20.h5 362 test: 363 suffix: b 364 TODO: broken 365 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/barycentricallyrefinedcube.h5 366 test: 367 suffix: c 368 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/blockcylinder-50.h5 369 test: 370 suffix: d 371 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/cube-hexahedra-refined.h5 372 test: 373 suffix: e 374 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/hybrid_hexwedge.h5 375 test: 376 suffix: f 377 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/square.h5 378 379 # load old format, save in new format, reload topology, distribute, load geometry and labels 380 testset: 381 suffix: 2 382 requires: !complex datafilespath 383 args: -dm_plex_name plex 384 args: -dm_plex_check_all -dm_plex_view_hdf5_storage_version 2.0.0 385 args: -dm_plex_interpolate 386 args: -load_dm_plex_check_all 387 args: -use_low_level_functions -load_dm_distribute 0 -distribute_after_topo_load -compare_boundary 388 args: -outfile ex56_2.h5 389 nsize: 3 390 test: 391 suffix: a 392 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/annulus-20.h5 393 test: 394 suffix: b 395 TODO: broken 396 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/barycentricallyrefinedcube.h5 397 test: 398 suffix: c 399 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/blockcylinder-50.h5 400 test: 401 suffix: d 402 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/cube-hexahedra-refined.h5 403 test: 404 suffix: e 405 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/hybrid_hexwedge.h5 406 test: 407 suffix: f 408 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/square.h5 409 410 # load old format, save in new format, reload topology, distribute, load geometry and labels 411 testset: 412 suffix: 3 413 requires: !complex datafilespath 414 args: -dm_plex_name plex 415 args: -dm_plex_view_hdf5_storage_version 2.0.0 416 args: -dm_plex_interpolate -load_dm_distribute 0 417 args: -use_low_level_functions -compare_pre_post 418 args: -outfile ex56_3.h5 419 nsize: 3 420 test: 421 suffix: a 422 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/annulus-20.h5 423 test: 424 suffix: b 425 TODO: broken 426 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/barycentricallyrefinedcube.h5 427 test: 428 suffix: c 429 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/blockcylinder-50.h5 430 test: 431 suffix: d 432 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/cube-hexahedra-refined.h5 433 test: 434 suffix: e 435 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/hybrid_hexwedge.h5 436 test: 437 suffix: f 438 args: -dm_plex_filename ${DATAFILESPATH}/meshes/hdf5-petsc/petsc-v3.16.0/v1.0.0/square.h5 439 TEST*/ 440