1 static char help[] = "Demonstrates HDF5 vector input/ouput\n\n"; 2 3 #include <petscsys.h> 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 #include <petscviewerhdf5.h> 7 8 int main(int argc, char **argv) { 9 PetscViewer viewer; 10 DM da; 11 Vec global, local, global2; 12 PetscMPIInt rank; 13 PetscBool flg; 14 PetscInt ndof; 15 16 /* 17 Every PETSc routine should begin with the PetscInitialize() routine. 18 argc, argv - These command line arguments are taken to extract the options 19 supplied to PETSc and options supplied to MPI. 20 help - When PETSc executable is invoked with the option -help, 21 it prints the various options that can be applied at 22 runtime. The user can use the "help" variable place 23 additional help messages in this printout. 24 */ 25 PetscFunctionBeginUser; 26 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 27 /* Get number of DOF's from command line */ 28 PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "DMDA VecView/VecLoad example", ""); 29 { 30 ndof = 1; 31 PetscOptionsBoundedInt("-ndof", "Number of DOF's in DMDA", "", ndof, &ndof, NULL, 1); 32 } 33 PetscOptionsEnd(); 34 35 /* Create a DMDA and an associated vector */ 36 PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, 100, 90, PETSC_DECIDE, PETSC_DECIDE, ndof, 1, NULL, NULL, &da)); 37 PetscCall(DMSetFromOptions(da)); 38 PetscCall(DMSetUp(da)); 39 PetscCall(DMCreateGlobalVector(da, &global)); 40 PetscCall(DMCreateLocalVector(da, &local)); 41 PetscCall(VecSet(global, -1.0)); 42 PetscCall(DMGlobalToLocalBegin(da, global, INSERT_VALUES, local)); 43 PetscCall(DMGlobalToLocalEnd(da, global, INSERT_VALUES, local)); 44 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 45 PetscCall(VecScale(local, rank + 1)); 46 PetscCall(DMLocalToGlobalBegin(da, local, ADD_VALUES, global)); 47 PetscCall(DMLocalToGlobalEnd(da, local, ADD_VALUES, global)); 48 49 /* Create the HDF5 viewer for writing */ 50 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "hdf5output.h5", FILE_MODE_WRITE, &viewer)); 51 PetscCall(PetscViewerSetFromOptions(viewer)); 52 53 /* Write the Vec without one extra dimension for BS */ 54 PetscCall(PetscViewerHDF5SetBaseDimension2(viewer, PETSC_FALSE)); 55 PetscCall(PetscObjectSetName((PetscObject)global, "noBsDim")); 56 PetscCall(VecView(global, viewer)); 57 58 /* Write the Vec with one extra, 1-sized, dimension for BS */ 59 PetscCall(PetscViewerHDF5SetBaseDimension2(viewer, PETSC_TRUE)); 60 PetscCall(PetscObjectSetName((PetscObject)global, "bsDim")); 61 PetscCall(VecView(global, viewer)); 62 63 PetscCall(PetscViewerDestroy(&viewer)); 64 PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD)); 65 PetscCall(VecDuplicate(global, &global2)); 66 67 /* Create the HDF5 viewer for reading */ 68 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, "hdf5output.h5", FILE_MODE_READ, &viewer)); 69 PetscCall(PetscViewerSetFromOptions(viewer)); 70 71 /* Load the Vec without the BS dim and compare */ 72 PetscCall(PetscObjectSetName((PetscObject)global2, "noBsDim")); 73 PetscCall(VecLoad(global2, viewer)); 74 75 PetscCall(VecEqual(global, global2, &flg)); 76 if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error: Vectors are not equal\n")); 77 78 /* Load the Vec with one extra, 1-sized, BS dim and compare */ 79 PetscCall(PetscObjectSetName((PetscObject)global2, "bsDim")); 80 PetscCall(VecLoad(global2, viewer)); 81 PetscCall(PetscViewerDestroy(&viewer)); 82 83 PetscCall(VecEqual(global, global2, &flg)); 84 if (!flg) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Error: Vectors are not equal\n")); 85 86 /* clean up and exit */ 87 PetscCall(VecDestroy(&local)); 88 PetscCall(VecDestroy(&global)); 89 PetscCall(VecDestroy(&global2)); 90 PetscCall(DMDestroy(&da)); 91 92 PetscCall(PetscFinalize()); 93 return 0; 94 } 95 96 /*TEST 97 98 build: 99 requires: hdf5 100 101 test: 102 nsize: 4 103 104 test: 105 nsize: 4 106 suffix: 2 107 args: -ndof 2 108 output_file: output/ex9_1.out 109 110 TEST*/ 111