Actual source code: ex45.c
2: static char help[] = "Demonstrates VecStrideSubSetScatter() and VecStrideSubSetGather().\n\n";
4: /*
5: Allows one to easily pull out some components of a multi-component vector and put them in another vector.
7: Note that these are special cases of VecScatter
8: */
10: /*
11: Include "petscvec.h" so that we can use vectors. Note that this file
12: automatically includes:
13: petscsys.h - base PETSc routines petscis.h - index sets
14: petscviewer.h - viewers
15: */
17: #include <petscvec.h>
19: int main(int argc, char **argv)
20: {
21: Vec v, s;
22: PetscInt i, start, end, n = 8;
23: PetscScalar value;
24: const PetscInt vidx[] = {1, 2}, sidx[] = {1, 0};
25: PetscInt miidx[2];
26: PetscReal mvidx[2];
29: PetscInitialize(&argc, &argv, (char *)0, help);
30: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
32: /*
33: Create multi-component vector with 4 components
34: */
35: VecCreate(PETSC_COMM_WORLD, &v);
36: VecSetSizes(v, PETSC_DECIDE, n);
37: VecSetBlockSize(v, 4);
38: VecSetFromOptions(v);
40: /*
41: Create double-component vectors
42: */
43: VecCreate(PETSC_COMM_WORLD, &s);
44: VecSetSizes(s, PETSC_DECIDE, n / 2);
45: VecSetBlockSize(s, 2);
46: VecSetFromOptions(s);
48: /*
49: Set the vector values
50: */
51: VecGetOwnershipRange(v, &start, &end);
52: for (i = start; i < end; i++) {
53: value = i;
54: VecSetValues(v, 1, &i, &value, INSERT_VALUES);
55: }
57: /*
58: Get the components from the large multi-component vector to the small multi-component vector,
59: scale the smaller vector and then move values back to the large vector
60: */
61: VecStrideSubSetGather(v, PETSC_DETERMINE, vidx, NULL, s, INSERT_VALUES);
62: VecView(s, PETSC_VIEWER_STDOUT_WORLD);
63: VecScale(s, 100.0);
65: VecStrideSubSetScatter(s, PETSC_DETERMINE, NULL, vidx, v, ADD_VALUES);
66: VecView(v, PETSC_VIEWER_STDOUT_WORLD);
68: /*
69: Get the components from the large multi-component vector to the small multi-component vector,
70: scale the smaller vector and then move values back to the large vector
71: */
72: VecStrideSubSetGather(v, 2, vidx, sidx, s, INSERT_VALUES);
73: VecView(s, PETSC_VIEWER_STDOUT_WORLD);
74: VecScale(s, 100.0);
76: VecStrideSubSetScatter(s, 2, sidx, vidx, v, ADD_VALUES);
77: VecView(v, PETSC_VIEWER_STDOUT_WORLD);
79: VecStrideMax(v, 1, &miidx[0], &mvidx[0]);
80: VecStrideMin(v, 1, &miidx[1], &mvidx[1]);
81: PetscPrintf(PETSC_COMM_WORLD, "Min/Max: %" PetscInt_FMT " %g, %" PetscInt_FMT " %g\n", miidx[0], (double)mvidx[0], miidx[1], (double)mvidx[1]);
82: /*
83: Free work space. All PETSc objects should be destroyed when they
84: are no longer needed.
85: */
86: VecDestroy(&v);
87: VecDestroy(&s);
88: PetscFinalize();
89: return 0;
90: }
92: /*TEST
94: test:
95: filter: grep -v type | grep -v " MPI process" | grep -v Process
96: diff_args: -j
97: nsize: 2
99: test:
100: filter: grep -v type | grep -v " MPI process" | grep -v Process
101: output_file: output/ex45_1.out
102: diff_args: -j
103: suffix: 2
104: nsize: 1
105: args: -vec_type {{seq mpi}}
107: TEST*/