Actual source code: ex16.c
2: static char help[] = "Demonstrates VecStrideScatter() and VecStrideGather() with subvectors that are also strided.\n\n";
4: /*
5: Include "petscvec.h" so that we can use vectors. Note that this file
6: automatically includes:
7: petscsys.h - base PETSc routines petscis.h - index sets
8: petscviewer.h - viewers
9: */
11: #include <petscvec.h>
13: int main(int argc, char **argv)
14: {
15: Vec v, s, r, vecs[2]; /* vectors */
16: PetscInt i, start, end, n = 20;
17: PetscScalar value;
19: PetscFunctionBeginUser;
20: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
21: PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
23: /*
24: Create multi-component vector with 2 components
25: */
26: PetscCall(VecCreate(PETSC_COMM_WORLD, &v));
27: PetscCall(VecSetSizes(v, PETSC_DECIDE, n));
28: PetscCall(VecSetBlockSize(v, 4));
29: PetscCall(VecSetFromOptions(v));
31: /*
32: Create double-component vectors
33: */
34: PetscCall(VecCreate(PETSC_COMM_WORLD, &s));
35: PetscCall(VecSetSizes(s, PETSC_DECIDE, n / 2));
36: PetscCall(VecSetBlockSize(s, 2));
37: PetscCall(VecSetFromOptions(s));
38: PetscCall(VecDuplicate(s, &r));
40: vecs[0] = s;
41: vecs[1] = r;
42: /*
43: Set the vector values
44: */
45: PetscCall(VecGetOwnershipRange(v, &start, &end));
46: for (i = start; i < end; i++) {
47: value = i;
48: PetscCall(VecSetValues(v, 1, &i, &value, INSERT_VALUES));
49: }
50: PetscCall(VecAssemblyBegin(v));
51: PetscCall(VecAssemblyEnd(v));
53: /*
54: Get the components from the multi-component vector to the other vectors
55: */
56: PetscCall(VecStrideGatherAll(v, vecs, INSERT_VALUES));
58: PetscCall(VecView(s, PETSC_VIEWER_STDOUT_WORLD));
59: PetscCall(VecView(r, PETSC_VIEWER_STDOUT_WORLD));
61: PetscCall(VecStrideScatterAll(vecs, v, ADD_VALUES));
63: PetscCall(VecView(v, PETSC_VIEWER_STDOUT_WORLD));
65: /*
66: Free work space. All PETSc objects should be destroyed when they
67: are no longer needed.
68: */
69: PetscCall(VecDestroy(&v));
70: PetscCall(VecDestroy(&s));
71: PetscCall(VecDestroy(&r));
72: PetscCall(PetscFinalize());
73: return 0;
74: }
76: /*TEST
78: test:
79: nsize: 2
81: TEST*/