Actual source code: ex51.c
1: static char help[] = "Test integrity of subvector data, use \n\
2: use -hdf5 to specify HDF5 viewer format for subvector I/O \n\n";
4: /*
5: Tests for transfer of data from subvectors to parent vectors after
6: loading data into subvector. This routine does the following : creates
7: a vector of size 50, sets it to 2 and saves it to disk. Creates a
8: vector of size 100, set it to 1 and extracts the last 50 elements
9: as a subvector. Loads the saved vector from disk into the subvector
10: and restores the subvector. To verify that the data has been loaded
11: into the parent vector, the sum of it's elements is calculated.
12: */
14: #include <petscvec.h>
15: #include <petscviewerhdf5.h>
17: int main(int argc,char **argv)
18: {
19: Vec testvec; /* parent vector of size 100 */
20: Vec loadvec; /* subvector extracted from the parent vector */
21: Vec writevec; /* vector used to save data to be loaded by loadvec */
22: IS loadis; /* index set to extract last 50 elements of testvec */
23: PetscInt low,high; /* used to store vecownership output */
24: PetscInt issize, isstart; /* index set params */
25: PetscInt skipuntil = 50; /* parameter to slice the last N elements of parent vec */
26: PetscViewer viewer; /* viewer for I/O */
27: PetscScalar sum; /* used to test sum of parent vector elements */
28: PetscBool usehdf5 = PETSC_FALSE;
31: PetscInitialize(&argc, &argv, (char*) 0, help);if (ierr) return ierr;
33: /* parse input options to determine I/O format */
34: PetscOptionsGetBool(NULL,NULL,"-hdf5",&usehdf5,NULL);
36: /* Create parent vector with 100 elements, set it to 1 */
37: VecCreate(PETSC_COMM_WORLD, &testvec);
38: VecSetSizes(testvec, PETSC_DECIDE,100);
39: VecSetUp(testvec);
40: VecSet(testvec, (PetscScalar) 1);
42: /* Create a vector with 50 elements, set it to 2. */
43: VecCreate(PETSC_COMM_WORLD, &writevec);
44: VecSetSizes(writevec, PETSC_DECIDE,50);
45: VecSetUp(writevec);
46: VecSet(writevec, (PetscScalar) 2);
47: PetscObjectSetName((PetscObject)writevec,"temp");
49: /* Save to disk in specified format, destroy vector & viewer */
50: if (usehdf5) {
51: PetscPrintf(PETSC_COMM_WORLD,"writing vector in hdf5 to vector.dat ...\n");
52: PetscViewerHDF5Open(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_WRITE,&viewer);
53: } else {
54: PetscPrintf(PETSC_COMM_WORLD,"writing vector in binary to vector.dat ...\n");
55: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_WRITE,&viewer);
56: }
57: VecView(writevec,viewer);
58: VecDestroy(&writevec);
59: PetscViewerDestroy(&viewer);
61: /* Create index sets on each mpi rank to select the last 50 elements of parent vec */
62: VecGetOwnershipRange(testvec, &low, &high);
63: if (low>=skipuntil) {
64: isstart = low;
65: issize = high - low;
66: } else if (low<=skipuntil && high>=skipuntil) {
67: isstart = skipuntil;
68: issize = high - skipuntil;
69: } else {
70: isstart = low;
71: issize = 0;
72: }
73: ISCreateStride(PETSC_COMM_WORLD, issize, isstart, 1, &loadis);
75: /* Create subvector using the index set created above */
76: VecGetSubVector(testvec, loadis, &loadvec);
77: PetscObjectSetName((PetscObject)loadvec,"temp");
79: /* Load the previously saved vector into the subvector, destroy viewer */
80: if (usehdf5) {
81: PetscPrintf(PETSC_COMM_WORLD,"reading vector in hdf5 from vector.dat ...\n");
82: PetscViewerHDF5Open(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_READ,&viewer);
83: } else {
84: PetscPrintf(PETSC_COMM_WORLD,"reading vector in binary from vector.dat ...\n");
85: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"vector.dat",FILE_MODE_READ,&viewer);
86: }
87: VecLoad(loadvec, viewer);
88: PetscViewerDestroy(&viewer);
90: /* Restore subvector to transfer loaded data into parent vector */
91: VecRestoreSubVector(testvec, loadis, &loadvec);
93: /* Compute sum of parent vector elements */
94: VecSum(testvec, &sum);
96: /* to verify that the loaded data has been transferred */
97: if (sum != 150) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_PLIB,"Data has not been transferred from subvector to parent vector");
98: PetscPrintf(PETSC_COMM_WORLD,"VecSum on parent vec is : %e\n",sum);
100: /* destroy parent vector, index set and exit */
101: VecDestroy(&testvec);
102: ISDestroy(&loadis);
103: PetscFinalize();
104: return ierr;
105: }
107: /*TEST
109: build:
110: requires: hdf5
112: test:
113: nsize: 4
115: test:
116: suffix: 2
117: nsize: 4
118: args: -hdf5
120: TEST*/