TABLE OF CONTENTS
ABINIT/ftgkk [ Functions ]
NAME
ftgkk
FUNCTION
If qtor=1 (q->r): Generates the Fourier transform of the recip space gkk matrices to obtain the real space ones. If qtor=0 (r->q): Generates the Fourier transform of the real space gkk matrices to obtain the reciprocal space ones.
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (MVer) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
gkqwrite = flag to write recip space matrix elements to disk gkrwrite = flag to write real space matrix elements to disk gprim(3,3)= Normalized coordinates in reciprocal space ikpt_phon0 = starting kpt number for forward FT. natom= Number of atoms in the unit cell nkpt_phon= Number of kpoints used for the FS ngkkband = number of bands kept in gkq and gkr matrix elements (=1 or nband) nkpt_used= number of FS kpoints used, starting at ikpt_phon0 nqpt= Number of q points in the Brillouin zone if qtor=0 this number is read in the input file nrpt= Number of R points in the Big Box qtor= ( q to r : see above ) rpt(3,nprt)= Canonical coordinates of the R points in the unit cell These coordinates are normalized (=> * acell(3)!!) qpt_full(3,nqpt)= Reduced coordinates of the q vectors in reciprocal space if qtor=0 these vectors are read in the input file unit_gkk_rpt = fortran unit for writing real-space matrix elements unitgkq = fortran unit for writing reciprocal-space matrix elements wghatm(natom,natom,nrpt) = Weights associated to a pair of atoms and to a R vector
OUTPUT
(see side effects)
SIDE EFFECTS
Input/output gkk_qpt(2,3*natom,nFSband,nFSband,nkpt_used,nqpt) = gkk matrices in recip space coming from the Derivative Data Base gkk_rpt(2,3*natom,nFSband,nFSband,nkpt_phon,nqpt) = gkk matrices in real space stored in file unit_gkk_rpt
PARENTS
get_all_gkr,interpolate_gkk,test_ftgkk
CHILDREN
NOTES
copied from ftiaf9.f recip to real space: real space is forced to disk file unit_gkk_rpt recip space depends on gkqwrite and unitgkq real to recip space: real space is forced to disk file unit_gkk_rpt recip space is necessarily in memory in gkk_qpt real space elements are complex, but could be reduced, as (-r) = (+r)*
SOURCE
70 #if defined HAVE_CONFIG_H 71 #include "config.h" 72 #endif 73 74 #include "abi_common.h" 75 76 77 subroutine ftgkk (wghatm,gkk_qpt,gkk_rpt,gkqwrite,gkrwrite,gprim,ikpt_phon0,& 78 & natom,nkpt_phon,ngkkband,nkpt_used,nqpt,nrpt,nsppol,& 79 & qtor,rpt,qpt_full,unit_gkk_rpt,unitgkq) 80 81 use defs_basis 82 use m_errors 83 use m_profiling_abi 84 use m_errors 85 86 !This section has been created automatically by the script Abilint (TD). 87 !Do not modify the following lines by hand. 88 #undef ABI_FUNC 89 #define ABI_FUNC 'ftgkk' 90 !End of the abilint section 91 92 implicit none 93 94 !Arguments ------------------------------- 95 !scalars 96 integer,intent(in) :: gkqwrite,gkrwrite,ikpt_phon0,nkpt_phon,natom,ngkkband 97 integer,intent(in) :: nkpt_used,nqpt,nrpt,nsppol,qtor,unit_gkk_rpt,unitgkq 98 !arrays 99 real(dp),intent(in) :: gprim(3,3),rpt(3,nrpt),qpt_full(3,nqpt) 100 real(dp),intent(in) :: wghatm(natom,natom,nrpt) 101 real(dp),intent(inout) :: gkk_qpt(2,ngkkband*ngkkband,3*natom*3*natom,nkpt_used,nsppol,nqpt) 102 real(dp),intent(inout) :: gkk_rpt(2,ngkkband*ngkkband,3*natom*3*natom,nkpt_used,nsppol,nrpt) 103 104 !Local variables ------------------------- 105 !scalars 106 integer :: ikpt_phon,iatom,ib1,ieffkpt_phon,ip,iqpt,irpt,isppol 107 integer :: jatom 108 real(dp) :: im,kr,re 109 character(len=500) :: message 110 !arrays 111 real(dp) :: coskr(nqpt,nrpt),ftwght(2,3*natom*3*natom) 112 real(dp) :: gkk_qpt_tmp(2,ngkkband*ngkkband,3*natom*3*natom,nkpt_used,nsppol) 113 real(dp) :: gkk_rpt_tmp(2,ngkkband*ngkkband,3*natom*3*natom,nkpt_phon,nsppol) 114 real(dp) :: kk(3),sinkr(nqpt,nrpt) 115 116 ! ********************************************************************* 117 118 !rewind (unit_gkk_rpt) 119 120 !prepare the phase factors 121 do iqpt=1,nqpt 122 ! Calculation of the k coordinates in Normalized Reciprocal 123 ! coordinates 124 kk(1)= qpt_full(1,iqpt)*gprim(1,1)+& 125 & qpt_full(2,iqpt)*gprim(1,2)+& 126 & qpt_full(3,iqpt)*gprim(1,3) 127 kk(2)= qpt_full(1,iqpt)*gprim(2,1)+& 128 & qpt_full(2,iqpt)*gprim(2,2)+& 129 & qpt_full(3,iqpt)*gprim(2,3) 130 kk(3)= qpt_full(1,iqpt)*gprim(3,1)+& 131 & qpt_full(2,iqpt)*gprim(3,2)+& 132 & qpt_full(3,iqpt)*gprim(3,3) 133 do irpt=1,nrpt 134 ! Product of k and r 135 kr = kk(1)*rpt(1,irpt)+& 136 & kk(2)*rpt(2,irpt)+& 137 & kk(3)*rpt(3,irpt) 138 coskr(iqpt,irpt)=cos(two_pi*kr) 139 sinkr(iqpt,irpt)=sin(two_pi*kr) 140 ! DEBUG 141 ! if (iqpt < 1000 .and. (irpt == 101 .or. irpt == 901)) then 142 ! write(std_out,*) iqpt,irpt,kk,rpt(:,irpt),coskr(iqpt,irpt), sinkr(iqpt,irpt) 143 ! end if 144 ! ENDDEBUG 145 end do 146 end do 147 148 149 150 !Recip to real space 151 if (qtor==1) then 152 ! 153 if (nkpt_used /= nkpt_phon) write(std_out,*) 'ftgkk: strange usage of nkpt_used for back FT!' 154 do irpt=1,nrpt 155 ! DEBUG 156 ! write(std_out,*) ' ftgkk : G->R irpt = ',irpt,' / ',nrpt 157 ! ENDDEBUG 158 gkk_rpt_tmp(:,:,:,:,:) = zero 159 160 do iqpt=1,nqpt 161 162 ! write(std_out,*) iqpt 163 164 if (gkqwrite == 0) then 165 gkk_qpt_tmp(:,:,:,:,:) = gkk_qpt(:,:,:,:,:,iqpt) 166 else 167 do ikpt_phon=1, nkpt_phon 168 read(unitgkq,REC=((iqpt-1)*nkpt_phon+ikpt_phon)) gkk_qpt_tmp(:,:,:,ikpt_phon,:) 169 end do 170 end if 171 ! Get the phase factor with normalization! 172 re=coskr(iqpt,irpt)/nqpt 173 im=sinkr(iqpt,irpt)/nqpt 174 do isppol=1,nsppol 175 do ikpt_phon=1,nkpt_used 176 ! DEBUG 177 ! write(std_out,*) ' ftgkk : G->R ikpt_phon = ',ikpt_phon,' / ',nkpt_used 178 ! ENDDEBUG 179 do ip=1,3*natom*3*natom 180 ! Real and imaginary part of the real-space gkk matrices -> exp(-i k.r) 181 do ib1=1,ngkkband*ngkkband 182 gkk_rpt_tmp(1,ib1,ip,ikpt_phon,isppol) = gkk_rpt_tmp(1,ib1,ip,ikpt_phon,isppol)& 183 & +re*gkk_qpt_tmp(1,ib1,ip,ikpt_phon,isppol) & 184 & +im*gkk_qpt_tmp(2,ib1,ip,ikpt_phon,isppol) 185 gkk_rpt_tmp(2,ib1,ip,ikpt_phon,isppol) = gkk_rpt_tmp(2,ib1,ip,ikpt_phon,isppol)& 186 & +re*gkk_qpt_tmp(2,ib1,ip,ikpt_phon,isppol) & 187 & -im*gkk_qpt_tmp(1,ib1,ip,ikpt_phon,isppol) 188 end do 189 end do 190 end do 191 end do 192 end do 193 if (gkrwrite == 0) then 194 gkk_rpt(:,:,:,:,:,irpt) = gkk_rpt_tmp(:,:,:,:,:) 195 else 196 write (unit_gkk_rpt,REC=irpt) gkk_rpt_tmp 197 end if 198 end do 199 200 ! Real space to recip space 201 else if (qtor==0) then 202 203 ! write(std_out,*) 'ftgkk : shape(gkk_qpt) = ', shape(gkk_qpt) 204 gkk_qpt(:,:,:,:,:,:)=zero 205 206 ! rewind (unit_gkk_rpt) 207 do irpt=1,nrpt 208 if (gkrwrite == 0) then 209 gkk_rpt_tmp(:,:,:,:,:) = gkk_rpt(:,:,:,:,:,irpt) 210 else 211 read(unit_gkk_rpt,REC=irpt) gkk_rpt_tmp 212 end if 213 214 215 do iqpt=1,nqpt 216 217 ! Avoid recalculating weights nkpt_used*9 times 218 do iatom=1,natom 219 do jatom=1,natom 220 ip = 3*((iatom-1)*natom+jatom-1) 221 ! copy same weight for all 3 directions 222 ftwght(1,ip+1:ip+3)=coskr(iqpt,irpt)*wghatm(iatom,jatom,irpt) 223 ftwght(2,ip+1:ip+3)=sinkr(iqpt,irpt)*wghatm(iatom,jatom,irpt) 224 end do 225 end do 226 227 228 229 do ip=1,3*natom*3*natom 230 ! Get phase factor 231 re = ftwght(1,ip) 232 im = ftwght(2,ip) 233 234 do isppol=1,nsppol 235 do ikpt_phon=1,nkpt_used 236 237 238 ! DEBUG 239 ! write(std_out,*) ' ftgkk : R->G ikpt_phon = ',ikpt_phon,' / ',nkpt_used 240 ! ENDDEBUG 241 ! effective FS kpt in real space array is ikpt_phon+ikpt_phon0-1 to allow for offset 242 ieffkpt_phon = ikpt_phon+ikpt_phon0-1 243 ! write(std_out,*) 'ftgkk :ikpt_phon,iqpt,ieffkpt_phon ', ikpt_phon,iqpt,ieffkpt_phon 244 245 do ib1=1,ngkkband*ngkkband 246 ! Real and imaginary part of the gamma matrices 247 gkk_qpt(1,ib1,ip,ikpt_phon,isppol,iqpt)=& 248 & gkk_qpt(1,ib1,ip,ikpt_phon,isppol,iqpt)& 249 & +re*gkk_rpt_tmp(1,ib1,ip,ieffkpt_phon,isppol)& 250 & -im*gkk_rpt_tmp(2,ib1,ip,ieffkpt_phon,isppol) 251 ! !DEBUG 252 gkk_qpt(2,ib1,ip,ikpt_phon,isppol,iqpt)=& 253 & gkk_qpt(2,ib1,ip,ikpt_phon,isppol,iqpt)& 254 & +im*gkk_rpt_tmp(1,ib1,ip,ieffkpt_phon,isppol)& 255 & +re*gkk_rpt_tmp(2,ib1,ip,ieffkpt_phon,isppol) 256 ! !ENDDEBUG 257 258 ! if (iqpt < 100 .and. irpt < 100 .and. & 259 ! & tmpgkkrim(irpt)**2+tmpgkkrre(irpt)**2 > tol6) then 260 ! write(std_out,'(2I4,2E16.8,x,2E16.8)') & 261 ! & iqpt,irpt,re,im,tmpgkkrre(irpt),tmpgkkrim(irpt) 262 ! end if 263 264 end do 265 end do 266 ! end ikpt_phon 267 end do 268 ! end isppol 269 ! write(std_out,'(a)') ' ftgkk :gkk_qpt :' 270 ! write(std_out,'(4E16.5)') gkk_qpt(:,1,1,,ikpt_phon,1:nqpt) 271 end do 272 ! end ip 273 end do 274 ! end iqpt 275 end do 276 ! end irpt 277 278 279 ! There is no other space to Fourier transform from ?? 280 else 281 write(message,'(a,a,a,i0,a)' )& 282 & 'The only allowed values for qtor are 0 or 1, while',ch10,& 283 & 'qtor=',qtor,' has been required.' 284 MSG_BUG(message) 285 end if 286 287 end subroutine ftgkk