Actual source code: subcomm.c


  2: /*
  3:      Provides utility routines for split MPI communicator.
  4: */
  5: #include <petscsys.h>
  6: #include <petscviewer.h>

  8: const char *const PetscSubcommTypes[] = {"GENERAL","CONTIGUOUS","INTERLACED","PetscSubcommType","PETSC_SUBCOMM_",NULL};

 10: static PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm);
 11: static PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm);

 13: /*@
 14:    PetscSubcommSetFromOptions - Allows setting options from a PetscSubcomm

 16:    Collective on PetscSubcomm

 18:    Input Parameter:
 19: .  psubcomm - PetscSubcomm context

 21:    Level: beginner

 23: @*/
 24: PetscErrorCode PetscSubcommSetFromOptions(PetscSubcomm psubcomm)
 25: {
 26:   PetscErrorCode   ierr;
 27:   PetscSubcommType type;
 28:   PetscBool        flg;


 32:   PetscOptionsBegin(psubcomm->parent,psubcomm->subcommprefix,"Options for PetscSubcomm",NULL);
 33:   PetscOptionsEnum("-psubcomm_type",NULL,NULL,PetscSubcommTypes,(PetscEnum)psubcomm->type,(PetscEnum*)&type,&flg);
 34:   if (flg && psubcomm->type != type) {
 35:     /* free old structures */
 36:     PetscCommDestroy(&(psubcomm)->dupparent);
 37:     PetscCommDestroy(&(psubcomm)->child);
 38:     PetscFree((psubcomm)->subsize);
 39:     switch (type) {
 40:     case PETSC_SUBCOMM_GENERAL:
 41:       SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Runtime option PETSC_SUBCOMM_GENERAL is not supported, use PetscSubcommSetTypeGeneral()");
 42:     case PETSC_SUBCOMM_CONTIGUOUS:
 43:       PetscSubcommCreate_contiguous(psubcomm);
 44:       break;
 45:     case PETSC_SUBCOMM_INTERLACED:
 46:       PetscSubcommCreate_interlaced(psubcomm);
 47:       break;
 48:     default:
 49:       SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"PetscSubcommType %s is not supported yet",PetscSubcommTypes[type]);
 50:     }
 51:   }

 53:   PetscOptionsName("-psubcomm_view","Triggers display of PetscSubcomm context","PetscSubcommView",&flg);
 54:   if (flg) {
 55:     PetscSubcommView(psubcomm,PETSC_VIEWER_STDOUT_(psubcomm->parent));
 56:   }
 57:   PetscOptionsEnd();
 58:   return 0;
 59: }

 61: /*@C
 62:   PetscSubcommSetOptionsPrefix - Sets the prefix used for searching for all
 63:   PetscSubcomm items in the options database.

 65:   Logically collective on PetscSubcomm.

 67:   Level: Intermediate

 69:   Input Parameters:
 70: +   psubcomm - PetscSubcomm context
 71: -   prefix - the prefix to prepend all PetscSubcomm item names with.

 73: @*/
 74: PetscErrorCode PetscSubcommSetOptionsPrefix(PetscSubcomm psubcomm,const char pre[])
 75: {
 76:    if (!pre) {
 77:     PetscFree(psubcomm->subcommprefix);
 78:   } else {
 80:     PetscFree(psubcomm->subcommprefix);
 81:     PetscStrallocpy(pre,&(psubcomm->subcommprefix));
 82:   }
 83:   return 0;
 84: }

 86: /*@C
 87:    PetscSubcommView - Views a PetscSubcomm of values as either ASCII text or a binary file

 89:    Collective on PetscSubcomm

 91:    Input Parameters:
 92: +  psubcomm - PetscSubcomm context
 93: -  viewer - location to view the values

 95:    Level: beginner
 96: @*/
 97: PetscErrorCode PetscSubcommView(PetscSubcomm psubcomm,PetscViewer viewer)
 98: {
 99:   PetscBool         iascii;
100:   PetscViewerFormat format;

102:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
103:   if (iascii) {
104:     PetscViewerGetFormat(viewer,&format);
105:     if (format == PETSC_VIEWER_DEFAULT) {
106:       MPI_Comm    comm=psubcomm->parent;
107:       PetscMPIInt rank,size,subsize,subrank,duprank;

109:       MPI_Comm_size(comm,&size);
110:       PetscViewerASCIIPrintf(viewer,"PetscSubcomm type %s with total %d MPI processes:\n",PetscSubcommTypes[psubcomm->type],size);
111:       MPI_Comm_rank(comm,&rank);
112:       MPI_Comm_size(psubcomm->child,&subsize);
113:       MPI_Comm_rank(psubcomm->child,&subrank);
114:       MPI_Comm_rank(psubcomm->dupparent,&duprank);
115:       PetscViewerASCIIPushSynchronized(viewer);
116:       PetscViewerASCIISynchronizedPrintf(viewer,"  [%d], color %d, sub-size %d, sub-rank %d, duprank %d\n",rank,psubcomm->color,subsize,subrank,duprank);
117:       PetscViewerFlush(viewer);
118:       PetscViewerASCIIPopSynchronized(viewer);
119:     }
120:   } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not supported yet");
121:   return 0;
122: }

124: /*@
125:   PetscSubcommSetNumber - Set total number of subcommunicators.

127:    Collective

129:    Input Parameters:
130: +  psubcomm - PetscSubcomm context
131: -  nsubcomm - the total number of subcommunicators in psubcomm

133:    Level: advanced

135: .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetType(),PetscSubcommSetTypeGeneral()
136: @*/
137: PetscErrorCode  PetscSubcommSetNumber(PetscSubcomm psubcomm,PetscInt nsubcomm)
138: {
139:   MPI_Comm       comm=psubcomm->parent;
140:   PetscMPIInt    msub,size;

143:   MPI_Comm_size(comm,&size);
144:   PetscMPIIntCast(nsubcomm,&msub);

147:   psubcomm->n = msub;
148:   return 0;
149: }

151: /*@
152:   PetscSubcommSetType - Set type of subcommunicators.

154:    Collective

156:    Input Parameters:
157: +  psubcomm - PetscSubcomm context
158: -  subcommtype - subcommunicator type, PETSC_SUBCOMM_CONTIGUOUS,PETSC_SUBCOMM_INTERLACED

160:    Level: advanced

162: .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetTypeGeneral(), PetscSubcommType
163: @*/
164: PetscErrorCode  PetscSubcommSetType(PetscSubcomm psubcomm,PetscSubcommType subcommtype)
165: {

169:   if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS) {
170:     PetscSubcommCreate_contiguous(psubcomm);
171:   } else if (subcommtype == PETSC_SUBCOMM_INTERLACED) {
172:     PetscSubcommCreate_interlaced(psubcomm);
173:   } else SETERRQ(psubcomm->parent,PETSC_ERR_SUP,"PetscSubcommType %s is not supported yet",PetscSubcommTypes[subcommtype]);
174:   return 0;
175: }

177: /*@
178:   PetscSubcommSetTypeGeneral - Set a PetscSubcomm from user's specifications

180:    Collective

182:    Input Parameters:
183: +  psubcomm - PetscSubcomm context
184: .  color   - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator.
185: -  subrank - rank in the subcommunicator

187:    Level: advanced

189: .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetType()
190: @*/
191: PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm,PetscMPIInt color,PetscMPIInt subrank)
192: {
193:   MPI_Comm       subcomm=0,dupcomm=0,comm=psubcomm->parent;
194:   PetscMPIInt    size,icolor,duprank,*recvbuf,sendbuf[3],mysubsize,rank,*subsize;
195:   PetscMPIInt    i,nsubcomm=psubcomm->n;


200:   MPI_Comm_split(comm,color,subrank,&subcomm);

202:   /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
203:   /* TODO: this can be done in an ostensibly scalale way (i.e., without allocating an array of size 'size') as is done in PetscObjectsCreateGlobalOrdering(). */
204:   MPI_Comm_size(comm,&size);
205:   PetscMalloc1(2*size,&recvbuf);

207:   MPI_Comm_rank(comm,&rank);
208:   MPI_Comm_size(subcomm,&mysubsize);

210:   sendbuf[0] = color;
211:   sendbuf[1] = mysubsize;
212:   MPI_Allgather(sendbuf,2,MPI_INT,recvbuf,2,MPI_INT,comm);

214:   PetscCalloc1(nsubcomm,&subsize);
215:   for (i=0; i<2*size; i+=2) {
216:     subsize[recvbuf[i]] = recvbuf[i+1];
217:   }
218:   PetscFree(recvbuf);

220:   duprank = 0;
221:   for (icolor=0; icolor<nsubcomm; icolor++) {
222:     if (icolor != color) { /* not color of this process */
223:       duprank += subsize[icolor];
224:     } else {
225:       duprank += subrank;
226:       break;
227:     }
228:   }
229:   MPI_Comm_split(comm,0,duprank,&dupcomm);

231:   PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);
232:   PetscCommDuplicate(subcomm,&psubcomm->child,NULL);
233:   MPI_Comm_free(&dupcomm);
234:   MPI_Comm_free(&subcomm);

236:   psubcomm->color   = color;
237:   psubcomm->subsize = subsize;
238:   psubcomm->type    = PETSC_SUBCOMM_GENERAL;
239:   return 0;
240: }

242: /*@
243:   PetscSubcommDestroy - Destroys a PetscSubcomm object

245:    Collective on PetscSubcomm

247:    Input Parameter:
248:    .  psubcomm - the PetscSubcomm context

250:    Level: advanced

252: .seealso: PetscSubcommCreate(),PetscSubcommSetType()
253: @*/
254: PetscErrorCode  PetscSubcommDestroy(PetscSubcomm *psubcomm)
255: {
256:   if (!*psubcomm) return 0;
257:   PetscCommDestroy(&(*psubcomm)->dupparent);
258:   PetscCommDestroy(&(*psubcomm)->child);
259:   PetscFree((*psubcomm)->subsize);
260:   if ((*psubcomm)->subcommprefix) PetscFree((*psubcomm)->subcommprefix);
261:   PetscFree((*psubcomm));
262:   return 0;
263: }

265: /*@
266:   PetscSubcommCreate - Create a PetscSubcomm context.

268:    Collective

270:    Input Parameter:
271: .  comm - MPI communicator

273:    Output Parameter:
274: .  psubcomm - location to store the PetscSubcomm context

276:    Level: advanced

278: .seealso: PetscSubcommDestroy(), PetscSubcommSetTypeGeneral(), PetscSubcommSetFromOptions(), PetscSubcommSetType(),
279:           PetscSubcommSetNumber()
280: @*/
281: PetscErrorCode  PetscSubcommCreate(MPI_Comm comm,PetscSubcomm *psubcomm)
282: {
283:   PetscMPIInt    rank,size;

285:   PetscNew(psubcomm);

287:   /* set defaults */
288:   MPI_Comm_rank(comm,&rank);
289:   MPI_Comm_size(comm,&size);

291:   (*psubcomm)->parent    = comm;
292:   (*psubcomm)->dupparent = comm;
293:   (*psubcomm)->child     = PETSC_COMM_SELF;
294:   (*psubcomm)->n         = size;
295:   (*psubcomm)->color     = rank;
296:   (*psubcomm)->subsize   = NULL;
297:   (*psubcomm)->type      = PETSC_SUBCOMM_INTERLACED;
298:   return 0;
299: }

301: /*@C
302:   PetscSubcommGetParent - Gets the communicator that was used to create the PetscSubcomm

304:    Collective

306:    Input Parameter:
307: .  scomm - the PetscSubcomm

309:    Output Parameter:
310: .  pcomm - location to store the parent communicator

312:    Level: intermediate

314: .seealso: PetscSubcommDestroy(), PetscSubcommSetTypeGeneral(), PetscSubcommSetFromOptions(), PetscSubcommSetType(),
315:           PetscSubcommSetNumber(), PetscSubcommGetChild(), PetscSubcommContiguousParent()
316: @*/
317: PetscErrorCode  PetscSubcommGetParent(PetscSubcomm scomm,MPI_Comm *pcomm)
318: {
319:   *pcomm = PetscSubcommParent(scomm);
320:   return 0;
321: }

323: /*@C
324:   PetscSubcommGetContiguousParent - Gets a communicator that that is a duplicate of the parent but has the ranks
325:                                     reordered by the order they are in the children

327:    Collective

329:    Input Parameter:
330: .  scomm - the PetscSubcomm

332:    Output Parameter:
333: .  pcomm - location to store the parent communicator

335:    Level: intermediate

337: .seealso: PetscSubcommDestroy(), PetscSubcommSetTypeGeneral(), PetscSubcommSetFromOptions(), PetscSubcommSetType(),
338:           PetscSubcommSetNumber(), PetscSubcommGetChild(), PetscSubcommContiguousParent()
339: @*/
340: PetscErrorCode  PetscSubcommGetContiguousParent(PetscSubcomm scomm,MPI_Comm *pcomm)
341: {
342:   *pcomm = PetscSubcommContiguousParent(scomm);
343:   return 0;
344: }

346: /*@C
347:   PetscSubcommGetChild - Gets the communicator created by the PetscSubcomm

349:    Collective

351:    Input Parameter:
352: .  scomm - the PetscSubcomm

354:    Output Parameter:
355: .  ccomm - location to store the child communicator

357:    Level: intermediate

359: .seealso: PetscSubcommDestroy(), PetscSubcommSetTypeGeneral(), PetscSubcommSetFromOptions(), PetscSubcommSetType(),
360:           PetscSubcommSetNumber(), PetscSubcommGetParent(), PetscSubcommContiguousParent()
361: @*/
362: PetscErrorCode  PetscSubcommGetChild(PetscSubcomm scomm,MPI_Comm *ccomm)
363: {
364:   *ccomm = PetscSubcommChild(scomm);
365:   return 0;
366: }

368: static PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm)
369: {
370:   PetscMPIInt    rank,size,*subsize,duprank=-1,subrank=-1;
371:   PetscMPIInt    np_subcomm,nleftover,i,color=-1,rankstart,nsubcomm=psubcomm->n;
372:   MPI_Comm       subcomm=0,dupcomm=0,comm=psubcomm->parent;

374:   MPI_Comm_rank(comm,&rank);
375:   MPI_Comm_size(comm,&size);

377:   /* get size of each subcommunicator */
378:   PetscMalloc1(1+nsubcomm,&subsize);

380:   np_subcomm = size/nsubcomm;
381:   nleftover  = size - nsubcomm*np_subcomm;
382:   for (i=0; i<nsubcomm; i++) {
383:     subsize[i] = np_subcomm;
384:     if (i<nleftover) subsize[i]++;
385:   }

387:   /* get color and subrank of this proc */
388:   rankstart = 0;
389:   for (i=0; i<nsubcomm; i++) {
390:     if (rank >= rankstart && rank < rankstart+subsize[i]) {
391:       color   = i;
392:       subrank = rank - rankstart;
393:       duprank = rank;
394:       break;
395:     } else rankstart += subsize[i];
396:   }

398:   MPI_Comm_split(comm,color,subrank,&subcomm);

400:   /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
401:   MPI_Comm_split(comm,0,duprank,&dupcomm);
402:   PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);
403:   PetscCommDuplicate(subcomm,&psubcomm->child,NULL);
404:   MPI_Comm_free(&dupcomm);
405:   MPI_Comm_free(&subcomm);

407:   psubcomm->color   = color;
408:   psubcomm->subsize = subsize;
409:   psubcomm->type    = PETSC_SUBCOMM_CONTIGUOUS;
410:   return 0;
411: }

413: /*
414:    Note:
415:    In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators
416:    by iteratively taking a process into a subcommunicator.
417:    Example: size=4, nsubcomm=(*psubcomm)->n=3
418:      comm=(*psubcomm)->parent:
419:       rank:     [0]  [1]  [2]  [3]
420:       color:     0    1    2    0

422:      subcomm=(*psubcomm)->comm:
423:       subrank:  [0]  [0]  [0]  [1]

425:      dupcomm=(*psubcomm)->dupparent:
426:       duprank:  [0]  [2]  [3]  [1]

428:      Here, subcomm[color = 0] has subsize=2, owns process [0] and [3]
429:            subcomm[color = 1] has subsize=1, owns process [1]
430:            subcomm[color = 2] has subsize=1, owns process [2]
431:            dupcomm has same number of processes as comm, and its duprank maps
432:            processes in subcomm contiguously into a 1d array:
433:             duprank: [0] [1]      [2]         [3]
434:             rank:    [0] [3]      [1]         [2]
435:                     subcomm[0] subcomm[1]  subcomm[2]
436: */

438: static PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm)
439: {
440:   PetscMPIInt    rank,size,*subsize,duprank,subrank;
441:   PetscMPIInt    np_subcomm,nleftover,i,j,color,nsubcomm=psubcomm->n;
442:   MPI_Comm       subcomm=0,dupcomm=0,comm=psubcomm->parent;

444:   MPI_Comm_rank(comm,&rank);
445:   MPI_Comm_size(comm,&size);

447:   /* get size of each subcommunicator */
448:   PetscMalloc1(1+nsubcomm,&subsize);

450:   np_subcomm = size/nsubcomm;
451:   nleftover  = size - nsubcomm*np_subcomm;
452:   for (i=0; i<nsubcomm; i++) {
453:     subsize[i] = np_subcomm;
454:     if (i<nleftover) subsize[i]++;
455:   }

457:   /* find color for this proc */
458:   color   = rank%nsubcomm;
459:   subrank = rank/nsubcomm;

461:   MPI_Comm_split(comm,color,subrank,&subcomm);

463:   j = 0; duprank = 0;
464:   for (i=0; i<nsubcomm; i++) {
465:     if (j == color) {
466:       duprank += subrank;
467:       break;
468:     }
469:     duprank += subsize[i]; j++;
470:   }

472:   /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
473:   MPI_Comm_split(comm,0,duprank,&dupcomm);
474:   PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);
475:   PetscCommDuplicate(subcomm,&psubcomm->child,NULL);
476:   MPI_Comm_free(&dupcomm);
477:   MPI_Comm_free(&subcomm);

479:   psubcomm->color   = color;
480:   psubcomm->subsize = subsize;
481:   psubcomm->type    = PETSC_SUBCOMM_INTERLACED;
482:   return 0;
483: }