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;
17: PetscInitialize(&argc, &args, (char *)0, help);
18: PetscOptionsGetString(NULL, NULL, "-filename", filename, sizeof(filename), &flg);
19: if (!flg) PetscStrcpy(filename, "x.h5");
20: VecCreate(PETSC_COMM_WORLD, &x);
21: VecSetFromOptions(x);
22: VecSetSizes(x, 11, PETSC_DETERMINE);
23: VecSet(x, 22.3);
25: PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_WRITE, &H5viewer);
26: PetscViewerSetFromOptions(H5viewer);
28: /* Write the Vec without one extra dimension for BS */
29: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_FALSE);
30: PetscObjectSetName((PetscObject)x, "noBsDim");
31: VecView(x, H5viewer);
33: /* Write the Vec with one extra, 1-sized, dimension for BS */
34: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_TRUE);
35: PetscObjectSetName((PetscObject)x, "bsDim");
36: VecView(x, H5viewer);
38: PetscViewerDestroy(&H5viewer);
39: VecDuplicate(x, &y);
41: /* Create the HDF5 viewer for reading */
42: PetscViewerHDF5Open(PETSC_COMM_WORLD, filename, FILE_MODE_READ, &H5viewer);
43: PetscViewerSetFromOptions(H5viewer);
45: /* Load the Vec without the BS dim and compare */
46: PetscObjectSetName((PetscObject)y, "noBsDim");
47: VecLoad(y, H5viewer);
48: VecAXPY(y, -1.0, x);
49: VecNorm(y, NORM_2, &norm);
50: VecNorm(x, NORM_2, &dnorm);
53: /* Load the Vec with one extra, 1-sized, BS dim and compare */
54: PetscObjectSetName((PetscObject)y, "bsDim");
55: VecLoad(y, H5viewer);
56: VecAXPY(y, -1.0, x);
57: VecNorm(y, NORM_2, &norm);
58: VecNorm(x, NORM_2, &dnorm);
61: PetscViewerDestroy(&H5viewer);
62: VecDestroy(&y);
63: VecDestroy(&x);
64: 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*/