Actual source code: ex47.c
petsc-3.7.5 2017-01-01
2: static char help[] = "Tests PetscViewerHDF5 VecView()/VecLoad() function.\n\n";
4: #include <petscviewer.h>
5: #include <petscviewerhdf5.h>
6: #include <petscvec.h>
10: int main(int argc,char **args)
11: {
13: Vec x,y;
14: PetscReal norm,dnorm;
15: PetscViewer H5viewer;
17: PetscInitialize(&argc,&args,(char*)0,help);
18: VecCreate(PETSC_COMM_WORLD,&x);
19: VecSetFromOptions(x);
20: VecSetSizes(x,11,PETSC_DETERMINE);
21: VecSet(x,22.3);
23: PetscViewerHDF5Open(PETSC_COMM_WORLD,"x.h5",FILE_MODE_WRITE,&H5viewer);
24: PetscViewerSetFromOptions(H5viewer);
26: /* Write the Vec without one extra dimension for BS */
27: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_FALSE);
28: PetscObjectSetName((PetscObject) x, "noBsDim");
29: VecView(x,H5viewer);
31: /* Write the Vec with one extra, 1-sized, dimension for BS */
32: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_TRUE);
33: PetscObjectSetName((PetscObject) x, "bsDim");
34: VecView(x,H5viewer);
36: PetscViewerDestroy(&H5viewer);
37: VecDuplicate(x,&y);
39: /* Create the HDF5 viewer for reading */
40: PetscViewerHDF5Open(PETSC_COMM_WORLD,"x.h5",FILE_MODE_READ,&H5viewer);
41: PetscViewerSetFromOptions(H5viewer);
43: /* Load the Vec without the BS dim and compare */
44: PetscObjectSetName((PetscObject) y, "noBsDim");
45: VecLoad(y,H5viewer);
46: VecAXPY(y,-1.0,x);
47: VecNorm(y,NORM_2,&norm);
48: VecNorm(x,NORM_2,&dnorm);
49: if (norm/dnorm > 1.e-6) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Vec read in 'noBsDim' does not match vector written out %g",(double)(norm/dnorm));
51: /* Load the Vec with one extra, 1-sized, BS dim and compare */
52: PetscObjectSetName((PetscObject) y, "bsDim");
53: VecLoad(y,H5viewer);
54: VecAXPY(y,-1.0,x);
55: VecNorm(y,NORM_2,&norm);
56: VecNorm(x,NORM_2,&dnorm);
57: if (norm/dnorm > 1.e-6) SETERRQ1(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Vec read in 'bsDim' does not match vector written out %g",(double)(norm/dnorm));
59: PetscViewerDestroy(&H5viewer);
60: VecDestroy(&y);
61: VecDestroy(&x);
62: PetscFinalize();
63: return 0;
64: }