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 { 26 PetscSubcommType type; 27 PetscBool flg; 28 29 PetscFunctionBegin; 30 PetscCheck(psubcomm,PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Must call PetscSubcommCreate firt"); 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(0); 57 } 58 59 /*@C 60 PetscSubcommSetOptionsPrefix - Sets the prefix used for searching for all 61 PetscSubcomm items in the options database. 62 63 Logically collective on PetscSubcomm. 64 65 Level: Intermediate 66 67 Input Parameters: 68 + psubcomm - PetscSubcomm context 69 - prefix - the prefix to prepend all PetscSubcomm item names with. 70 71 @*/ 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(0); 83 } 84 85 /*@C 86 PetscSubcommView - Views a PetscSubcomm of values as either ASCII text or a binary file 87 88 Collective on PetscSubcomm 89 90 Input Parameters: 91 + psubcomm - PetscSubcomm context 92 - viewer - location to view the values 93 94 Level: beginner 95 @*/ 96 PetscErrorCode PetscSubcommView(PetscSubcomm psubcomm,PetscViewer viewer) 97 { 98 PetscBool iascii; 99 PetscViewerFormat format; 100 101 PetscFunctionBegin; 102 PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii)); 103 if (iascii) { 104 PetscCall(PetscViewerGetFormat(viewer,&format)); 105 if (format == PETSC_VIEWER_DEFAULT) { 106 MPI_Comm comm=psubcomm->parent; 107 PetscMPIInt rank,size,subsize,subrank,duprank; 108 109 PetscCallMPI(MPI_Comm_size(comm,&size)); 110 PetscCall(PetscViewerASCIIPrintf(viewer,"PetscSubcomm type %s with total %d MPI processes:\n",PetscSubcommTypes[psubcomm->type],size)); 111 PetscCallMPI(MPI_Comm_rank(comm,&rank)); 112 PetscCallMPI(MPI_Comm_size(psubcomm->child,&subsize)); 113 PetscCallMPI(MPI_Comm_rank(psubcomm->child,&subrank)); 114 PetscCallMPI(MPI_Comm_rank(psubcomm->dupparent,&duprank)); 115 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 116 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer," [%d], color %d, sub-size %d, sub-rank %d, duprank %d\n",rank,psubcomm->color,subsize,subrank,duprank)); 117 PetscCall(PetscViewerFlush(viewer)); 118 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 119 } 120 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not supported yet"); 121 PetscFunctionReturn(0); 122 } 123 124 /*@ 125 PetscSubcommSetNumber - Set total number of subcommunicators. 126 127 Collective 128 129 Input Parameters: 130 + psubcomm - PetscSubcomm context 131 - nsubcomm - the total number of subcommunicators in psubcomm 132 133 Level: advanced 134 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; 141 142 PetscFunctionBegin; 143 PetscCheck(psubcomm,PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate() first"); 144 PetscCallMPI(MPI_Comm_size(comm,&size)); 145 PetscCall(PetscMPIIntCast(nsubcomm,&msub)); 146 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); 147 148 psubcomm->n = msub; 149 PetscFunctionReturn(0); 150 } 151 152 /*@ 153 PetscSubcommSetType - Set type of subcommunicators. 154 155 Collective 156 157 Input Parameters: 158 + psubcomm - PetscSubcomm context 159 - subcommtype - subcommunicator type, PETSC_SUBCOMM_CONTIGUOUS,PETSC_SUBCOMM_INTERLACED 160 161 Level: advanced 162 163 .seealso: `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommType` 164 @*/ 165 PetscErrorCode PetscSubcommSetType(PetscSubcomm psubcomm,PetscSubcommType subcommtype) 166 { 167 PetscFunctionBegin; 168 PetscCheck(psubcomm,PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 169 PetscCheck(psubcomm->n >= 1,PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %d is incorrect. Call PetscSubcommSetNumber()",psubcomm->n); 170 171 if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS) { 172 PetscCall(PetscSubcommCreate_contiguous(psubcomm)); 173 } else if (subcommtype == PETSC_SUBCOMM_INTERLACED) { 174 PetscCall(PetscSubcommCreate_interlaced(psubcomm)); 175 } else SETERRQ(psubcomm->parent,PETSC_ERR_SUP,"PetscSubcommType %s is not supported yet",PetscSubcommTypes[subcommtype]); 176 PetscFunctionReturn(0); 177 } 178 179 /*@ 180 PetscSubcommSetTypeGeneral - Set a PetscSubcomm from user's specifications 181 182 Collective 183 184 Input Parameters: 185 + psubcomm - PetscSubcomm context 186 . color - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator. 187 - subrank - rank in the subcommunicator 188 189 Level: advanced 190 191 .seealso: `PetscSubcommCreate()`, `PetscSubcommDestroy()`, `PetscSubcommSetNumber()`, `PetscSubcommSetType()` 192 @*/ 193 PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm,PetscMPIInt color,PetscMPIInt subrank) 194 { 195 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 196 PetscMPIInt size,icolor,duprank,*recvbuf,sendbuf[3],mysubsize,rank,*subsize; 197 PetscMPIInt i,nsubcomm=psubcomm->n; 198 199 PetscFunctionBegin; 200 PetscCheck(psubcomm,PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 201 PetscCheck(nsubcomm >= 1,PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %d is incorrect. Call PetscSubcommSetNumber()",nsubcomm); 202 203 PetscCallMPI(MPI_Comm_split(comm,color,subrank,&subcomm)); 204 205 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 206 /* TODO: this can be done in an ostensibly scalale way (i.e., without allocating an array of size 'size') as is done in PetscObjectsCreateGlobalOrdering(). */ 207 PetscCallMPI(MPI_Comm_size(comm,&size)); 208 PetscCall(PetscMalloc1(2*size,&recvbuf)); 209 210 PetscCallMPI(MPI_Comm_rank(comm,&rank)); 211 PetscCallMPI(MPI_Comm_size(subcomm,&mysubsize)); 212 213 sendbuf[0] = color; 214 sendbuf[1] = mysubsize; 215 PetscCallMPI(MPI_Allgather(sendbuf,2,MPI_INT,recvbuf,2,MPI_INT,comm)); 216 217 PetscCall(PetscCalloc1(nsubcomm,&subsize)); 218 for (i=0; i<2*size; i+=2) { 219 subsize[recvbuf[i]] = recvbuf[i+1]; 220 } 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(0); 243 } 244 245 /*@ 246 PetscSubcommDestroy - Destroys a PetscSubcomm object 247 248 Collective on PetscSubcomm 249 250 Input Parameter: 251 . psubcomm - the PetscSubcomm context 252 253 Level: advanced 254 255 .seealso: `PetscSubcommCreate()`, `PetscSubcommSetType()` 256 @*/ 257 PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm) 258 { 259 PetscFunctionBegin; 260 if (!*psubcomm) PetscFunctionReturn(0); 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(0); 267 } 268 269 /*@ 270 PetscSubcommCreate - Create a PetscSubcomm context. 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: `PetscSubcommDestroy()`, `PetscSubcommSetTypeGeneral()`, `PetscSubcommSetFromOptions()`, `PetscSubcommSetType()`, 283 `PetscSubcommSetNumber()` 284 @*/ 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(0); 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 @*/ 322 PetscErrorCode PetscSubcommGetParent(PetscSubcomm scomm,MPI_Comm *pcomm) 323 { 324 *pcomm = PetscSubcommParent(scomm); 325 return 0; 326 } 327 328 /*@C 329 PetscSubcommGetContiguousParent - Gets a communicator that 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 @*/ 345 PetscErrorCode PetscSubcommGetContiguousParent(PetscSubcomm scomm,MPI_Comm *pcomm) 346 { 347 *pcomm = PetscSubcommContiguousParent(scomm); 348 return 0; 349 } 350 351 /*@C 352 PetscSubcommGetChild - Gets the communicator 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 @*/ 367 PetscErrorCode PetscSubcommGetChild(PetscSubcomm scomm,MPI_Comm *ccomm) 368 { 369 *ccomm = PetscSubcommChild(scomm); 370 return 0; 371 } 372 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(0); 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 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; duprank = 0; 471 for (i=0; i<nsubcomm; i++) { 472 if (j == color) { 473 duprank += subrank; 474 break; 475 } 476 duprank += subsize[i]; j++; 477 } 478 479 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 480 PetscCallMPI(MPI_Comm_split(comm,0,duprank,&dupcomm)); 481 PetscCall(PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL)); 482 PetscCall(PetscCommDuplicate(subcomm,&psubcomm->child,NULL)); 483 PetscCallMPI(MPI_Comm_free(&dupcomm)); 484 PetscCallMPI(MPI_Comm_free(&subcomm)); 485 486 psubcomm->color = color; 487 psubcomm->subsize = subsize; 488 psubcomm->type = PETSC_SUBCOMM_INTERLACED; 489 PetscFunctionReturn(0); 490 } 491