Actual source code: ex47.c
2: static char help[] = "Tests PetscViewerHDF5 VecView()/VecLoad() function.\n\n";
4: #include <petscviewer.h>
5: #include <petscviewerhdf5.h>
6: #include <petscvec.h>
8: int main(int argc, char **args)
9: {
10: Vec x, y;
11: PetscReal norm, dnorm;
12: PetscViewer H5viewer;
13: char filename[PETSC_MAX_PATH_LEN];
14: PetscBool flg;
16: PetscFunctionBeginUser;
17: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
18: PetscCall(PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg));
19: if (!flg) PetscCall(PetscStrncpy(filename, "x.h5", sizeof(filename)));
20: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
21: PetscCall(VecSetFromOptions(x));
22: PetscCall(VecSetSizes(x, 11, PETSC_DETERMINE));
23: PetscCall(VecSet(x, 22.3));
25: PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_WRITE, &H5viewer));
26: PetscCall(PetscViewerSetFromOptions(H5viewer));
28: /* Write the Vec without one extra dimension for BS */
29: PetscCall(PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_FALSE));
30: PetscCall(PetscObjectSetName((PetscObject)x, "noBsDim"));
31: PetscCall(VecView(x, H5viewer));
33: /* Write the Vec with one extra, 1-sized, dimension for BS */
34: PetscCall(PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_TRUE));
35: PetscCall(PetscObjectSetName((PetscObject)x, "bsDim"));
36: PetscCall(VecView(x, H5viewer));
38: PetscCall(PetscViewerDestroy(&H5viewer));
39: PetscCall(VecDuplicate(x, &y));
41: /* Create the HDF5 viewer for reading */
42: PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_READ, &H5viewer));
43: PetscCall(PetscViewerSetFromOptions(H5viewer));
45: /* Load the Vec without the BS dim and compare */
46: PetscCall(PetscObjectSetName((PetscObject)y, "noBsDim"));
47: PetscCall(VecLoad(y, H5viewer));
48: PetscCall(VecAXPY(y, -1.0, x));
49: PetscCall(VecNorm(y, NORM_2, &norm));
50: PetscCall(VecNorm(x, NORM_2, &dnorm));
51: PetscCheck(norm / dnorm <= 1.e-6, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Vec read in 'noBsDim' does not match vector written out %g", (double)(norm / dnorm));
53: /* Load the Vec with one extra, 1-sized, BS dim and compare */
54: PetscCall(PetscObjectSetName((PetscObject)y, "bsDim"));
55: PetscCall(VecLoad(y, H5viewer));
56: PetscCall(VecAXPY(y, -1.0, x));
57: PetscCall(VecNorm(y, NORM_2, &norm));
58: PetscCall(VecNorm(x, NORM_2, &dnorm));
59: PetscCheck(norm / dnorm <= 1.e-6, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Vec read in 'bsDim' does not match vector written out %g", (double)(norm / dnorm));
61: PetscCall(PetscViewerDestroy(&H5viewer));
62: PetscCall(VecDestroy(&y));
63: PetscCall(VecDestroy(&x));
64: PetscCall(PetscFinalize());
65: return 0;
66: }
68: /*TEST
70: build:
71: requires: hdf5
73: test:
74: requires: hdf5
76: test:
77: suffix: 2
78: nsize: 4
80: test:
81: suffix: 3
82: nsize: 4
83: args: -viewer_hdf5_base_dimension2
85: test:
86: suffix: 4
87: nsize: 4
88: args: -viewer_hdf5_sp_output
90: TEST*/