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