TABLE OF CONTENTS
ABINIT/distrb2_hf [ Functions ]
NAME
distrb2_hf
FUNCTION
This routine creates the tabs of repartition of processors for sharing the jobs on occupied states (labeled by k-points, bands and spin indices) for an Hartree-Fock calculation.
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (CMartins) 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
nbandhf = maximum number of occupied bands nkpthf = number of k-points in full BZ 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(nkpthf,nbandhf,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
mpi_setup
CHILDREN
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 54 subroutine distrb2_hf(nbandhf,nkpthf, nproc, nsppol, mpi_enreg) 55 56 use defs_basis 57 use defs_abitypes 58 use m_errors 59 use m_profiling_abi 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'distrb2_hf' 65 !End of the abilint section 66 67 implicit none 68 69 !Arguments ------------------------------------ 70 integer,intent(in) :: nbandhf,nkpthf,nproc,nsppol 71 type(MPI_type),intent(inout) :: mpi_enreg 72 73 !Local variables------------------------------- 74 integer :: ind,iiband,iikpt,iistep,nproc_hf 75 character(len=500) :: message 76 77 !****************************************************************** 78 79 nproc_hf=mpi_enreg%nproc_hf 80 81 !* Initialize distrb_hf (the array exists necessarily) 82 do iiband=1,nbandhf 83 do iikpt=1,nkpthf 84 mpi_enreg%distrb_hf(iikpt,iiband,1)=nproc_hf-1 85 end do 86 end do 87 88 !* End of the routine for an empty communication space 89 if (nproc==0) return 90 91 !*** Testing section *** 92 93 if (nsppol==2) then 94 !* Check that the distribution over (spin,k point) allow to consider spin up and spin dn independently. 95 if (mpi_enreg%nproc_kpt/=1.and.mod(mpi_enreg%nproc_kpt,2)/=0) then 96 MSG_ERROR( 'The variable nproc_kpt is not even but nssppol= 2') 97 !* In this case, one processor will carry both spin. (this will create pbms for the calculation) 98 end if 99 !* Check that the number of band is the same for each spin, at each k-point. (It should be) 100 !* do iikpt=1,nkpthf 101 !* if (nband(iikpt)/=nband(iikpt+nkpthf)) then 102 !* message = ' WARNING - the number of bands is different for spin up or spin down. ' 103 !* MSG_ERROR(message) 104 !* end if 105 !* end do 106 !* If one of this test is not good for one proc then other procs fall in deadlock, according to distrb2. 107 !* What does it mean ??? 108 end if 109 110 111 !* Check if nkpthf and nproc_hf match 112 if (nproc_hf>nkpthf*nbandhf) then 113 !* There are too many processors with respect to nkpthf*nbandhf 114 write(message, '(a,a,i4,a,i4,a,i4,a,a)' ) ch10,& 115 & 'nproc_hf=',nproc_hf,' >= nkpthf=',nkpthf,'* nbandhf=',nbandhf,ch10,& 116 & 'The number of processors is larger than nkpthf*nbandhf. This is a waste.' 117 MSG_WARNING(message) 118 119 else if(mod(nkpthf*nbandhf,nproc_hf)/=0) then 120 !* nkpthf*nbandhf is not a multiple of nproc_hf 121 write(message, '(2a,i5,a,i5,3a)' ) ch10,& 122 & 'nkpthf*nbandhf (', nkpthf*nbandhf, ') is not a multiple of nproc_hf (',nproc_hf, ')', ch10,& 123 & 'The parallelisation may not be efficient.' 124 MSG_WARNING(message) 125 end if 126 127 !*** End of testing section *** 128 129 !*** Initialize the processor distribution from a simple algorithm *** 130 131 if (nproc_hf<nkpthf) then 132 !* In this case, a parallelization over kpts only. 133 iistep=nkpthf/nproc_hf 134 if (mod(nkpthf,nproc_hf) /=0) iistep=iistep+1 135 ind=0 136 do iikpt=1,nkpthf 137 !*** Only the first "nbandhf" bands are considered (they are assumed to be the only occupied ones) 138 do iiband=1,nbandhf 139 mpi_enreg%distrb_hf(iikpt,iiband,1)=ind/iistep 140 end do 141 ind=ind+1 142 end do 143 144 else 145 !* In this case, a parallelization over all the occupied states is possible. 146 if (nproc_hf < nbandhf*nkpthf) then 147 iistep=(nbandhf*nkpthf)/nproc_hf; 148 if (mod((nbandhf*nkpthf),nproc_hf) /=0) iistep=iistep+1 149 else 150 iistep=1 151 end if 152 ind=0 153 do iikpt=1,nkpthf 154 !*** Only the first "nbandhf" bands are considered (they are assumed to be the only occupied ones) 155 do iiband=1,nbandhf 156 mpi_enreg%distrb_hf(iikpt,iiband,1)=ind/iistep 157 ind=ind+1 158 end do 159 end do 160 end if 161 162 !*** Initialization of processor distribution from a file (simple copy from distrb2, not yet implemented) *** 163 164 ! !* Inquire whether there exist a file containing the processor distribution 165 ! if (first) then 166 ! ! Case first time : test file to do 167 ! ! Open the file containing the (k-points,bands) distribution 168 ! open(unit=temp_unit,file='kpt_distrb_hf',form='formatted',status='old',iostat=ios) 169 ! if(ios==0) then 170 ! ! 'kpt_distrb_hf' file exists 171 ! file_exist=1 172 ! close(temp_unit) 173 ! else 174 ! file_exist=0 175 ! end if 176 ! first=.false. 177 ! end if 178 ! 179 ! !* Initialize the processor distribution, either from a file, or from an algorithm 180 ! if (file_exist == 1) then 181 ! !* Read (k-points,bands) distribution out of the file 182 ! open(unit=temp_unit,file='kpt_distrb_hf',form='formatted',status='old',iostat=ios) 183 ! rewind(unit=temp_unit) 184 ! read(temp_unit,*) mpi_enreg%distrb_hf 185 ! close(temp_unit) 186 ! !* Determine the range of processors requested 187 ! proc_max=0 188 ! proc_min=nproc_hf 189 ! do iikpt=1,nkpthf 190 ! mband_occ_k = mband_occ(iikpt+(iisppol-1)*nkpthf) 191 ! proc_max=maxval(mpi_enreg%distrb_hf(iikpt,1:mband_occ_k,1)) 192 ! proc_min=minval(mpi_enreg%distrb_hf(iikpt,1:mband_occ_k,1)) 193 ! end do 194 ! 195 ! if(proc_max>(nproc_hf-1)) then 196 ! !* Too much proc. requested 197 ! write(message, '(a,a,a,i4,a,a,a)' )& 198 ! & ' The number of processors mentioned in the kpt_distrb file',ch10,& 199 ! & ' must be lower or equal to the actual number of processors =',& 200 ! & nproc_hf-1,ch10,& 201 ! & ' Action : change the kpt_distrb file, or increase the',& 202 ! & ' number of processors.' 203 ! MSG_ERROR(message) 204 ! end if 205 ! 206 ! if(proc_max/=(nproc_hf-1)) then 207 ! !* Too few proc. used 208 ! write(message, '(a,i4,a,a,a,i4,a,a,a)' )& 209 ! & ' Only ',proc_max+1,' processors are used (from kpt_distrb file),',ch10,& 210 ! & ' when',nproc_hf,' processors are available.',ch10,& 211 ! & ' Action : adjust number of processors and kpt_distrb file.' 212 ! MSG_ERROR(message) 213 ! end if 214 ! 215 ! if(proc_min<0) then 216 ! write(message, '(a,a,a)' )& 217 ! & ' The number of processors must be bigger than 0 in kpt_distrb file.',ch10,& 218 ! & ' Action : modify kpt_distrb file.' 219 ! MSG_ERROR(message) 220 ! end if 221 ! else 222 ! !* The file does not exist... 223 ! end if ! file_exist 224 225 !DEBUG 226 !write(std_out,*)' distrb2_hf: exit ' 227 !ENDDEBUG 228 229 end subroutine distrb2_hf