xref: /petsc/src/vec/is/sf/tests/ex22.c (revision 66af8762ec03dbef0e079729eb2a1734a35ed7ff)
1 static const char help[] = "Test PetscSFFetchAndOp \n\n";
2 
3 #include <petscvec.h>
4 #include <petscsf.h>
5 
6 int main(int argc, char *argv[])
7 {
8   PetscInt           n, N = 12;
9   PetscInt          *indices;
10   IS                 ix, iy;
11   VecScatter         vscat;
12   Vec                x, y, z;
13   PetscInt           rstart, rend;
14   const PetscScalar *xarray;
15   PetscScalar       *yarray, *zarray;
16   PetscMemType       xmtype, ymtype, zmtype;
17 
18   PetscFunctionBeginUser;
19   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
20 
21   PetscCall(VecCreateFromOptions(PETSC_COMM_WORLD, NULL, 1, PETSC_DECIDE, N, &x));
22   PetscCall(VecDuplicate(x, &y));
23   PetscCall(VecDuplicate(x, &z));
24   PetscCall(VecGetLocalSize(x, &n));
25 
26   PetscCall(VecGetOwnershipRange(x, &rstart, &rend));
27   PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, rstart, 1, &ix));
28   PetscCall(PetscMalloc1(n, &indices));
29   for (PetscInt i = rstart; i < rend; i++) indices[i - rstart] = i / 4;
30   PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, n, indices, PETSC_OWN_POINTER, &iy));
31 
32   // connect y[0] to x[0..3], y[1] to x[4..7], etc
33   PetscCall(VecScatterCreate(y, iy, x, ix, &vscat)); // y has roots, x has leaves
34 
35   PetscCall(VecSet(x, 1.0));
36   PetscCall(VecSet(y, 2.0));
37 
38   PetscCall(VecGetArrayReadAndMemType(x, &xarray, &xmtype));
39   PetscCall(VecGetArrayAndMemType(y, &yarray, &ymtype));
40   PetscCall(VecGetArrayWriteAndMemType(z, &zarray, &zmtype));
41 
42   PetscCall(PetscSFFetchAndOpWithMemTypeBegin(vscat, MPIU_SCALAR, ymtype, yarray, xmtype, xarray, zmtype, zarray, MPI_SUM));
43   PetscCall(PetscSFFetchAndOpEnd(vscat, MPIU_SCALAR, yarray, xarray, zarray, MPI_SUM));
44 
45   PetscCall(VecRestoreArrayReadAndMemType(x, &xarray));
46   PetscCall(VecRestoreArrayAndMemType(y, &yarray));
47   PetscCall(VecRestoreArrayWriteAndMemType(z, &zarray));
48 
49   PetscCall(VecView(y, PETSC_VIEWER_STDOUT_WORLD));
50   PetscCall(VecView(z, PETSC_VIEWER_STDOUT_WORLD));
51   PetscCall(ISDestroy(&ix));
52   PetscCall(ISDestroy(&iy));
53   PetscCall(VecDestroy(&x));
54   PetscCall(VecDestroy(&y));
55   PetscCall(VecDestroy(&z));
56   PetscCall(VecScatterDestroy(&vscat));
57   PetscCall(PetscFinalize());
58 }
59 
60 /*TEST
61   testset:
62     nsize: {{1 4}}
63     # since FetchAndOp on complex would need to be atomic in this test
64     requires: !complex
65     output_file: output/ex22.out
66     filter: grep -v "type" | grep -v "Process" |grep -v "Vec Object"
67 
68     test:
69       suffix: cuda
70       requires: cuda
71       args: -vec_type cuda
72 
73     test:
74       suffix: hip
75       requires: hip
76       args: -vec_type hip
77 
78     test:
79       suffix: kok
80       requires: kokkos_kernels
81       args: -vec_type kokkos
82 
83 TEST*/
84