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