1 2 /* 3 Provides utility routines for split MPI communicator. 4 */ 5 #include <petscsys.h> /*I "petscsys.h" I*/ 6 7 extern PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm); 8 extern PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm); 9 10 #undef __FUNCT__ 11 #define __FUNCT__ "PetscSubcommSetNumber" 12 /*@C 13 PetscSubcommSetNumber - Set total number of subcommunicators. 14 15 Collective on MPI_Comm 16 17 Input Parameter: 18 + psubcomm - PetscSubcomm context 19 - nsubcomm - the total number of subcommunicators in psubcomm 20 21 Level: advanced 22 23 .keywords: communicator 24 25 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetType(),PetscSubcommSetTypeGeneral() 26 @*/ 27 PetscErrorCode PetscSubcommSetNumber(PetscSubcomm psubcomm,PetscInt nsubcomm) 28 { 29 PetscErrorCode ierr; 30 MPI_Comm comm=psubcomm->parent; 31 PetscMPIInt rank,size; 32 33 PetscFunctionBegin; 34 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate() first"); 35 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 36 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 37 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); 38 39 psubcomm->n = nsubcomm; 40 PetscFunctionReturn(0); 41 } 42 43 #undef __FUNCT__ 44 #define __FUNCT__ "PetscSubcommSetType" 45 /*@C 46 PetscSubcommSetType - Set type of subcommunicators. 47 48 Collective on MPI_Comm 49 50 Input Parameter: 51 + psubcomm - PetscSubcomm context 52 - subcommtype - subcommunicator type, PETSC_SUBCOMM_CONTIGUOUS,PETSC_SUBCOMM_INTERLACED 53 54 Level: advanced 55 56 .keywords: communicator 57 58 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetTypeGeneral() 59 @*/ 60 PetscErrorCode PetscSubcommSetType(PetscSubcomm psubcomm,const PetscSubcommType subcommtype) 61 { 62 PetscErrorCode ierr; 63 64 PetscFunctionBegin; 65 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 66 if (psubcomm->n < 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %D is incorrect. Call PetscSubcommSetNumber()",psubcomm->n); 67 68 if (subcommtype == PETSC_SUBCOMM_CONTIGUOUS){ 69 ierr = PetscSubcommCreate_contiguous(psubcomm);CHKERRQ(ierr); 70 } else if (subcommtype == PETSC_SUBCOMM_INTERLACED){ 71 ierr = PetscSubcommCreate_interlaced(psubcomm);CHKERRQ(ierr); 72 } else { 73 SETERRQ1(psubcomm->parent,PETSC_ERR_SUP,"PetscSubcommType %D is not supported yet",subcommtype); 74 } 75 PetscFunctionReturn(0); 76 } 77 78 #undef __FUNCT__ 79 #define __FUNCT__ "PetscSubcommSetTypeGeneral" 80 /*@C 81 PetscSubcommSetTypeGeneral - Set type of subcommunicators from user's specifications 82 83 Collective on MPI_Comm 84 85 Input Parameter: 86 + psubcomm - PetscSubcomm context 87 . color - control of subset assignment (nonnegative integer). Processes with the same color are in the same subcommunicator. 88 . subrank - rank in the subcommunicator 89 - duprank - rank in the dupparent (see PetscSubcomm) 90 91 Level: advanced 92 93 .keywords: communicator, create 94 95 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetType() 96 @*/ 97 PetscErrorCode PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm,PetscMPIInt color,PetscMPIInt subrank,PetscMPIInt duprank) 98 { 99 PetscErrorCode ierr; 100 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 101 PetscMPIInt size; 102 103 PetscFunctionBegin; 104 if (!psubcomm) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NULL,"PetscSubcomm is not created. Call PetscSubcommCreate()"); 105 if (psubcomm->n < 1) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG,"number of subcommunicators %D is incorrect. Call PetscSubcommSetNumber()",psubcomm->n); 106 107 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 108 109 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm 110 if duprank is not a valid number, then dupcomm is not created - not all applications require dupcomm! */ 111 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 112 if (duprank == PETSC_DECIDE){ 113 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"duprank==PETSC_DECIDE is not supported yet"); 114 } else if (duprank >= 0 && duprank < size){ 115 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 116 } 117 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,PETSC_NULL);CHKERRQ(ierr); 118 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,PETSC_NULL);CHKERRQ(ierr); 119 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 120 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 121 psubcomm->color = color; 122 PetscFunctionReturn(0); 123 } 124 125 #undef __FUNCT__ 126 #define __FUNCT__ "PetscSubcommDestroy" 127 PetscErrorCode PetscSubcommDestroy(PetscSubcomm *psubcomm) 128 { 129 PetscErrorCode ierr; 130 131 PetscFunctionBegin; 132 if (!*psubcomm) PetscFunctionReturn(0); 133 ierr = PetscCommDestroy(&(*psubcomm)->dupparent);CHKERRQ(ierr); 134 ierr = PetscCommDestroy(&(*psubcomm)->comm);CHKERRQ(ierr); 135 ierr = PetscFree((*psubcomm));CHKERRQ(ierr); 136 PetscFunctionReturn(0); 137 } 138 139 #undef __FUNCT__ 140 #define __FUNCT__ "PetscSubcommCreate" 141 /*@C 142 PetscSubcommCreate - Create a PetscSubcomm context. 143 144 Collective on MPI_Comm 145 146 Input Parameter: 147 + comm - MPI communicator 148 . nsubcomm - the number of subcommunicators to be created 149 - subcommtype - subcommunicator type 150 151 Output Parameter: 152 . psubcomm - location to store the PetscSubcomm context 153 154 Level: advanced 155 156 .keywords: communicator, create 157 158 .seealso: PetscSubcommDestroy() 159 @*/ 160 PetscErrorCode PetscSubcommCreate(MPI_Comm comm,PetscSubcomm *psubcomm) 161 { 162 PetscErrorCode ierr; 163 164 PetscFunctionBegin; 165 166 ierr = PetscNew(struct _n_PetscSubcomm,psubcomm);CHKERRQ(ierr); 167 (*psubcomm)->parent = comm; 168 PetscFunctionReturn(0); 169 } 170 171 #undef __FUNCT__ 172 #define __FUNCT__ "PetscSubcommCreate_contiguous" 173 PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm) 174 { 175 PetscErrorCode ierr; 176 PetscMPIInt rank,size,*subsize,duprank=-1,subrank=-1; 177 PetscInt np_subcomm,nleftover,i,color=-1,rankstart,nsubcomm=psubcomm->n; 178 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 179 180 PetscFunctionBegin; 181 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 182 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 183 184 /* get size of each subcommunicator */ 185 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 186 np_subcomm = size/nsubcomm; 187 nleftover = size - nsubcomm*np_subcomm; 188 for (i=0; i<nsubcomm; i++){ 189 subsize[i] = np_subcomm; 190 if (i<nleftover) subsize[i]++; 191 } 192 193 /* get color and subrank of this proc */ 194 rankstart = 0; 195 for (i=0; i<nsubcomm; i++){ 196 if ( rank >= rankstart && rank < rankstart+subsize[i]) { 197 color = i; 198 subrank = rank - rankstart; 199 duprank = rank; 200 break; 201 } else { 202 rankstart += subsize[i]; 203 } 204 } 205 ierr = PetscFree(subsize);CHKERRQ(ierr); 206 207 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 208 209 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 210 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 211 212 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,PETSC_NULL);CHKERRQ(ierr); 213 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,PETSC_NULL);CHKERRQ(ierr); 214 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 215 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 216 psubcomm->color = color; 217 PetscFunctionReturn(0); 218 } 219 220 #undef __FUNCT__ 221 #define __FUNCT__ "PetscSubcommCreate_interlaced" 222 /* 223 Note: 224 In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators 225 by iteratively taking a process into a subcommunicator. 226 Example: size=4, nsubcomm=(*psubcomm)->n=3 227 comm=(*psubcomm)->parent: 228 rank: [0] [1] [2] [3] 229 color: 0 1 2 0 230 231 subcomm=(*psubcomm)->comm: 232 subrank: [0] [0] [0] [1] 233 234 dupcomm=(*psubcomm)->dupparent: 235 duprank: [0] [2] [3] [1] 236 237 Here, subcomm[color = 0] has subsize=2, owns process [0] and [3] 238 subcomm[color = 1] has subsize=1, owns process [1] 239 subcomm[color = 2] has subsize=1, owns process [2] 240 dupcomm has same number of processes as comm, and its duprank maps 241 processes in subcomm contiguously into a 1d array: 242 duprank: [0] [1] [2] [3] 243 rank: [0] [3] [1] [2] 244 subcomm[0] subcomm[1] subcomm[2] 245 */ 246 247 PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm) 248 { 249 PetscErrorCode ierr; 250 PetscMPIInt rank,size,*subsize,duprank,subrank; 251 PetscInt np_subcomm,nleftover,i,j,color,nsubcomm=psubcomm->n; 252 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 253 254 PetscFunctionBegin; 255 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 256 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 257 258 /* get size of each subcommunicator */ 259 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 260 np_subcomm = size/nsubcomm; 261 nleftover = size - nsubcomm*np_subcomm; 262 for (i=0; i<nsubcomm; i++){ 263 subsize[i] = np_subcomm; 264 if (i<nleftover) subsize[i]++; 265 } 266 267 /* find color for this proc */ 268 color = rank%nsubcomm; 269 subrank = rank/nsubcomm; 270 271 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 272 273 j = 0; duprank = 0; 274 for (i=0; i<nsubcomm; i++){ 275 if (j == color){ 276 duprank += subrank; 277 break; 278 } 279 duprank += subsize[i]; j++; 280 } 281 ierr = PetscFree(subsize);CHKERRQ(ierr); 282 283 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 284 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 285 286 ierr = PetscCommDuplicate(dupcomm,&psubcomm->dupparent,PETSC_NULL);CHKERRQ(ierr); 287 ierr = PetscCommDuplicate(subcomm,&psubcomm->comm,PETSC_NULL);CHKERRQ(ierr); 288 ierr = MPI_Comm_free(&dupcomm);CHKERRQ(ierr); 289 ierr = MPI_Comm_free(&subcomm);CHKERRQ(ierr); 290 psubcomm->color = color; 291 PetscFunctionReturn(0); 292 } 293 294 295