1 static char help[] = "Demonstrates a scatter with a stride and general index set.\n\n";
2
3 #include <petscvec.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 PetscInt n = 6, idx1[3] = {0, 1, 2}, loc[6] = {0, 1, 2, 3, 4, 5};
8 PetscScalar two = 2.0, vals[6] = {10, 11, 12, 13, 14, 15};
9 Vec x, y;
10 IS is1, is2;
11 VecScatter ctx = 0;
12
13 PetscFunctionBeginUser;
14 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
15
16 /* create two vectors */
17 PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &x));
18 PetscCall(VecDuplicate(x, &y));
19
20 /* create two index sets */
21 PetscCall(ISCreateStride(PETSC_COMM_SELF, 3, 0, 2, &is1));
22 PetscCall(ISCreateGeneral(PETSC_COMM_SELF, 3, idx1, PETSC_COPY_VALUES, &is2));
23
24 PetscCall(VecSetValues(x, 6, loc, vals, INSERT_VALUES));
25 PetscCall(VecView(x, PETSC_VIEWER_STDOUT_SELF));
26 PetscCall(PetscPrintf(PETSC_COMM_SELF, "----\n"));
27 PetscCall(VecSet(y, two));
28 PetscCall(VecScatterCreate(x, is1, y, is2, &ctx));
29 PetscCall(VecScatterBegin(ctx, x, y, INSERT_VALUES, SCATTER_FORWARD));
30 PetscCall(VecScatterEnd(ctx, x, y, INSERT_VALUES, SCATTER_FORWARD));
31 PetscCall(VecScatterDestroy(&ctx));
32
33 PetscCall(VecView(y, PETSC_VIEWER_STDOUT_SELF));
34
35 PetscCall(ISDestroy(&is1));
36 PetscCall(ISDestroy(&is2));
37 PetscCall(VecDestroy(&x));
38 PetscCall(VecDestroy(&y));
39
40 PetscCall(PetscFinalize());
41 return 0;
42 }
43
44 /*TEST
45
46 test:
47
48 TEST*/
49