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