1 static char help[] = "Scatters from a parallel vector to a parallel vector.\n\n";
2
3 #include <petscvec.h>
4
main(int argc,char ** argv)5 int main(int argc, char **argv)
6 {
7 PetscInt n = 5, N, i;
8 PetscMPIInt size, rank;
9 PetscScalar value, zero = 0.0;
10 Vec x, y;
11 IS is1, is2;
12 VecScatter ctx = 0;
13
14 PetscFunctionBeginUser;
15 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
16 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
17 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
18
19 /* create two vectors */
20 N = size * n;
21 PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
22 PetscCall(VecSetSizes(y, n, PETSC_DECIDE));
23 PetscCall(VecSetFromOptions(y));
24
25 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
26 PetscCall(VecSetSizes(x, n, PETSC_DECIDE));
27 PetscCall(VecSetFromOptions(x));
28
29 /* create two index sets */
30 PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, n * rank, 1, &is1));
31 PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, (n * (rank + 1)) % N, 1, &is2));
32
33 /* fill local part of parallel vector x */
34 value = (PetscScalar)(rank + 1);
35 for (i = n * rank; i < n * (rank + 1); i++) PetscCall(VecSetValues(x, 1, &i, &value, INSERT_VALUES));
36 PetscCall(VecAssemblyBegin(x));
37 PetscCall(VecAssemblyEnd(x));
38
39 PetscCall(VecSet(y, zero));
40
41 PetscCall(VecScatterCreate(x, is1, y, is2, &ctx));
42 PetscCall(VecScatterBegin(ctx, x, y, ADD_VALUES, SCATTER_FORWARD));
43 PetscCall(VecScatterEnd(ctx, x, y, ADD_VALUES, SCATTER_FORWARD));
44 PetscCall(VecScatterDestroy(&ctx));
45
46 PetscCall(VecView(y, PETSC_VIEWER_STDOUT_WORLD));
47
48 PetscCall(VecDestroy(&x));
49 PetscCall(VecDestroy(&y));
50 PetscCall(ISDestroy(&is1));
51 PetscCall(ISDestroy(&is2));
52
53 PetscCall(PetscFinalize());
54 return 0;
55 }
56
57 /*TEST
58
59 testset:
60 nsize: 4
61 output_file: output/ex22_1.out
62 filter: grep -v " type:"
63 diff_args: -j
64 test:
65 suffix: standard
66 args: -vec_type standard
67 test:
68 requires: cuda
69 suffix: cuda
70 args: -vec_type cuda
71 test:
72 requires: viennacl
73 suffix: viennacl
74 args: -vec_type viennacl
75 test:
76 requires: kokkos_kernels
77 suffix: kokkos
78 args: -vec_type kokkos
79 test:
80 requires: hip
81 suffix: hip
82 args: -vec_type hip
83
84 TEST*/
85