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