Actual source code: ex19.c
petsc-3.8.4 2018-03-24
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, y1, y2, y3, y4;
15: PetscViewer viewer;
16: PetscMPIInt rank;
17: PetscInt i, nlocal, n = 6;
18: PetscScalar *array;
19: PetscBool equal;
23: PetscInitialize(&argc, &argv, (char*) 0, help);if (ierr) return ierr;
24: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
25: PetscOptionsGetInt(NULL,NULL, "-n", &n, NULL);
27: VecCreate(PETSC_COMM_WORLD, &x1);
28: PetscObjectSetName((PetscObject) x1, "TestVec");
29: VecSetSizes(x1, PETSC_DECIDE, n);
30: VecSetFromOptions(x1);
32: VecGetLocalSize(x1, &nlocal);
33: VecGetArray(x1, &array);
34: for (i = 0; i < nlocal; i++) array[i] = rank + 1;
35: VecRestoreArray(x1, &array);
36: VecAssemblyBegin(x1);
37: VecAssemblyEnd(x1);
39: VecCreate(PETSC_COMM_WORLD, &x2);
40: PetscObjectSetName((PetscObject) x2, "TestVec2");
41: VecSetSizes(x2, PETSC_DECIDE, n);
42: VecSetBlockSize(x2, 2);
43: VecSetFromOptions(x2);
44: VecCopy(x1, x2);
46: PetscViewerHDF5Open(PETSC_COMM_WORLD, "ex19.h5", FILE_MODE_WRITE, &viewer);
47: PetscViewerHDF5PushGroup(viewer, "/");
48: VecView(x1, viewer);
49: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
50: VecView(x2, viewer);
51: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
52: PetscViewerHDF5SetTimestep(viewer, 0);
53: VecView(x2, viewer);
54: PetscViewerHDF5SetTimestep(viewer, 1);
55: VecView(x2, viewer);
56: PetscViewerHDF5PopGroup(viewer);
57: PetscViewerHDF5PopGroup(viewer);
58: PetscViewerHDF5PopGroup(viewer);
59: PetscViewerDestroy(&viewer);
61: VecCreate(PETSC_COMM_WORLD, &y1);
62: PetscObjectSetName((PetscObject) y1, "TestVec");
63: VecSetSizes(y1, PETSC_DECIDE, n);
64: VecSetFromOptions(y1);
66: VecCreate(PETSC_COMM_WORLD,&y2);
67: VecSetBlockSize(y2, 2);
68: PetscObjectSetName((PetscObject) y2, "TestVec2");
70: VecCreate(PETSC_COMM_WORLD,&y3);
71: VecSetBlockSize(y3, 2);
72: PetscObjectSetName((PetscObject) y3, "TestVec2");
74: VecCreate(PETSC_COMM_WORLD,&y4);
75: VecSetBlockSize(y4, 2);
76: PetscObjectSetName((PetscObject) y4, "TestVec2");
78: PetscViewerHDF5Open(PETSC_COMM_WORLD, "ex19.h5", FILE_MODE_READ, &viewer);
79: PetscViewerHDF5PushGroup(viewer, "/");
80: VecLoad(y1, viewer);
82: PetscViewerHDF5PushGroup(viewer, "/testBlockSize");
83: VecLoad(y2, viewer);
84: PetscViewerHDF5PushGroup(viewer, "/testTimestep");
85: PetscViewerHDF5SetTimestep(viewer, 0);
86: VecLoad(y3, viewer);
87: PetscViewerHDF5SetTimestep(viewer, 1);
88: VecLoad(y4, viewer);
89: PetscViewerHDF5PopGroup(viewer);
90: PetscViewerHDF5PopGroup(viewer);
91: PetscViewerHDF5PopGroup(viewer);
92: PetscViewerDestroy(&viewer);
94: VecEqual(x1, y1, &equal);
95: if (!equal) {
96: VecView(x1, PETSC_VIEWER_STDOUT_WORLD);
97: VecView(y1, PETSC_VIEWER_STDOUT_WORLD);
98: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB, "Error in HDF5 viewer");
99: }
100: VecEqual(x2, y2, &equal);
101: if (!equal) {
102: VecView(x2, PETSC_VIEWER_STDOUT_WORLD);
103: VecView(y2, PETSC_VIEWER_STDOUT_WORLD);
104: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB, "Error in HDF5 viewer");
105: }
107: VecDestroy(&x1);
108: VecDestroy(&x2);
109: VecDestroy(&y1);
110: VecDestroy(&y2);
111: VecDestroy(&y3);
112: VecDestroy(&y4);
113: PetscFinalize();
114: return ierr;
115: }