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>
13: int main(int argc,char **argv)
14: {
15: Vec x1, x2, y1, y2, y3, y4;
16: PetscViewer viewer;
17: PetscMPIInt rank;
18: PetscInt i, nlocal, n = 6;
19: PetscScalar *array;
20: PetscBool equal;
24: PetscInitialize(&argc, &argv, (char *) 0, help);
25: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
26: PetscOptionsGetInt(PETSC_NULL, "-n", &n, PETSC_NULL);
28: VecCreate(PETSC_COMM_WORLD, &x1);
29: PetscObjectSetName((PetscObject) x1, "TestVec");
30: VecSetSizes(x1, PETSC_DECIDE, n);
31: VecSetFromOptions(x1);
33: VecGetLocalSize(x1, &nlocal);
34: VecGetArray(x1, &array);
35: for(i = 0; i < nlocal; i++) {
36: array[i] = rank + 1;
37: }
38: VecRestoreArray(x1, &array);
39: VecAssemblyBegin(x1);
40: VecAssemblyEnd(x1);
42: VecDuplicate(x1, &x2);
43: PetscObjectSetName((PetscObject) x2, "TestVec");
44: VecCopy(x1, x2);
45: VecSetBlockSize(x2, 2);
47: PetscViewerHDF5Open(PETSC_COMM_WORLD, "ex19.h5", FILE_MODE_WRITE, &viewer);
48: PetscViewerHDF5PushGroup(viewer, "/");
49: VecView(x1, viewer);
50: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
51: VecView(x2, viewer);
52: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
53: PetscViewerHDF5SetTimestep(viewer, 0);
54: VecView(x2, viewer);
55: PetscViewerHDF5SetTimestep(viewer, 1);
56: VecView(x2, viewer);
57: PetscViewerHDF5PopGroup(viewer);
58: PetscViewerHDF5PopGroup(viewer);
59: PetscViewerHDF5PopGroup(viewer);
60: PetscViewerDestroy(&viewer);
62: VecCreate(PETSC_COMM_WORLD, &y1);
63: PetscObjectSetName((PetscObject) y1, "TestVec");
64: VecSetSizes(y1, PETSC_DECIDE, n);
65: VecSetFromOptions(y1);
66: VecDuplicate(y1, &y2);
67: PetscObjectSetName((PetscObject) y2, "TestVec");
68: VecDuplicate(y1, &y3);
69: PetscObjectSetName((PetscObject) y3, "TestVec");
70: VecDuplicate(y1, &y4);
71: PetscObjectSetName((PetscObject) y4, "TestVec");
73: PetscViewerHDF5Open(PETSC_COMM_WORLD, "ex19.h5", FILE_MODE_READ, &viewer);
74: PetscViewerHDF5PushGroup(viewer, "/");
75: VecLoad(y1, viewer);
76: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
77: VecLoad(y2, viewer);
78: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
79: PetscViewerHDF5SetTimestep(viewer, 0);
80: VecLoad(y3, viewer);
81: PetscViewerHDF5SetTimestep(viewer, 1);
82: VecLoad(y4, viewer);
83: PetscViewerHDF5PopGroup(viewer);
84: PetscViewerHDF5PopGroup(viewer);
85: PetscViewerHDF5PopGroup(viewer);
86: PetscViewerDestroy(&viewer);
88: VecEqual(x1, y1, &equal);
89: if (!equal) {
90: VecView(x1, PETSC_VIEWER_STDOUT_WORLD);
91: VecView(y1, PETSC_VIEWER_STDOUT_WORLD);
92: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB, "Error in HDF5 viewer");
93: }
94: VecEqual(x2, y2, &equal);
95: if (!equal) {
96: VecView(x2, PETSC_VIEWER_STDOUT_WORLD);
97: VecView(y2, PETSC_VIEWER_STDOUT_WORLD);
98: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB, "Error in HDF5 viewer");
99: }
101: VecDestroy(&x1);
102: VecDestroy(&x2);
103: VecDestroy(&y1);
104: VecDestroy(&y2);
105: VecDestroy(&y3);
106: VecDestroy(&y4);
107: PetscFinalize();
108: return(0);
109: }
110: