Actual source code: sfbasic.c

  1: #include <../src/vec/is/sf/impls/basic/sfbasic.h>
  2: #include <../src/vec/is/sf/impls/basic/sfpack.h>
  3: #include <petsc/private/viewerimpl.h>

  5: /*===================================================================================*/
  6: /*              SF public interface implementations                                  */
  7: /*===================================================================================*/
  8: PETSC_INTERN PetscErrorCode PetscSFSetUp_Basic(PetscSF sf)
  9: {
 11:   PetscSF_Basic  *bas = (PetscSF_Basic*)sf->data;
 12:   PetscInt       *rlengths,*ilengths,i,nRemoteRootRanks,nRemoteLeafRanks;
 13:   PetscMPIInt    rank,niranks,*iranks,tag;
 14:   MPI_Comm       comm;
 15:   MPI_Group      group;
 16:   MPI_Request    *rootreqs,*leafreqs;

 19:   MPI_Comm_group(PETSC_COMM_SELF,&group);
 20:   PetscSFSetUpRanks(sf,group);
 21:   MPI_Group_free(&group);
 22:   PetscObjectGetComm((PetscObject)sf,&comm);
 23:   PetscObjectGetNewTag((PetscObject)sf,&tag);
 24:   MPI_Comm_rank(comm,&rank);
 25:   /*
 26:    * Inform roots about how many leaves and from which ranks
 27:    */
 28:   PetscMalloc1(sf->nranks,&rlengths);
 29:   /* Determine number, sending ranks and length of incoming */
 30:   for (i=0; i<sf->nranks; i++) {
 31:     rlengths[i] = sf->roffset[i+1] - sf->roffset[i]; /* Number of roots referenced by my leaves; for rank sf->ranks[i] */
 32:   }
 33:   nRemoteRootRanks = sf->nranks-sf->ndranks;
 34:   PetscCommBuildTwoSided(comm,1,MPIU_INT,nRemoteRootRanks,sf->ranks+sf->ndranks,rlengths+sf->ndranks,&niranks,&iranks,(void**)&ilengths);

 36:   /* Sort iranks. See use of VecScatterGetRemoteOrdered_Private() in MatGetBrowsOfAoCols_MPIAIJ() on why.
 37:      We could sort ranks there at the price of allocating extra working arrays. Presumably, niranks is
 38:      small and the sorting is cheap.
 39:    */
 40:   PetscSortMPIIntWithIntArray(niranks,iranks,ilengths);

 42:   /* Partition into distinguished and non-distinguished incoming ranks */
 43:   bas->ndiranks = sf->ndranks;
 44:   bas->niranks = bas->ndiranks + niranks;
 45:   PetscMalloc2(bas->niranks,&bas->iranks,bas->niranks+1,&bas->ioffset);
 46:   bas->ioffset[0] = 0;
 47:   for (i=0; i<bas->ndiranks; i++) {
 48:     bas->iranks[i] = sf->ranks[i];
 49:     bas->ioffset[i+1] = bas->ioffset[i] + rlengths[i];
 50:   }
 51:   if (bas->ndiranks > 1 || (bas->ndiranks == 1 && bas->iranks[0] != rank)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Broken setup for shared ranks");
 52:   for (; i<bas->niranks; i++) {
 53:     bas->iranks[i] = iranks[i-bas->ndiranks];
 54:     bas->ioffset[i+1] = bas->ioffset[i] + ilengths[i-bas->ndiranks];
 55:   }
 56:   bas->itotal = bas->ioffset[i];
 57:   PetscFree(rlengths);
 58:   PetscFree(iranks);
 59:   PetscFree(ilengths);

 61:   /* Send leaf identities to roots */
 62:   nRemoteLeafRanks = bas->niranks-bas->ndiranks;
 63:   PetscMalloc1(bas->itotal,&bas->irootloc);
 64:   PetscMalloc2(nRemoteLeafRanks,&rootreqs,nRemoteRootRanks,&leafreqs);
 65:   for (i=bas->ndiranks; i<bas->niranks; i++) {
 66:     MPI_Irecv(bas->irootloc+bas->ioffset[i],bas->ioffset[i+1]-bas->ioffset[i],MPIU_INT,bas->iranks[i],tag,comm,&rootreqs[i-bas->ndiranks]);
 67:   }
 68:   for (i=0; i<sf->nranks; i++) {
 69:     PetscMPIInt npoints;
 70:     PetscMPIIntCast(sf->roffset[i+1] - sf->roffset[i],&npoints);
 71:     if (i < sf->ndranks) {
 72:       if (sf->ranks[i] != rank) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Cannot interpret distinguished leaf rank");
 73:       if (bas->iranks[0] != rank) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Cannot interpret distinguished root rank");
 74:       if (npoints != bas->ioffset[1]-bas->ioffset[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Distinguished rank exchange has mismatched lengths");
 75:       PetscArraycpy(bas->irootloc+bas->ioffset[0],sf->rremote+sf->roffset[i],npoints);
 76:       continue;
 77:     }
 78:     MPI_Isend(sf->rremote+sf->roffset[i],npoints,MPIU_INT,sf->ranks[i],tag,comm,&leafreqs[i-sf->ndranks]);
 79:   }
 80:   MPI_Waitall(nRemoteLeafRanks,rootreqs,MPI_STATUSES_IGNORE);
 81:   MPI_Waitall(nRemoteRootRanks,leafreqs,MPI_STATUSES_IGNORE);

 83:   sf->nleafreqs  = nRemoteRootRanks;
 84:   bas->nrootreqs = nRemoteLeafRanks;
 85:   sf->persistent = PETSC_TRUE;

 87:   /* Setup fields related to packing, such as rootbuflen[] */
 88:   PetscSFSetUpPackFields(sf);
 89:   PetscFree2(rootreqs,leafreqs);
 90:   return(0);
 91: }

 93: PETSC_INTERN PetscErrorCode PetscSFReset_Basic(PetscSF sf)
 94: {
 95:   PetscErrorCode    ierr;
 96:   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
 97:   PetscSFLink       link = bas->avail,next;

100:   if (bas->inuse) SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Outstanding operation has not been completed");
101:   PetscFree2(bas->iranks,bas->ioffset);
102:   PetscFree(bas->irootloc);

104:  #if defined(PETSC_HAVE_DEVICE)
105:   for (PetscInt i=0; i<2; i++) {PetscSFFree(sf,PETSC_MEMTYPE_DEVICE,bas->irootloc_d[i]);}
106:  #endif

108:  #if defined(PETSC_HAVE_NVSHMEM)
109:   PetscSFReset_Basic_NVSHMEM(sf);
110:  #endif

112:   for (; link; link=next) {next = link->next; PetscSFLinkDestroy(sf,link);}
113:   bas->avail = NULL;
114:   PetscSFResetPackFields(sf);
115:   return(0);
116: }

118: PETSC_INTERN PetscErrorCode PetscSFDestroy_Basic(PetscSF sf)
119: {

123:   PetscSFReset_Basic(sf);
124:   PetscFree(sf->data);
125:   return(0);
126: }

128: #if defined(PETSC_USE_SINGLE_LIBRARY)
129: #include <petscmat.h>

131: PETSC_INTERN PetscErrorCode PetscSFView_Basic_PatternAndSizes(PetscSF sf,PetscViewer viewer)
132: {
133:   PetscErrorCode       ierr;
134:   PetscSF_Basic        *bas = (PetscSF_Basic*)sf->data;
135:   PetscInt             i,nrootranks,ndrootranks;
136:   const PetscInt       *rootoffset;
137:   PetscMPIInt          rank,size;
138:   const PetscMPIInt    *rootranks;
139:   MPI_Comm             comm = PetscObjectComm((PetscObject)sf);
140:   PetscScalar          unitbytes;
141:   Mat                  A;

144:   MPI_Comm_size(comm,&size);
145:   MPI_Comm_rank(comm,&rank);
146:   /* PetscSFView is most useful for the SF used in VecScatterBegin/End in MatMult etc, where we do
147:     PetscSFBcast, i.e., roots send data to leaves.  We dump the communication pattern into a matrix
148:     in senders' view point: how many bytes I will send to my neighbors.

150:     Looking at a column of the matrix, one can also know how many bytes the rank will receive from others.

152:     If PetscSFLink bas->inuse is available, we can use that to get tree vertex size. But that would give
153:     different interpretations for the same SF for different data types. Since we most care about VecScatter,
154:     we uniformly treat each vertex as a PetscScalar.
155:   */
156:   unitbytes = (PetscScalar)sizeof(PetscScalar);

158:   PetscSFGetRootInfo_Basic(sf,&nrootranks,&ndrootranks,&rootranks,&rootoffset,NULL);
159:   MatCreateAIJ(comm,1,1,size,size,1,NULL,nrootranks-ndrootranks,NULL,&A);
160:   MatSetOptionsPrefix(A,"__petsc_internal__"); /* To prevent the internal A from taking any command line options */
161:   for (i=0; i<nrootranks; i++) {
162:     MatSetValue(A,(PetscInt)rank,bas->iranks[i],(rootoffset[i+1]-rootoffset[i])*unitbytes,INSERT_VALUES);
163:   }
164:   MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
165:   MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
166:   MatView(A,viewer);
167:   MatDestroy(&A);
168:   return(0);
169: }
170: #endif

172: PETSC_INTERN PetscErrorCode PetscSFView_Basic(PetscSF sf,PetscViewer viewer)
173: {
175:   PetscBool      isascii;

178:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
179:   if (isascii && viewer->format != PETSC_VIEWER_ASCII_MATLAB) {PetscViewerASCIIPrintf(viewer,"  MultiSF sort=%s\n",sf->rankorder ? "rank-order" : "unordered");}
180: #if defined(PETSC_USE_SINGLE_LIBRARY)
181:   else {
182:     PetscBool  isdraw,isbinary;
183:     PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
184:     PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
185:     if ((isascii && viewer->format == PETSC_VIEWER_ASCII_MATLAB) || isdraw || isbinary) {
186:       PetscSFView_Basic_PatternAndSizes(sf,viewer);
187:     }
188:   }
189: #endif
190:   return(0);
191: }

193: static PetscErrorCode PetscSFBcastBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,const void *rootdata,PetscMemType leafmtype,void *leafdata,MPI_Op op)
194: {
195:   PetscErrorCode    ierr;
196:   PetscSFLink       link = NULL;

199:   /* Create a communication link, which provides buffers, MPI requests etc (if MPI is used) */
200:   PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,PETSCSF_BCAST,&link);
201:   /* Pack rootdata to rootbuf for remote communication */
202:   PetscSFLinkPackRootData(sf,link,PETSCSF_REMOTE,rootdata);
203:   /* Start communcation, e.g., post MPI_Isend */
204:   PetscSFLinkStartCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
205:   /* Do local scatter (i.e., self to self communication), which overlaps with the remote communication above */
206:   PetscSFLinkScatterLocal(sf,link,PETSCSF_../../../../../..2LEAF,(void*)rootdata,leafdata,op);
207:   return(0);
208: }

210: PETSC_INTERN PetscErrorCode PetscSFBcastEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata,MPI_Op op)
211: {
212:   PetscErrorCode    ierr;
213:   PetscSFLink       link = NULL;

216:   /* Retrieve the link used in XxxBegin() with root/leafdata as key */
217:   PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
218:   /* Finish remote communication, e.g., post MPI_Waitall */
219:   PetscSFLinkFinishCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
220:   /* Unpack data in leafbuf to leafdata for remote communication */
221:   PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafdata,op);
222:   /* Recycle the link */
223:   PetscSFLinkReclaim(sf,&link);
224:   return(0);
225: }

227: /* Shared by ReduceBegin and FetchAndOpBegin */
228: PETSC_STATIC_INLINE PetscErrorCode PetscSFLeafToRootBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType leafmtype,const void *leafdata,PetscMemType rootmtype,void *rootdata,MPI_Op op,PetscSFOperation sfop,PetscSFLink *out)
229: {
230:   PetscErrorCode    ierr;
231:   PetscSFLink       link = NULL;

234:   PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,sfop,&link);
235:   PetscSFLinkPackLeafData(sf,link,PETSCSF_REMOTE,leafdata);
236:   PetscSFLinkStartCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
237:   *out = link;
238:   return(0);
239: }

241: /* leaf -> root with reduction */
242: static PetscErrorCode PetscSFReduceBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType leafmtype,const void *leafdata,PetscMemType rootmtype,void *rootdata,MPI_Op op)
243: {
244:   PetscErrorCode    ierr;
245:   PetscSFLink       link = NULL;

248:   PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_REDUCE,&link);
249:   PetscSFLinkScatterLocal(sf,link,PETSCSF_LEAF2../../../../../..,rootdata,(void*)leafdata,op);
250:   return(0);
251: }

253: PETSC_INTERN PetscErrorCode PetscSFReduceEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
254: {
255:   PetscErrorCode    ierr;
256:   PetscSFLink       link = NULL;

259:   PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
260:   PetscSFLinkFinishCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
261:   PetscSFLinkUnpackRootData(sf,link,PETSCSF_REMOTE,rootdata,op);
262:   PetscSFLinkReclaim(sf,&link);
263:   return(0);
264: }

266: PETSC_INTERN PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,void *rootdata,PetscMemType leafmtype,const void *leafdata,void *leafupdate,MPI_Op op)
267: {
268:   PetscErrorCode    ierr;
269:   PetscSFLink       link = NULL;

272:   PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_FETCH,&link);
273:   PetscSFLinkFetchAndOpLocal(sf,link,rootdata,leafdata,leafupdate,op);
274:   return(0);
275: }

277: static PetscErrorCode PetscSFFetchAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
278: {
279:   PetscErrorCode    ierr;
280:   PetscSFLink       link = NULL;

283:   PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
284:   /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
285:   PetscSFLinkFinishCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
286:   /* Do fetch-and-op, the (remote) update results are in rootbuf */
287:   PetscSFLinkFetchAndOpRemote(sf,link,rootdata,op);
288:   /* Bcast rootbuf to leafupdate */
289:   PetscSFLinkStartCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
290:   PetscSFLinkFinishCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
291:   /* Unpack and insert fetched data into leaves */
292:   PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafupdate,MPI_REPLACE);
293:   PetscSFLinkReclaim(sf,&link);
294:   return(0);
295: }

297: PETSC_INTERN PetscErrorCode PetscSFGetLeafRanks_Basic(PetscSF sf,PetscInt *niranks,const PetscMPIInt **iranks,const PetscInt **ioffset,const PetscInt **irootloc)
298: {
299:   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;

302:   if (niranks)  *niranks  = bas->niranks;
303:   if (iranks)   *iranks   = bas->iranks;
304:   if (ioffset)  *ioffset  = bas->ioffset;
305:   if (irootloc) *irootloc = bas->irootloc;
306:   return(0);
307: }

309: /* An optimized PetscSFCreateEmbeddedRootSF. We aggresively make use of the established communication on sf.
310:    We need one bcast on sf, and no communication anymore to build the embedded sf. Note that selected[]
311:    was sorted before calling the routine.
312:  */
313: PETSC_INTERN PetscErrorCode PetscSFCreateEmbeddedRootSF_Basic(PetscSF sf,PetscInt nselected,const PetscInt *selected,PetscSF *newsf)
314: {
315:   PetscSF           esf;
316:   PetscInt          esf_nranks,esf_ndranks,*esf_roffset,*esf_rmine,*esf_rremote;
317:   PetscInt          i,j,p,q,nroots,esf_nleaves,*new_ilocal,nranks,ndranks,niranks,ndiranks,minleaf,maxleaf,maxlocal;
318:   char              *rootdata,*leafdata,*leafmem; /* Only stores 0 or 1, so we can save memory with char */
319:   PetscMPIInt       *esf_ranks;
320:   const PetscMPIInt *ranks,*iranks;
321:   const PetscInt    *roffset,*rmine,*rremote,*ioffset,*irootloc;
322:   PetscBool         connected;
323:   PetscSFNode       *new_iremote;
324:   PetscSF_Basic     *bas;
325:   PetscErrorCode    ierr;

328:   PetscSFCreate(PetscObjectComm((PetscObject)sf),&esf);
329:   PetscSFSetFromOptions(esf);
330:   PetscSFSetType(esf,PETSCSFBASIC); /* This optimized routine can only create a basic sf */

332:   /* Find out which leaves are still connected to roots in the embedded sf by doing a Bcast */
333:   PetscSFGetGraph(sf,&nroots,NULL,NULL,NULL);
334:   PetscSFGetLeafRange(sf,&minleaf,&maxleaf);
335:   maxlocal = maxleaf - minleaf + 1;
336:   PetscCalloc2(nroots,&rootdata,maxlocal,&leafmem);
337:   leafdata = leafmem - minleaf;
338:   /* Tag selected roots */
339:   for (i=0; i<nselected; ++i) rootdata[selected[i]] = 1;

341:   PetscSFBcastBegin(sf,MPI_CHAR,rootdata,leafdata,MPI_REPLACE);
342:   PetscSFBcastEnd(sf,MPI_CHAR,rootdata,leafdata,MPI_REPLACE);
343:   PetscSFGetLeafInfo_Basic(sf,&nranks,&ndranks,&ranks,&roffset,&rmine,&rremote); /* Get send info */
344:   esf_nranks = esf_ndranks = esf_nleaves = 0;
345:   for (i=0; i<nranks; i++) {
346:     connected = PETSC_FALSE; /* Is this process still connected to this remote root rank? */
347:     for (j=roffset[i]; j<roffset[i+1]; j++) {if (leafdata[rmine[j]]) {esf_nleaves++; connected = PETSC_TRUE;}}
348:     if (connected) {esf_nranks++; if (i < ndranks) esf_ndranks++;}
349:   }

351:   /* Set graph of esf and also set up its outgoing communication (i.e., send info), which is usually done by PetscSFSetUpRanks */
352:   PetscMalloc1(esf_nleaves,&new_ilocal);
353:   PetscMalloc1(esf_nleaves,&new_iremote);
354:   PetscMalloc4(esf_nranks,&esf_ranks,esf_nranks+1,&esf_roffset,esf_nleaves,&esf_rmine,esf_nleaves,&esf_rremote);
355:   p    = 0; /* Counter for connected root ranks */
356:   q    = 0; /* Counter for connected leaves */
357:   esf_roffset[0] = 0;
358:   for (i=0; i<nranks; i++) { /* Scan leaf data again to fill esf arrays */
359:     connected = PETSC_FALSE;
360:     for (j=roffset[i]; j<roffset[i+1]; j++) {
361:       if (leafdata[rmine[j]]) {
362:         esf_rmine[q]         = new_ilocal[q] = rmine[j];
363:         esf_rremote[q]       = rremote[j];
364:         new_iremote[q].index = rremote[j];
365:         new_iremote[q].rank  = ranks[i];
366:         connected            = PETSC_TRUE;
367:         q++;
368:       }
369:     }
370:     if (connected) {
371:       esf_ranks[p]     = ranks[i];
372:       esf_roffset[p+1] = q;
373:       p++;
374:     }
375:   }

377:   /* SetGraph internally resets the SF, so we only set its fields after the call */
378:   PetscSFSetGraph(esf,nroots,esf_nleaves,new_ilocal,PETSC_OWN_POINTER,new_iremote,PETSC_OWN_POINTER);
379:   esf->nranks    = esf_nranks;
380:   esf->ndranks   = esf_ndranks;
381:   esf->ranks     = esf_ranks;
382:   esf->roffset   = esf_roffset;
383:   esf->rmine     = esf_rmine;
384:   esf->rremote   = esf_rremote;
385:   esf->nleafreqs = esf_nranks - esf_ndranks;

387:   /* Set up the incoming communication (i.e., recv info) stored in esf->data, which is usually done by PetscSFSetUp_Basic */
388:   bas  = (PetscSF_Basic*)esf->data;
389:   PetscSFGetRootInfo_Basic(sf,&niranks,&ndiranks,&iranks,&ioffset,&irootloc); /* Get recv info */
390:   /* Embedded sf always has simpler communication than the original one. We might allocate longer arrays than needed here. But we
391:      we do not care since these arrays are usually short. The benefit is we can fill these arrays by just parsing irootloc once.
392:    */
393:   PetscMalloc2(niranks,&bas->iranks,niranks+1,&bas->ioffset);
394:   PetscMalloc1(ioffset[niranks],&bas->irootloc);
395:   bas->niranks = bas->ndiranks = bas->ioffset[0] = 0;
396:   p = 0; /* Counter for connected leaf ranks */
397:   q = 0; /* Counter for connected roots */
398:   for (i=0; i<niranks; i++) {
399:     connected = PETSC_FALSE; /* Is the current process still connected to this remote leaf rank? */
400:     for (j=ioffset[i]; j<ioffset[i+1]; j++) {
401:       if (rootdata[irootloc[j]]) {
402:         bas->irootloc[q++] = irootloc[j];
403:         connected = PETSC_TRUE;
404:       }
405:     }
406:     if (connected) {
407:       bas->niranks++;
408:       if (i<ndiranks) bas->ndiranks++; /* Note that order of ranks (including distinguished ranks) is kept */
409:       bas->iranks[p]    = iranks[i];
410:       bas->ioffset[p+1] = q;
411:       p++;
412:     }
413:   }
414:   bas->itotal     = q;
415:   bas->nrootreqs  = bas->niranks - bas->ndiranks;
416:   esf->persistent = PETSC_TRUE;
417:   /* Setup packing related fields */
418:   PetscSFSetUpPackFields(esf);

420:   /* Copy from PetscSFSetUp(), since this method wants to skip PetscSFSetUp(). */
421: #if defined(PETSC_HAVE_CUDA)
422:   if (esf->backend == PETSCSF_BACKEND_CUDA) {
423:     esf->ops->Malloc = PetscSFMalloc_CUDA;
424:     esf->ops->Free   = PetscSFFree_CUDA;
425:   }
426: #endif

428: #if defined(PETSC_HAVE_HIP)
429:   /* TODO: Needs debugging */
430:   if (esf->backend == PETSCSF_BACKEND_HIP) {
431:     esf->ops->Malloc = PetscSFMalloc_HIP;
432:     esf->ops->Free   = PetscSFFree_HIP;
433:   }
434: #endif

436: #if defined(PETSC_HAVE_KOKKOS)
437:   if (esf->backend == PETSCSF_BACKEND_KOKKOS) {
438:     esf->ops->Malloc = PetscSFMalloc_Kokkos;
439:     esf->ops->Free   = PetscSFFree_Kokkos;
440:   }
441: #endif
442:   esf->setupcalled = PETSC_TRUE; /* We have done setup ourselves! */
443:   PetscFree2(rootdata,leafmem);
444:   *newsf = esf;
445:   return(0);
446: }

448: PETSC_EXTERN PetscErrorCode PetscSFCreate_Basic(PetscSF sf)
449: {
450:   PetscSF_Basic  *dat;

454:   sf->ops->SetUp                = PetscSFSetUp_Basic;
455:   sf->ops->Reset                = PetscSFReset_Basic;
456:   sf->ops->Destroy              = PetscSFDestroy_Basic;
457:   sf->ops->View                 = PetscSFView_Basic;
458:   sf->ops->BcastBegin           = PetscSFBcastBegin_Basic;
459:   sf->ops->BcastEnd             = PetscSFBcastEnd_Basic;
460:   sf->ops->ReduceBegin          = PetscSFReduceBegin_Basic;
461:   sf->ops->ReduceEnd            = PetscSFReduceEnd_Basic;
462:   sf->ops->FetchAndOpBegin      = PetscSFFetchAndOpBegin_Basic;
463:   sf->ops->FetchAndOpEnd        = PetscSFFetchAndOpEnd_Basic;
464:   sf->ops->GetLeafRanks         = PetscSFGetLeafRanks_Basic;
465:   sf->ops->CreateEmbeddedRootSF = PetscSFCreateEmbeddedRootSF_Basic;

467:   PetscNewLog(sf,&dat);
468:   sf->data = (void*)dat;
469:   return(0);
470: }