TABLE OF CONTENTS
ABINIT/accumulate_chi0sumrule [ Functions ]
NAME
accumulate_chi0sumrule
FUNCTION
Accumulate the contribution to the sum rule for Im chi0 arising from a single transition. Eventually symmetrize it using the symmetry operations of the little group of q.
INPUTS
ik_bz=Index of the k-point in the full BZ whose contribution has to be added and symmetrized. symchi=1 if we are summing over IBZ_q and symmetrization has to be performed. npwe=Number of G vectors in chi0. npwepG0=Number of G vectors in the "enlarged" sphere to treat umklapp. factor=factor entering the expression. delta_ene=Transition energy. Ltg_q=<littlegroup_t>= Structure gathering information on the little group of the external q. Gsph_epsG0=<gsphere_t>=Info on the G-sphere for chi0. rhotwg(npwepG0)=Fouriet transform of u_{b1 k-q} u_{b2 k} in the "enlarged" sphere.
OUTPUT
See SIDES EFFECTS SIDES EFFECTS chi0sumrule(npwe)= In input the sum rule calculated so far, In output the contribution of this transition is accounted for, and, eventually, symmetrized. using the symmetry operations of the little group of the external q.
PARENTS
cchi0,cchi0q0
CHILDREN
SOURCE
286 subroutine accumulate_chi0sumrule(ik_bz,symchi,npwe,factor,delta_ene,& 287 & Ltg_q,Gsph_epsG0,npwepG0,rhotwg,chi0sumrule) 288 289 use defs_basis 290 use m_errors 291 use m_profiling_abi 292 293 use m_gsphere, only : gsphere_t 294 use m_bz_mesh, only : littlegroup_t 295 296 !This section has been created automatically by the script Abilint (TD). 297 !Do not modify the following lines by hand. 298 #undef ABI_FUNC 299 #define ABI_FUNC 'accumulate_chi0sumrule' 300 !End of the abilint section 301 302 implicit none 303 304 !Arguments ------------------------------------ 305 !scalars 306 integer,intent(in) :: ik_bz,npwe,npwepG0,symchi 307 real(dp),intent(in) :: delta_ene,factor 308 type(gsphere_t),intent(in) :: Gsph_epsG0 309 type(littlegroup_t),target,intent(in) :: Ltg_q 310 !arrays 311 real(dp),intent(inout) :: chi0sumrule(npwe) 312 complex(gwpc),intent(in) :: rhotwg(npwepG0) 313 314 !Local variables------------------------------- 315 !scalars 316 integer :: isym,itim 317 character(len=500) :: msg 318 !arrays 319 integer,allocatable :: Sm1_gmG0(:) 320 integer, ABI_CONTIGUOUS pointer :: gmG0(:) 321 complex(gwpc),allocatable :: rhotwg_sym(:) 322 323 !************************************************************************ 324 325 ! Accumulating the sum rule on chi0. 326 ! Eq.(5.284) in G. D. Mahan Many-Particle Physics 3rd edition 327 328 SELECT CASE (symchi) 329 330 CASE (0) ! Do not use symmetries, sum is performed in the full BZ. 331 chi0sumrule(:)=chi0sumrule(:) + factor*delta_ene*ABS(rhotwg(1:npwe))**2 332 333 CASE (1) ! Symmetrize the contribution in the full BZ. 334 ABI_ALLOCATE(rhotwg_sym,(npwe)) 335 ABI_ALLOCATE(Sm1_gmG0,(npwe)) 336 337 do itim=1,Ltg_q%timrev 338 do isym=1,Ltg_q%nsym_sg 339 if (Ltg_q%wtksym(itim,isym,ik_bz)==1) then 340 ! === This operation belongs to the little group and has to be used to reconstruct the BZ === 341 ! * In the following 2 lines mind the slicing (1:npwe) 342 ! 343 gmG0 => Ltg_q%igmG0(1:npwe,itim,isym) 344 Sm1_gmG0(1:npwe)=Gsph_epsG0%rottbm1(gmG0(1:npwe),itim,isym) 345 rhotwg_sym(1:npwe)=rhotwg(Sm1_gmG0) 346 347 chi0sumrule(:)=chi0sumrule(:) + factor*delta_ene*ABS(rhotwg_sym(1:npwe))**2 348 end if 349 end do !isym 350 end do !itim 351 352 ABI_DEALLOCATE(rhotwg_sym) 353 ABI_DEALLOCATE(Sm1_gmG0) 354 355 CASE DEFAULT 356 write(msg,'(a,i3)')'Wrong value for symchi= ',symchi 357 MSG_BUG(msg) 358 END SELECT 359 360 end subroutine accumulate_chi0sumrule
ABINIT/completechi0_deltapart [ Functions ]
NAME
completechi0_deltapart
FUNCTION
Apply the delta part of the completeness correction to chi0
COPYRIGHT
Copyright (C) 1999-2017 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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
ik_bz=Index of the k-point in the full BZ whose contribution has to be added and symmetrized. qzero=.TRUE. is long wave-length limit. symchi=1 if we are summing over IBZ_q and symmetrization has to be performed. npwe=Number of G vectors in chi0. npwvec=MAX number of G. nomega=Number of frequencies. nspinor=Number of spinorial components. nfftot=Total Number of points in the FFT ngfft(18)=Info on the FFT. igfft0(npwvec)=Index of each G in the FFT array. Gsph_FFT=<gsphere_t>=Info on the largest G-sphere contained in the FFT box used for wavefunctions. Ltg_q=<littlegroup_t>= Structure gathering information on the little group of the external q. green_enhigh_w=Approximated frequency dependent part of the Green function entering equation (TODO put reference) wfwfg=Fourier components of u_{kb1}.u_{kb2}
OUTPUT
See SIDES EFFECTS SIDES EFFECTS chi0(npwe,npwe,nomega)= In input chi0 calculated so far, In output the "delta part" of the completeness correction is added.
PARENTS
cchi0,cchi0q0
CHILDREN
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 subroutine completechi0_deltapart(ik_bz,qzero,symchi,npwe,npwvec,nomega,nspinor,& 54 & nfftot,ngfft,igfft0,Gsph_FFT,Ltg_q,green_enhigh_w,wfwfg,chi0) 55 56 use defs_basis 57 use m_errors 58 use m_profiling_abi 59 60 use m_gsphere, only : gsphere_t, gsph_gmg_idx, gsph_gmg_fftidx 61 use m_bz_mesh, only : littlegroup_t 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 'completechi0_deltapart' 67 use interfaces_14_hidewrite 68 !End of the abilint section 69 70 implicit none 71 72 !Arguments ------------------------------------ 73 !scalars 74 integer,intent(in) :: ik_bz,nfftot,nomega,npwe,npwvec,nspinor,symchi 75 logical,intent(in) :: qzero 76 type(gsphere_t),intent(in) :: Gsph_FFT 77 type(littlegroup_t),intent(in) :: Ltg_q 78 !arrays 79 integer,intent(in) :: igfft0(npwvec),ngfft(18) 80 complex(dpc),intent(in) :: green_enhigh_w(nomega) 81 complex(gwpc),intent(in) :: wfwfg(nfftot*nspinor**2) 82 complex(gwpc),intent(inout) :: chi0(npwe,npwe,nomega) 83 84 !Local variables ------------------------------ 85 !scalars 86 integer,save :: enough=0 87 integer :: iSm1_g1mg2,iSm1_g1mg2_fft,ig,gmg_sph,gmg_fft 88 integer :: igp,igstart,isym,itim,outofbox_wfn 89 complex(gwpc) :: phmGt 90 character(len=500) :: msg 91 92 !************************************************************************ 93 94 igstart=1; if (qzero) igstart=2 95 outofbox_wfn=0 96 97 SELECT CASE (symchi) 98 99 CASE (0) ! Do not use symmetries. 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 do igp=igstart,npwe 103 do ig=igstart,npwe 104 105 gmg_fft = gsph_gmg_fftidx(Gsph_FFT,ig,igp,ngfft) 106 if (gmg_fft==0) then 107 outofbox_wfn=outofbox_wfn+1; CYCLE 108 end if 109 110 chi0(ig,igp,:) = chi0(ig,igp,:) + wfwfg(gmg_fft)*green_enhigh_w(:) 111 end do 112 end do 113 114 CASE (1) 115 ! Symmetrize the integrand in the full BZ. 116 ! * <Sk b|e^{-i(G1-G2}.r}|b Sk> = e^{-i(G1-G2).\tau} <k b|e^{-i(S^{-1}(G1-G2).r)|b k> 117 ! * green_enhigh_w in invariant under symmetry 118 ! * We symmetrize using the operations of the little group of q since this routine 119 ! is called inside a sum over IBZ_q, it would be possible to symmetrize 120 ! this term by just summing over the IBZ and rotating the matrix elements. 121 ! * Time-reversal does not lead to a complex conjugated since bra and ket are the same. 122 ! 123 do igp=igstart,npwe 124 do ig=igstart,npwe 125 126 ! Get the index of G1-G2. 127 gmg_sph = gsph_gmg_idx(Gsph_FFT,ig,igp) 128 if (gmg_sph==0) then 129 outofbox_wfn=outofbox_wfn+1; CYCLE 130 end if 131 132 do itim=1,Ltg_q%timrev 133 do isym=1,Ltg_q%nsym_sg 134 if (Ltg_q%wtksym(itim,isym,ik_bz)==1) then 135 ! * This operation belongs to the little group and has to be used to reconstruct the BZ. 136 ! * Time-reversal in not used to rotate (G1-G2) see comment above. 137 phmGt = Gsph_FFT%phmGt (gmg_sph,isym) 138 iSm1_g1mg2 = Gsph_FFT%rottbm1(gmg_sph,1,isym) 139 iSm1_g1mg2_fft = igfft0(iSm1_g1mg2) 140 141 chi0(ig,igp,:) = chi0(ig,igp,:) + phmGt*wfwfg(iSm1_g1mg2_fft)*green_enhigh_w(:) 142 end if 143 end do !isym 144 end do !itim 145 146 end do !igp 147 end do !ig 148 149 CASE DEFAULT 150 MSG_BUG("Wrong value of symchi") 151 END SELECT 152 153 if (outofbox_wfn/=0) then 154 enough=enough+1 155 if (enough<=50) then 156 write(msg,'(a,i0)')' Number of G1-G2 pairs outside the G-sphere for Wfns = ',outofbox_wfn 157 MSG_WARNING(msg) 158 if (enough==50) then 159 call wrtout(std_out,' ========== Stop writing Warnings ==========','COLL') 160 end if 161 end if 162 end if 163 164 end subroutine completechi0_deltapart
ABINIT/output_chi0sumrule [ Functions ]
NAME
output_chi0sumrule
FUNCTION
Calculate and output the value of the sum rule for the non-interacting polarizability chi0
INPUTS
OUTPUT
(for writing routines, no output) otherwise, should be described
PARENTS
screening
CHILDREN
SOURCE
190 subroutine output_chi0sumrule(qeq0,iq,npwe,omegaplasma,chi0sumrule,epsm1_w0,vc_sqrt) 191 192 use defs_basis 193 use m_profiling_abi 194 195 !This section has been created automatically by the script Abilint (TD). 196 !Do not modify the following lines by hand. 197 #undef ABI_FUNC 198 #define ABI_FUNC 'output_chi0sumrule' 199 use interfaces_14_hidewrite 200 !End of the abilint section 201 202 implicit none 203 204 !Arguments ------------------------------------ 205 !scalars 206 integer,intent(in) :: iq,npwe 207 real(dp),intent(in) :: omegaplasma 208 logical,intent(in) :: qeq0 209 !arrays 210 real(dp),intent(inout) :: chi0sumrule(npwe) 211 complex(gwpc),intent(in) :: epsm1_w0(npwe,npwe),vc_sqrt(npwe) 212 213 !Local variables ------------------------------ 214 !scalars 215 integer :: ig,igstart 216 real(dp) :: average,norm 217 character(len=500) :: msg 218 219 !************************************************************************ 220 221 igstart=1; if (qeq0) igstart=2 222 ! 223 ! The sumrule reads: 224 ! $ \int d\omega \omega v * Im[ \chi_0(\omega) ] = \pi/2 * w_p^2 $. 225 chi0sumrule(igstart:npwe) = chi0sumrule(igstart:npwe) * vc_sqrt(igstart:npwe)**2 226 ! 227 ! Calculate a weighted average of the fulfilment of the sumrule on epsilon 228 ! The weight is given according to the significance of each q+G in the 229 ! subsequent GW calculation: It is proportional to v * (epsm1 -1 ) 230 average = zero; norm = zero 231 do ig=igstart,npwe 232 average = average + chi0sumrule(ig) * real( vc_sqrt(ig)**2 * (epsm1_w0(ig,ig) - 1.0_dp ) ) 233 norm = norm + real( vc_sqrt(ig)**2 * (epsm1_w0(ig,ig) - 1.0_dp ) ) 234 !average = average + chi0sumrule(ig) * real( (epsm1_w0(ig,ig) - 1.0_dp ) ) 235 !norm = norm + real( (epsm1_w0(ig,ig) - 1.0_dp ) ) 236 !write(203,'(i4,8(2x,e12.6))') ig,1.0_dp/vc_sqrt(ig),chi0sumrule(ig)/ (0.5d0*omegaplasma**2*pi) 237 end do 238 239 if (abs(norm)>tol8) then 240 write(msg,'(1x,a,i4,a,f10.2,2x,a)')& 241 ' Average fulfillment of the sum rule on Im[epsilon] for q-point ',& 242 iq,' :',average/norm/(0.5_dp*omegaplasma**2*pi)*100.0_dp,'[%]' 243 call wrtout(std_out,msg,'COLL'); call wrtout(ab_out, msg,'COLL') 244 end if 245 246 end subroutine output_chi0sumrule