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;
25: PetscInitialize(&argc, &argv, (char*) 0, help);if (ierr) return ierr;
26: comm = PETSC_COMM_WORLD;
27: MPI_Comm_rank(comm, &rank);
28: PetscOptionsGetInt(NULL,NULL, "-n", &n, NULL);
29: PetscOptionsGetInt(NULL,NULL, "-n_timesteps", &n_timesteps, NULL);
30: if (n_timesteps < 0) SETERRQ(comm, PETSC_ERR_USER_INPUT, "-n_timesteps must be nonnegative");
32: /* create, initialize and write vectors */
33: PetscRandomCreate(comm, &rand);
34: PetscRandomSetFromOptions(rand);
35: PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_WRITE, &viewer);
37: VecCreate(comm, &x1);
38: PetscObjectSetName((PetscObject) x1, "x1");
39: VecSetSizes(x1, PETSC_DECIDE, n);
40: VecSetFromOptions(x1);
41: VecSetRandom(x1, rand);
42: VecView(x1, viewer);
44: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
45: VecCreate(comm, &x2);
46: PetscObjectSetName((PetscObject) x2, "x2");
47: VecSetSizes(x2, PETSC_DECIDE, n);
48: VecSetBlockSize(x2, 2);
49: VecSetFromOptions(x2);
50: VecSetRandom(x2, rand);
51: VecView(x2, viewer);
52: PetscViewerHDF5PopGroup(viewer);
54: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
55: PetscViewerHDF5PushTimestepping(viewer);
57: VecDuplicateVecs(x1, n_timesteps, &x3ts);
58: for (i=0; i<n_timesteps; i++) {
59: PetscObjectSetName((PetscObject) x3ts[i], "x3");
60: VecSetRandom(x3ts[i], rand);
61: VecView(x3ts[i], viewer);
62: PetscViewerHDF5IncrementTimestep(viewer);
63: }
65: PetscViewerHDF5PushGroup(viewer, "testBlockSize");
66: VecDuplicateVecs(x2, n_timesteps, &x4ts);
67: for (i=0; i<n_timesteps; i++) {
68: PetscObjectSetName((PetscObject) x4ts[i], "x4");
69: VecSetRandom(x4ts[i], rand);
70: PetscViewerHDF5SetTimestep(viewer, i);
71: VecView(x4ts[i], viewer);
72: }
73: PetscViewerHDF5PopGroup(viewer);
75: PetscViewerHDF5PopTimestepping(viewer);
76: PetscViewerHDF5PopGroup(viewer);
78: PetscViewerDestroy(&viewer);
79: PetscRandomDestroy(&rand);
81: /* read and compare */
82: PetscViewerHDF5Open(comm, "ex19.h5", FILE_MODE_READ, &viewer);
84: VecDuplicate(x1, &x1r);
85: PetscObjectSetName((PetscObject) x1r, "x1");
86: VecLoad(x1r, viewer);
87: VecEqual(x1, x1r, &equal);
88: if (!equal) {
89: VecView(x1, PETSC_VIEWER_STDOUT_WORLD);
90: VecView(x1r, PETSC_VIEWER_STDOUT_WORLD);
91: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x1 != x1r");
92: }
94: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
95: VecDuplicate(x2, &x2r);
96: PetscObjectSetName((PetscObject) x2r, "x2");
97: VecLoad(x2r, viewer);
98: VecEqual(x2, x2r, &equal);
99: if (!equal) {
100: VecView(x2, PETSC_VIEWER_STDOUT_WORLD);
101: VecView(x2r, PETSC_VIEWER_STDOUT_WORLD);
102: SETERRQ(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x2 != x2r");
103: }
104: PetscViewerHDF5PopGroup(viewer);
106: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
107: PetscViewerHDF5PushTimestepping(viewer);
109: VecDuplicate(x1, &x3r);
110: PetscObjectSetName((PetscObject) x3r, "x3");
111: for (i=0; i<n_timesteps; i++) {
112: PetscViewerHDF5SetTimestep(viewer, i);
113: VecLoad(x3r, viewer);
114: VecEqual(x3r, x3ts[i], &equal);
115: if (!equal) {
116: VecView(x3r, PETSC_VIEWER_STDOUT_WORLD);
117: VecView(x3ts[i], PETSC_VIEWER_STDOUT_WORLD);
118: SETERRQ1(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x3ts[%D] != x3r", i);
119: }
120: }
122: PetscViewerHDF5PushGroup(viewer, "testBlockSize");
123: VecDuplicate(x2, &x4r);
124: PetscObjectSetName((PetscObject) x4r, "x4");
125: PetscViewerHDF5SetTimestep(viewer, 0);
126: for (i=0; i<n_timesteps; i++) {
127: VecLoad(x4r, viewer);
128: VecEqual(x4r, x4ts[i], &equal);
129: if (!equal) {
130: VecView(x4r, PETSC_VIEWER_STDOUT_WORLD);
131: VecView(x4ts[i], PETSC_VIEWER_STDOUT_WORLD);
132: SETERRQ1(comm, PETSC_ERR_PLIB, "Error in HDF5 viewer: x4ts[%D] != x4r", i);
133: }
134: PetscViewerHDF5IncrementTimestep(viewer);
135: }
136: PetscViewerHDF5PopGroup(viewer);
138: PetscViewerHDF5PopTimestepping(viewer);
139: PetscViewerHDF5PopGroup(viewer);
141: /* cleanup */
142: PetscViewerDestroy(&viewer);
143: VecDestroy(&x1);
144: VecDestroy(&x2);
145: VecDestroyVecs(n_timesteps, &x3ts);
146: VecDestroyVecs(n_timesteps, &x4ts);
147: VecDestroy(&x1r);
148: VecDestroy(&x2r);
149: VecDestroy(&x3r);
150: VecDestroy(&x4r);
151: PetscFinalize();
152: return ierr;
153: }
155: /*TEST
157: build:
158: requires: hdf5
160: test:
161: nsize: {{1 2 3 4}}
163: TEST*/