1 2 /* 3 Provides utility routines for split MPI communicator. 4 */ 5 #include <petscsys.h> /*I "petscsys.h" I*/ 6 #include <petsc-private/threadcommimpl.h> /* Petsc_ThreadComm_keyval */ 7 #include <petscviewer.h> 8 9 const char *const PetscSubcommTypes[] = {"GENERAL","CONTIGUOUS","INTERLACED","PetscSubcommType","PETSC_SUBCOMM_",0}; 10 11 extern PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm); 12 extern PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm); 13 #undef __FUNCT__ 14 #define __FUNCT__ "PetscSubcommSetFromOptions" 15 PetscErrorCode PetscSubcommSetFromOptions(PetscSubcomm psubcomm) 16 { 17 PetscErrorCode ierr; 18 PetscInt type; 19 const char *psubcommTypes[3] = {"general","contiguous","interlaced"}; 20 PetscBool flg; 21 22 PetscFunctionBegin; 23 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"Must call PetscSubcommCreate firt"); 24 ierr = PetscOptionsEList("-psubcomm_type","PETSc subcommunicator","PetscSubcommSetType",psubcommTypes,3,psubcommTypes[2],&type,&flg);CHKERRQ(ierr); 25 if (flg && psubcomm->type != type) { 26 /* free old structures */ 27 ierr = PetscCommDestroy(&(psubcomm)->dupparent);CHKERRQ(ierr); 28 ierr = PetscCommDestroy(&(psubcomm)->comm);CHKERRQ(ierr); 29 ierr = PetscFree((psubcomm)->subsize);CHKERRQ(ierr); 30 switch (type) { 31 case 1: 32 ierr = PetscSubcommCreate_contiguous(psubcomm);CHKERRQ(ierr); 33 break; 34 case 2: 35 ierr = PetscSubcommCreate_interlaced(psubcomm);CHKERRQ(ierr); 36 break; 37 default: 38 SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"PetscSubcommType %D is not supported yet",type); 39 } 40 } 41 42 ierr = PetscOptionsHasName(NULL, "-psubcomm_view", &flg);CHKERRQ(ierr); 43 if (flg) { 44 ierr = PetscSubcommView(psubcomm,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 45 } 46 PetscFunctionReturn(0); 47 } 48 49 #undef __FUNCT__ 50 #define __FUNCT__ "PetscSubcommView" 51 PetscErrorCode PetscSubcommView(PetscSubcomm psubcomm,PetscViewer viewer) 52 { 53 PetscErrorCode ierr; 54 PetscBool iascii; 55 PetscViewerFormat format; 56 57 PetscFunctionBegin; 58 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 59 if (iascii) { 60 ierr = PetscViewerGetFormat(viewer,&format);CHKERRQ(ierr); 61 if (format == PETSC_VIEWER_DEFAULT) { 62 MPI_Comm comm=psubcomm->parent; 63 PetscMPIInt rank,size,subsize,subrank,duprank; 64 65 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 66 ierr = PetscViewerASCIIPrintf(viewer,"PetscSubcomm type %s with total %D MPI processes:\n",PetscSubcommTypes[psubcomm->type],size);CHKERRQ(ierr); 67 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 68 ierr = MPI_Comm_size(psubcomm->comm,&subsize);CHKERRQ(ierr); 69 ierr = MPI_Comm_rank(psubcomm->comm,&subrank);CHKERRQ(ierr); 70 ierr = MPI_Comm_rank(psubcomm->dupparent,&duprank);CHKERRQ(ierr); 71 ierr = PetscSynchronizedPrintf(comm," [%D], color %D, sub-size %D, sub-rank %D, duprank %D\n",rank,psubcomm->color,subsize,subrank,duprank); 72 ierr = PetscSynchronizedFlush(comm);CHKERRQ(ierr); 73 } 74 } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Not supported yet"); 75 PetscFunctionReturn(0); 76 } 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "PetscSubcommSetNumber" 80 /*@C 81 PetscSubcommSetNumber - Set total number of subcommunicators. 82 83 Collective on MPI_Comm 84 85 Input Parameter: 86 + psubcomm - PetscSubcomm context 87 - nsubcomm - the total number of subcommunicators in psubcomm 88 89 Level: advanced 90 91 .keywords: communicator 92 93 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetType(),PetscSubcommSetTypeGeneral() 94 @*/ 95 PetscErrorCode PetscSubcommSetNumber(PetscSubcomm psubcomm,PetscInt nsubcomm) 96 { 97 PetscErrorCode ierr; 98 MPI_Comm comm=psubcomm->parent; 99 PetscMPIInt rank,size; 100 101 PetscFunctionBegin; 102 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate() first"); 103 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 104 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 105 if (nsubcomm < 1 || nsubcomm > size) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Num of subcommunicators %D cannot be < 1 or > input comm size %D",nsubcomm,size); 106 107 psubcomm->n = nsubcomm; 108 PetscFunctionReturn(0); 109 } 110 111 #undef __FUNCT__ 112 #define __FUNCT__ "PetscSubcommSetType" 113 /*@C 114 PetscSubcommSetType - Set type of subcommunicators. 115 116 Collective on MPI_Comm 117 118 Input Parameter: 119 + psubcomm - PetscSubcomm context 120 - subcommtype - subcommunicator type, PETSC_SUBCOMM_CONTIGUOUS,PETSC_SUBCOMM_INTERLACED 121 122 Level: advanced 123 124 .keywords: communicator 125 126 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetTypeGeneral() 127 @*/ 128 PetscErrorCode PetscSubcommSetType(PetscSubcomm psubcomm,PetscSubcommType subcommtype) 129 { 130 PetscErrorCode ierr; 131 132 PetscFunctionBegin; 133 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 134 if (psubcomm->n < 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %D is incorrect. Call PetscSubcommSetNumber()",psubcomm->n); 135 136 if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS) { 137 ierr = PetscSubcommCreate_contiguous(psubcomm);CHKERRQ(ierr); 138 } else if (subcommtype == PETSC_SUBCOMM_INTERLACED) { 139 ierr = PetscSubcommCreate_interlaced(psubcomm);CHKERRQ(ierr); 140 } else SETERRQ1(psubcomm->parent,PETSC_ERR_SUP,"PetscSubcommType %D is not supported yet",subcommtype); 141 PetscFunctionReturn(0); 142 } 143 144 #undef __FUNCT__ 145 #define __FUNCT__ "PetscSubcommSetTypeGeneral" 146 /*@C 147 PetscSubcommSetTypeGeneral - Set a PetscSubcomm from user's specifications 148 149 Collective on MPI_Comm 150 151 Input Parameter: 152 + psubcomm - PetscSubcomm context 153 . color - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator. 154 - subrank - rank in the subcommunicator 155 156 Level: advanced 157 158 .keywords: communicator, create 159 160 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetType() 161 @*/ 162 PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm,PetscMPIInt color,PetscMPIInt subrank) 163 { 164 PetscErrorCode ierr; 165 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 166 PetscMPIInt size,icolor,duprank,*recvbuf,sendbuf[3],subsize,rank,subsizes[psubcomm->n]; 167 PetscInt i; 168 169 PetscFunctionBegin; 170 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 171 if (psubcomm->n < 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %D is incorrect. Call PetscSubcommSetNumber()",psubcomm->n); 172 173 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 174 175 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 176 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 177 ierr = PetscMalloc(3*size*sizeof(PetscMPIInt),&recvbuf);CHKERRQ(ierr); 178 179 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 180 ierr = MPI_Comm_size(subcomm,&subsize);CHKERRQ(ierr); 181 182 sendbuf[0] = color; 183 sendbuf[1] = subrank; 184 sendbuf[2] = subsize; 185 ierr = MPI_Allgather(sendbuf,3,MPI_INT,recvbuf,3,MPI_INT,comm);CHKERRQ(ierr); 186 187 /* printf("[%d] (color,subrank,subsize): ",rank); */ 188 for (i=0; i<3*size; i+=3) { 189 /* printf(" (%d %d %d);",recvbuf[i],recvbuf[i+1],recvbuf[i+2]); */ 190 subsizes[recvbuf[i]] = recvbuf[i+2]; 191 } 192 ierr = PetscFree(recvbuf);CHKERRQ(ierr); 193 194 duprank = 0; 195 for (icolor=0; icolor<psubcomm->n; icolor++) { 196 if (icolor != color) { /* not color of this process */ 197 duprank += subsizes[icolor]; 198 } else { 199 duprank += subrank; 200 break; 201 } 202 } 203 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 204 205 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);CHKERRQ(ierr); 206 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,NULL);CHKERRQ(ierr); 207 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 208 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 209 210 psubcomm->color = color; 211 PetscFunctionReturn(0); 212 } 213 214 #undef __FUNCT__ 215 #define __FUNCT__ "PetscSubcommDestroy" 216 PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm) 217 { 218 PetscErrorCode ierr; 219 220 PetscFunctionBegin; 221 if (!*psubcomm) PetscFunctionReturn(0); 222 ierr = PetscCommDestroy(&(*psubcomm)->dupparent);CHKERRQ(ierr); 223 ierr = PetscCommDestroy(&(*psubcomm)->comm);CHKERRQ(ierr); 224 ierr = PetscFree((*psubcomm)->subsize);CHKERRQ(ierr); 225 ierr = PetscFree((*psubcomm));CHKERRQ(ierr); 226 PetscFunctionReturn(0); 227 } 228 229 #undef __FUNCT__ 230 #define __FUNCT__ "PetscSubcommCreate" 231 /*@C 232 PetscSubcommCreate - Create a PetscSubcomm context. 233 234 Collective on MPI_Comm 235 236 Input Parameter: 237 . comm - MPI communicator 238 239 Output Parameter: 240 . psubcomm - location to store the PetscSubcomm context 241 242 Level: advanced 243 244 .keywords: communicator, create 245 246 .seealso: PetscSubcommDestroy() 247 @*/ 248 PetscErrorCode PetscSubcommCreate(MPI_Comm comm,PetscSubcomm *psubcomm) 249 { 250 PetscErrorCode ierr; 251 PetscMPIInt rank,size; 252 253 PetscFunctionBegin; 254 ierr = PetscNew(struct _n_PetscSubcomm,psubcomm);CHKERRQ(ierr); 255 256 /* set defaults */ 257 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 258 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 259 260 (*psubcomm)->parent = comm; 261 (*psubcomm)->dupparent = comm; 262 (*psubcomm)->comm = PETSC_COMM_SELF; 263 (*psubcomm)->n = size; 264 (*psubcomm)->color = rank; 265 (*psubcomm)->subsize = NULL; 266 (*psubcomm)->type = PETSC_SUBCOMM_INTERLACED; 267 PetscFunctionReturn(0); 268 } 269 270 #undef __FUNCT__ 271 #define __FUNCT__ "PetscSubcommCreate_contiguous" 272 PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm) 273 { 274 PetscErrorCode ierr; 275 PetscMPIInt rank,size,*subsize,duprank=-1,subrank=-1; 276 PetscInt np_subcomm,nleftover,i,color=-1,rankstart,nsubcomm=psubcomm->n; 277 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 278 279 PetscFunctionBegin; 280 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 281 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 282 283 /* get size of each subcommunicator */ 284 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 285 286 np_subcomm = size/nsubcomm; 287 nleftover = size - nsubcomm*np_subcomm; 288 for (i=0; i<nsubcomm; i++) { 289 subsize[i] = np_subcomm; 290 if (i<nleftover) subsize[i]++; 291 } 292 293 /* get color and subrank of this proc */ 294 rankstart = 0; 295 for (i=0; i<nsubcomm; i++) { 296 if (rank >= rankstart && rank < rankstart+subsize[i]) { 297 color = i; 298 subrank = rank - rankstart; 299 duprank = rank; 300 break; 301 } else rankstart += subsize[i]; 302 } 303 304 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 305 306 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 307 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 308 { 309 PetscThreadComm tcomm; 310 ierr = PetscCommGetThreadComm(comm,&tcomm);CHKERRQ(ierr); 311 ierr = MPI_Attr_put(dupcomm,Petsc_ThreadComm_keyval,tcomm);CHKERRQ(ierr); 312 tcomm->refct++; 313 ierr = MPI_Attr_put(subcomm,Petsc_ThreadComm_keyval,tcomm);CHKERRQ(ierr); 314 tcomm->refct++; 315 } 316 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);CHKERRQ(ierr); 317 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,NULL);CHKERRQ(ierr); 318 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 319 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 320 321 psubcomm->color = color; 322 psubcomm->subsize = subsize; 323 psubcomm->type = PETSC_SUBCOMM_CONTIGUOUS; 324 PetscFunctionReturn(0); 325 } 326 327 #undef __FUNCT__ 328 #define __FUNCT__ "PetscSubcommCreate_interlaced" 329 /* 330 Note: 331 In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators 332 by iteratively taking a process into a subcommunicator. 333 Example: size=4, nsubcomm=(*psubcomm)->n=3 334 comm=(*psubcomm)->parent: 335 rank: [0] [1] [2] [3] 336 color: 0 1 2 0 337 338 subcomm=(*psubcomm)->comm: 339 subrank: [0] [0] [0] [1] 340 341 dupcomm=(*psubcomm)->dupparent: 342 duprank: [0] [2] [3] [1] 343 344 Here, subcomm[color = 0] has subsize=2, owns process [0] and [3] 345 subcomm[color = 1] has subsize=1, owns process [1] 346 subcomm[color = 2] has subsize=1, owns process [2] 347 dupcomm has same number of processes as comm, and its duprank maps 348 processes in subcomm contiguously into a 1d array: 349 duprank: [0] [1] [2] [3] 350 rank: [0] [3] [1] [2] 351 subcomm[0] subcomm[1] subcomm[2] 352 */ 353 354 PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm) 355 { 356 PetscErrorCode ierr; 357 PetscMPIInt rank,size,*subsize,duprank,subrank; 358 PetscInt np_subcomm,nleftover,i,j,color,nsubcomm=psubcomm->n; 359 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 360 361 PetscFunctionBegin; 362 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 363 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 364 365 /* get size of each subcommunicator */ 366 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 367 368 np_subcomm = size/nsubcomm; 369 nleftover = size - nsubcomm*np_subcomm; 370 for (i=0; i<nsubcomm; i++) { 371 subsize[i] = np_subcomm; 372 if (i<nleftover) subsize[i]++; 373 } 374 375 /* find color for this proc */ 376 color = rank%nsubcomm; 377 subrank = rank/nsubcomm; 378 379 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 380 381 j = 0; duprank = 0; 382 for (i=0; i<nsubcomm; i++) { 383 if (j == color) { 384 duprank += subrank; 385 break; 386 } 387 duprank += subsize[i]; j++; 388 } 389 390 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 391 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 392 { 393 PetscThreadComm tcomm; 394 ierr = PetscCommGetThreadComm(comm,&tcomm);CHKERRQ(ierr); 395 ierr = MPI_Attr_put(dupcomm,Petsc_ThreadComm_keyval,tcomm);CHKERRQ(ierr); 396 tcomm->refct++; 397 ierr = MPI_Attr_put(subcomm,Petsc_ThreadComm_keyval,tcomm);CHKERRQ(ierr); 398 tcomm->refct++; 399 } 400 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,NULL);CHKERRQ(ierr); 401 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,NULL);CHKERRQ(ierr); 402 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 403 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 404 405 psubcomm->color = color; 406 psubcomm->subsize = subsize; 407 psubcomm->type = PETSC_SUBCOMM_INTERLACED; 408 PetscFunctionReturn(0); 409 } 410 411 412