TABLE OF CONTENTS
ABINIT/calc_coh [ Functions ]
NAME
calc_coh
FUNCTION
Calculates the partial contribution to the COH part of the COHSEX self-energy for a given q-point.
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (FB,MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
iqibz=index of the irreducible q-point in the array qibz, point which is related by a symmetry operation to the point q summed over (see csigme). This index is also used to treat the integrable coulombian singularity at q=0 ngfft(18)=contain all needed information about 3D FFT for GW wavefuntions, see ~abinit/doc/variables/vargs.htm#ngfft nsig_ab=Number of components in the self-energy operator (1 for collinear magnetism) npwc=number of plane waves in $\tilde epsilon^{-1}$ nspinor=Number of spinorial components. nfftot=number of points in real space i_sz=contribution arising from the integrable coulombian singularity at q==0 (see csigme for the method used), note that in case of 3-D systems the factor 4pi in the coulombian potential is included in the definition of i_sz gvec(3,npwc)=G vectors in reduced coordinates vc_sqrt(npwc)= square root of the coulombian matrix elements for this q-point epsm1q_o(npwc,npwc)= contains $\tilde epsilon^{-1}(q,w=0) - \delta_{G Gp}$ for the particular q-point considered in the sum wfg2_jk(nsig_ab*nfftot)= Fourier Transform of $\u_{jb k}^*(r) u_{kb k}$ jb,kb=left and righ band indeces definining the left and right states where the partial contribution to the matrix element of $\Sigma_{COH}$ is evaluated
OUTPUT
sigcohme=partial contribution to the matrix element of $<jb k \sigma|\Sigma_{COH} | kb k \sigma>$ coming from this single q-point
PARENTS
cohsex_me
CHILDREN
wrtout
SOURCE
48 #if defined HAVE_CONFIG_H 49 #include "config.h" 50 #endif 51 52 #include "abi_common.h" 53 54 subroutine calc_coh(nspinor,nsig_ab,nfftot,ngfft,npwc,gvec,wfg2_jk,epsm1q_o,vc_sqrt,i_sz,iqibz,same_band,sigcohme) 55 56 use defs_basis 57 use m_profiling_abi 58 use m_errors 59 60 use m_fstrings, only : sjoin, itoa 61 use m_gwdefs, only : czero_gw 62 63 !This section has been created automatically by the script Abilint (TD). 64 !Do not modify the following lines by hand. 65 #undef ABI_FUNC 66 #define ABI_FUNC 'calc_coh' 67 use interfaces_14_hidewrite 68 !End of the abilint section 69 70 implicit none 71 72 !Arguments ------------------------------------ 73 !scalars 74 integer,intent(in) :: iqibz,nfftot,npwc,nsig_ab,nspinor 75 real(dp),intent(in) :: i_sz 76 logical,intent(in) :: same_band 77 !arrays 78 integer,intent(in) :: gvec(3,npwc),ngfft(18) 79 complex(gwpc),intent(in) :: epsm1q_o(npwc,npwc),vc_sqrt(npwc) 80 complex(gwpc),intent(in) :: wfg2_jk(nfftot*nsig_ab) 81 complex(gwpc),intent(out) :: sigcohme(nsig_ab) 82 83 !Local variables------------------------------- 84 !scalars 85 integer,save :: enough=0 86 integer :: ig,ig4,ig4x,ig4y,ig4z,igp,igmin,ispinor,spad,outofbox 87 !arrays 88 integer :: g2mg1(3) 89 90 ! ************************************************************************* 91 92 DBG_ENTER("COLL") 93 94 ! === Partial contribution to the matrix element of Sigma_c === 95 ! * For nspinor==2, the closure relation reads: 96 ! $\sum_s \psi_a^*(1)\psi_b(2) = \delta_{ab} \delta(1-2)$ 97 ! where a,b are the spinor components. As a consequence, Sigma_{COH} is always 98 ! diagonal in spin-space and only diagonal matrix elements have to be calculated. 99 ! MG TODO wfg2_jk should be calculated on an augmented FFT box to avoid spurious wrapping of G1-G2. 100 ! MG: One has to make sure G1-G2 is still in the FFT mesh for each G1 and G2 in chi0 (not always true) 101 ! MODULO wraps G1-G2 in the FFT box but the Fourier components are not periodic! 102 103 ! * Treat the case q --> 0 adequately. 104 ! TODO Better treatment of wings, check cutoff in the coulombian interaction. 105 igmin=1; if (iqibz==1) igmin=2 106 107 sigcohme(:)=czero_gw 108 109 do ispinor=1,nspinor 110 spad=(ispinor-1)*nfftot 111 outofbox=0 112 113 do igp=igmin,npwc 114 do ig=igmin,npwc 115 116 g2mg1 = gvec(:,igp)-gvec(:,ig) 117 if (ANY(g2mg1(:)>ngfft(1:3)/2) .or. ANY(g2mg1(:)<-(ngfft(1:3)-1)/2)) then 118 outofbox = outofbox+1; CYCLE 119 end if 120 121 ig4x=MODULO(g2mg1(1),ngfft(1)) 122 ig4y=MODULO(g2mg1(2),ngfft(2)) 123 ig4z=MODULO(g2mg1(3),ngfft(3)) 124 ig4= 1+ig4x+ig4y*ngfft(1)+ig4z*ngfft(1)*ngfft(2) 125 126 sigcohme(ispinor) = sigcohme(ispinor) + & 127 & half*wfg2_jk(spad+ig4)*epsm1q_o(ig,igp)*vc_sqrt(ig)*vc_sqrt(igp) 128 end do !ig 129 end do !igp 130 131 if (iqibz==1.and.same_band) then 132 sigcohme(ispinor) = sigcohme(ispinor) + half*wfg2_jk(spad+ig4)*epsm1q_o(1,1)*i_sz 133 end if 134 end do !ispinor 135 136 if (outofbox/=0) then 137 enough=enough+1 138 if (enough<=50) then 139 MSG_WARNING(sjoin(' Number of G1-G2 pairs outside the G-sphere for Wfns:', itoa(outofbox))) 140 if (enough==50) then 141 call wrtout(std_out,' ========== Stop writing Warnings ==========') 142 end if 143 end if 144 end if 145 146 DBG_EXIT("COLL") 147 148 end subroutine calc_coh
ABINIT/calc_coh_comp [ Functions ]
NAME
calc_coh_comp
FUNCTION
Calculates the COH-like contribution to the self-energy when the extrapolar technique and the closure relation is used to reduce the number of empty states to be summed over in the Green function entering the definition of the GW self-energy.
COPYRIGHT
Copyright (C) 2005-2018 ABINIT group (FB,MG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
iqibz=index of the irreducible q-point in the array qibz, point which is related by a symmetry operation to the point q summed over (see csigme). This index is also used to treat the integrable coulombian singularity at q=0 ngfft(18)=contain all needed information about 3D FFT for GW wavefuntions, see ~abinit/doc/variables/vargs.htm#ngfft nsig_ab=Number of components in the self-energy operator (1 for collinear magnetism) npwc=number of plane waves in $\tilde epsilon^{-1}$ nspinor=Number of spinorial components. i_sz=contribution arising from the integrable coulombian singularity at q==0 (see csigme for the method used), note that in case of 3-D systems the factor 4pi in the coulombian potential is included in the definition of i_sz gvec(3,npwc)=G vectors in reduced coordinates vc_sqrt(npwc)= square root of the coulombian matrix elements for this q-point botsq = Plasmon-pole parameters otq = PPm parameters
OUTPUT
sigcohme=partial contribution to the matrix element of $<jb k|\Sigma_{COH}| kb k>$ coming from this single q-point for completeness trick
SIDE EFFECTS
PARENTS
calc_sigc_me
CHILDREN
wrtout
SOURCE
200 subroutine calc_coh_comp(iqibz,i_sz,same_band,nspinor,nsig_ab,ediff,npwc,gvec,& 201 & ngfft,nfftot,wfg2_jk,vc_sqrt,botsq,otq,sigcohme) 202 203 use defs_basis 204 use m_profiling_abi 205 use m_errors 206 207 use m_fstrings, only : sjoin, itoa 208 use m_gwdefs, only : czero_gw 209 210 !This section has been created automatically by the script Abilint (TD). 211 !Do not modify the following lines by hand. 212 #undef ABI_FUNC 213 #define ABI_FUNC 'calc_coh_comp' 214 use interfaces_14_hidewrite 215 !End of the abilint section 216 217 implicit none 218 219 !Arguments ------------------------------------ 220 !scalars 221 integer,intent(in) :: iqibz,npwc,nsig_ab,nspinor,nfftot 222 real(dp),intent(in) :: i_sz,ediff 223 logical,intent(in) :: same_band 224 !arrays 225 integer,intent(in) :: gvec(3,npwc),ngfft(18) 226 complex(gwpc),intent(in) :: botsq(npwc,npwc),otq(npwc,npwc) 227 complex(gwpc),intent(in) :: vc_sqrt(npwc) 228 complex(gwpc),intent(in) :: wfg2_jk(nfftot*nsig_ab) 229 complex(gwpc),intent(out) :: sigcohme(nsig_ab) 230 231 !Local variables------------------------------- 232 !scalars 233 integer,save :: enough=0 234 integer :: ig,ig4,ig4x,ig4y,ig4z,igp,igmin,ispinor,ngfft1,ngfft2,ngfft3,spad,outofbox 235 !arrays 236 integer :: g2mg1(3) 237 238 ! ************************************************************************* 239 240 DBG_ENTER("COLL") 241 242 ! === Treat the case q --> 0 adequately === 243 ! TODO Better treatment of wings 244 igmin=1 ; if (iqibz==1) igmin=2 245 ! 246 ! === Partial contribution to the matrix element of Sigma_c === 247 ! * For nspinor==2, the closure relation reads: 248 ! $\sum_s \psi_a^*(1)\psi_b(2) = \delta_{ab} \delta(1-2)$ 249 ! where a,b are the spinor components. As a consequence, Sigma_{COH} is always 250 ! diagonal in spin-space and only diagonal matrix elements have to be calculated. 251 ! MG TODO wfg2_jk should be calculated on an augmented FFT box to avoid spurious wrapping of G1-G2. 252 ! 253 ngfft1 = ngfft(1); ngfft2 = ngfft(2); ngfft3 = ngfft(3) 254 sigcohme(:) = czero_gw 255 256 do ispinor=1,nspinor 257 spad=(ispinor-1)*nfftot 258 outofbox=0 259 260 do igp=igmin,npwc 261 do ig=igmin,npwc 262 263 g2mg1 = gvec(:,igp)-gvec(:,ig) 264 if (ANY(g2mg1(:)>ngfft(1:3)/2) .or. ANY(g2mg1(:)<-(ngfft(1:3)-1)/2)) then 265 outofbox = outofbox+1; CYCLE 266 end if 267 268 ig4x=MODULO(g2mg1(1),ngfft1) 269 ig4y=MODULO(g2mg1(2),ngfft2) 270 ig4z=MODULO(g2mg1(3),ngfft3) 271 ig4= 1+ig4x+ig4y*ngfft1+ig4z*ngfft1*ngfft2 272 273 !MG where is neta here, ediff, otq might be close to zero depending on gwecomp 274 sigcohme(ispinor) = sigcohme(ispinor) + & 275 & half*wfg2_jk(spad+ig4)*vc_sqrt(ig)*vc_sqrt(igp) * botsq(ig,igp) / ( otq(ig,igp) * ( ediff -otq(ig,igp) ) ) 276 end do 277 end do 278 279 if (iqibz==1.and.same_band) then 280 sigcohme(ispinor) = sigcohme(ispinor) + half*wfg2_jk(spad+ig4)*i_sz*botsq(1,1) / ( otq(1,1) * (ediff -otq(1,1)) ) 281 end if 282 end do !ispinor 283 284 if (outofbox/=0) then 285 enough=enough+1 286 if (enough<=50) then 287 MSG_WARNING(sjoin('Number of G1-G2 pairs outside the G-sphere for Wfns: ',itoa(outofbox))) 288 if (enough==50) then 289 call wrtout(std_out,' ========== Stop writing Warnings ==========') 290 end if 291 end if 292 end if 293 294 DBG_EXIT("COLL") 295 296 end subroutine calc_coh_comp