TABLE OF CONTENTS
ABINIT/initmv [ Functions ]
NAME
initmv
FUNCTION
Initialize finite difference calculation of the ddk im dfptnl_mv.f
COPYRIGHT
Copyright (C) 1999-2018 ABINIT group (MVeithen) 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
dtset <type(dataset_type)> = all input variables in this dataset gmet(3,3) = reciprocal space metric tensor in bohr**-2 kg(3,mpw*mkmem) = reduced (integer) coordinates of G vecs in basis sphere kneigh(30,nkpt2) = index of the neighbours of each k-point kg_neigh(30,nkpt2,3) = necessary to construct the vector joining a k-point to its nearest neighbour in case of a single k-point, a line of k-points or a plane of k-points. See getshell.F90 for details kptindex(2,nkpt3)= index of the k-points in the reduced BZ related to a k-point in the full BZ kpt3(3,nkpt3) = reduced coordinates of k-points in the full BZ mband = maximum number of bands mkmem = number of k points which can fit in memory mpi_enreg = informations about MPI parallelization mpw = maximum number of plane waves nband(nkpt*nsppol)=number of bands at each k point, for each polarization nkpt2 = number of k-points in the reduced BZ nkpt3 = number of k-points in the full BZ nneigh = total number of neighbours required to evaluate the finite difference formula npwarr(nkpt2)=number of planewaves at each k point nsppol = number of spin polarizations occ(mband*nkpt*nsppol) = occupation number for each band for each k
OUTPUT
cgindex(nkpt2,nsppol) = for each k-point, cgindex tores the location of the WF in the cg array me = index of the current processor ineigh = index of a neighbour ikpt_loc = index of the iteration on ikpt on the current processor ikpt_rbz = index of a k-point in the reduced BZ pwind(mpw,nneigh,mkmem) = array used to compute the overlap matrix smat between k-points (see initberry.f for more explanations) COMMENTS
PARENTS
nonlinear
CHILDREN
kpgio,xmpi_sum
SOURCE
62 #if defined HAVE_CONFIG_H 63 #include "config.h" 64 #endif 65 66 #include "abi_common.h" 67 68 subroutine initmv(cgindex,dtset,gmet,kg,kneigh,kg_neigh,kptindex,& 69 & kpt3,mband,mkmem,mpi_enreg,mpw,nband,nkpt2,& 70 & nkpt3,nneigh,npwarr,nsppol,occ,pwind) 71 72 use defs_basis 73 use defs_abitypes 74 use m_xmpi 75 use m_errors 76 use m_profiling_abi 77 78 use m_kg, only : kpgio 79 80 !This section has been created automatically by the script Abilint (TD). 81 !Do not modify the following lines by hand. 82 #undef ABI_FUNC 83 #define ABI_FUNC 'initmv' 84 use interfaces_32_util 85 !End of the abilint section 86 87 implicit none 88 89 !Arguments ------------------------------------ 90 !scalars 91 integer,intent(in) :: mband,mkmem,mpw,nkpt2,nkpt3,nneigh,nsppol 92 type(MPI_type),intent(inout) :: mpi_enreg 93 type(dataset_type),intent(in) :: dtset 94 !arrays 95 integer,intent(in) :: kg(3,mpw*mkmem),kneigh(30,nkpt2),kg_neigh(30,nkpt2,3) 96 integer,intent(in) :: nband(nkpt2*nsppol),npwarr(nkpt2),kptindex(2,nkpt3) 97 integer,intent(out) :: cgindex(nkpt2,nsppol),pwind(mpw,nneigh,mkmem) 98 real(dp),intent(in) :: gmet(3,3),kpt3(3,nkpt3),occ(mband*nkpt2*nsppol) 99 100 !Local variables------------------------------- 101 !scalars 102 integer :: flag,iband,icg,ierr,ikg,ikg1,ikpt,ikpt2,ikpt_loc,ikpt_rbz 103 integer :: index,ineigh,ipw,isppol,jpw,nband_k,mband_occ,mband_occ_k,npw_k 104 integer :: npw_k1,orig,spaceComm 105 real(dp) :: ecut_eff,sdeg 106 character(len=500) :: message 107 !arrays 108 integer :: dg(3) 109 integer,allocatable :: kg1(:,:),kg1_k(:,:),npwar1(:),npwtot(:) 110 real(dp) :: dk(3),dk_(3) 111 real(dp),allocatable :: kpt1(:,:) 112 113 !************************************************************************ 114 115 !DEBUG 116 !write(std_out,*)' initmv : enter ' 117 !ENDDEBUG 118 119 if (xmpi_paral== 1) then 120 ! BEGIN TF_CHANGES 121 spaceComm=mpi_enreg%comm_cell 122 ! END TF_CHANGES 123 mpi_enreg%kpt_loc2ibz_sp(:,:,:) = 0 124 mpi_enreg%mkmem(:) = 0 125 end if 126 127 ecut_eff = dtset%ecut*(dtset%dilatmx)**2 128 ABI_ALLOCATE(kg1_k,(3,mpw)) 129 ABI_ALLOCATE(kg1,(3,mkmem*mpw)) 130 ABI_ALLOCATE(kpt1,(3,nkpt2)) 131 ABI_ALLOCATE(npwar1,(nkpt2)) 132 ABI_ALLOCATE(npwtot,(nkpt2)) 133 kg1_k(:,:) = 0 134 pwind(:,:,:) = 0 135 cgindex(:,:) = 0 136 137 !Compute the number of occupied bands. 138 !Check that it is the same for every k-point and that 139 !nband(ikpt) is equal to this value 140 141 if (nsppol == 1) then 142 sdeg = two 143 else if (nsppol == 2) then 144 sdeg = one 145 end if 146 147 !DEBUG 148 !write(std_out,*)' list of nband ' 149 !do isppol = 1, nsppol 150 !do ikpt = 1, nkpt2 151 ! 152 !nband_k = nband(ikpt + (isppol - 1)*nkpt2) 153 !write(std_out,*)' isppol, ikpt, nband_k=',isppol, ikpt, nband_k 154 !end do 155 !end do 156 !ENDDEBUG 157 158 159 160 161 162 163 164 index = 0 165 do isppol = 1, nsppol 166 do ikpt = 1, nkpt2 167 168 mband_occ_k = 0 169 nband_k = nband(ikpt + (isppol - 1)*nkpt2) 170 171 do iband = 1, nband_k 172 index = index + 1 173 if (abs(occ(index) - sdeg) < tol8) mband_occ_k = mband_occ_k + 1 174 end do 175 176 if (nband_k /= mband_occ_k) then 177 write(message,'(a,a,a)')& 178 & ' In a non-linear response calculation, nband must be equal ',ch10,& 179 & ' to the number of valence bands.' 180 MSG_ERROR(message) 181 end if 182 183 ! Note that the number of bands can be different for spin up and spin down 184 if (ikpt > 1) then 185 if (mband_occ /= mband_occ_k) then 186 message = 'The number of valence bands is not the same for every k-point' 187 MSG_ERROR(message) 188 end if 189 else 190 mband_occ = mband_occ_k 191 end if 192 193 end do ! close loop over ikpt 194 end do ! close loop over isppol 195 196 !Find the location of each wavefunction 197 198 icg = 0 199 do isppol = 1, nsppol 200 do ikpt = 1, nkpt2 201 202 203 ! fab: inserted the shift due to the spin... 204 205 nband_k = dtset%nband(ikpt+(isppol - 1)*nkpt2) 206 npw_k = npwarr(ikpt) 207 208 if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,isppol,mpi_enreg%me)) cycle 209 210 cgindex(ikpt,isppol) = icg 211 icg = icg + dtset%nspinor*npw_k*nband_k 212 213 end do 214 end do 215 216 217 !Build pwind 218 219 do ineigh = 1, nneigh 220 221 do ikpt = 1, nkpt2 222 ikpt2 = kneigh(ineigh,ikpt) 223 ikpt_rbz = kptindex(1,ikpt2) ! index of the k-point in the reduced BZ 224 kpt1(:,ikpt) = dtset%kptns(:,ikpt_rbz) 225 end do 226 227 ! Set up the basis sphere of plane waves at kpt1 228 kg1(:,:) = 0 229 call kpgio(ecut_eff,dtset%exchn2n3d,gmet,dtset%istwfk,kg1,& 230 & kpt1,mkmem,dtset%nband,nkpt2,'PERS',mpi_enreg,mpw,& 231 & npwar1,npwtot,dtset%nsppol) 232 233 ikg = 0 ; ikg1 = 0 ; ikpt_loc = 0 234 235 if(dtset%nsppol/=1)then 236 if(mpi_enreg%nproc/=1)then 237 message = ' At present, non-linear response calculations for spin-polarized system cannot be done in parallel.' 238 MSG_ERROR(message) 239 else 240 isppol=1 241 end if 242 else 243 isppol=1 244 end if 245 246 do ikpt = 1, nkpt2 247 248 nband_k = dtset%nband(ikpt+(isppol - 1)*nkpt2) 249 ikpt2 = kneigh(ineigh,ikpt) 250 ikpt_rbz = kptindex(1,ikpt2) ! index of the k-point in the reduced BZ 251 252 if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,-1,mpi_enreg%me)) cycle 253 254 ikpt_loc = ikpt_loc + 1 255 256 mpi_enreg%kpt_loc2ibz_sp(mpi_enreg%me, ikpt_loc, 1) = ikpt 257 258 flag = 0 259 npw_k = npwarr(ikpt) 260 npw_k1 = npwarr(ikpt_rbz) 261 dk_(:) = kpt3(:,ikpt2) - dtset%kptns(:,ikpt) 262 dk(:) = dk_(:) - nint(dk_(:)) + real(kg_neigh(ineigh,ikpt,:),dp) 263 dg(:) = nint(dk(:) - dk_(:)) 264 265 266 if (kptindex(2,ikpt2) == 0) then 267 kg1_k(:,1:npw_k1) = kg1(:,ikg1+1:ikg1+npw_k1) 268 if (dg(1)==0.and.dg(2)==0.and.dg(3)==0) flag = 1 269 else 270 kg1_k(:,1:npw_k1) = -1*kg1(:,ikg1+1:ikg1+npw_k1) 271 end if 272 273 orig = 1 274 do ipw = 1, npw_k 275 do jpw = orig, npw_k1 276 277 if ((kg(1,ikg + ipw) == kg1_k(1,jpw) - dg(1)).and. & 278 & (kg(2,ikg + ipw) == kg1_k(2,jpw) - dg(2)).and. & 279 & (kg(3,ikg + ipw) == kg1_k(3,jpw) - dg(3))) then 280 281 pwind(ipw,ineigh,ikpt_loc) = jpw 282 if (flag == 1) orig = jpw + 1 283 exit 284 285 end if 286 287 end do 288 end do 289 290 ikg = ikg + npw_k 291 ikg1 = ikg1 + npw_k1 292 293 end do ! close loop over k-points 294 295 end do ! close loop over ineigh 296 mpi_enreg%mkmem(mpi_enreg%me) = mkmem 297 298 call xmpi_sum(mpi_enreg%kpt_loc2ibz_sp,spaceComm,ierr) 299 call xmpi_sum(mpi_enreg%mkmem,spaceComm,ierr) 300 301 302 !---------------------------------------------------------------------------- 303 304 ABI_DEALLOCATE(kg1) 305 ABI_DEALLOCATE(kg1_k) 306 ABI_DEALLOCATE(kpt1) 307 ABI_DEALLOCATE(npwar1) 308 ABI_DEALLOCATE(npwtot) 309 310 end subroutine initmv