Actual source code: ex19.c
1: static char help[] = "Parallel HDF5 Vec Viewing.\n\n";
3: /*T
4: Concepts: vectors^viewing
5: Concepts: viewers^hdf5
6: Processors: n
7: T*/
9: #include <petscvec.h>
10: #include <petscviewerhdf5.h>
12: int main(int argc,char **argv)
13: {
14: Vec x1, x2, *x3ts, *x4ts;
15: Vec x1r, x2r, x3r, x4r;
16: PetscViewer viewer;
17: PetscRandom rand;
18: PetscMPIInt rank;
19: PetscInt i, n = 6, n_timesteps = 5;
20: PetscBool equal;
21: MPI_Comm comm;
23: PetscInitialize(&argc, &argv, (char*) 0, help);
24: comm = PETSC_COMM_WORLD;
25: MPI_Comm_rank(comm, &rank);
26: PetscOptionsGetInt(NULL,NULL, "-n", &n, NULL);
27: PetscOptionsGetInt(NULL,NULL, "-n_timesteps", &n_timesteps, NULL);
30: /* create, initialize and write vectors */
31: PetscRandomCreate(comm, &rand);
32: PetscRandomSetFromOptions(rand);
33: PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_WRITE, &viewer);
35: VecCreate(comm, &x1);
36: PetscObjectSetName((PetscObject) x1, "x1");
37: VecSetSizes(x1, PETSC_DECIDE, n);
38: VecSetFromOptions(x1);
39: VecSetRandom(x1, rand);
40: VecView(x1, viewer);
42: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
43: VecCreate(comm, &x2);
44: PetscObjectSetName((PetscObject) x2, "x2");
45: VecSetSizes(x2, PETSC_DECIDE, n);
46: VecSetBlockSize(x2, 2);
47: VecSetFromOptions(x2);
48: VecSetRandom(x2, rand);
49: VecView(x2, viewer);
50: PetscViewerHDF5PopGroup(viewer);
52: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
53: PetscViewerHDF5PushTimestepping(viewer);
55: VecDuplicateVecs(x1, n_timesteps, &x3ts);
56: for (i=0; i<n_timesteps; i++) {
57: PetscObjectSetName((PetscObject) x3ts[i], "x3");
58: VecSetRandom(x3ts[i], rand);
59: VecView(x3ts[i], viewer);
60: PetscViewerHDF5IncrementTimestep(viewer);
61: }
63: PetscViewerHDF5PushGroup(viewer, "testBlockSize");
64: VecDuplicateVecs(x2, n_timesteps, &x4ts);
65: for (i=0; i<n_timesteps; i++) {
66: PetscObjectSetName((PetscObject) x4ts[i], "x4");
67: VecSetRandom(x4ts[i], rand);
68: PetscViewerHDF5SetTimestep(viewer, i);
69: VecView(x4ts[i], viewer);
70: }
71: PetscViewerHDF5PopGroup(viewer);
73: PetscViewerHDF5PopTimestepping(viewer);
74: PetscViewerHDF5PopGroup(viewer);
76: PetscViewerDestroy(&viewer);
77: PetscRandomDestroy(&rand);
79: /* read and compare */
80: PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_READ, &viewer);
82: VecDuplicate(x1, &x1r);
83: PetscObjectSetName((PetscObject) x1r, "x1");
84: VecLoad(x1r, viewer);
85: VecEqual(x1, x1r, &equal);
86: if (!equal) {
87: VecView(x1, PETSC_VIEWER_STDOUT_WORLD);
88: VecView(x1r, PETSC_VIEWER_STDOUT_WORLD);
89: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x1 != x1r");
90: }
92: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
93: VecDuplicate(x2, &x2r);
94: PetscObjectSetName((PetscObject) x2r, "x2");
95: VecLoad(x2r, viewer);
96: VecEqual(x2, x2r, &equal);
97: if (!equal) {
98: VecView(x2, PETSC_VIEWER_STDOUT_WORLD);
99: VecView(x2r, PETSC_VIEWER_STDOUT_WORLD);
100: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x2 != x2r");
101: }
102: PetscViewerHDF5PopGroup(viewer);
104: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
105: PetscViewerHDF5PushTimestepping(viewer);
107: VecDuplicate(x1, &x3r);
108: PetscObjectSetName((PetscObject) x3r, "x3");
109: for (i=0; i<n_timesteps; i++) {
110: PetscViewerHDF5SetTimestep(viewer, i);
111: VecLoad(x3r, viewer);
112: VecEqual(x3r, x3ts[i], &equal);
113: if (!equal) {
114: VecView(x3r, PETSC_VIEWER_STDOUT_WORLD);
115: VecView(x3ts[i], PETSC_VIEWER_STDOUT_WORLD);
116: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x3ts[%" PetscInt_FMT "] != x3r", i);
117: }
118: }
120: PetscViewerHDF5PushGroup(viewer, "testBlockSize");
121: VecDuplicate(x2, &x4r);
122: PetscObjectSetName((PetscObject) x4r, "x4");
123: PetscViewerHDF5SetTimestep(viewer, 0);
124: for (i=0; i<n_timesteps; i++) {
125: VecLoad(x4r, viewer);
126: VecEqual(x4r, x4ts[i], &equal);
127: if (!equal) {
128: VecView(x4r, PETSC_VIEWER_STDOUT_WORLD);
129: VecView(x4ts[i], PETSC_VIEWER_STDOUT_WORLD);
130: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x4ts[%" PetscInt_FMT "] != x4r", i);
131: }
132: PetscViewerHDF5IncrementTimestep(viewer);
133: }
134: PetscViewerHDF5PopGroup(viewer);
136: PetscViewerHDF5PopTimestepping(viewer);
137: PetscViewerHDF5PopGroup(viewer);
139: /* cleanup */
140: PetscViewerDestroy(&viewer);
141: VecDestroy(&x1);
142: VecDestroy(&x2);
143: VecDestroyVecs(n_timesteps, &x3ts);
144: VecDestroyVecs(n_timesteps, &x4ts);
145: VecDestroy(&x1r);
146: VecDestroy(&x2r);
147: VecDestroy(&x3r);
148: VecDestroy(&x4r);
149: PetscFinalize();
150: return 0;
151: }
153: /*TEST
155: build:
156: requires: hdf5
158: test:
159: nsize: {{1 2 3 4}}
161: TEST*/