1 /*
2 Provides utility routines for split MPI communicator.
3 */
4 #include <petscsys.h> /*I "petscsys.h" I*/
5 #include <petscviewer.h>
6
7 const char *const PetscSubcommTypes[] = {"GENERAL", "CONTIGUOUS", "INTERLACED", "PetscSubcommType", "PETSC_SUBCOMM_", NULL};
8
9 static PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm);
10 static PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm);
11
12 /*@
13 PetscSubcommSetFromOptions - Allows setting options for a `PetscSubcomm`
14
15 Collective
16
17 Input Parameter:
18 . psubcomm - `PetscSubcomm` context
19
20 Level: beginner
21
22 .seealso: `PetscSubcomm`, `PetscSubcommCreate()`
23 @*/
PetscSubcommSetFromOptions(PetscSubcomm psubcomm)24 PetscErrorCode PetscSubcommSetFromOptions(PetscSubcomm psubcomm)
25 {
26 PetscSubcommType type;
27 PetscBool flg;
28
29 PetscFunctionBegin;
30 PetscCheck(psubcomm, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "Must call PetscSubcommCreate first");
31
32 PetscOptionsBegin(psubcomm->parent, psubcomm->subcommprefix, "Options for PetscSubcomm", NULL);
33 PetscCall(PetscOptionsEnum("-psubcomm_type", NULL, NULL, PetscSubcommTypes, (PetscEnum)psubcomm->type, (PetscEnum *)&type, &flg));
34 if (flg && psubcomm->type != type) {
35 /* free old structures */
36 PetscCall(PetscCommDestroy(&(psubcomm)->dupparent));
37 PetscCall(PetscCommDestroy(&(psubcomm)->child));
38 PetscCall(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 PetscCall(PetscSubcommCreate_contiguous(psubcomm));
44 break;
45 case PETSC_SUBCOMM_INTERLACED:
46 PetscCall(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 }
52
53 PetscCall(PetscOptionsName("-psubcomm_view", "Triggers display of PetscSubcomm context", "PetscSubcommView", &flg));
54 if (flg) PetscCall(PetscSubcommView(psubcomm, PETSC_VIEWER_STDOUT_(psubcomm->parent)));
55 PetscOptionsEnd();
56 PetscFunctionReturn(PETSC_SUCCESS);
57 }
58
59 /*@
60 PetscSubcommSetOptionsPrefix - Sets the prefix used for searching for options in the options database for this object
61
62 Logically Collective
63
64 Level: intermediate
65
66 Input Parameters:
67 + psubcomm - `PetscSubcomm` context
68 - pre - the prefix to prepend all `PetscSubcomm` item names with.
69
70 .seealso: `PetscSubcomm`, `PetscSubcommCreate()`
71 @*/
PetscSubcommSetOptionsPrefix(PetscSubcomm psubcomm,const char pre[])72 PetscErrorCode PetscSubcommSetOptionsPrefix(PetscSubcomm psubcomm, const char pre[])
73 {
74 PetscFunctionBegin;
75 if (!pre) {
76 PetscCall(PetscFree(psubcomm->subcommprefix));
77 } else {
78 PetscCheck(pre[0] != '-', PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Options prefix should not begin with a hyphen");
79 PetscCall(PetscFree(psubcomm->subcommprefix));
80 PetscCall(PetscStrallocpy(pre, &psubcomm->subcommprefix));
81 }
82 PetscFunctionReturn(PETSC_SUCCESS);
83 }
84
85 /*@
86 PetscSubcommView - Views a `PetscSubcomm`
87
88 Collective
89
90 Input Parameters:
91 + psubcomm - `PetscSubcomm` context
92 - viewer - `PetscViewer` to display the information
93
94 Level: beginner
95
96 .seealso: `PetscSubcomm`, `PetscSubcommCreate()`, `PetscViewer`
97 @*/
PetscSubcommView(PetscSubcomm psubcomm,PetscViewer viewer)98 PetscErrorCode PetscSubcommView(PetscSubcomm psubcomm, PetscViewer viewer)
99 {
100 PetscBool isascii;
101 PetscViewerFormat format;
102
103 PetscFunctionBegin;
104 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
105 PetscCheck(isascii, PETSC_COMM_SELF, PETSC_ERR_SUP, "Not supported yet");
106 PetscCall(PetscViewerGetFormat(viewer, &format));
107 if (format == PETSC_VIEWER_DEFAULT) {
108 MPI_Comm comm = psubcomm->parent;
109 PetscMPIInt rank, size, subsize, subrank, duprank;
110
111 PetscCallMPI(MPI_Comm_size(comm, &size));
112 PetscCall(PetscViewerASCIIPrintf(viewer, "PetscSubcomm type %s with total %d MPI processes:\n", PetscSubcommTypes[psubcomm->type], size));
113 PetscCallMPI(MPI_Comm_rank(comm, &rank));
114 PetscCallMPI(MPI_Comm_size(psubcomm->child, &subsize));
115 PetscCallMPI(MPI_Comm_rank(psubcomm->child, &subrank));
116 PetscCallMPI(MPI_Comm_rank(psubcomm->dupparent, &duprank));
117 PetscCall(PetscViewerASCIIPushSynchronized(viewer));
118 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, " [%d], color %d, sub-size %d, sub-rank %d, duprank %d\n", rank, psubcomm->color, subsize, subrank, duprank));
119 PetscCall(PetscViewerFlush(viewer));
120 PetscCall(PetscViewerASCIIPopSynchronized(viewer));
121 }
122 PetscFunctionReturn(PETSC_SUCCESS);
123 }
124
125 /*@
126 PetscSubcommSetNumber - Set total number of subcommunicators desired in the given `PetscSubcomm`
127
128 Collective
129
130 Input Parameters:
131 + psubcomm - `PetscSubcomm` context
132 - nsubcomm - the total number of subcommunicators in psubcomm
133
134 Level: advanced
135
136 .seealso: `PetscSubcomm`, `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetType()`, `PetscSubcommSetTypeGeneral()`
137 @*/
PetscSubcommSetNumber(PetscSubcomm psubcomm,PetscInt nsubcomm)138 PetscErrorCode PetscSubcommSetNumber(PetscSubcomm psubcomm, PetscInt nsubcomm)
139 {
140 MPI_Comm comm = psubcomm->parent;
141 PetscMPIInt msub, size;
142
143 PetscFunctionBegin;
144 PetscCheck(psubcomm, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "PetscSubcomm is not created. Call PetscSubcommCreate() first");
145 PetscCallMPI(MPI_Comm_size(comm, &size));
146 PetscCall(PetscMPIIntCast(nsubcomm, &msub));
147 PetscCheck(msub >= 1 && msub <= size, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Num of subcommunicators %d cannot be < 1 or > input comm size %d", msub, size);
148
149 psubcomm->n = msub;
150 PetscFunctionReturn(PETSC_SUCCESS);
151 }
152
153 /*@
154 PetscSubcommSetType - Set the way the original MPI communicator is divided up in the `PetscSubcomm`
155
156 Collective
157
158 Input Parameters:
159 + psubcomm - `PetscSubcomm` context
160 - subcommtype - `PetscSubcommType` `PETSC_SUBCOMM_CONTIGUOUS` or `PETSC_SUBCOMM_INTERLACED`
161
162 Level: advanced
163
164 .seealso: `PetscSubcommType`, `PETSC_SUBCOMM_CONTIGUOUS`, `PETSC_SUBCOMM_INTERLACED`,
165 `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetTypeGeneral()`
166 @*/
PetscSubcommSetType(PetscSubcomm psubcomm,PetscSubcommType subcommtype)167 PetscErrorCode PetscSubcommSetType(PetscSubcomm psubcomm, PetscSubcommType subcommtype)
168 {
169 PetscFunctionBegin;
170 PetscCheck(psubcomm, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "PetscSubcomm is not created. Call PetscSubcommCreate()");
171 PetscCheck(psubcomm->n >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "number of subcommunicators %d is incorrect. Call PetscSubcommSetNumber()", psubcomm->n);
172
173 if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS) {
174 PetscCall(PetscSubcommCreate_contiguous(psubcomm));
175 } else if (subcommtype == PETSC_SUBCOMM_INTERLACED) {
176 PetscCall(PetscSubcommCreate_interlaced(psubcomm));
177 } else SETERRQ(psubcomm->parent, PETSC_ERR_SUP, "PetscSubcommType %s is not supported yet", PetscSubcommTypes[subcommtype]);
178 PetscFunctionReturn(PETSC_SUCCESS);
179 }
180
181 /*@
182 PetscSubcommSetTypeGeneral - Divides up a communicator based on a specific user's specification
183
184 Collective
185
186 Input Parameters:
187 + psubcomm - `PetscSubcomm` context
188 . color - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator.
189 - subrank - rank in the subcommunicator
190
191 Level: advanced
192
193 .seealso: `PetscSubcommType`, `PETSC_SUBCOMM_CONTIGUOUS`, `PETSC_SUBCOMM_INTERLACED`, `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetType()`
194 @*/
PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm,PetscMPIInt color,PetscMPIInt subrank)195 PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm, PetscMPIInt color, PetscMPIInt subrank)
196 {
197 MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
198 PetscMPIInt size, icolor, duprank, *recvbuf, sendbuf[3], mysubsize, rank, *subsize;
199 PetscMPIInt i, nsubcomm = psubcomm->n;
200
201 PetscFunctionBegin;
202 PetscCheck(psubcomm, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "PetscSubcomm is not created. Call PetscSubcommCreate()");
203 PetscCheck(nsubcomm >= 1, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "number of subcommunicators %d is incorrect. Call PetscSubcommSetNumber()", nsubcomm);
204
205 PetscCallMPI(MPI_Comm_split(comm, color, subrank, &subcomm));
206
207 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
208 /* TODO: this can be done in an ostensibly scalale way (i.e., without allocating an array of size 'size') as is done in PetscObjectsCreateGlobalOrdering(). */
209 PetscCallMPI(MPI_Comm_size(comm, &size));
210 PetscCall(PetscMalloc1(2 * size, &recvbuf));
211
212 PetscCallMPI(MPI_Comm_rank(comm, &rank));
213 PetscCallMPI(MPI_Comm_size(subcomm, &mysubsize));
214
215 sendbuf[0] = color;
216 sendbuf[1] = mysubsize;
217 PetscCallMPI(MPI_Allgather(sendbuf, 2, MPI_INT, recvbuf, 2, MPI_INT, comm));
218
219 PetscCall(PetscCalloc1(nsubcomm, &subsize));
220 for (i = 0; i < 2 * size; i += 2) subsize[recvbuf[i]] = recvbuf[i + 1];
221 PetscCall(PetscFree(recvbuf));
222
223 duprank = 0;
224 for (icolor = 0; icolor < nsubcomm; icolor++) {
225 if (icolor != color) { /* not color of this process */
226 duprank += subsize[icolor];
227 } else {
228 duprank += subrank;
229 break;
230 }
231 }
232 PetscCallMPI(MPI_Comm_split(comm, 0, duprank, &dupcomm));
233
234 PetscCall(PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL));
235 PetscCall(PetscCommDuplicate(subcomm, &psubcomm->child, NULL));
236 PetscCallMPI(MPI_Comm_free(&dupcomm));
237 PetscCallMPI(MPI_Comm_free(&subcomm));
238
239 psubcomm->color = color;
240 psubcomm->subsize = subsize;
241 psubcomm->type = PETSC_SUBCOMM_GENERAL;
242 PetscFunctionReturn(PETSC_SUCCESS);
243 }
244
245 /*@
246 PetscSubcommDestroy - Destroys a `PetscSubcomm` object
247
248 Collective
249
250 Input Parameter:
251 . psubcomm - the `PetscSubcomm` context
252
253 Level: advanced
254
255 .seealso: `PetscSubcommCreate()`, `PetscSubcommSetType()`
256 @*/
PetscSubcommDestroy(PetscSubcomm * psubcomm)257 PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm)
258 {
259 PetscFunctionBegin;
260 if (!*psubcomm) PetscFunctionReturn(PETSC_SUCCESS);
261 PetscCall(PetscCommDestroy(&(*psubcomm)->dupparent));
262 PetscCall(PetscCommDestroy(&(*psubcomm)->child));
263 PetscCall(PetscFree((*psubcomm)->subsize));
264 if ((*psubcomm)->subcommprefix) PetscCall(PetscFree((*psubcomm)->subcommprefix));
265 PetscCall(PetscFree(*psubcomm));
266 PetscFunctionReturn(PETSC_SUCCESS);
267 }
268
269 /*@
270 PetscSubcommCreate - Create a `PetscSubcomm` context. This object is used to manage the division of a `MPI_Comm` into subcommunicators
271
272 Collective
273
274 Input Parameter:
275 . comm - MPI communicator
276
277 Output Parameter:
278 . psubcomm - location to store the `PetscSubcomm` context
279
280 Level: advanced
281
282 .seealso: `PetscSubcomm`, `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
283 `PetscSubcommSetNumber()`
284 @*/
PetscSubcommCreate(MPI_Comm comm,PetscSubcomm * psubcomm)285 PetscErrorCode PetscSubcommCreate(MPI_Comm comm, PetscSubcomm *psubcomm)
286 {
287 PetscMPIInt rank, size;
288
289 PetscFunctionBegin;
290 PetscCall(PetscNew(psubcomm));
291
292 /* set defaults */
293 PetscCallMPI(MPI_Comm_rank(comm, &rank));
294 PetscCallMPI(MPI_Comm_size(comm, &size));
295
296 (*psubcomm)->parent = comm;
297 (*psubcomm)->dupparent = comm;
298 (*psubcomm)->child = PETSC_COMM_SELF;
299 (*psubcomm)->n = size;
300 (*psubcomm)->color = rank;
301 (*psubcomm)->subsize = NULL;
302 (*psubcomm)->type = PETSC_SUBCOMM_INTERLACED;
303 PetscFunctionReturn(PETSC_SUCCESS);
304 }
305
306 /*@C
307 PetscSubcommGetParent - Gets the communicator that was used to create the `PetscSubcomm`
308
309 Collective
310
311 Input Parameter:
312 . scomm - the `PetscSubcomm`
313
314 Output Parameter:
315 . pcomm - location to store the parent communicator
316
317 Level: intermediate
318
319 .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
320 `PetscSubcommSetNumber()`, `PetscSubcommGetChild()`, `PetscSubcommContiguousParent()`
321 @*/
PetscSubcommGetParent(PetscSubcomm scomm,MPI_Comm * pcomm)322 PetscErrorCode PetscSubcommGetParent(PetscSubcomm scomm, MPI_Comm *pcomm)
323 {
324 *pcomm = PetscSubcommParent(scomm);
325 return PETSC_SUCCESS;
326 }
327
328 /*@C
329 PetscSubcommGetContiguousParent - Gets a communicator that is a duplicate of the parent but has the ranks
330 reordered by the order they are in the children
331
332 Collective
333
334 Input Parameter:
335 . scomm - the `PetscSubcomm`
336
337 Output Parameter:
338 . pcomm - location to store the parent communicator
339
340 Level: intermediate
341
342 .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
343 `PetscSubcommSetNumber()`, `PetscSubcommGetChild()`, `PetscSubcommContiguousParent()`
344 @*/
PetscSubcommGetContiguousParent(PetscSubcomm scomm,MPI_Comm * pcomm)345 PetscErrorCode PetscSubcommGetContiguousParent(PetscSubcomm scomm, MPI_Comm *pcomm)
346 {
347 *pcomm = PetscSubcommContiguousParent(scomm);
348 return PETSC_SUCCESS;
349 }
350
351 /*@C
352 PetscSubcommGetChild - Gets the communicator created by the `PetscSubcomm`. This is part of one of the subcommunicators created by the `PetscSubcomm`
353
354 Collective
355
356 Input Parameter:
357 . scomm - the `PetscSubcomm`
358
359 Output Parameter:
360 . ccomm - location to store the child communicator
361
362 Level: intermediate
363
364 .seealso: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`,
365 `PetscSubcommSetNumber()`, `PetscSubcommGetParent()`, `PetscSubcommContiguousParent()`
366 @*/
PetscSubcommGetChild(PetscSubcomm scomm,MPI_Comm * ccomm)367 PetscErrorCode PetscSubcommGetChild(PetscSubcomm scomm, MPI_Comm *ccomm)
368 {
369 *ccomm = PetscSubcommChild(scomm);
370 return PETSC_SUCCESS;
371 }
372
PetscSubcommCreate_contiguous(PetscSubcomm psubcomm)373 static PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm)
374 {
375 PetscMPIInt rank, size, *subsize, duprank = -1, subrank = -1;
376 PetscMPIInt np_subcomm, nleftover, i, color = -1, rankstart, nsubcomm = psubcomm->n;
377 MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
378
379 PetscFunctionBegin;
380 PetscCallMPI(MPI_Comm_rank(comm, &rank));
381 PetscCallMPI(MPI_Comm_size(comm, &size));
382
383 /* get size of each subcommunicator */
384 PetscCall(PetscMalloc1(1 + nsubcomm, &subsize));
385
386 np_subcomm = size / nsubcomm;
387 nleftover = size - nsubcomm * np_subcomm;
388 for (i = 0; i < nsubcomm; i++) {
389 subsize[i] = np_subcomm;
390 if (i < nleftover) subsize[i]++;
391 }
392
393 /* get color and subrank of this proc */
394 rankstart = 0;
395 for (i = 0; i < nsubcomm; i++) {
396 if (rank >= rankstart && rank < rankstart + subsize[i]) {
397 color = i;
398 subrank = rank - rankstart;
399 duprank = rank;
400 break;
401 } else rankstart += subsize[i];
402 }
403
404 PetscCallMPI(MPI_Comm_split(comm, color, subrank, &subcomm));
405
406 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
407 PetscCallMPI(MPI_Comm_split(comm, 0, duprank, &dupcomm));
408 PetscCall(PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL));
409 PetscCall(PetscCommDuplicate(subcomm, &psubcomm->child, NULL));
410 PetscCallMPI(MPI_Comm_free(&dupcomm));
411 PetscCallMPI(MPI_Comm_free(&subcomm));
412
413 psubcomm->color = color;
414 psubcomm->subsize = subsize;
415 psubcomm->type = PETSC_SUBCOMM_CONTIGUOUS;
416 PetscFunctionReturn(PETSC_SUCCESS);
417 }
418
419 /*
420 Note:
421 In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators
422 by iteratively taking a process into a subcommunicator.
423 Example: size=4, nsubcomm=(*psubcomm)->n=3
424 comm=(*psubcomm)->parent:
425 rank: [0] [1] [2] [3]
426 color: 0 1 2 0
427
428 subcomm=(*psubcomm)->comm:
429 subrank: [0] [0] [0] [1]
430
431 dupcomm=(*psubcomm)->dupparent:
432 duprank: [0] [2] [3] [1]
433
434 Here, subcomm[color = 0] has subsize=2, owns process [0] and [3]
435 subcomm[color = 1] has subsize=1, owns process [1]
436 subcomm[color = 2] has subsize=1, owns process [2]
437 dupcomm has same number of processes as comm, and its duprank maps
438 processes in subcomm contiguously into a 1d array:
439 duprank: [0] [1] [2] [3]
440 rank: [0] [3] [1] [2]
441 subcomm[0] subcomm[1] subcomm[2]
442 */
443
PetscSubcommCreate_interlaced(PetscSubcomm psubcomm)444 static PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm)
445 {
446 PetscMPIInt rank, size, *subsize, duprank, subrank;
447 PetscMPIInt np_subcomm, nleftover, i, j, color, nsubcomm = psubcomm->n;
448 MPI_Comm subcomm = 0, dupcomm = 0, comm = psubcomm->parent;
449
450 PetscFunctionBegin;
451 PetscCallMPI(MPI_Comm_rank(comm, &rank));
452 PetscCallMPI(MPI_Comm_size(comm, &size));
453
454 /* get size of each subcommunicator */
455 PetscCall(PetscMalloc1(1 + nsubcomm, &subsize));
456
457 np_subcomm = size / nsubcomm;
458 nleftover = size - nsubcomm * np_subcomm;
459 for (i = 0; i < nsubcomm; i++) {
460 subsize[i] = np_subcomm;
461 if (i < nleftover) subsize[i]++;
462 }
463
464 /* find color for this proc */
465 color = rank % nsubcomm;
466 subrank = rank / nsubcomm;
467
468 PetscCallMPI(MPI_Comm_split(comm, color, subrank, &subcomm));
469
470 j = 0;
471 duprank = 0;
472 for (i = 0; i < nsubcomm; i++) {
473 if (j == color) {
474 duprank += subrank;
475 break;
476 }
477 duprank += subsize[i];
478 j++;
479 }
480
481 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */
482 PetscCallMPI(MPI_Comm_split(comm, 0, duprank, &dupcomm));
483 PetscCall(PetscCommDuplicate(dupcomm, &psubcomm->dupparent, NULL));
484 PetscCall(PetscCommDuplicate(subcomm, &psubcomm->child, NULL));
485 PetscCallMPI(MPI_Comm_free(&dupcomm));
486 PetscCallMPI(MPI_Comm_free(&subcomm));
487
488 psubcomm->color = color;
489 psubcomm->subsize = subsize;
490 psubcomm->type = PETSC_SUBCOMM_INTERLACED;
491 PetscFunctionReturn(PETSC_SUCCESS);
492 }
493