1 #if !defined(__SFBASIC_H) 2 #define __SFBASIC_H 3 4 #include <../src/vec/is/sf/impls/basic/sfpack.h> 5 6 #define SFBASICHEADER \ 7 PetscMPIInt niranks; /* Number of incoming ranks (ranks accessing my roots) */ \ 8 PetscMPIInt ndiranks; /* Number of incoming ranks (ranks accessing my roots) in distinguished set */ \ 9 PetscMPIInt *iranks; /* Array of ranks that reference my roots */ \ 10 PetscInt itotal; /* Total number of graph edges referencing my roots */ \ 11 PetscInt *ioffset; /* Array of length niranks+1 holding offset in irootloc[] for each rank */ \ 12 PetscInt *irootloc; /* Incoming roots referenced by ranks starting at ioffset[rank] */ \ 13 PetscInt *irootloc_d; /* irootloc in device memory when needed */ \ 14 PetscSFPackOpt rootpackopt; /* Optimization plans to (un)pack roots based on patterns in irootloc[]. NULL for no plans */ \ 15 PetscSFPackOpt selfrootpackopt; /* Optimization plans to (un)pack roots connected to local leaves */ \ 16 PetscSFPack avail; /* One or more entries per MPI Datatype, lazily constructed */ \ 17 PetscSFPack inuse; /* Buffers being used for transactions that have not yet completed */ \ 18 PetscBool selfrootdups; /* Indices of roots in irootloc[0,ioffset[ndiranks]) have dups, implying theads working ... */ \ 19 /* ... on these roots in parallel may have data race. */ \ 20 PetscBool remoterootdups /* Indices of roots in irootloc[ioffset[ndiranks],ioffset[niranks]) have dups */ 21 22 typedef struct { 23 SFBASICHEADER; 24 } PetscSF_Basic; 25 26 PETSC_STATIC_INLINE PetscErrorCode PetscSFGetRootInfo_Basic(PetscSF sf,PetscInt *nrootranks,PetscInt *ndrootranks,const PetscMPIInt **rootranks,const PetscInt **rootoffset,const PetscInt **rootloc) 27 { 28 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 29 30 PetscFunctionBegin; 31 if (nrootranks) *nrootranks = bas->niranks; 32 if (ndrootranks) *ndrootranks = bas->ndiranks; 33 if (rootranks) *rootranks = bas->iranks; 34 if (rootoffset) *rootoffset = bas->ioffset; 35 if (rootloc) *rootloc = bas->irootloc; 36 PetscFunctionReturn(0); 37 } 38 39 PETSC_STATIC_INLINE PetscErrorCode PetscSFGetLeafInfo_Basic(PetscSF sf,PetscInt *nleafranks,PetscInt *ndleafranks,const PetscMPIInt **leafranks,const PetscInt **leafoffset,const PetscInt **leafloc,const PetscInt **leafrremote) 40 { 41 PetscFunctionBegin; 42 if (nleafranks) *nleafranks = sf->nranks; 43 if (ndleafranks) *ndleafranks = sf->ndranks; 44 if (leafranks) *leafranks = sf->ranks; 45 if (leafoffset) *leafoffset = sf->roffset; 46 if (leafloc) *leafloc = sf->rmine; 47 if (leafrremote) *leafrremote = sf->rremote; 48 PetscFunctionReturn(0); 49 } 50 51 /* Get root locations either on Host or Device */ 52 PETSC_STATIC_INLINE PetscErrorCode PetscSFGetRootIndicesWithMemType_Basic(PetscSF sf,PetscMemType mtype, const PetscInt **rootloc) 53 { 54 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 55 PetscFunctionBegin; 56 if (rootloc) { 57 if (mtype == PETSC_MEMTYPE_HOST) *rootloc = bas->irootloc; 58 #if defined(PETSC_HAVE_CUDA) 59 else if (mtype == PETSC_MEMTYPE_DEVICE) { 60 if (!bas->irootloc_d) { 61 cudaError_t err; 62 PetscErrorCode ierr; 63 size_t size = bas->ioffset[bas->niranks]*sizeof(PetscInt); 64 err = cudaMalloc((void **)&bas->irootloc_d,size);CHKERRCUDA(err); 65 ierr = PetscMemcpyWithMemType(PETSC_MEMTYPE_DEVICE,PETSC_MEMTYPE_HOST,bas->irootloc_d,bas->irootloc,size);CHKERRQ(ierr); 66 } 67 *rootloc = bas->irootloc_d; 68 } 69 #endif 70 else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Wrong PetscMemType in getting rootloc %d",(int)mtype); 71 } 72 PetscFunctionReturn(0); 73 } 74 75 /* Get leaf locations either on Host (CPU) or Device (GPU) */ 76 PETSC_STATIC_INLINE PetscErrorCode PetscSFGetLeafIndicesWithMemType_Basic(PetscSF sf,PetscMemType mtype, const PetscInt **leafloc) 77 { 78 PetscFunctionBegin; 79 if (leafloc) { 80 if (mtype == PETSC_MEMTYPE_HOST) *leafloc = sf->rmine; 81 #if defined(PETSC_HAVE_CUDA) 82 else if (mtype == PETSC_MEMTYPE_DEVICE) { 83 if (!sf->rmine_d) { 84 cudaError_t err; 85 PetscErrorCode ierr; 86 size_t size = sf->roffset[sf->nranks]*sizeof(PetscInt); 87 err = cudaMalloc((void **)&sf->rmine_d,size);CHKERRCUDA(err); 88 ierr = PetscMemcpyWithMemType(PETSC_MEMTYPE_DEVICE,PETSC_MEMTYPE_HOST,sf->rmine_d,sf->rmine,size);CHKERRQ(ierr); 89 } 90 *leafloc = sf->rmine_d; 91 } 92 #endif 93 else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"Wrong PetscMemType in getting leafloc %d",(int)mtype); 94 } 95 PetscFunctionReturn(0); 96 } 97 98 /* A convenience struct to provide info to the following (un)packing routines so that we can treat packings to self and to remote in one loop. */ 99 typedef struct _n_PackInfo { 100 PetscInt count; /* Number of entries to pack, unpack etc. */ 101 const PetscInt *idx; /* Indices of the entries. NULL for contiguous indices [0,count-1) */ 102 PetscSFPackOpt opt; /* Pack optimizations */ 103 char *buf; /* The contiguous buffer where we pack to or unpack from */ 104 PetscBool atomic; /* Whether the unpack routine needs to use atomics */ 105 } PackInfo; 106 107 /* Utility routine to pack selected entries of rootdata into root buffer. 108 Input Arguments: 109 + sf - The SF this packing works on. 110 . link - The PetscSFPack, which gives the memtype of the roots and also provides root buffer. 111 . rootloc - Indices of the roots, only meaningful if the root space is sparse 112 . rootdata - Where to read the roots. 113 - sparse - Is the root space sparse (for SFBasic, SFNeighbor) or dense (for SFAllgatherv etc) 114 */ 115 PETSC_STATIC_INLINE PetscErrorCode PetscSFPackRootData(PetscSF sf,PetscSFPack link,const PetscInt *rootloc,const void *rootdata,PetscBool sparse) 116 { 117 PetscErrorCode ierr; 118 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 119 PetscInt i; 120 PetscErrorCode (*Pack)(PetscInt,const PetscInt*,PetscSFPack,PetscSFPackOpt,const void*,void*); 121 PackInfo p[2]; 122 123 PetscFunctionBegin; 124 if (sparse) { 125 /* For SFBasic and SFNeighbor, whose root space is sparse and have separate buffers for self and remote. */ 126 p[0].count = link->selfbuflen; p[1].count = link->rootbuflen; 127 p[0].idx = rootloc; p[1].idx = rootloc+bas->ioffset[bas->ndiranks]; 128 p[0].opt = bas->selfrootpackopt; p[1].opt = bas->rootpackopt; 129 p[0].buf = link->selfbuf[link->rootmtype]; p[1].buf = link->rootbuf[link->rootmtype]; 130 p[0].atomic = PETSC_FALSE; p[1].atomic = PETSC_FALSE; 131 } else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_PLIB,"SFAllgather etc should directly use rootdata instead of packing it"); 132 133 ierr = PetscSFPackGetPack(link,link->rootmtype,&Pack);CHKERRQ(ierr); 134 /* Only do packing when count != 0 so that we can avoid invoking empty CUDA kernels */ 135 for (i=0; i<2; i++) {if (p[i].count) {ierr = (*Pack)(p[i].count,p[i].idx,link,p[i].opt,rootdata,p[i].buf);CHKERRQ(ierr);}} 136 137 #if defined(PETSC_HAVE_CUDA) 138 if (link->rootmtype == PETSC_MEMTYPE_DEVICE && p[1].count) { /* We only care p[1], which is for remote and involves MPI */ 139 /* Without use_gpu_aware_mpi, we need to copy rootbuf on device to rootbuf on host. The cudaStreamSynchronize() 140 is to make sure rootbuf is ready before MPI communicatioin starts. 141 */ 142 cudaError_t err; 143 if (!use_gpu_aware_mpi) { 144 err = cudaMemcpyAsync(link->rootbuf[PETSC_MEMTYPE_HOST],link->rootbuf[PETSC_MEMTYPE_DEVICE],link->rootbuflen*link->unitbytes,cudaMemcpyDeviceToHost,link->stream);CHKERRCUDA(err); 145 ierr = PetscLogGpuToCpu(link->rootbuflen*link->unitbytes);CHKERRQ(ierr); 146 } 147 err = cudaStreamSynchronize(link->stream);CHKERRCUDA(err); 148 } 149 #endif 150 PetscFunctionReturn(0); 151 } 152 153 /* Utility routine to pack selected entries of leafdata into leaf buffer */ 154 PETSC_STATIC_INLINE PetscErrorCode PetscSFPackLeafData(PetscSF sf,PetscSFPack link,const PetscInt *leafloc,const void *leafdata,PetscBool sparse) 155 { 156 PetscErrorCode ierr; 157 PetscInt i; 158 PetscErrorCode (*Pack)(PetscInt,const PetscInt*,PetscSFPack,PetscSFPackOpt,const void*,void*); 159 PackInfo p[2]; 160 161 PetscFunctionBegin; 162 if (sparse) { 163 p[0].count = link->selfbuflen; p[1].count = link->leafbuflen; 164 p[0].idx = leafloc; p[1].idx = leafloc+sf->roffset[sf->ndranks]; 165 p[0].opt = sf->selfleafpackopt; p[1].opt = sf->leafpackopt; 166 p[0].buf = link->selfbuf[link->leafmtype]; p[1].buf = link->leafbuf[link->leafmtype]; 167 p[0].atomic = PETSC_FALSE; p[1].atomic = PETSC_FALSE; 168 } else SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_PLIB,"SFAllgather etc should directly use leafdata instead of packing it"); 169 170 ierr = PetscSFPackGetPack(link,link->leafmtype,&Pack);CHKERRQ(ierr); 171 for (i=0; i<2; i++) {if (p[i].count) {ierr = (*Pack)(p[i].count,p[i].idx,link,p[i].opt,leafdata,p[i].buf);CHKERRQ(ierr);}} 172 #if defined(PETSC_HAVE_CUDA) 173 if (link->leafmtype == PETSC_MEMTYPE_DEVICE && p[1].count) { 174 cudaError_t err; 175 if (!use_gpu_aware_mpi) { 176 err = cudaMemcpyAsync(link->leafbuf[PETSC_MEMTYPE_HOST],link->leafbuf[PETSC_MEMTYPE_DEVICE],link->leafbuflen*link->unitbytes,cudaMemcpyDeviceToHost,link->stream);CHKERRCUDA(err); 177 ierr = PetscLogGpuToCpu(link->leafbuflen*link->unitbytes);CHKERRQ(ierr); 178 } 179 err = cudaStreamSynchronize(link->stream);CHKERRCUDA(err); 180 } 181 #endif 182 PetscFunctionReturn(0); 183 } 184 185 /* Utility routine to unpack data from root buffer and Op it into selected entries of rootdata */ 186 PETSC_STATIC_INLINE PetscErrorCode PetscSFUnpackAndOpRootData(PetscSF sf,PetscSFPack link,const PetscInt *rootloc,void *rootdata,MPI_Op op,PetscBool sparse) 187 { 188 PetscErrorCode ierr; 189 PetscInt i; 190 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 191 PetscErrorCode (*UnpackAndOp)(PetscInt,const PetscInt*,PetscSFPack,PetscSFPackOpt,void*,const void*); 192 PackInfo p[2]; 193 194 PetscFunctionBegin; 195 #if defined(PETSC_HAVE_CUDA) 196 if (!use_gpu_aware_mpi && link->rootmtype == PETSC_MEMTYPE_DEVICE) { /* Copy roots from host to device when needed */ 197 cudaError_t err; 198 if (!link->rootbuf[PETSC_MEMTYPE_DEVICE]) {ierr = PetscMallocWithMemType(PETSC_MEMTYPE_DEVICE,link->rootbuflen*link->unitbytes,(void**)&link->rootbuf[PETSC_MEMTYPE_DEVICE]);CHKERRQ(ierr);} 199 err = cudaMemcpyAsync(link->rootbuf[PETSC_MEMTYPE_DEVICE],link->rootbuf[PETSC_MEMTYPE_HOST],link->rootbuflen*link->unitbytes,cudaMemcpyHostToDevice,link->stream);CHKERRCUDA(err); 200 ierr = PetscLogCpuToGpu(link->rootbuflen*link->unitbytes);CHKERRQ(ierr); 201 } 202 #endif 203 204 if (sparse) { 205 p[0].count = link->selfbuflen; p[1].count = link->rootbuflen; 206 p[0].idx = rootloc; p[1].idx = rootloc+bas->ioffset[bas->ndiranks]; 207 p[0].opt = bas->selfrootpackopt; p[1].opt = bas->rootpackopt; 208 p[0].buf = link->selfbuf[link->rootmtype]; p[1].buf = link->rootbuf[link->rootmtype]; 209 p[0].atomic = bas->selfrootdups; p[1].atomic = bas->remoterootdups; 210 } else { 211 p[0].count = 0; 212 p[1].count = sf->nroots; 213 p[1].idx = NULL; 214 p[1].opt = NULL; 215 p[1].buf = link->rootbuf[link->rootmtype]; /* Might just allocated by PetscMallocWithMemType() above */ 216 p[1].atomic = PETSC_FALSE; 217 } 218 219 for (i=0; i<2; i++) { 220 if (p[i].count) { 221 ierr = PetscSFPackGetUnpackAndOp(link,link->rootmtype,op,p[i].atomic,&UnpackAndOp);CHKERRQ(ierr); 222 if (UnpackAndOp) {ierr = (*UnpackAndOp)(p[i].count,p[i].idx,link,p[i].opt,rootdata,p[i].buf);CHKERRQ(ierr);} 223 #if defined(PETSC_HAVE_MPI_REDUCE_LOCAL) 224 else { 225 PetscInt j; 226 PetscMPIInt n; 227 if (p[i].idx) { 228 /* Note if done on GPU, this can be very slow due to the huge number of kernel calls. The op is likely user-defined. We must 229 use link->unit (instead of link->basicunit) as the datatype and 1 (instead of link->bs) as the count in MPI_Reduce_local. 230 */ 231 for (j=0; j<p[i].count; j++) {ierr = MPI_Reduce_local(p[i].buf+j*link->unitbytes,(char *)rootdata+p[i].idx[j]*link->unitbytes,1,link->unit,op);CHKERRQ(ierr);} 232 } else { 233 ierr = PetscMPIIntCast(p[i].count,&n);CHKERRQ(ierr); 234 ierr = MPI_Reduce_local(p[i].buf,rootdata,n,link->unit,op);CHKERRQ(ierr); 235 } 236 } 237 #else 238 else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No unpacking reduction operation for this MPI_Op"); 239 #endif 240 } 241 } 242 #if defined(PETSC_HAVE_CUDA) 243 /* Make sure rootdata is ready to use by SF client */ 244 if (link->rootmtype == PETSC_MEMTYPE_DEVICE && p[1].count) {cudaError_t err = cudaStreamSynchronize(link->stream);CHKERRCUDA(err);} 245 #endif 246 PetscFunctionReturn(0); 247 } 248 249 /* Utility routine to unpack data from leaf buffer and Op it into selected entries of leafdata */ 250 PETSC_STATIC_INLINE PetscErrorCode PetscSFUnpackAndOpLeafData(PetscSF sf,PetscSFPack link,const PetscInt *leafloc,void *leafdata,MPI_Op op,PetscBool sparse) 251 { 252 PetscErrorCode ierr; 253 PetscInt i; 254 PetscErrorCode (*UnpackAndOp)(PetscInt,const PetscInt*,PetscSFPack,PetscSFPackOpt,void*,const void*); 255 PackInfo p[2]; 256 257 PetscFunctionBegin; 258 #if defined(PETSC_HAVE_CUDA) 259 if (!use_gpu_aware_mpi && link->leafmtype == PETSC_MEMTYPE_DEVICE) { 260 cudaError_t err; 261 if (!link->leafbuf[PETSC_MEMTYPE_DEVICE]) {ierr = PetscMallocWithMemType(PETSC_MEMTYPE_DEVICE,link->leafbuflen*link->unitbytes,(void**)&link->leafbuf[PETSC_MEMTYPE_DEVICE]);CHKERRQ(ierr);} 262 err = cudaMemcpyAsync(link->leafbuf[PETSC_MEMTYPE_DEVICE],link->leafbuf[PETSC_MEMTYPE_HOST],link->leafbuflen*link->unitbytes,cudaMemcpyHostToDevice,link->stream);CHKERRCUDA(err); 263 ierr = PetscLogCpuToGpu(link->leafbuflen*link->unitbytes);CHKERRQ(ierr); 264 } 265 #endif 266 267 if (sparse) { 268 p[0].count = link->selfbuflen; p[1].count = link->leafbuflen; 269 p[0].idx = leafloc; p[1].idx = leafloc+sf->roffset[sf->ndranks]; 270 p[0].opt = sf->selfleafpackopt; p[1].opt = sf->leafpackopt; 271 p[0].buf = link->selfbuf[link->leafmtype]; p[1].buf = link->leafbuf[link->leafmtype]; 272 p[0].atomic = sf->selfleafdups; p[1].atomic = sf->remoteleafdups; 273 } else { 274 p[0].count = 0; 275 p[1].count = sf->nleaves; 276 p[1].idx = NULL; 277 p[1].opt = NULL; 278 p[1].buf = link->leafbuf[link->leafmtype]; 279 p[1].atomic = PETSC_FALSE; 280 } 281 282 for (i=0; i<2; i++) { 283 if (p[i].count) { 284 ierr = PetscSFPackGetUnpackAndOp(link,link->leafmtype,op,p[i].atomic,&UnpackAndOp);CHKERRQ(ierr); 285 if (UnpackAndOp) {ierr = (*UnpackAndOp)(p[i].count,p[i].idx,link,p[i].opt,leafdata,p[i].buf);CHKERRQ(ierr);} 286 #if defined(PETSC_HAVE_MPI_REDUCE_LOCAL) 287 else { 288 PetscInt j; 289 PetscMPIInt n; 290 if (p[i].idx) { 291 for (j=0; j<p[i].count; j++) {ierr = MPI_Reduce_local(p[i].buf+j*link->unitbytes,(char *)leafdata+p[i].idx[j]*link->unitbytes,1,link->unit,op);CHKERRQ(ierr);} 292 } else { 293 ierr = PetscMPIIntCast(p[i].count,&n);CHKERRQ(ierr); 294 ierr = MPI_Reduce_local(p[i].buf,leafdata,n,link->unit,op);CHKERRQ(ierr); 295 } 296 } 297 #else 298 else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"No unpacking reduction operation for this MPI_Op"); 299 #endif 300 } 301 } 302 #if defined(PETSC_HAVE_CUDA) 303 if (link->leafmtype == PETSC_MEMTYPE_DEVICE && p[1].count) {cudaError_t err = cudaStreamSynchronize(link->stream);CHKERRCUDA(err);} 304 #endif 305 PetscFunctionReturn(0); 306 } 307 308 /* Utility routine to fetch and Op selected entries of rootdata */ 309 PETSC_STATIC_INLINE PetscErrorCode PetscSFFetchAndOpRootData(PetscSF sf,PetscSFPack link,const PetscInt *rootloc,void *rootdata,MPI_Op op,PetscBool sparse) 310 { 311 PetscErrorCode ierr; 312 PetscInt i; 313 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 314 PetscErrorCode (*FetchAndOp)(PetscInt,const PetscInt*,PetscSFPack,PetscSFPackOpt,void*,void*); 315 PackInfo p[2]; 316 317 PetscFunctionBegin; 318 #if defined(PETSC_HAVE_CUDA) 319 if (!use_gpu_aware_mpi && link->rootmtype == PETSC_MEMTYPE_DEVICE) { 320 cudaError_t err; 321 if (!link->rootbuf[PETSC_MEMTYPE_DEVICE]) {ierr = PetscMallocWithMemType(PETSC_MEMTYPE_DEVICE,link->rootbuflen*link->unitbytes,(void**)&link->rootbuf[PETSC_MEMTYPE_DEVICE]);CHKERRQ(ierr);} 322 err = cudaMemcpyAsync(link->rootbuf[PETSC_MEMTYPE_DEVICE],link->rootbuf[PETSC_MEMTYPE_HOST],link->rootbuflen*link->unitbytes,cudaMemcpyHostToDevice,link->stream);CHKERRCUDA(err); 323 ierr = PetscLogCpuToGpu(link->rootbuflen*link->unitbytes);CHKERRQ(ierr); 324 } 325 #endif 326 327 if (sparse) { 328 /* For SFBasic and SFNeighbor, whose root space is sparse and have separate buffers for self and remote. */ 329 p[0].count = link->selfbuflen; p[1].count = link->rootbuflen; 330 p[0].idx = rootloc; p[1].idx = rootloc+bas->ioffset[bas->ndiranks]; 331 p[0].opt = bas->selfrootpackopt; p[1].opt = bas->rootpackopt; 332 p[0].buf = link->selfbuf[link->rootmtype]; p[1].buf = link->rootbuf[link->rootmtype]; 333 p[0].atomic = bas->selfrootdups; p[1].atomic = bas->remoterootdups; 334 } else { 335 p[0].count = 0; 336 p[1].count = sf->nroots; 337 p[1].idx = NULL; 338 p[1].opt = NULL; 339 p[1].buf = link->rootbuf[link->rootmtype]; 340 p[1].atomic = PETSC_FALSE; 341 } 342 343 for (i=0; i<2; i++) { 344 if (p[i].count) { 345 ierr = PetscSFPackGetFetchAndOp(link,link->rootmtype,op,p[i].atomic,&FetchAndOp);CHKERRQ(ierr); 346 ierr = (*FetchAndOp)(p[i].count,p[i].idx,link,p[i].opt,rootdata,p[i].buf);CHKERRQ(ierr); 347 } 348 } 349 #if defined(PETSC_HAVE_CUDA) 350 if (link->rootmtype == PETSC_MEMTYPE_DEVICE && p[1].count) {cudaError_t err = cudaStreamSynchronize(link->stream);CHKERRCUDA(err);} 351 #endif 352 PetscFunctionReturn(0); 353 } 354 355 PETSC_STATIC_INLINE PetscErrorCode PetscSFPackSetupOptimizations_Basic(PetscSF sf) 356 { 357 PetscErrorCode ierr; 358 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 359 360 PetscFunctionBegin; 361 ierr = PetscSFPackOptCreate(sf->ndranks, sf->roffset, sf->rmine, &sf->selfleafpackopt);CHKERRQ(ierr); 362 ierr = PetscSFPackOptCreate(sf->nranks-sf->ndranks, sf->roffset+sf->ndranks, sf->rmine, &sf->leafpackopt);CHKERRQ(ierr); 363 ierr = PetscSFPackOptCreate(bas->ndiranks, bas->ioffset, bas->irootloc,&bas->selfrootpackopt);CHKERRQ(ierr); 364 ierr = PetscSFPackOptCreate(bas->niranks-bas->ndiranks,bas->ioffset+bas->ndiranks,bas->irootloc,&bas->rootpackopt);CHKERRQ(ierr); 365 366 #if defined(PETSC_HAVE_CUDA) 367 /* Check duplicates in irootloc[] so CUDA packing kernels can use cheaper regular operations 368 instead of atomics to unpack data on leaves/roots, when they know there is not data race. 369 */ 370 ierr = PetscCheckDupsInt(sf->roffset[sf->ndranks], sf->rmine, &sf->selfleafdups);CHKERRQ(ierr); 371 ierr = PetscCheckDupsInt(sf->roffset[sf->nranks]-sf->roffset[sf->ndranks], sf->rmine+sf->roffset[sf->ndranks], &sf->remoteleafdups);CHKERRQ(ierr); 372 ierr = PetscCheckDupsInt(bas->ioffset[bas->ndiranks], bas->irootloc, &bas->selfrootdups);CHKERRQ(ierr); 373 ierr = PetscCheckDupsInt(bas->ioffset[bas->niranks]-bas->ioffset[bas->ndiranks],bas->irootloc+bas->ioffset[bas->ndiranks],&bas->remoterootdups);CHKERRQ(ierr); 374 #endif 375 PetscFunctionReturn(0); 376 } 377 378 PETSC_STATIC_INLINE PetscErrorCode PetscSFPackDestroyOptimizations_Basic(PetscSF sf) 379 { 380 PetscErrorCode ierr; 381 PetscSF_Basic *bas = (PetscSF_Basic*)sf->data; 382 383 PetscFunctionBegin; 384 ierr = PetscSFPackOptDestroy(&sf->leafpackopt);CHKERRQ(ierr); 385 ierr = PetscSFPackOptDestroy(&sf->selfleafpackopt);CHKERRQ(ierr); 386 ierr = PetscSFPackOptDestroy(&bas->rootpackopt);CHKERRQ(ierr); 387 ierr = PetscSFPackOptDestroy(&bas->selfrootpackopt);CHKERRQ(ierr); 388 #if defined(PETSC_HAVE_CUDA) 389 sf->selfleafdups = PETSC_TRUE; 390 sf->remoteleafdups = PETSC_TRUE; 391 bas->selfrootdups = PETSC_TRUE; /* The default is assuming there are dups so that atomics are used. */ 392 bas->remoterootdups = PETSC_TRUE; 393 #endif 394 PetscFunctionReturn(0); 395 } 396 397 PETSC_INTERN PetscErrorCode PetscSFSetUp_Basic(PetscSF); 398 PETSC_INTERN PetscErrorCode PetscSFView_Basic(PetscSF,PetscViewer); 399 PETSC_INTERN PetscErrorCode PetscSFReset_Basic(PetscSF); 400 PETSC_INTERN PetscErrorCode PetscSFDestroy_Basic(PetscSF); 401 PETSC_INTERN PetscErrorCode PetscSFBcastAndOpEnd_Basic (PetscSF,MPI_Datatype,PetscMemType,const void*,PetscMemType, void*, MPI_Op); 402 PETSC_INTERN PetscErrorCode PetscSFReduceEnd_Basic (PetscSF,MPI_Datatype,PetscMemType,const void*,PetscMemType, void*, MPI_Op); 403 PETSC_INTERN PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF,MPI_Datatype,PetscMemType, void*,PetscMemType,const void*,void*,MPI_Op); 404 PETSC_INTERN PetscErrorCode PetscSFCreateEmbeddedSF_Basic(PetscSF,PetscInt,const PetscInt*,PetscSF*); 405 PETSC_INTERN PetscErrorCode PetscSFCreateEmbeddedLeafSF_Basic(PetscSF,PetscInt,const PetscInt*,PetscSF*); 406 PETSC_INTERN PetscErrorCode PetscSFGetLeafRanks_Basic(PetscSF,PetscInt*,const PetscMPIInt**,const PetscInt**,const PetscInt**); 407 PETSC_INTERN PetscErrorCode PetscSFPackGet_Basic_Common(PetscSF,MPI_Datatype,PetscMemType,const void*,PetscMemType,const void*,PetscInt,PetscInt,PetscSFPack*); 408 PETSC_INTERN PetscErrorCode PetscSFSetFromOptions_Basic(PetscOptionItems*,PetscSF); 409 #endif 410