TABLE OF CONTENTS
ABINIT/distrb2 [ Functions ]
NAME
distrb2
FUNCTION
This routine creates the tabs of repartition of processors for sharing the jobs on k-points, spins and bands.
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (AR,XG,MB) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
mband = maximum number of bands nband(nkpt*nsppol) = number of bands per k point, for each spin nkpt = number of k-points nproc= number of processors available for this distribution nsppol = 1 for unpolarized, 2 for polarized
SIDE EFFECTS
mpi_enreg = informations about MPI parallelization mpi_enreg%proc_distrb(nkpt,mband,nsppol)=number of the processor that will treat each band in each k point. mpi_enreg%nproc_kpt is set
NOTES
For the time being, the band parallelisation works only when the number of bands is identical for spin up and spin down at the same k point. The problem is the most clearly seen in the kpgio routine, where a different parallel repartition of k points for spin up and spin down would conflict with the present computation of k+G sphere, independent of the spin.
PARENTS
dfpt_looppert,eig2stern,eig2tot,mpi_setup
CHILDREN
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 subroutine distrb2(mband,nband,nkpt,nproc,nsppol,mpi_enreg) 54 55 use defs_basis 56 use defs_abitypes 57 use m_profiling_abi 58 use m_errors 59 60 use m_io_tools, only : file_exists, open_file 61 62 !This section has been created automatically by the script Abilint (TD). 63 !Do not modify the following lines by hand. 64 #undef ABI_FUNC 65 #define ABI_FUNC 'distrb2' 66 use interfaces_32_util 67 !End of the abilint section 68 69 implicit none 70 71 !Arguments ------------------------------------ 72 integer,intent(in) :: mband,nkpt,nproc,nsppol 73 integer,intent(in) :: nband(nkpt*nsppol) 74 type(MPI_type),intent(inout) :: mpi_enreg 75 76 !Local variables------------------------------- 77 integer :: inb,inb1,ind,ind0,nband_k,proc_max,proc_min 78 integer :: iiband,iikpt,iisppol,ikpt_this_proc,nbsteps,nproc_kpt,temp_unit 79 integer :: kpt_distrb(nkpt) 80 logical,save :: first=.true.,has_file 81 character(len=500) :: message 82 83 !****************************************************************** 84 85 nproc_kpt=mpi_enreg%nproc_kpt 86 if (mpi_enreg%paral_pert==1) nproc_kpt=nproc 87 88 !Initialization of proc_distrb 89 do iisppol=1,nsppol 90 do iiband=1,mband 91 do iikpt=1,nkpt 92 mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=nproc_kpt-1 93 end do 94 end do 95 end do 96 !That s all for an empty communication space 97 if (nproc==0) return 98 99 !Some checks 100 if (mpi_enreg%paralbd==0) then 101 ! Check if nkpt and nproc_kpt match 102 if(nproc_kpt>nkpt*nsppol) then 103 ! Too much proc. with respect to nkpt 104 write(message,'(a,i0,a,i0,a,i0,2a)')& 105 & 'nproc_kpt=',nproc_kpt,' >= nkpt=',nkpt,'* nsppol=',nsppol,ch10,& 106 & 'The number of processors is larger than nkpt*nsppol. This is a waste.' 107 MSG_WARNING(message) 108 else if(mod(nkpt*nsppol,nproc_kpt)/=0) then 109 ! nkpt not a multiple of nproc_kpt 110 write(message,'(a,i0,a,i0,3a)')& 111 & 'nkpt*nsppol (', nkpt*nsppol, ') is not a multiple of nproc_kpt (',nproc_kpt, ')', ch10,& 112 & 'The k-point parallelisation is not efficient.' 113 MSG_WARNING(message) 114 end if 115 end if 116 117 !Inquire whether there exist a file containing the processor distribution 118 if (first) then 119 ! Case first time : test file to do 120 ! Open the file containing the k-point distribution 121 first=.false.; has_file = file_exists("kpt_distrb") 122 end if 123 124 !Initialize the processor distribution, either from a file, or from an algorithm 125 if (has_file) then 126 if (open_file('kpt_distrb',message,newunit=temp_unit,form='formatted',status='old') /= 0) then 127 MSG_ERROR(message) 128 end if 129 rewind(unit=temp_unit) 130 if (mpi_enreg%paralbd == 1) then 131 ! -> read bands distribution 132 read(temp_unit,*) mpi_enreg%proc_distrb 133 else 134 read(temp_unit,*) kpt_distrb 135 end if 136 close(temp_unit) 137 proc_max=0 138 proc_min=nproc_kpt 139 ! -> determine the range of proc. requested 140 if (mpi_enreg%paralbd == 1) then 141 do iisppol=1,nsppol 142 do iikpt=1,nkpt 143 nband_k = nband(iikpt+(iisppol-1)*nkpt) 144 proc_max=maxval(mpi_enreg%proc_distrb(iikpt,1:nband_k,iisppol)) 145 proc_min=minval(mpi_enreg%proc_distrb(iikpt,1:nband_k,iisppol)) 146 end do 147 end do 148 else 149 proc_max=maxval(kpt_distrb(1:nkpt)) 150 proc_min=minval(kpt_distrb(1:nkpt)) 151 ! -> fill the tab proc_distrb with kpt_distrb 152 do iisppol=1,nsppol 153 do iikpt=1,nkpt 154 nband_k = nband(iikpt+(iisppol-1)*nkpt) 155 do iiband=1,nband_k 156 mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=kpt_distrb(iikpt) 157 end do 158 end do 159 end do 160 end if ! mpi_enreg%paralbd 161 162 if(proc_max>(nproc_kpt-1)) then 163 ! Too much proc. requested 164 write(message, '(a,a,a,i0,a,a,a)' )& 165 & 'The number of processors mentioned in the kpt_distrb file',ch10,& 166 & 'must be lower or equal to the actual number of processors =',nproc_kpt-1,ch10,& 167 & 'Action: change the kpt_distrb file, or increase the',' number of processors.' 168 MSG_ERROR(message) 169 end if 170 171 if(proc_max/=(nproc_kpt-1)) then 172 ! Too few proc. used 173 write(message, '(a,i0,a,a,a,i0,a,a,a)' )& 174 & 'Only ',proc_max+1,' processors are used (from kpt_distrb file),',ch10,& 175 & 'when',nproc_kpt,' processors are available.',ch10,& 176 & 'Action : adjust number of processors and kpt_distrb file.' 177 MSG_ERROR(message) 178 end if 179 180 if(proc_min<0) then 181 write(message, '(a,a,a)' )& 182 & 'The number of processors must be bigger than 0 in kpt_distrb file.',ch10,& 183 & 'Action: modify kpt_distrb file.' 184 MSG_ERROR(message) 185 end if 186 187 else 188 ! 'kpt_distrb' file does not exist 189 190 if (mpi_enreg%paralbd==1) then 191 192 ! No possible band parallelization 193 if (nproc<(nkpt*nsppol)) then 194 195 ! Does not allow a processor to treat different spins 196 ind0=0 197 inb1=(nkpt*nsppol)/nproc;if (mod((nkpt*nsppol),nproc)/=0) inb1=inb1+1 198 do iikpt=1,nkpt 199 nband_k=nband(iikpt) 200 ind=ind0/inb1 201 do iiband=1,nband_k 202 mpi_enreg%proc_distrb(iikpt,iiband,1)=ind 203 if (nsppol==2) mpi_enreg%proc_distrb(iikpt,iiband,2)=nproc-ind-1 204 end do 205 ind0=ind0+1 206 end do 207 208 ! MT130831 : OLD CODING 209 ! do iisppol=1,nsppol;do iikpt=1,nkpt 210 ! ind=(iikpt+(iisppol-1)*nkpt-1)/inb1 211 ! nband_k=nband(iikpt+(iisppol-1)*nkpt) 212 ! do iiband=1,nband_k 213 ! mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=ind 214 ! end do;end do;end do 215 ! MT130831 : END OF OLD CODING 216 217 ! Possible band parallelization 218 else 219 ! Does not allow a processor to treat different spins 220 ind0=0;inb=nproc/(nkpt*nsppol) 221 do iikpt=1,nkpt 222 nband_k=nband(iikpt) 223 inb1=nband_k/inb;if (mod(nband_k,inb)/=0) inb1=inb1+1 224 do iiband=1,nband_k 225 ind=(iiband-1)/inb1+ind0 226 mpi_enreg%proc_distrb(iikpt,iiband,1)=ind 227 if (nsppol==2) mpi_enreg%proc_distrb(iikpt,iiband,2)=nproc-ind-1 228 end do 229 ind0=ind+1 230 end do 231 232 ! MT130831 : OLD CODING 233 ! ind0=0;inb=nproc/(nkpt*nsppol) 234 ! do iisppol=1,nsppol;do iikpt=1,nkpt 235 ! nband_k=nband(iikpt+(iisppol-1)*nkpt) 236 ! inb1=nband_k/inb;if (mod(nband_k,inb)/=0) inb1=inb1+1 237 ! do iiband=1,nband_k 238 ! ind=(iiband-1)/inb1+ind0 239 ! mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=ind 240 ! end do 241 ! ind0=ind+1 242 ! end do;end do 243 ! MT130831 : END OF OLD CODING 244 245 end if 246 247 ! XG060807 : OLD CODING 248 ! ind=0 249 ! do iisppol=1,nsppol;do iikpt=1,nkpt 250 ! nband_k=nband(iikpt+(iisppol-1)*nkpt) 251 ! do iiband=1,nband_k 252 ! mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=ind/nbsteps 253 ! ind=ind+1 254 ! end do;end do;end do 255 ! XG060807 : END OF OLD CODING 256 257 elseif (mpi_enreg%paralbd==0) then 258 259 ! Does not allow a processor to treat different spins 260 ind0=0 261 nbsteps=(nsppol*nkpt)/nproc_kpt 262 if (mod((nsppol*nkpt),nproc_kpt)/=0) nbsteps=nbsteps+1 263 do iikpt=1,nkpt 264 nband_k=nband(iikpt) 265 ind=ind0/nbsteps 266 do iiband=1,nband_k 267 mpi_enreg%proc_distrb(iikpt,iiband,1)=ind 268 if (nsppol==2) mpi_enreg%proc_distrb(iikpt,iiband,2)=nproc_kpt-ind-1 269 end do 270 ind0=ind0+1 271 end do 272 273 ! XG060807 : OLD CODING 274 ! ind=0 275 ! do iisppol=1,nsppol;do iikpt=1,nkpt 276 ! nband_k = nband(iikpt+(iisppol-1)*nkpt) 277 ! do iiband=1,nband_k 278 ! Distribute k-points homogeneously 279 ! proc_distrb(iikpt,iiband,iisppol)=mod(iikpt-1,nproc_kpt) 280 ! mpi_enreg%proc_distrb(iikpt,iiband,iisppol)=ind/nbsteps 281 ! end do 282 ! ind=ind + 1 283 ! end do;end do 284 ! XG060807 : END OF OLD CODING 285 286 end if ! mpi_enreg%paralbd 287 288 end if ! has_file 289 290 mpi_enreg%my_kpttab(:)=0 291 mpi_enreg%my_isppoltab(:)=0 292 do iisppol=1,nsppol 293 ikpt_this_proc=0 294 do iikpt=1,nkpt 295 nband_k=nband(iikpt+(iisppol-1)*nkpt) 296 if(proc_distrb_cycle(mpi_enreg%proc_distrb,iikpt,1,nband_k,iisppol,mpi_enreg%me_kpt)) cycle 297 ikpt_this_proc=ikpt_this_proc+1 298 ! This test should be done when dataset are read and slipt of work do between processor 299 ! If this test is not good for one proc then other procs fall in deadlock->so PERS and MPI_ABORT 300 ! if (ikpt_this_proc > mkmem) then 301 ! message = ' this bandfft tab cannot be allocated !' 302 ! MSG_BUG(message) 303 ! end if 304 mpi_enreg%my_kpttab(iikpt)=ikpt_this_proc 305 mpi_enreg%my_isppoltab(iisppol)=1 306 end do 307 end do 308 309 end subroutine distrb2