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 { 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(0); 40 } 41 42 struct FCtx { 43 PetscMPIInt rank; 44 PetscMPIInt nto; 45 PetscMPIInt *toranks; 46 Unit *todata; 47 PetscSegBuffer seg; 48 }; 49 50 static PetscErrorCode FSend(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt tonum,PetscMPIInt rank,void *todata,MPI_Request req[],void *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(0); 60 } 61 62 static PetscErrorCode FRecv(MPI_Comm comm,const PetscMPIInt tag[],PetscMPIInt rank,void *fromdata,MPI_Request req[],void *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(0); 76 } 77 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 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 87 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 88 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 89 90 verbose = PETSC_FALSE; 91 PetscCall(PetscOptionsGetBool(NULL,NULL,"-verbose",&verbose,NULL)); 92 build_twosided_f = PETSC_FALSE; 93 PetscCall(PetscOptionsGetBool(NULL,NULL,"-build_twosided_f",&build_twosided_f,NULL)); 94 95 for (i=1,nto=0; i<size; i*=2) nto++; 96 PetscCall(PetscMalloc2(nto,&todata,nto,&toranks)); 97 for (n=0,i=1; i<size; n++,i*=2) { 98 toranks[n] = (rank+i) % size; 99 todata[n].rank = (rank+i) % size; 100 todata[n].value = (PetscScalar)rank; 101 todata[n].ok[0] = 'o'; 102 todata[n].ok[1] = 'k'; 103 todata[n].ok[2] = 0; 104 } 105 if (verbose) { 106 for (i=0; i<nto; i++) { 107 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 } 109 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT)); 110 } 111 112 PetscCall(MakeDatatype(&dtype)); 113 114 if (build_twosided_f) { 115 struct FCtx fctx; 116 PetscMPIInt *todummy,*fromdummy; 117 fctx.rank = rank; 118 fctx.nto = nto; 119 fctx.toranks = toranks; 120 fctx.todata = todata; 121 PetscCall(PetscSegBufferCreate(sizeof(Unit),1,&fctx.seg)); 122 PetscCall(PetscMalloc1(nto,&todummy)); 123 for (i=0; i<nto; i++) todummy[i] = rank; 124 PetscCall(PetscCommBuildTwoSidedF(PETSC_COMM_WORLD,1,MPI_INT,nto,toranks,todummy,&nfrom,&fromranks,&fromdummy,2,FSend,FRecv,&fctx)); 125 PetscCall(PetscFree(todummy)); 126 PetscCall(PetscFree(fromdummy)); 127 PetscCall(PetscSegBufferExtractAlloc(fctx.seg,&fromdata)); 128 PetscCall(PetscSegBufferDestroy(&fctx.seg)); 129 } else { 130 PetscCall(PetscCommBuildTwoSided(PETSC_COMM_WORLD,1,dtype,nto,toranks,todata,&nfrom,&fromranks,&fromdata)); 131 } 132 PetscCallMPI(MPI_Type_free(&dtype)); 133 134 if (verbose) { 135 PetscInt *iranks,*iperm; 136 PetscCall(PetscMalloc2(nfrom,&iranks,nfrom,&iperm)); 137 for (i=0; i<nfrom; i++) { 138 iranks[i] = fromranks[i]; 139 iperm[i] = i; 140 } 141 /* Receive ordering is non-deterministic in general, so sort to make verbose output deterministic. */ 142 PetscCall(PetscSortIntWithPermutation(nfrom,iranks,iperm)); 143 for (i=0; i<nfrom; i++) { 144 PetscInt ip = iperm[i]; 145 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)); 146 } 147 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT)); 148 PetscCall(PetscFree2(iranks,iperm)); 149 } 150 151 PetscCheck(nto == nfrom,PETSC_COMM_SELF,PETSC_ERR_PLIB,"[%d] From ranks %d does not match To ranks %d",rank,nto,nfrom); 152 for (i=1; i<size; i*=2) { 153 PetscMPIInt expected_rank = (rank-i+size)%size; 154 PetscBool flg; 155 for (n=0; n<nfrom; n++) { 156 if (expected_rank == fromranks[n]) goto found; 157 } 158 SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_PLIB,"[%d] Could not find expected from rank %d",rank,expected_rank); 159 found: 160 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); 161 PetscCall(PetscStrcmp(fromdata[n].ok,"ok",&flg)); 162 PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"[%d] Got string %s from rank %d",rank,fromdata[n].ok,expected_rank); 163 } 164 PetscCall(PetscFree2(todata,toranks)); 165 PetscCall(PetscFree(fromdata)); 166 PetscCall(PetscFree(fromranks)); 167 PetscCall(PetscFinalize()); 168 return 0; 169 } 170 171 /*TEST 172 173 test: 174 nsize: 4 175 args: -verbose -build_twosided allreduce 176 177 test: 178 suffix: f 179 nsize: 4 180 args: -verbose -build_twosided_f -build_twosided allreduce 181 output_file: output/ex8_1.out 182 183 test: 184 suffix: f_ibarrier 185 nsize: 4 186 args: -verbose -build_twosided_f -build_twosided ibarrier 187 output_file: output/ex8_1.out 188 requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES) 189 190 test: 191 suffix: ibarrier 192 nsize: 4 193 args: -verbose -build_twosided ibarrier 194 output_file: output/ex8_1.out 195 requires: defined(PETSC_HAVE_MPI_NONBLOCKING_COLLECTIVES) 196 197 test: 198 suffix: redscatter 199 requires: mpi_reduce_scatter_block 200 nsize: 4 201 args: -verbose -build_twosided redscatter 202 output_file: output/ex8_1.out 203 204 TEST*/ 205