Actual source code: ex47.c
petsc-3.12.0 2019-09-29
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: {
11: Vec x,y;
12: PetscReal norm,dnorm;
13: PetscViewer H5viewer;
15: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
16: VecCreate(PETSC_COMM_WORLD,&x);
17: VecSetFromOptions(x);
18: VecSetSizes(x,11,PETSC_DETERMINE);
19: VecSet(x,22.3);
21: PetscViewerHDF5Open(PETSC_COMM_WORLD,"x.h5",FILE_MODE_WRITE,&H5viewer);
22: PetscViewerSetFromOptions(H5viewer);
24: /* Write the Vec without one extra dimension for BS */
25: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_FALSE);
26: PetscObjectSetName((PetscObject) x, "noBsDim");
27: VecView(x,H5viewer);
29: /* Write the Vec with one extra, 1-sized, dimension for BS */
30: PetscViewerHDF5SetBaseDimension2(H5viewer, PETSC_TRUE);
31: PetscObjectSetName((PetscObject) x, "bsDim");
32: VecView(x,H5viewer);
34: PetscViewerDestroy(&H5viewer);
35: VecDuplicate(x,&y);
37: /* Create the HDF5 viewer for reading */
38: PetscViewerHDF5Open(PETSC_COMM_WORLD,"x.h5",FILE_MODE_READ,&H5viewer);
39: PetscViewerSetFromOptions(H5viewer);
41: /* Load the Vec without the BS dim and compare */
42: PetscObjectSetName((PetscObject) y, "noBsDim");
43: VecLoad(y,H5viewer);
44: VecAXPY(y,-1.0,x);
45: VecNorm(y,NORM_2,&norm);
46: VecNorm(x,NORM_2,&dnorm);
47: 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));
49: /* Load the Vec with one extra, 1-sized, BS dim and compare */
50: PetscObjectSetName((PetscObject) y, "bsDim");
51: VecLoad(y,H5viewer);
52: VecAXPY(y,-1.0,x);
53: VecNorm(y,NORM_2,&norm);
54: VecNorm(x,NORM_2,&dnorm);
55: 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));
57: PetscViewerDestroy(&H5viewer);
58: VecDestroy(&y);
59: VecDestroy(&x);
60: PetscFinalize();
61: return ierr;
62: }
65: /*TEST
67: build:
68: requires: hdf5
70: test:
71: requires: hdf5
73: test:
74: suffix: 2
75: nsize: 4
77: test:
78: suffix: 3
79: nsize: 4
80: args: -viewer_hdf5_base_dimension2
82: test:
83: suffix: 4
84: nsize: 4
85: args: -viewer_hdf5_sp_output
87: TEST*/