1 static const char help[] = "Test overlapped communication on a single star forest (PetscSF)\n\n";
2
3 #include <petscvec.h>
4 #include <petscsf.h>
5 #include <petscviewer.h>
6
main(int argc,char ** argv)7 int main(int argc, char **argv)
8 {
9 PetscSF sf;
10 Vec A, Aout;
11 Vec B, Bout;
12 PetscScalar *bufA;
13 PetscScalar *bufAout;
14 PetscScalar *bufB;
15 PetscScalar *bufBout;
16 PetscMPIInt rank, size;
17 PetscInt nroots, nleaves;
18 PetscInt i;
19 PetscInt *ilocal;
20 PetscSFNode *iremote;
21
22 PetscFunctionBeginUser;
23 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
24 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
25 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
26
27 PetscCheck(size == 2, PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for two MPI processes");
28
29 PetscCall(PetscSFCreate(PETSC_COMM_WORLD, &sf));
30 PetscCall(PetscSFSetFromOptions(sf));
31
32 nleaves = 2;
33 nroots = 1;
34 PetscCall(PetscMalloc1(nleaves, &ilocal));
35
36 for (i = 0; i < nleaves; i++) ilocal[i] = i;
37
38 PetscCall(PetscMalloc1(nleaves, &iremote));
39 if (rank == 0) {
40 iremote[0].rank = 0;
41 iremote[0].index = 0;
42 iremote[1].rank = 1;
43 iremote[1].index = 0;
44 } else {
45 iremote[0].rank = 1;
46 iremote[0].index = 0;
47 iremote[1].rank = 0;
48 iremote[1].index = 0;
49 }
50 PetscCall(PetscSFSetGraph(sf, nroots, nleaves, ilocal, PETSC_OWN_POINTER, iremote, PETSC_OWN_POINTER));
51 PetscCall(PetscSFSetUp(sf));
52 PetscCall(PetscSFView(sf, PETSC_VIEWER_STDOUT_WORLD));
53 PetscCall(VecCreate(PETSC_COMM_WORLD, &A));
54 PetscCall(VecSetSizes(A, 2, PETSC_DETERMINE));
55 PetscCall(VecSetFromOptions(A));
56 PetscCall(VecSetUp(A));
57
58 PetscCall(VecDuplicate(A, &B));
59 PetscCall(VecDuplicate(A, &Aout));
60 PetscCall(VecDuplicate(A, &Bout));
61 PetscCall(VecGetArray(A, &bufA));
62 PetscCall(VecGetArray(B, &bufB));
63 for (i = 0; i < 2; i++) {
64 bufA[i] = (PetscScalar)rank;
65 bufB[i] = (PetscScalar)rank + 10.0;
66 }
67 PetscCall(VecRestoreArray(A, &bufA));
68 PetscCall(VecRestoreArray(B, &bufB));
69
70 PetscCall(VecGetArrayRead(A, (const PetscScalar **)&bufA));
71 PetscCall(VecGetArrayRead(B, (const PetscScalar **)&bufB));
72 PetscCall(VecGetArray(Aout, &bufAout));
73 PetscCall(VecGetArray(Bout, &bufBout));
74 PetscCall(PetscSFBcastBegin(sf, MPIU_SCALAR, (const void *)bufA, (void *)bufAout, MPI_REPLACE));
75 PetscCall(PetscSFBcastBegin(sf, MPIU_SCALAR, (const void *)bufB, (void *)bufBout, MPI_REPLACE));
76 PetscCall(PetscSFBcastEnd(sf, MPIU_SCALAR, (const void *)bufA, (void *)bufAout, MPI_REPLACE));
77 PetscCall(PetscSFBcastEnd(sf, MPIU_SCALAR, (const void *)bufB, (void *)bufBout, MPI_REPLACE));
78 PetscCall(VecRestoreArrayRead(A, (const PetscScalar **)&bufA));
79 PetscCall(VecRestoreArrayRead(B, (const PetscScalar **)&bufB));
80 PetscCall(VecRestoreArray(Aout, &bufAout));
81 PetscCall(VecRestoreArray(Bout, &bufBout));
82
83 PetscCall(VecView(Aout, PETSC_VIEWER_STDOUT_WORLD));
84 PetscCall(VecView(Bout, PETSC_VIEWER_STDOUT_WORLD));
85 PetscCall(VecDestroy(&A));
86 PetscCall(VecDestroy(&B));
87 PetscCall(VecDestroy(&Aout));
88 PetscCall(VecDestroy(&Bout));
89 PetscCall(PetscSFDestroy(&sf));
90
91 PetscCall(PetscFinalize());
92 return 0;
93 }
94
95 /*TEST
96
97 test:
98 suffix: basic
99 nsize: 2
100 filter: grep -v "type" | grep -v "sort"
101 args: -sf_type basic
102
103 test:
104 suffix: window
105 nsize: 2
106 filter: grep -v "type" | grep -v "sort"
107 output_file: output/ex2_basic.out
108 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create dynamic allocate}}
109 requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW)
110
111 # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes
112 test:
113 suffix: window_shared
114 nsize: 2
115 filter: grep -v "type" | grep -v "sort"
116 output_file: output/ex2_basic.out
117 args: -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared
118 requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH) defined(PETSC_HAVE_MPI_ONE_SIDED)
119
120 TEST*/
121