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