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