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