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: {
10: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
11: PetscInt *rlengths,*ilengths,i,nRemoteRootRanks,nRemoteLeafRanks;
12: PetscMPIInt rank,niranks,*iranks,tag;
13: MPI_Comm comm;
14: MPI_Group group;
15: MPI_Request *rootreqs,*leafreqs;
17: MPI_Comm_group(PETSC_COMM_SELF,&group);
18: PetscSFSetUpRanks(sf,group);
19: MPI_Group_free(&group);
20: PetscObjectGetComm((PetscObject)sf,&comm);
21: PetscObjectGetNewTag((PetscObject)sf,&tag);
22: MPI_Comm_rank(comm,&rank);
23: /*
24: * Inform roots about how many leaves and from which ranks
25: */
26: PetscMalloc1(sf->nranks,&rlengths);
27: /* Determine number, sending ranks and length of incoming */
28: for (i=0; i<sf->nranks; i++) {
29: rlengths[i] = sf->roffset[i+1] - sf->roffset[i]; /* Number of roots referenced by my leaves; for rank sf->ranks[i] */
30: }
31: nRemoteRootRanks = sf->nranks-sf->ndranks;
32: PetscCommBuildTwoSided(comm,1,MPIU_INT,nRemoteRootRanks,sf->ranks+sf->ndranks,rlengths+sf->ndranks,&niranks,&iranks,(void**)&ilengths);
34: /* Sort iranks. See use of VecScatterGetRemoteOrdered_Private() in MatGetBrowsOfAoCols_MPIAIJ() on why.
35: We could sort ranks there at the price of allocating extra working arrays. Presumably, niranks is
36: small and the sorting is cheap.
37: */
38: PetscSortMPIIntWithIntArray(niranks,iranks,ilengths);
40: /* Partition into distinguished and non-distinguished incoming ranks */
41: bas->ndiranks = sf->ndranks;
42: bas->niranks = bas->ndiranks + niranks;
43: PetscMalloc2(bas->niranks,&bas->iranks,bas->niranks+1,&bas->ioffset);
44: bas->ioffset[0] = 0;
45: for (i=0; i<bas->ndiranks; i++) {
46: bas->iranks[i] = sf->ranks[i];
47: bas->ioffset[i+1] = bas->ioffset[i] + rlengths[i];
48: }
50: for (; i<bas->niranks; i++) {
51: bas->iranks[i] = iranks[i-bas->ndiranks];
52: bas->ioffset[i+1] = bas->ioffset[i] + ilengths[i-bas->ndiranks];
53: }
54: bas->itotal = bas->ioffset[i];
55: PetscFree(rlengths);
56: PetscFree(iranks);
57: PetscFree(ilengths);
59: /* Send leaf identities to roots */
60: nRemoteLeafRanks = bas->niranks-bas->ndiranks;
61: PetscMalloc1(bas->itotal,&bas->irootloc);
62: PetscMalloc2(nRemoteLeafRanks,&rootreqs,nRemoteRootRanks,&leafreqs);
63: for (i=bas->ndiranks; i<bas->niranks; i++) {
64: MPIU_Irecv(bas->irootloc+bas->ioffset[i],bas->ioffset[i+1]-bas->ioffset[i],MPIU_INT,bas->iranks[i],tag,comm,&rootreqs[i-bas->ndiranks]);
65: }
66: for (i=0; i<sf->nranks; i++) {
67: PetscInt npoints = sf->roffset[i+1] - sf->roffset[i];
68: if (i < sf->ndranks) {
72: PetscArraycpy(bas->irootloc+bas->ioffset[0],sf->rremote+sf->roffset[i],npoints);
73: continue;
74: }
75: MPIU_Isend(sf->rremote+sf->roffset[i],npoints,MPIU_INT,sf->ranks[i],tag,comm,&leafreqs[i-sf->ndranks]);
76: }
77: MPI_Waitall(nRemoteLeafRanks,rootreqs,MPI_STATUSES_IGNORE);
78: MPI_Waitall(nRemoteRootRanks,leafreqs,MPI_STATUSES_IGNORE);
80: sf->nleafreqs = nRemoteRootRanks;
81: bas->nrootreqs = nRemoteLeafRanks;
82: sf->persistent = PETSC_TRUE;
84: /* Setup fields related to packing, such as rootbuflen[] */
85: PetscSFSetUpPackFields(sf);
86: PetscFree2(rootreqs,leafreqs);
87: return 0;
88: }
90: PETSC_INTERN PetscErrorCode PetscSFReset_Basic(PetscSF sf)
91: {
92: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
93: PetscSFLink link = bas->avail,next;
96: PetscFree2(bas->iranks,bas->ioffset);
97: PetscFree(bas->irootloc);
99: #if defined(PETSC_HAVE_DEVICE)
100: for (PetscInt i=0; i<2; i++) PetscSFFree(sf,PETSC_MEMTYPE_DEVICE,bas->irootloc_d[i]);
101: #endif
103: #if defined(PETSC_HAVE_NVSHMEM)
104: PetscSFReset_Basic_NVSHMEM(sf);
105: #endif
107: for (; link; link=next) {next = link->next; PetscSFLinkDestroy(sf,link);}
108: bas->avail = NULL;
109: PetscSFResetPackFields(sf);
110: return 0;
111: }
113: PETSC_INTERN PetscErrorCode PetscSFDestroy_Basic(PetscSF sf)
114: {
115: PetscSFReset_Basic(sf);
116: PetscFree(sf->data);
117: return 0;
118: }
120: #if defined(PETSC_USE_SINGLE_LIBRARY)
121: #include <petscmat.h>
123: PETSC_INTERN PetscErrorCode PetscSFView_Basic_PatternAndSizes(PetscSF sf,PetscViewer viewer)
124: {
125: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
126: PetscInt i,nrootranks,ndrootranks;
127: const PetscInt *rootoffset;
128: PetscMPIInt rank,size;
129: const PetscMPIInt *rootranks;
130: MPI_Comm comm = PetscObjectComm((PetscObject)sf);
131: PetscScalar unitbytes;
132: Mat A;
134: MPI_Comm_size(comm,&size);
135: MPI_Comm_rank(comm,&rank);
136: /* PetscSFView is most useful for the SF used in VecScatterBegin/End in MatMult etc, where we do
137: PetscSFBcast, i.e., roots send data to leaves. We dump the communication pattern into a matrix
138: in senders' view point: how many bytes I will send to my neighbors.
140: Looking at a column of the matrix, one can also know how many bytes the rank will receive from others.
142: If PetscSFLink bas->inuse is available, we can use that to get tree vertex size. But that would give
143: different interpretations for the same SF for different data types. Since we most care about VecScatter,
144: we uniformly treat each vertex as a PetscScalar.
145: */
146: unitbytes = (PetscScalar)sizeof(PetscScalar);
148: PetscSFGetRootInfo_Basic(sf,&nrootranks,&ndrootranks,&rootranks,&rootoffset,NULL);
149: MatCreateAIJ(comm,1,1,size,size,1,NULL,nrootranks-ndrootranks,NULL,&A);
150: MatSetOptionsPrefix(A,"__petsc_internal__"); /* To prevent the internal A from taking any command line options */
151: for (i=0; i<nrootranks; i++) {
152: MatSetValue(A,(PetscInt)rank,bas->iranks[i],(rootoffset[i+1]-rootoffset[i])*unitbytes,INSERT_VALUES);
153: }
154: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
155: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
156: MatView(A,viewer);
157: MatDestroy(&A);
158: return 0;
159: }
160: #endif
162: PETSC_INTERN PetscErrorCode PetscSFView_Basic(PetscSF sf,PetscViewer viewer)
163: {
164: PetscBool isascii;
166: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
167: if (isascii && viewer->format != PETSC_VIEWER_ASCII_MATLAB) PetscViewerASCIIPrintf(viewer," MultiSF sort=%s\n",sf->rankorder ? "rank-order" : "unordered");
168: #if defined(PETSC_USE_SINGLE_LIBRARY)
169: else {
170: PetscBool isdraw,isbinary;
171: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
172: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
173: if ((isascii && viewer->format == PETSC_VIEWER_ASCII_MATLAB) || isdraw || isbinary) {
174: PetscSFView_Basic_PatternAndSizes(sf,viewer);
175: }
176: }
177: #endif
178: return 0;
179: }
181: static PetscErrorCode PetscSFBcastBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,const void *rootdata,PetscMemType leafmtype,void *leafdata,MPI_Op op)
182: {
183: PetscSFLink link = NULL;
185: /* Create a communication link, which provides buffers, MPI requests etc (if MPI is used) */
186: PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,PETSCSF_BCAST,&link);
187: /* Pack rootdata to rootbuf for remote communication */
188: PetscSFLinkPackRootData(sf,link,PETSCSF_REMOTE,rootdata);
189: /* Start communcation, e.g., post MPI_Isend */
190: PetscSFLinkStartCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
191: /* Do local scatter (i.e., self to self communication), which overlaps with the remote communication above */
192: PetscSFLinkScatterLocal(sf,link,PETSCSF_../../../../../..2LEAF,(void*)rootdata,leafdata,op);
193: return 0;
194: }
196: PETSC_INTERN PetscErrorCode PetscSFBcastEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *rootdata,void *leafdata,MPI_Op op)
197: {
198: PetscSFLink link = NULL;
200: /* Retrieve the link used in XxxBegin() with root/leafdata as key */
201: PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
202: /* Finish remote communication, e.g., post MPI_Waitall */
203: PetscSFLinkFinishCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
204: /* Unpack data in leafbuf to leafdata for remote communication */
205: PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafdata,op);
206: /* Recycle the link */
207: PetscSFLinkReclaim(sf,&link);
208: return 0;
209: }
211: /* Shared by ReduceBegin and FetchAndOpBegin */
212: 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)
213: {
214: PetscSFLink link = NULL;
216: PetscSFLinkCreate(sf,unit,rootmtype,rootdata,leafmtype,leafdata,op,sfop,&link);
217: PetscSFLinkPackLeafData(sf,link,PETSCSF_REMOTE,leafdata);
218: PetscSFLinkStartCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
219: *out = link;
220: return 0;
221: }
223: /* leaf -> root with reduction */
224: static PetscErrorCode PetscSFReduceBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType leafmtype,const void *leafdata,PetscMemType rootmtype,void *rootdata,MPI_Op op)
225: {
226: PetscSFLink link = NULL;
228: PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_REDUCE,&link);
229: PetscSFLinkScatterLocal(sf,link,PETSCSF_LEAF2../../../../../..,rootdata,(void*)leafdata,op);
230: return 0;
231: }
233: PETSC_INTERN PetscErrorCode PetscSFReduceEnd_Basic(PetscSF sf,MPI_Datatype unit,const void *leafdata,void *rootdata,MPI_Op op)
234: {
235: PetscSFLink link = NULL;
237: PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
238: PetscSFLinkFinishCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
239: PetscSFLinkUnpackRootData(sf,link,PETSCSF_REMOTE,rootdata,op);
240: PetscSFLinkReclaim(sf,&link);
241: return 0;
242: }
244: PETSC_INTERN PetscErrorCode PetscSFFetchAndOpBegin_Basic(PetscSF sf,MPI_Datatype unit,PetscMemType rootmtype,void *rootdata,PetscMemType leafmtype,const void *leafdata,void *leafupdate,MPI_Op op)
245: {
246: PetscSFLink link = NULL;
248: PetscSFLeafToRootBegin_Basic(sf,unit,leafmtype,leafdata,rootmtype,rootdata,op,PETSCSF_FETCH,&link);
249: PetscSFLinkFetchAndOpLocal(sf,link,rootdata,leafdata,leafupdate,op);
250: return 0;
251: }
253: static PetscErrorCode PetscSFFetchAndOpEnd_Basic(PetscSF sf,MPI_Datatype unit,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op op)
254: {
255: PetscSFLink link = NULL;
257: PetscSFLinkGetInUse(sf,unit,rootdata,leafdata,PETSC_OWN_POINTER,&link);
258: /* This implementation could be changed to unpack as receives arrive, at the cost of non-determinism */
259: PetscSFLinkFinishCommunication(sf,link,PETSCSF_LEAF2../../../../../..);
260: /* Do fetch-and-op, the (remote) update results are in rootbuf */
261: PetscSFLinkFetchAndOpRemote(sf,link,rootdata,op);
262: /* Bcast rootbuf to leafupdate */
263: PetscSFLinkStartCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
264: PetscSFLinkFinishCommunication(sf,link,PETSCSF_../../../../../..2LEAF);
265: /* Unpack and insert fetched data into leaves */
266: PetscSFLinkUnpackLeafData(sf,link,PETSCSF_REMOTE,leafupdate,MPI_REPLACE);
267: PetscSFLinkReclaim(sf,&link);
268: return 0;
269: }
271: PETSC_INTERN PetscErrorCode PetscSFGetLeafRanks_Basic(PetscSF sf,PetscInt *niranks,const PetscMPIInt **iranks,const PetscInt **ioffset,const PetscInt **irootloc)
272: {
273: PetscSF_Basic *bas = (PetscSF_Basic*)sf->data;
275: if (niranks) *niranks = bas->niranks;
276: if (iranks) *iranks = bas->iranks;
277: if (ioffset) *ioffset = bas->ioffset;
278: if (irootloc) *irootloc = bas->irootloc;
279: return 0;
280: }
282: /* An optimized PetscSFCreateEmbeddedRootSF. We aggresively make use of the established communication on sf.
283: We need one bcast on sf, and no communication anymore to build the embedded sf. Note that selected[]
284: was sorted before calling the routine.
285: */
286: PETSC_INTERN PetscErrorCode PetscSFCreateEmbeddedRootSF_Basic(PetscSF sf,PetscInt nselected,const PetscInt *selected,PetscSF *newsf)
287: {
288: PetscSF esf;
289: PetscInt esf_nranks,esf_ndranks,*esf_roffset,*esf_rmine,*esf_rremote;
290: PetscInt i,j,p,q,nroots,esf_nleaves,*new_ilocal,nranks,ndranks,niranks,ndiranks,minleaf,maxleaf,maxlocal;
291: char *rootdata,*leafdata,*leafmem; /* Only stores 0 or 1, so we can save memory with char */
292: PetscMPIInt *esf_ranks;
293: const PetscMPIInt *ranks,*iranks;
294: const PetscInt *roffset,*rmine,*rremote,*ioffset,*irootloc;
295: PetscBool connected;
296: PetscSFNode *new_iremote;
297: PetscSF_Basic *bas;
299: PetscSFCreate(PetscObjectComm((PetscObject)sf),&esf);
300: PetscSFSetFromOptions(esf);
301: PetscSFSetType(esf,PETSCSFBASIC); /* This optimized routine can only create a basic sf */
303: /* Find out which leaves are still connected to roots in the embedded sf by doing a Bcast */
304: PetscSFGetGraph(sf,&nroots,NULL,NULL,NULL);
305: PetscSFGetLeafRange(sf,&minleaf,&maxleaf);
306: maxlocal = maxleaf - minleaf + 1;
307: PetscCalloc2(nroots,&rootdata,maxlocal,&leafmem);
308: leafdata = leafmem - minleaf;
309: /* Tag selected roots */
310: for (i=0; i<nselected; ++i) rootdata[selected[i]] = 1;
312: PetscSFBcastBegin(sf,MPI_CHAR,rootdata,leafdata,MPI_REPLACE);
313: PetscSFBcastEnd(sf,MPI_CHAR,rootdata,leafdata,MPI_REPLACE);
314: PetscSFGetLeafInfo_Basic(sf,&nranks,&ndranks,&ranks,&roffset,&rmine,&rremote); /* Get send info */
315: esf_nranks = esf_ndranks = esf_nleaves = 0;
316: for (i=0; i<nranks; i++) {
317: connected = PETSC_FALSE; /* Is this process still connected to this remote root rank? */
318: for (j=roffset[i]; j<roffset[i+1]; j++) {if (leafdata[rmine[j]]) {esf_nleaves++; connected = PETSC_TRUE;}}
319: if (connected) {esf_nranks++; if (i < ndranks) esf_ndranks++;}
320: }
322: /* Set graph of esf and also set up its outgoing communication (i.e., send info), which is usually done by PetscSFSetUpRanks */
323: PetscMalloc1(esf_nleaves,&new_ilocal);
324: PetscMalloc1(esf_nleaves,&new_iremote);
325: PetscMalloc4(esf_nranks,&esf_ranks,esf_nranks+1,&esf_roffset,esf_nleaves,&esf_rmine,esf_nleaves,&esf_rremote);
326: p = 0; /* Counter for connected root ranks */
327: q = 0; /* Counter for connected leaves */
328: esf_roffset[0] = 0;
329: for (i=0; i<nranks; i++) { /* Scan leaf data again to fill esf arrays */
330: connected = PETSC_FALSE;
331: for (j=roffset[i]; j<roffset[i+1]; j++) {
332: if (leafdata[rmine[j]]) {
333: esf_rmine[q] = new_ilocal[q] = rmine[j];
334: esf_rremote[q] = rremote[j];
335: new_iremote[q].index = rremote[j];
336: new_iremote[q].rank = ranks[i];
337: connected = PETSC_TRUE;
338: q++;
339: }
340: }
341: if (connected) {
342: esf_ranks[p] = ranks[i];
343: esf_roffset[p+1] = q;
344: p++;
345: }
346: }
348: /* SetGraph internally resets the SF, so we only set its fields after the call */
349: PetscSFSetGraph(esf,nroots,esf_nleaves,new_ilocal,PETSC_OWN_POINTER,new_iremote,PETSC_OWN_POINTER);
350: esf->nranks = esf_nranks;
351: esf->ndranks = esf_ndranks;
352: esf->ranks = esf_ranks;
353: esf->roffset = esf_roffset;
354: esf->rmine = esf_rmine;
355: esf->rremote = esf_rremote;
356: esf->nleafreqs = esf_nranks - esf_ndranks;
358: /* Set up the incoming communication (i.e., recv info) stored in esf->data, which is usually done by PetscSFSetUp_Basic */
359: bas = (PetscSF_Basic*)esf->data;
360: PetscSFGetRootInfo_Basic(sf,&niranks,&ndiranks,&iranks,&ioffset,&irootloc); /* Get recv info */
361: /* Embedded sf always has simpler communication than the original one. We might allocate longer arrays than needed here. But we
362: we do not care since these arrays are usually short. The benefit is we can fill these arrays by just parsing irootloc once.
363: */
364: PetscMalloc2(niranks,&bas->iranks,niranks+1,&bas->ioffset);
365: PetscMalloc1(ioffset[niranks],&bas->irootloc);
366: bas->niranks = bas->ndiranks = bas->ioffset[0] = 0;
367: p = 0; /* Counter for connected leaf ranks */
368: q = 0; /* Counter for connected roots */
369: for (i=0; i<niranks; i++) {
370: connected = PETSC_FALSE; /* Is the current process still connected to this remote leaf rank? */
371: for (j=ioffset[i]; j<ioffset[i+1]; j++) {
372: if (rootdata[irootloc[j]]) {
373: bas->irootloc[q++] = irootloc[j];
374: connected = PETSC_TRUE;
375: }
376: }
377: if (connected) {
378: bas->niranks++;
379: if (i<ndiranks) bas->ndiranks++; /* Note that order of ranks (including distinguished ranks) is kept */
380: bas->iranks[p] = iranks[i];
381: bas->ioffset[p+1] = q;
382: p++;
383: }
384: }
385: bas->itotal = q;
386: bas->nrootreqs = bas->niranks - bas->ndiranks;
387: esf->persistent = PETSC_TRUE;
388: /* Setup packing related fields */
389: PetscSFSetUpPackFields(esf);
391: /* Copy from PetscSFSetUp(), since this method wants to skip PetscSFSetUp(). */
392: #if defined(PETSC_HAVE_CUDA)
393: if (esf->backend == PETSCSF_BACKEND_CUDA) {
394: esf->ops->Malloc = PetscSFMalloc_CUDA;
395: esf->ops->Free = PetscSFFree_CUDA;
396: }
397: #endif
399: #if defined(PETSC_HAVE_HIP)
400: /* TODO: Needs debugging */
401: if (esf->backend == PETSCSF_BACKEND_HIP) {
402: esf->ops->Malloc = PetscSFMalloc_HIP;
403: esf->ops->Free = PetscSFFree_HIP;
404: }
405: #endif
407: #if defined(PETSC_HAVE_KOKKOS)
408: if (esf->backend == PETSCSF_BACKEND_KOKKOS) {
409: esf->ops->Malloc = PetscSFMalloc_Kokkos;
410: esf->ops->Free = PetscSFFree_Kokkos;
411: }
412: #endif
413: esf->setupcalled = PETSC_TRUE; /* We have done setup ourselves! */
414: PetscFree2(rootdata,leafmem);
415: *newsf = esf;
416: return 0;
417: }
419: PETSC_EXTERN PetscErrorCode PetscSFCreate_Basic(PetscSF sf)
420: {
421: PetscSF_Basic *dat;
423: sf->ops->SetUp = PetscSFSetUp_Basic;
424: sf->ops->Reset = PetscSFReset_Basic;
425: sf->ops->Destroy = PetscSFDestroy_Basic;
426: sf->ops->View = PetscSFView_Basic;
427: sf->ops->BcastBegin = PetscSFBcastBegin_Basic;
428: sf->ops->BcastEnd = PetscSFBcastEnd_Basic;
429: sf->ops->ReduceBegin = PetscSFReduceBegin_Basic;
430: sf->ops->ReduceEnd = PetscSFReduceEnd_Basic;
431: sf->ops->FetchAndOpBegin = PetscSFFetchAndOpBegin_Basic;
432: sf->ops->FetchAndOpEnd = PetscSFFetchAndOpEnd_Basic;
433: sf->ops->GetLeafRanks = PetscSFGetLeafRanks_Basic;
434: sf->ops->CreateEmbeddedRootSF = PetscSFCreateEmbeddedRootSF_Basic;
436: PetscNewLog(sf,&dat);
437: sf->data = (void*)dat;
438: return 0;
439: }