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 for a `PetscSubcomm`
16: Collective
18: Input Parameter:
19: . psubcomm - `PetscSubcomm` context
21: Level: beginner
23: .seealso: `PetscSubcomm`, `PetscSubcommCreate()`
24: @*/
25: PetscErrorCode PetscSubcommSetFromOptions(PetscSubcomm psubcomm)
26: {
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) PetscSubcommView(psubcomm, PETSC_VIEWER_STDOUT_(psubcomm->parent));
55: PetscOptionsEnd();
56: return 0;
57: }
59: /*@C
60: PetscSubcommSetOptionsPrefix - Sets the prefix used for searching for options in the options database for this object
62: Logically collective
64: Level: Intermediate
66: Input Parameters:
67: + psubcomm - `PetscSubcomm` context
68: - prefix - the prefix to prepend all PetscSubcomm item names with.
70: .seealso: `PetscSubcomm`, `PetscSubcommCreate()`
71: @*/
72: PetscErrorCode PetscSubcommSetOptionsPrefix(PetscSubcomm psubcomm, const char pre[])
73: {
74: if (!pre) {
75: PetscFree(psubcomm->subcommprefix);
76: } else {
78: PetscFree(psubcomm->subcommprefix);
79: PetscStrallocpy(pre, &(psubcomm->subcommprefix));
80: }
81: return 0;
82: }
84: /*@C
85: PetscSubcommView - Views a `PetscSubcomm`
87: Collective
89: Input Parameters:
90: + psubcomm - `PetscSubcomm` context
91: - viewer - `PetscViewer` to display the information
93: Level: beginner
95: .seealso: `PetscSubcomm`, `PetscSubcommCreate()`, `PetscViewer`
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 desired in the given `PetscSubcomm`
127: Collective
129: Input Parameters:
130: + psubcomm - `PetscSubcomm` context
131: - nsubcomm - the total number of subcommunicators in psubcomm
133: Level: advanced
135: .seealso: `PetscSubcomm`, `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 the way the original MPI communicator is divided up in the `PetscSubcomm`
154: Collective
156: Input Parameters:
157: + psubcomm - `PetscSubcomm` context
158: - subcommtype - `PetscSubcommType` `PETSC_SUBCOMM_CONTIGUOUS` or `PETSC_SUBCOMM_INTERLACED`
160: Level: advanced
162: .seealso: `PetscSubcommType`, `PETSC_SUBCOMM_CONTIGUOUS`, `PETSC_SUBCOMM_INTERLACED`,
163: `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommType`
164: @*/
165: PetscErrorCode PetscSubcommSetType(PetscSubcomm psubcomm, PetscSubcommType subcommtype)
166: {
170: if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS) {
171: PetscSubcommCreate_contiguous(psubcomm);
172: } else if (subcommtype == PETSC_SUBCOMM_INTERLACED) {
173: PetscSubcommCreate_interlaced(psubcomm);
174: } else SETERRQ(psubcomm->parent, PETSC_ERR_SUP, "PetscSubcommType %s is not supported yet", PetscSubcommTypes[subcommtype]);
175: return 0;
176: }
178: /*@
179: PetscSubcommSetTypeGeneral - Divides up a communicator based on a specific user's specification
181: Collective
183: Input Parameters:
184: + psubcomm - `PetscSubcomm` context
185: . color - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator.
186: - subrank - rank in the subcommunicator
188: Level: advanced
190: .seealso: `PetscSubcommType`, `PETSC_SUBCOMM_CONTIGUOUS`, `PETSC_SUBCOMM_INTERLACED`, `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetType()`
191: @*/
192: PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm, PetscMPIInt color, PetscMPIInt subrank)
193: {
194: MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
195: PetscMPIInt size, icolor, duprank, *recvbuf, sendbuf[3], mysubsize, rank, *subsize;
196: PetscMPIInt i, nsubcomm = psubcomm->n;
201: MPI_Comm_split(comm, color, subrank, &subcomm);
203: /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
204: /* TODO: this can be done in an ostensibly scalale way (i.e., without allocating an array of size 'size') as is done in PetscObjectsCreateGlobalOrdering(). */
205: MPI_Comm_size(comm, &size);
206: PetscMalloc1(2 * size, &recvbuf);
208: MPI_Comm_rank(comm, &rank);
209: MPI_Comm_size(subcomm, &mysubsize);
211: sendbuf[0] = color;
212: sendbuf[1] = mysubsize;
213: MPI_Allgather(sendbuf, 2, MPI_INT, recvbuf, 2, MPI_INT, comm);
215: PetscCalloc1(nsubcomm, &subsize);
216: for (i = 0; i < 2 * size; i += 2) subsize[recvbuf[i]] = recvbuf[i + 1];
217: PetscFree(recvbuf);
219: duprank = 0;
220: for (icolor = 0; icolor < nsubcomm; icolor++) {
221: if (icolor != color) { /* not color of this process */
222: duprank += subsize[icolor];
223: } else {
224: duprank += subrank;
225: break;
226: }
227: }
228: MPI_Comm_split(comm, 0, duprank, &dupcomm);
230: PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL);
231: PetscCommDuplicate(subcomm, &psubcomm->child, NULL);
232: MPI_Comm_free(&dupcomm);
233: MPI_Comm_free(&subcomm);
235: psubcomm->color = color;
236: psubcomm->subsize = subsize;
237: psubcomm->type = PETSC_SUBCOMM_GENERAL;
238: return 0;
239: }
241: /*@
242: PetscSubcommDestroy - Destroys a `PetscSubcomm` object
244: Collective
246: Input Parameter:
247: . psubcomm - the `PetscSubcomm` context
249: Level: advanced
251: .seealso: `PetscSubcommCreate()`, `PetscSubcommSetType()`
252: @*/
253: PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm)
254: {
255: if (!*psubcomm) return 0;
256: PetscCommDestroy(&(*psubcomm)->dupparent);
257: PetscCommDestroy(&(*psubcomm)->child);
258: PetscFree((*psubcomm)->subsize);
259: if ((*psubcomm)->subcommprefix) PetscFree((*psubcomm)->subcommprefix);
260: PetscFree((*psubcomm));
261: return 0;
262: }
264: /*@
265: PetscSubcommCreate - Create a `PetscSubcomm` context. This object is used to manage the division of a `MPI_Comm` into subcommunicators
267: Collective
269: Input Parameter:
270: . comm - MPI communicator
272: Output Parameter:
273: . psubcomm - location to store the `PetscSubcomm` context
275: Level: advanced
277: .seealso: `PetscSubcomm`, `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
278: `PetscSubcommSetNumber()`
279: @*/
280: PetscErrorCode PetscSubcommCreate(MPI_Comm comm, PetscSubcomm *psubcomm)
281: {
282: PetscMPIInt rank, size;
284: PetscNew(psubcomm);
286: /* set defaults */
287: MPI_Comm_rank(comm, &rank);
288: MPI_Comm_size(comm, &size);
290: (*psubcomm)->parent = comm;
291: (*psubcomm)->dupparent = comm;
292: (*psubcomm)->child = PETSC_COMM_SELF;
293: (*psubcomm)->n = size;
294: (*psubcomm)->color = rank;
295: (*psubcomm)->subsize = NULL;
296: (*psubcomm)->type = PETSC_SUBCOMM_INTERLACED;
297: return 0;
298: }
300: /*@C
301: PetscSubcommGetParent - Gets the communicator that was used to create the `PetscSubcomm`
303: Collective
305: Input Parameter:
306: . scomm - the `PetscSubcomm`
308: Output Parameter:
309: . pcomm - location to store the parent communicator
311: Level: intermediate
313: .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
314: `PetscSubcommSetNumber()`, `PetscSubcommGetChild()`, `PetscSubcommContiguousParent()`
315: @*/
316: PetscErrorCode PetscSubcommGetParent(PetscSubcomm scomm, MPI_Comm *pcomm)
317: {
318: *pcomm = PetscSubcommParent(scomm);
319: return 0;
320: }
322: /*@C
323: PetscSubcommGetContiguousParent - Gets a communicator that that is a duplicate of the parent but has the ranks
324: reordered by the order they are in the children
326: Collective
328: Input Parameter:
329: . scomm - the `PetscSubcomm`
331: Output Parameter:
332: . pcomm - location to store the parent communicator
334: Level: intermediate
336: .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
337: `PetscSubcommSetNumber()`, `PetscSubcommGetChild()`, `PetscSubcommContiguousParent()`
338: @*/
339: PetscErrorCode PetscSubcommGetContiguousParent(PetscSubcomm scomm, MPI_Comm *pcomm)
340: {
341: *pcomm = PetscSubcommContiguousParent(scomm);
342: return 0;
343: }
345: /*@C
346: PetscSubcommGetChild - Gets the communicator created by the `PetscSubcomm`. This is part of one of the subcommunicators created by the `PetscSubcomm`
348: Collective
350: Input Parameter:
351: . scomm - the `PetscSubcomm`
353: Output Parameter:
354: . ccomm - location to store the child communicator
356: Level: intermediate
358: .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
359: `PetscSubcommSetNumber()`, `PetscSubcommGetParent()`, `PetscSubcommContiguousParent()`
360: @*/
361: PetscErrorCode PetscSubcommGetChild(PetscSubcomm scomm, MPI_Comm *ccomm)
362: {
363: *ccomm = PetscSubcommChild(scomm);
364: return 0;
365: }
367: static PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm)
368: {
369: PetscMPIInt rank, size, *subsize, duprank = -1, subrank = -1;
370: PetscMPIInt np_subcomm, nleftover, i, color = -1, rankstart, nsubcomm = psubcomm->n;
371: MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
373: MPI_Comm_rank(comm, &rank);
374: MPI_Comm_size(comm, &size);
376: /* get size of each subcommunicator */
377: PetscMalloc1(1 + nsubcomm, &subsize);
379: np_subcomm = size / nsubcomm;
380: nleftover = size - nsubcomm * np_subcomm;
381: for (i = 0; i < nsubcomm; i++) {
382: subsize[i] = np_subcomm;
383: if (i < nleftover) subsize[i]++;
384: }
386: /* get color and subrank of this proc */
387: rankstart = 0;
388: for (i = 0; i < nsubcomm; i++) {
389: if (rank >= rankstart && rank < rankstart + subsize[i]) {
390: color = i;
391: subrank = rank - rankstart;
392: duprank = rank;
393: break;
394: } else rankstart += subsize[i];
395: }
397: MPI_Comm_split(comm, color, subrank, &subcomm);
399: /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
400: MPI_Comm_split(comm, 0, duprank, &dupcomm);
401: PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL);
402: PetscCommDuplicate(subcomm, &psubcomm->child, NULL);
403: MPI_Comm_free(&dupcomm);
404: MPI_Comm_free(&subcomm);
406: psubcomm->color = color;
407: psubcomm->subsize = subsize;
408: psubcomm->type = PETSC_SUBCOMM_CONTIGUOUS;
409: return 0;
410: }
412: /*
413: Note:
414: In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators
415: by iteratively taking a process into a subcommunicator.
416: Example: size=4, nsubcomm=(*psubcomm)->n=3
417: comm=(*psubcomm)->parent:
418: rank: [0] [1] [2] [3]
419: color: 0 1 2 0
421: subcomm=(*psubcomm)->comm:
422: subrank: [0] [0] [0] [1]
424: dupcomm=(*psubcomm)->dupparent:
425: duprank: [0] [2] [3] [1]
427: Here, subcomm[color = 0] has subsize=2, owns process [0] and [3]
428: subcomm[color = 1] has subsize=1, owns process [1]
429: subcomm[color = 2] has subsize=1, owns process [2]
430: dupcomm has same number of processes as comm, and its duprank maps
431: processes in subcomm contiguously into a 1d array:
432: duprank: [0] [1] [2] [3]
433: rank: [0] [3] [1] [2]
434: subcomm[0] subcomm[1] subcomm[2]
435: */
437: static PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm)
438: {
439: PetscMPIInt rank, size, *subsize, duprank, subrank;
440: PetscMPIInt np_subcomm, nleftover, i, j, color, nsubcomm = psubcomm->n;
441: MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
443: MPI_Comm_rank(comm, &rank);
444: MPI_Comm_size(comm, &size);
446: /* get size of each subcommunicator */
447: PetscMalloc1(1 + nsubcomm, &subsize);
449: np_subcomm = size / nsubcomm;
450: nleftover = size - nsubcomm * np_subcomm;
451: for (i = 0; i < nsubcomm; i++) {
452: subsize[i] = np_subcomm;
453: if (i < nleftover) subsize[i]++;
454: }
456: /* find color for this proc */
457: color = rank % nsubcomm;
458: subrank = rank / nsubcomm;
460: MPI_Comm_split(comm, color, subrank, &subcomm);
462: j = 0;
463: duprank = 0;
464: for (i = 0; i < nsubcomm; i++) {
465: if (j == color) {
466: duprank += subrank;
467: break;
468: }
469: duprank += subsize[i];
470: j++;
471: }
473: /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
474: MPI_Comm_split(comm, 0, duprank, &dupcomm);
475: PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL);
476: PetscCommDuplicate(subcomm, &psubcomm->child, NULL);
477: MPI_Comm_free(&dupcomm);
478: MPI_Comm_free(&subcomm);
480: psubcomm->color = color;
481: psubcomm->subsize = subsize;
482: psubcomm->type = PETSC_SUBCOMM_INTERLACED;
483: return 0;
484: }