1 #define PETSC_DLL 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 PETSCMAT_DLLEXPORT 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 53 54 Level: advanced 55 56 .keywords: communicator 57 58 .seealso: PetscSubcommCreate(),PetscSubcommDestroy(),PetscSubcommSetNumber(),PetscSubcommSetTypeGeneral() 59 @*/ 60 PetscErrorCode PETSCMAT_DLLEXPORT 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 PetscErrorCode PETSCMAT_DLLEXPORT PetscSubcommSetTypeGeneral(PetscSubcomm psubcomm) 81 { 82 PetscFunctionBegin; 83 PetscFunctionReturn(0); 84 } 85 86 #undef __FUNCT__ 87 #define __FUNCT__ "PetscSubcommDestroy" 88 PetscErrorCode PETSCMAT_DLLEXPORT PetscSubcommDestroy(PetscSubcomm psubcomm) 89 { 90 PetscErrorCode ierr; 91 92 PetscFunctionBegin; 93 ierr = MPI_Comm_free(&psubcomm->dupparent);CHKERRQ(ierr); 94 ierr = MPI_Comm_free(&psubcomm->comm);CHKERRQ(ierr); 95 ierr = PetscFree(psubcomm);CHKERRQ(ierr); 96 PetscFunctionReturn(0); 97 } 98 99 #undef __FUNCT__ 100 #define __FUNCT__ "PetscSubcommCreate" 101 /*@C 102 PetscSubcommCreate - Create a PetscSubcomm context. 103 104 Collective on MPI_Comm 105 106 Input Parameter: 107 + comm - MPI communicator 108 . nsubcomm - the number of subcommunicators to be created 109 - subcommtype - subcommunicator type 110 111 Output Parameter: 112 . psubcomm - location to store the PetscSubcomm context 113 114 Level: advanced 115 116 .keywords: communicator, create 117 118 .seealso: PetscSubcommDestroy() 119 @*/ 120 PetscErrorCode PETSCMAT_DLLEXPORT PetscSubcommCreate(MPI_Comm comm,PetscSubcomm *psubcomm) 121 { 122 PetscErrorCode ierr; 123 124 PetscFunctionBegin; 125 126 ierr = PetscNew(struct _n_PetscSubcomm,psubcomm);CHKERRQ(ierr); 127 (*psubcomm)->parent = comm; 128 PetscFunctionReturn(0); 129 } 130 131 #undef __FUNCT__ 132 #define __FUNCT__ "PetscSubcommCreate_contiguous" 133 PetscErrorCode PetscSubcommCreate_contiguous(PetscSubcomm psubcomm) 134 { 135 PetscErrorCode ierr; 136 PetscMPIInt rank,size,*subsize,duprank=-1,subrank=-1; 137 PetscInt np_subcomm,nleftover,i,color=-1,rankstart,nsubcomm=psubcomm->n; 138 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 139 140 PetscFunctionBegin; 141 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 142 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 143 144 /* get size of each subcommunicator */ 145 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 146 np_subcomm = size/nsubcomm; 147 nleftover = size - nsubcomm*np_subcomm; 148 for (i=0; i<nsubcomm; i++){ 149 subsize[i] = np_subcomm; 150 if (i<nleftover) subsize[i]++; 151 } 152 153 /* get color and subrank of this proc */ 154 rankstart = 0; 155 for (i=0; i<nsubcomm; i++){ 156 if ( rank >= rankstart && rank < rankstart+subsize[i]) { 157 color = i; 158 subrank = rank - rankstart; 159 duprank = rank; 160 break; 161 } else { 162 rankstart += subsize[i]; 163 } 164 } 165 ierr = PetscFree(subsize);CHKERRQ(ierr); 166 167 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 168 169 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 170 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 171 172 psubcomm->dupparent = dupcomm; 173 psubcomm->comm = subcomm; 174 psubcomm->color = color; 175 PetscFunctionReturn(0); 176 } 177 178 #undef __FUNCT__ 179 #define __FUNCT__ "PetscSubcommCreate_interlaced" 180 /* 181 Note: 182 In PCREDUNDANT, to avoid data scattering from subcomm back to original comm, we create subcommunicators 183 by iteratively taking a process into a subcommunicator. 184 Example: size=4, nsubcomm=(*psubcomm)->n=3 185 comm=(*psubcomm)->parent: 186 rank: [0] [1] [2] [3] 187 color: 0 1 2 0 188 189 subcomm=(*psubcomm)->comm: 190 subrank: [0] [0] [0] [1] 191 192 dupcomm=(*psubcomm)->dupparent: 193 duprank: [0] [2] [3] [1] 194 195 Here, subcomm[color = 0] has subsize=2, owns process [0] and [3] 196 subcomm[color = 1] has subsize=1, owns process [1] 197 subcomm[color = 2] has subsize=1, owns process [2] 198 dupcomm has same number of processes as comm, and its duprank maps 199 processes in subcomm contiguously into a 1d array: 200 duprank: [0] [1] [2] [3] 201 rank: [0] [3] [1] [2] 202 subcomm[0] subcomm[1] subcomm[2] 203 */ 204 205 PetscErrorCode PetscSubcommCreate_interlaced(PetscSubcomm psubcomm) 206 { 207 PetscErrorCode ierr; 208 PetscMPIInt rank,size,*subsize,duprank,subrank; 209 PetscInt np_subcomm,nleftover,i,j,color,nsubcomm=psubcomm->n; 210 MPI_Comm subcomm=0,dupcomm=0,comm=psubcomm->parent; 211 212 PetscFunctionBegin; 213 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 214 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 215 216 /* get size of each subcommunicator */ 217 ierr = PetscMalloc((1+nsubcomm)*sizeof(PetscMPIInt),&subsize);CHKERRQ(ierr); 218 np_subcomm = size/nsubcomm; 219 nleftover = size - nsubcomm*np_subcomm; 220 for (i=0; i<nsubcomm; i++){ 221 subsize[i] = np_subcomm; 222 if (i<nleftover) subsize[i]++; 223 } 224 225 /* find color for this proc */ 226 color = rank%nsubcomm; 227 subrank = rank/nsubcomm; 228 229 ierr = MPI_Comm_split(comm,color,subrank,&subcomm);CHKERRQ(ierr); 230 231 j = 0; duprank = 0; 232 for (i=0; i<nsubcomm; i++){ 233 if (j == color){ 234 duprank += subrank; 235 break; 236 } 237 duprank += subsize[i]; j++; 238 } 239 ierr = PetscFree(subsize);CHKERRQ(ierr); 240 241 /* create dupcomm with same size as comm, but its rank, duprank, maps subcomm's contiguously into dupcomm */ 242 ierr = MPI_Comm_split(comm,0,duprank,&dupcomm);CHKERRQ(ierr); 243 244 psubcomm->dupparent = dupcomm; 245 psubcomm->comm = subcomm; 246 psubcomm->color = color; 247 PetscFunctionReturn(0); 248 } 249 250 251