Actual source code: ex23.c
2: static char help[] = "Scatters from a parallel vector to a sequential vector.\n\
3: Using a blocked send and a strided receive.\n\n";
5: /*
6: 0 1 2 3 | 4 5 6 7 || 8 9 10 11
8: Scatter first and third block to first processor and
9: second and third block to second processor
10: */
12: #include <petscvec.h>
14: int main(int argc, char **argv)
15: {
16: PetscInt i, blocks[2], nlocal;
17: PetscMPIInt size, rank;
18: PetscScalar value;
19: Vec x, y;
20: IS is1, is2;
21: VecScatter ctx = 0;
22: PetscViewer subviewer;
24: PetscFunctionBeginUser;
25: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
26: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
27: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
29: PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Must run with 2 processors");
31: /* create two vectors */
32: if (rank == 0) nlocal = 8;
33: else nlocal = 4;
34: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
35: PetscCall(VecSetSizes(x, nlocal, 12));
36: PetscCall(VecSetFromOptions(x));
37: PetscCall(VecCreate(PETSC_COMM_SELF, &y));
38: PetscCall(VecSetSizes(y, 8, PETSC_DECIDE));
39: PetscCall(VecSetFromOptions(y));
41: /* create two index sets */
42: if (rank == 0) {
43: blocks[0] = 0;
44: blocks[1] = 2;
45: } else {
46: blocks[0] = 1;
47: blocks[1] = 2;
48: }
49: PetscCall(ISCreateBlock(PETSC_COMM_SELF, 4, 2, blocks, PETSC_COPY_VALUES, &is1));
50: PetscCall(ISCreateStride(PETSC_COMM_SELF, 8, 0, 1, &is2));
52: for (i = 0; i < 12; i++) {
53: value = i;
54: PetscCall(VecSetValues(x, 1, &i, &value, INSERT_VALUES));
55: }
56: PetscCall(VecAssemblyBegin(x));
57: PetscCall(VecAssemblyEnd(x));
59: PetscCall(VecScatterCreate(x, is1, y, is2, &ctx));
60: PetscCall(VecScatterBegin(ctx, x, y, INSERT_VALUES, SCATTER_FORWARD));
61: PetscCall(VecScatterEnd(ctx, x, y, INSERT_VALUES, SCATTER_FORWARD));
62: PetscCall(VecScatterDestroy(&ctx));
64: PetscCall(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &subviewer));
65: PetscCall(VecView(y, subviewer));
66: PetscCall(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &subviewer));
68: PetscCall(VecDestroy(&x));
69: PetscCall(VecDestroy(&y));
70: PetscCall(ISDestroy(&is1));
71: PetscCall(ISDestroy(&is2));
73: PetscCall(PetscFinalize());
74: return 0;
75: }
77: /*TEST
79: testset:
80: nsize: 2
81: output_file: output/ex23_1.out
82: filter: grep -v " type:"
83: diff_args: -j
84: test:
85: suffix: standard
86: args: -vec_type standard
87: test:
88: requires: cuda
89: suffix: cuda
90: args: -vec_type cuda
91: test:
92: requires: viennacl
93: suffix: viennacl
94: args: -vec_type viennacl
95: test:
96: requires: kokkos_kernels
97: suffix: kokkos
98: args: -vec_type kokkos
99: test:
100: requires: hip
101: suffix: hip
102: args: -vec_type hip
104: TEST*/