xref: /petsc/src/vec/is/sf/impls/basic/sfbasic.c (revision 7fd2d3dbf8105d3f2e002a5ca11f019cd0ad7420)
1 
2 #include <../src/vec/is/sf/impls/basic/sfbasic.h>
3 #include <../src/vec/is/sf/impls/basic/sfpack.h>
4 
5 /*===================================================================================*/
6 /*              SF public interface implementations                                  */
7 /*===================================================================================*/
8 PETSC_INTERN PetscErrorCode PetscSFSetUp_Basic(PetscSF sf)
9 {
10   PetscErrorCode ierr;
11   PetscSF_Basic  *bas = (PetscSF_Basic*)sf->data;
12   PetscInt       *rlengths,*ilengths,i;
13   PetscMPIInt    rank,niranks,*iranks,tag;
14   MPI_Comm       comm;
15   MPI_Group      group;
16   MPI_Request    *rootreqs,*leafreqs;
17 
18   PetscFunctionBegin;
19   ierr = MPI_Comm_group(PETSC_COMM_SELF,&group);CHKERRQ(ierr);
20   ierr = PetscSFSetUpRanks(sf,group);CHKERRQ(ierr);
21   ierr = MPI_Group_free(&group);CHKERRQ(ierr);
22   ierr = PetscObjectGetComm((PetscObject)sf,&comm);CHKERRQ(ierr);
23   ierr = PetscObjectGetNewTag((PetscObject)sf,&tag);CHKERRQ(ierr);
24   ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr);
25   /*
26    * Inform roots about how many leaves and from which ranks
27    */
28   ierr = PetscMalloc1(sf->nranks,&rlengths);CHKERRQ(ierr);
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   ierr = PetscCommBuildTwoSided(comm,1,MPIU_INT,sf->nranks-sf->ndranks,sf->ranks+sf->ndranks,rlengths+sf->ndranks,&niranks,&iranks,(void**)&ilengths);CHKERRQ(ierr);
34 
35   /* Sort iranks. See use of VecScatterGetRemoteOrdered_Private() in MatGetBrowsOfAoCols_MPIAIJ() on why.
36      We could sort ranks there at the price of allocating extra working arrays. Presumably, niranks is
37      small and the sorting is cheap.
38    */
39   ierr = PetscSortMPIIntWithIntArray(niranks,iranks,ilengths);CHKERRQ(ierr);
40 
41   /* Partition into distinguished and non-distinguished incoming ranks */
42   bas->ndiranks = sf->ndranks;
43   bas->niranks = bas->ndiranks + niranks;
44   ierr = PetscMalloc2(bas->niranks,&bas->iranks,bas->niranks+1,&bas->ioffset);CHKERRQ(ierr);
45   bas->ioffset[0] = 0;
46   for (i=0; i<bas->ndiranks; i++) {
47     bas->iranks[i] = sf->ranks[i];
48     bas->ioffset[i+1] = bas->ioffset[i] + rlengths[i];
49   }
50   if (bas->ndiranks > 1 || (bas->ndiranks == 1 && bas->iranks[0] != rank)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Broken setup for shared ranks");
51   for (; i<bas->niranks; i++) {
52     bas->iranks[i] = iranks[i-bas->ndiranks];
53     bas->ioffset[i+1] = bas->ioffset[i] + ilengths[i-bas->ndiranks];
54   }
55   bas->itotal = bas->ioffset[i];
56   ierr = PetscFree(rlengths);CHKERRQ(ierr);
57   ierr = PetscFree(iranks);CHKERRQ(ierr);
58   ierr = PetscFree(ilengths);CHKERRQ(ierr);
59 
60   /* Send leaf identities to roots */
61   ierr = PetscMalloc1(bas->itotal,&bas->irootloc);CHKERRQ(ierr);
62   ierr = PetscMalloc2(bas->niranks-bas->ndiranks,&rootreqs,sf->nranks-sf->ndranks,&leafreqs);CHKERRQ(ierr);
63   for (i=bas->ndiranks; i<bas->niranks; i++) {
64     ierr = 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]);CHKERRQ(ierr);
65   }
66   for (i=0; i<sf->nranks; i++) {
67     PetscMPIInt npoints;
68     ierr = PetscMPIIntCast(sf->roffset[i+1] - sf->roffset[i],&npoints);CHKERRQ(ierr);
69     if (i < sf->ndranks) {
70       if (sf->ranks[i] != rank) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Cannot interpret distinguished leaf rank");
71       if (bas->iranks[0] != rank) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Cannot interpret distinguished root rank");
72       if (npoints != bas->ioffset[1]-bas->ioffset[0]) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Distinguished rank exchange has mismatched lengths");
73       ierr = PetscArraycpy(bas->irootloc+bas->ioffset[0],sf->rremote+sf->roffset[i],npoints);CHKERRQ(ierr);
74       continue;
75     }
76     ierr = MPI_Isend(sf->rremote+sf->roffset[i],npoints,MPIU_INT,sf->ranks[i],tag,comm,&leafreqs[i-sf->ndranks]);CHKERRQ(ierr);
77   }
78   ierr = MPI_Waitall(bas->niranks-bas->ndiranks,rootreqs,MPI_STATUSES_IGNORE);CHKERRQ(ierr);
79   ierr = MPI_Waitall(sf->nranks-sf->ndranks,leafreqs,MPI_STATUSES_IGNORE);CHKERRQ(ierr);
80   ierr = PetscFree2(rootreqs,leafreqs);CHKERRQ(ierr);
81 
82   sf->nleafreqs  = sf->nranks - sf->ndranks;
83   bas->nrootreqs = bas->niranks - bas->ndiranks;
84   sf->persistent = PETSC_TRUE;
85 
86   /* Setup fields related to packing */
87   ierr = PetscSFSetUpPackFields(sf);CHKERRQ(ierr);
88   PetscFunctionReturn(0);
89 }
90 
91 PETSC_INTERN PetscErrorCode PetscSFReset_Basic(PetscSF sf)
92 {
93   PetscErrorCode    ierr;
94   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
95 
96   PetscFunctionBegin;
97   if (bas->inuse) SETERRQ(PetscObjectComm((PetscObject)sf),PETSC_ERR_ARG_WRONGSTATE,"Outstanding operation has not been completed");
98   ierr = PetscFree2(bas->iranks,bas->ioffset);CHKERRQ(ierr);
99   ierr = PetscFree(bas->irootloc);CHKERRQ(ierr);
100 #if defined(PETSC_HAVE_DEVICE)
101   for (PetscInt i=0; i<2; i++) {ierr = PetscSFFree(PETSC_MEMTYPE_DEVICE,bas->irootloc_d[i]);CHKERRQ(ierr);}
102 #endif
103   ierr = PetscSFLinkDestroy(sf,&bas->avail);CHKERRQ(ierr);
104   ierr = PetscSFResetPackFields(sf);CHKERRQ(ierr);
105   PetscFunctionReturn(0);
106 }
107 
108 PETSC_INTERN PetscErrorCode PetscSFDestroy_Basic(PetscSF sf)
109 {
110   PetscErrorCode ierr;
111 
112   PetscFunctionBegin;
113   ierr = PetscSFReset_Basic(sf);CHKERRQ(ierr);
114   ierr = PetscFree(sf->data);CHKERRQ(ierr);
115   PetscFunctionReturn(0);
116 }
117 
118 PETSC_INTERN PetscErrorCode PetscSFView_Basic(PetscSF sf,PetscViewer viewer)
119 {
120   PetscErrorCode ierr;
121   PetscBool      iascii;
122 
123   PetscFunctionBegin;
124   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr);
125   if (iascii) {ierr = PetscViewerASCIIPrintf(viewer,"  sort=%s\n",sf->rankorder ? "rank-order" : "unordered");CHKERRQ(ierr);}
126   PetscFunctionReturn(0);
127 }
128 
129 static PetscErrorCode PetscSFBcastAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,const void *rootdata,PetscMemType leafmtype,void *leafdata,MPI_Op op)
130 {
131   PetscErrorCode    ierr;
132   PetscSFLink       link = NULL;
133   MPI_Request       *rootreqs = NULL,*leafreqs = NULL;
134   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
135 
136   PetscFunctionBegin;
137   /* Create a communication link, which provides buffers & MPI requests etc */
138   ierr = PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,PETSCSF_BCAST,&link);CHKERRQ(ierr);
139   /* Get MPI requests from the link. We do not need buffers explicitly since we use persistent MPI */
140   ierr = PetscSFLinkGetMPIBuffersAndRequests(sf,link,PETSCSF_ROOT2LEAF,NULL,NULL,&rootreqs,&leafreqs);CHKERRQ(ierr);
141   /* Post Irecv for remote */
142   ierr = MPI_Startall_irecv(sf->leafbuflen[PETSCSF_REMOTE],unit,sf->nleafreqs,leafreqs);CHKERRQ(ierr);
143   /* Pack rootdata and do Isend for remote */
144   ierr = PetscSFLinkPackRootData(sf,link,PETSCSF_REMOTE,rootdata);CHKERRQ(ierr);
145   ierr = MPI_Startall_isend(bas->rootbuflen[PETSCSF_REMOTE],unit,bas->nrootreqs,rootreqs);CHKERRQ(ierr);
146   /* Do local BcastAndOp, which overlaps with the irecv/isend above */
147   ierr = PetscSFLinkBcastAndOpLocal(sf,link,rootdata,leafdata,op);CHKERRQ(ierr);
148   PetscFunctionReturn(0);
149 }
150 
151 PETSC_INTERN PetscErrorCode PetscSFBcastAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata,MPI_Op op)
152 {
153   PetscErrorCode    ierr;
154   PetscSFLink       link = NULL;
155 
156   PetscFunctionBegin;
157   /* Retrieve the link used in XxxBegin() with root/leafdata as key */
158   ierr = PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);CHKERRQ(ierr);
159   /* Wait for the completion of mpi */
160   ierr = PetscSFLinkMPIWaitall(sf,link,PETSCSF_ROOT2LEAF);CHKERRQ(ierr);
161   /* Unpack leafdata and reclaim the link */
162   ierr = PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafdata,op);CHKERRQ(ierr);
163   ierr = PetscSFLinkReclaim(sf,&link);CHKERRQ(ierr);
164   PetscFunctionReturn(0);
165 }
166 
167 /* Shared by ReduceBegin and FetchAndOpBegin */
168 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)
169 {
170   PetscErrorCode    ierr;
171   PetscSFLink       link;
172   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
173   MPI_Request       *rootreqs = NULL,*leafreqs = NULL;
174 
175   PetscFunctionBegin;
176   ierr = PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,sfop,&link);CHKERRQ(ierr);
177   ierr = PetscSFLinkGetMPIBuffersAndRequests(sf,link,PETSCSF_LEAF2ROOT,NULL,NULL,&rootreqs,&leafreqs);CHKERRQ(ierr);
178   ierr = MPI_Startall_irecv(bas->rootbuflen[PETSCSF_REMOTE],unit,bas->nrootreqs,rootreqs);CHKERRQ(ierr);
179   ierr = PetscSFLinkPackLeafData(sf,link,PETSCSF_REMOTE,leafdata);CHKERRQ(ierr);
180   ierr = MPI_Startall_isend(sf->leafbuflen[PETSCSF_REMOTE],unit,sf->nleafreqs,leafreqs);CHKERRQ(ierr);
181   *out = link;
182   PetscFunctionReturn(0);
183 }
184 
185 /* leaf -> root with reduction */
186 static PetscErrorCode PetscSFReduceBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType leafmtype,const void *leafdata,PetscMemType rootmtype,void *rootdata,MPI_Op op)
187 {
188   PetscErrorCode    ierr;
189   PetscSFLink       link = NULL;
190 
191   PetscFunctionBegin;
192   ierr = PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_REDUCE,&link);CHKERRQ(ierr);
193   ierr = PetscSFLinkReduceLocal(sf,link,leafdata,rootdata,op);CHKERRQ(ierr);
194   PetscFunctionReturn(0);
195 }
196 
197 PETSC_INTERN PetscErrorCode PetscSFReduceEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
198 {
199   PetscErrorCode    ierr;
200   PetscSFLink       link = NULL;
201 
202   PetscFunctionBegin;
203   ierr = PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);CHKERRQ(ierr);
204   ierr = PetscSFLinkMPIWaitall(sf,link,PETSCSF_LEAF2ROOT);CHKERRQ(ierr);
205   ierr = PetscSFLinkUnpackRootData(sf,link,PETSCSF_REMOTE,rootdata,op);CHKERRQ(ierr);
206   ierr = PetscSFLinkReclaim(sf,&link);CHKERRQ(ierr);
207   PetscFunctionReturn(0);
208 }
209 
210 PETSC_INTERN PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,void *rootdata,PetscMemType leafmtype,const void *leafdata,void *leafupdate,MPI_Op op)
211 {
212   PetscErrorCode    ierr;
213   PetscSFLink       link = NULL;
214 
215   PetscFunctionBegin;
216   ierr = PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_FETCH,&link);CHKERRQ(ierr);
217   ierr = PetscSFLinkFetchAndOpLocal(sf,link,rootdata,leafdata,leafupdate,op);CHKERRQ(ierr);
218   PetscFunctionReturn(0);
219 }
220 
221 static PetscErrorCode PetscSFFetchAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
222 {
223   PetscErrorCode    ierr;
224   PetscSFLink       link = NULL;
225   MPI_Request       *rootreqs = NULL,*leafreqs = NULL;
226   PetscSF_Basic     *bas = (PetscSF_Basic*)sf->data;
227 
228   PetscFunctionBegin;
229   ierr = PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);CHKERRQ(ierr);
230   /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
231   ierr = PetscSFLinkMPIWaitall(sf,link,PETSCSF_LEAF2ROOT);CHKERRQ(ierr);
232   /* Do fetch-and-op, the (remote) update results are in rootbuf */
233   ierr = PetscSFLinkFetchRootData(sf,link,PETSCSF_REMOTE,rootdata,op);CHKERRQ(ierr);
234 
235   /* Bcast rootbuf to leafupdate */
236   ierr = PetscSFLinkGetMPIBuffersAndRequests(sf,link,PETSCSF_ROOT2LEAF,NULL,NULL,&rootreqs,&leafreqs);CHKERRQ(ierr);
237   /* Post leaf receives and root sends */
238   ierr = MPI_Startall_irecv(sf->leafbuflen[PETSCSF_REMOTE],unit,sf->nleafreqs,leafreqs);CHKERRQ(ierr);
239   ierr = MPI_Startall_isend(bas->rootbuflen[PETSCSF_REMOTE],unit,bas->nrootreqs,rootreqs);CHKERRQ(ierr);
240   /* Unpack and insert fetched data into leaves */
241   ierr = PetscSFLinkMPIWaitall(sf,link,PETSCSF_ROOT2LEAF);CHKERRQ(ierr);
242   ierr = PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafupdate,MPIU_REPLACE);CHKERRQ(ierr);
243   ierr = PetscSFLinkReclaim(sf,&link);CHKERRQ(ierr);
244   PetscFunctionReturn(0);
245 }
246 
247 PETSC_INTERN PetscErrorCode PetscSFGetLeafRanks_Basic(PetscSF sf,PetscInt *niranks,const PetscMPIInt **iranks,const PetscInt **ioffset,const PetscInt **irootloc)
248 {
249   PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
250 
251   PetscFunctionBegin;
252   if (niranks)  *niranks  = bas->niranks;
253   if (iranks)   *iranks   = bas->iranks;
254   if (ioffset)  *ioffset  = bas->ioffset;
255   if (irootloc) *irootloc = bas->irootloc;
256   PetscFunctionReturn(0);
257 }
258 
259 /* An optimized PetscSFCreateEmbeddedSF. We aggresively make use of the established communication on sf.
260    We need one bcast on sf, and no communication anymore to build the embedded sf. Note that selected[]
261    was sorted before calling the routine.
262  */
263 PETSC_INTERN PetscErrorCode PetscSFCreateEmbeddedSF_Basic(PetscSF sf,PetscInt nselected,const PetscInt *selected,PetscSF *newsf)
264 {
265   PetscSF           esf;
266   PetscInt          esf_nranks,esf_ndranks,*esf_roffset,*esf_rmine,*esf_rremote;
267   PetscInt          i,j,p,q,nroots,esf_nleaves,*new_ilocal,nranks,ndranks,niranks,ndiranks,minleaf,maxleaf,maxlocal;
268   char              *rootdata,*leafdata,*leafmem; /* Only stores 0 or 1, so we can save memory with char */
269   PetscMPIInt       *esf_ranks;
270   const PetscMPIInt *ranks,*iranks;
271   const PetscInt    *roffset,*rmine,*rremote,*ioffset,*irootloc;
272   PetscBool         connected;
273   PetscSFNode       *new_iremote;
274   PetscSF_Basic     *bas;
275   PetscErrorCode    ierr;
276 
277   PetscFunctionBegin;
278   ierr = PetscSFCreate(PetscObjectComm((PetscObject)sf),&esf);CHKERRQ(ierr);
279   ierr = PetscSFSetType(esf,PETSCSFBASIC);CHKERRQ(ierr); /* This optimized routine can only create a basic sf */
280 
281   /* Find out which leaves are still connected to roots in the embedded sf by doing a Bcast */
282   ierr = PetscSFGetGraph(sf,&nroots,NULL,NULL,NULL);CHKERRQ(ierr);
283   ierr = PetscSFGetLeafRange(sf,&minleaf,&maxleaf);CHKERRQ(ierr);
284   maxlocal = maxleaf - minleaf + 1;
285   ierr = PetscCalloc2(nroots,&rootdata,maxlocal,&leafmem);CHKERRQ(ierr);
286   leafdata = leafmem - minleaf;
287   /* Tag selected roots */
288   for (i=0; i<nselected; ++i) rootdata[selected[i]] = 1;
289 
290   ierr = PetscSFBcastBegin(sf,MPI_CHAR,rootdata,leafdata);CHKERRQ(ierr);
291   ierr = PetscSFBcastEnd(sf,MPI_CHAR,rootdata,leafdata);CHKERRQ(ierr);
292   ierr = PetscSFGetLeafInfo_Basic(sf,&nranks,&ndranks,&ranks,&roffset,&rmine,&rremote);CHKERRQ(ierr); /* Get send info */
293   esf_nranks = esf_ndranks = esf_nleaves = 0;
294   for (i=0; i<nranks; i++) {
295     connected = PETSC_FALSE; /* Is this process still connected to this remote root rank? */
296     for (j=roffset[i]; j<roffset[i+1]; j++) {if (leafdata[rmine[j]]) {esf_nleaves++; connected = PETSC_TRUE;}}
297     if (connected) {esf_nranks++; if (i < ndranks) esf_ndranks++;}
298   }
299 
300   /* Set graph of esf and also set up its outgoing communication (i.e., send info), which is usually done by PetscSFSetUpRanks */
301   ierr = PetscMalloc1(esf_nleaves,&new_ilocal);CHKERRQ(ierr);
302   ierr = PetscMalloc1(esf_nleaves,&new_iremote);CHKERRQ(ierr);
303   ierr = PetscMalloc4(esf_nranks,&esf_ranks,esf_nranks+1,&esf_roffset,esf_nleaves,&esf_rmine,esf_nleaves,&esf_rremote);CHKERRQ(ierr);
304   p    = 0; /* Counter for connected root ranks */
305   q    = 0; /* Counter for connected leaves */
306   esf_roffset[0] = 0;
307   for (i=0; i<nranks; i++) { /* Scan leaf data again to fill esf arrays */
308     connected = PETSC_FALSE;
309     for (j=roffset[i]; j<roffset[i+1]; j++) {
310       if (leafdata[rmine[j]]) {
311         esf_rmine[q]         = new_ilocal[q] = rmine[j];
312         esf_rremote[q]       = rremote[j];
313         new_iremote[q].index = rremote[j];
314         new_iremote[q].rank  = ranks[i];
315         connected            = PETSC_TRUE;
316         q++;
317       }
318     }
319     if (connected) {
320       esf_ranks[p]     = ranks[i];
321       esf_roffset[p+1] = q;
322       p++;
323     }
324   }
325 
326   /* SetGraph internally resets the SF, so we only set its fields after the call */
327   ierr           = PetscSFSetGraph(esf,nroots,esf_nleaves,new_ilocal,PETSC_OWN_POINTER,new_iremote,PETSC_OWN_POINTER);CHKERRQ(ierr);
328   esf->nranks    = esf_nranks;
329   esf->ndranks   = esf_ndranks;
330   esf->ranks     = esf_ranks;
331   esf->roffset   = esf_roffset;
332   esf->rmine     = esf_rmine;
333   esf->rremote   = esf_rremote;
334   esf->nleafreqs = esf_nranks - esf_ndranks;
335 
336   /* Set up the incoming communication (i.e., recv info) stored in esf->data, which is usually done by PetscSFSetUp_Basic */
337   bas  = (PetscSF_Basic*)esf->data;
338   ierr = PetscSFGetRootInfo_Basic(sf,&niranks,&ndiranks,&iranks,&ioffset,&irootloc);CHKERRQ(ierr); /* Get recv info */
339   /* Embedded sf always has simpler communication than the original one. We might allocate longer arrays than needed here. But we
340      we do not care since these arrays are usually short. The benefit is we can fill these arrays by just parsing irootloc once.
341    */
342   ierr = PetscMalloc2(niranks,&bas->iranks,niranks+1,&bas->ioffset);CHKERRQ(ierr);
343   ierr = PetscMalloc1(ioffset[niranks],&bas->irootloc);CHKERRQ(ierr);
344   bas->niranks = bas->ndiranks = bas->ioffset[0] = 0;
345   p = 0; /* Counter for connected leaf ranks */
346   q = 0; /* Counter for connected roots */
347   for (i=0; i<niranks; i++) {
348     connected = PETSC_FALSE; /* Is the current process still connected to this remote leaf rank? */
349     for (j=ioffset[i]; j<ioffset[i+1]; j++) {
350       if (rootdata[irootloc[j]]) {
351         bas->irootloc[q++] = irootloc[j];
352         connected = PETSC_TRUE;
353       }
354     }
355     if (connected) {
356       bas->niranks++;
357       if (i<ndiranks) bas->ndiranks++; /* Note that order of ranks (including distinguished ranks) is kept */
358       bas->iranks[p]    = iranks[i];
359       bas->ioffset[p+1] = q;
360       p++;
361     }
362   }
363   bas->itotal     = q;
364   bas->nrootreqs  = bas->niranks - bas->ndiranks;
365   esf->persistent = PETSC_TRUE;
366   /* Setup packing related fields */
367   ierr = PetscSFSetUpPackFields(esf);CHKERRQ(ierr);
368 
369   esf->setupcalled = PETSC_TRUE; /* We have done setup ourselves! */
370   ierr = PetscFree2(rootdata,leafmem);CHKERRQ(ierr);
371   *newsf = esf;
372   PetscFunctionReturn(0);
373 }
374 
375 PETSC_EXTERN PetscErrorCode PetscSFCreate_Basic(PetscSF sf)
376 {
377   PetscSF_Basic  *dat;
378   PetscErrorCode ierr;
379 
380   PetscFunctionBegin;
381   sf->ops->SetUp                = PetscSFSetUp_Basic;
382   sf->ops->Reset                = PetscSFReset_Basic;
383   sf->ops->Destroy              = PetscSFDestroy_Basic;
384   sf->ops->View                 = PetscSFView_Basic;
385   sf->ops->BcastAndOpBegin      = PetscSFBcastAndOpBegin_Basic;
386   sf->ops->BcastAndOpEnd        = PetscSFBcastAndOpEnd_Basic;
387   sf->ops->ReduceBegin          = PetscSFReduceBegin_Basic;
388   sf->ops->ReduceEnd            = PetscSFReduceEnd_Basic;
389   sf->ops->FetchAndOpBegin      = PetscSFFetchAndOpBegin_Basic;
390   sf->ops->FetchAndOpEnd        = PetscSFFetchAndOpEnd_Basic;
391   sf->ops->GetLeafRanks         = PetscSFGetLeafRanks_Basic;
392   sf->ops->CreateEmbeddedSF     = PetscSFCreateEmbeddedSF_Basic;
393 
394   ierr = PetscNewLog(sf,&dat);CHKERRQ(ierr);
395   sf->data = (void*)dat;
396   PetscFunctionReturn(0);
397 }
398