1 static char help[] = "Demonstrates BuildTwoSided functions.\n";
2
3 #include <petscsys.h>
4
5 typedef struct {
6 PetscInt rank;
7 PetscScalar value;
8 char ok[3];
9 } Unit;
10
MakeDatatype(MPI_Datatype * dtype)11 static PetscErrorCode MakeDatatype(MPI_Datatype *dtype)
12 {
13 MPI_Datatype dtypes[3], tmptype;
14 PetscMPIInt lengths[3];
15 MPI_Aint displs[3];
16 Unit dummy;
17
18 PetscFunctionBegin;
19 dtypes[0] = MPIU_INT;
20 dtypes[1] = MPIU_SCALAR;
21 dtypes[2] = MPI_CHAR;
22 lengths[0] = 1;
23 lengths[1] = 1;
24 lengths[2] = 3;
25 /* Curse the evil beings that made std::complex a non-POD type. */
26 displs[0] = (char *)&dummy.rank - (char *)&dummy; /* offsetof(Unit,rank); */
27 displs[1] = (char *)&dummy.value - (char *)&dummy; /* offsetof(Unit,value); */
28 displs[2] = (char *)&dummy.ok - (char *)&dummy; /* offsetof(Unit,ok); */
29 PetscCallMPI(MPI_Type_create_struct(3, lengths, displs, dtypes, &tmptype));
30 PetscCallMPI(MPI_Type_commit(&tmptype));
31 PetscCallMPI(MPI_Type_create_resized(tmptype, 0, sizeof(Unit), dtype));
32 PetscCallMPI(MPI_Type_commit(dtype));
33 PetscCallMPI(MPI_Type_free(&tmptype));
34 {
35 MPI_Aint lb, extent;
36 PetscCallMPI(MPI_Type_get_extent(*dtype, &lb, &extent));
37 PetscCheck(extent == sizeof(Unit), PETSC_COMM_WORLD, PETSC_ERR_LIB, "New type has extent %d != sizeof(Unit) %d", (int)extent, (int)sizeof(Unit));
38 }
39 PetscFunctionReturn(PETSC_SUCCESS);
40 }
41
42 struct FCtx {
43 PetscMPIInt rank;
44 PetscMPIInt nto;
45 PetscMPIInt *toranks;
46 Unit *todata;
47 PetscSegBuffer seg;
48 };
49
FSend(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt tonum,PetscMPIInt rank,void * todata,MPI_Request req[],PetscCtx ctx)50 static PetscErrorCode FSend(MPI_Comm comm, const PetscMPIInt tag[], PetscMPIInt tonum, PetscMPIInt rank, void *todata, MPI_Request req[], PetscCtx ctx)
51 {
52 struct FCtx *fctx = (struct FCtx *)ctx;
53
54 PetscFunctionBegin;
55 PetscCheck(rank == fctx->toranks[tonum], PETSC_COMM_SELF, PETSC_ERR_PLIB, "Rank %d does not match toranks[%d] %d", rank, tonum, fctx->toranks[tonum]);
56 PetscCheck(fctx->rank == *(PetscMPIInt *)todata, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Todata %d does not match rank %d", *(PetscMPIInt *)todata, fctx->rank);
57 PetscCallMPI(MPI_Isend(&fctx->todata[tonum].rank, 1, MPIU_INT, rank, tag[0], comm, &req[0]));
58 PetscCallMPI(MPI_Isend(&fctx->todata[tonum].value, 1, MPIU_SCALAR, rank, tag[1], comm, &req[1]));
59 PetscFunctionReturn(PETSC_SUCCESS);
60 }
61
FRecv(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt rank,void * fromdata,MPI_Request req[],PetscCtx ctx)62 static PetscErrorCode FRecv(MPI_Comm comm, const PetscMPIInt tag[], PetscMPIInt rank, void *fromdata, MPI_Request req[], PetscCtx ctx)
63 {
64 struct FCtx *fctx = (struct FCtx *)ctx;
65 Unit *buf;
66
67 PetscFunctionBegin;
68 PetscCheck(*(PetscMPIInt *)fromdata == rank, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Dummy data %d from rank %d corrupt", *(PetscMPIInt *)fromdata, rank);
69 PetscCall(PetscSegBufferGet(fctx->seg, 1, &buf));
70 PetscCallMPI(MPI_Irecv(&buf->rank, 1, MPIU_INT, rank, tag[0], comm, &req[0]));
71 PetscCallMPI(MPI_Irecv(&buf->value, 1, MPIU_SCALAR, rank, tag[1], comm, &req[1]));
72 buf->ok[0] = 'o';
73 buf->ok[1] = 'k';
74 buf->ok[2] = 0;
75 PetscFunctionReturn(PETSC_SUCCESS);
76 }
77
main(int argc,char ** argv)78 int main(int argc, char **argv)
79 {
80 PetscMPIInt rank, size, *toranks, *fromranks, nto, nfrom;
81 PetscInt i, n;
82 PetscBool verbose, build_twosided_f;
83 Unit *todata, *fromdata;
84 MPI_Datatype dtype;
85
86 PetscFunctionBeginUser;
87 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
88 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
89 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
90
91 verbose = PETSC_FALSE;
92 PetscCall(PetscOptionsGetBool(NULL, NULL, "-verbose", &verbose, NULL));
93 build_twosided_f = PETSC_FALSE;
94 PetscCall(PetscOptionsGetBool(NULL, NULL, "-build_twosided_f", &build_twosided_f, NULL));
95
96 for (i = 1, nto = 0; i < size; i *= 2) nto++;
97 PetscCall(PetscMalloc2(nto, &todata, nto, &toranks));
98 for (n = 0, i = 1; i < size; n++, i *= 2) {
99 toranks[n] = (rank + i) % size;
100 todata[n].rank = (rank + i) % size;
101 todata[n].value = (PetscScalar)rank;
102 todata[n].ok[0] = 'o';
103 todata[n].ok[1] = 'k';
104 todata[n].ok[2] = 0;
105 }
106 if (verbose) {
107 for (i = 0; i < nto; i++) PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] TO %d: {%" PetscInt_FMT ", %g, \"%s\"}\n", rank, toranks[i], todata[i].rank, (double)PetscRealPart(todata[i].value), todata[i].ok));
108 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT));
109 }
110
111 PetscCall(MakeDatatype(&dtype));
112
113 if (build_twosided_f) {
114 struct FCtx fctx;
115 PetscMPIInt *todummy, *fromdummy;
116 fctx.rank = rank;
117 fctx.nto = nto;
118 fctx.toranks = toranks;
119 fctx.todata = todata;
120 PetscCall(PetscSegBufferCreate(sizeof(Unit), 1, &fctx.seg));
121 PetscCall(PetscMalloc1(nto, &todummy));
122 for (i = 0; i < nto; i++) todummy[i] = rank;
123 PetscCall(PetscCommBuildTwoSidedF(PETSC_COMM_WORLD, 1, MPI_INT, nto, toranks, todummy, &nfrom, &fromranks, &fromdummy, 2, FSend, FRecv, &fctx));
124 PetscCall(PetscFree(todummy));
125 PetscCall(PetscFree(fromdummy));
126 PetscCall(PetscSegBufferExtractAlloc(fctx.seg, &fromdata));
127 PetscCall(PetscSegBufferDestroy(&fctx.seg));
128 } else {
129 PetscCall(PetscCommBuildTwoSided(PETSC_COMM_WORLD, 1, dtype, nto, toranks, todata, &nfrom, &fromranks, &fromdata));
130 }
131 PetscCallMPI(MPI_Type_free(&dtype));
132
133 if (verbose) {
134 PetscInt *iranks, *iperm;
135 PetscCall(PetscMalloc2(nfrom, &iranks, nfrom, &iperm));
136 for (i = 0; i < nfrom; i++) {
137 iranks[i] = fromranks[i];
138 iperm[i] = i;
139 }
140 /* Receive ordering is non-deterministic in general, so sort to make verbose output deterministic. */
141 PetscCall(PetscSortIntWithPermutation(nfrom, iranks, iperm));
142 for (i = 0; i < nfrom; i++) {
143 PetscInt ip = iperm[i];
144 PetscCall(PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] FROM %d: {%" PetscInt_FMT ", %g, \"%s\"}\n", rank, fromranks[ip], fromdata[ip].rank, (double)PetscRealPart(fromdata[ip].value), fromdata[ip].ok));
145 }
146 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD, PETSC_STDOUT));
147 PetscCall(PetscFree2(iranks, iperm));
148 }
149
150 PetscCheck(nto == nfrom, PETSC_COMM_SELF, PETSC_ERR_PLIB, "[%d] From ranks %d does not match To ranks %d", rank, nto, nfrom);
151 for (i = 1; i < size; i *= 2) {
152 PetscMPIInt expected_rank = (rank - i + size) % size;
153 PetscBool flg;
154 for (n = 0; n < nfrom; n++) {
155 if (expected_rank == fromranks[n]) goto found;
156 }
157 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "[%d] Could not find expected from rank %d", rank, expected_rank);
158 found:
159 PetscCheck(PetscRealPart(fromdata[n].value) == expected_rank, PETSC_COMM_SELF, PETSC_ERR_PLIB, "[%d] Got data %g from rank %d", rank, (double)PetscRealPart(fromdata[n].value), expected_rank);
160 PetscCall(PetscStrcmp(fromdata[n].ok, "ok", &flg));
161 PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "[%d] Got string %s from rank %d", rank, fromdata[n].ok, expected_rank);
162 }
163 PetscCall(PetscFree2(todata, toranks));
164 PetscCall(PetscFree(fromdata));
165 PetscCall(PetscFree(fromranks));
166 PetscCall(PetscFinalize());
167 return 0;
168 }
169
170 /*TEST
171
172 test:
173 nsize: 4
174 args: -verbose -build_twosided allreduce
175
176 test:
177 suffix: f
178 nsize: 4
179 args: -verbose -build_twosided_f -build_twosided allreduce
180 output_file: output/ex8_1.out
181
182 test:
183 suffix: f_ibarrier
184 nsize: 4
185 args: -verbose -build_twosided_f -build_twosided ibarrier
186 output_file: output/ex8_1.out
187 requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES)
188
189 test:
190 suffix: ibarrier
191 nsize: 4
192 args: -verbose -build_twosided ibarrier
193 output_file: output/ex8_1.out
194 requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES)
195
196 test:
197 suffix: redscatter
198 requires: mpi_reduce_scatter_block
199 nsize: 4
200 args: -verbose -build_twosided redscatter
201 output_file: output/ex8_1.out
202
203 TEST*/
204