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