TABLE OF CONTENTS
- ABINIT/m_gwdefs
- m_gwdefs/em1params_free
- m_gwdefs/em1params_t
- m_gwdefs/g0g0w
- m_gwdefs/sigijtab_free
- m_gwdefs/sigma_is_herm
- m_gwdefs/sigma_needs_ppm
- m_gwdefs/sigma_needs_w
- m_gwdefs/sigma_type_from_key
- m_gwdefs/sigparams_free
- m_gwdefs/sigparams_t
ABINIT/m_gwdefs [ Modules ]
NAME
m_gwdefs
FUNCTION
This module contains definitions for a number of named constants used in the GW part of abinit
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (MG, FB, GMR, VO, LR, RWG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_gwdefs 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 28 implicit none 29 30 private 31 32 ! Unit number for formatted files produced by GW calculations. 33 ! These files are not supposed to be read by abinit therefore 34 ! their names and unit numbers are not defined in the Dtfil% structure. 35 integer,public,parameter :: unt_gw = 21 ! GW corrections 36 integer,public,parameter :: unt_sig = 22 ! Self-energy as a function of frequency 37 integer,public,parameter :: unt_sgr = 23 ! Derivative wrt omega of the Self-energy 38 integer,public,parameter :: unt_sgm = 20 ! Sigma on the Matsubara axis 39 integer,public,parameter :: unt_sigc = 24 ! Sigma_c as a function of (epsilon_i) MRM 40 integer,public,parameter :: unt_gwdiag = 40 ! GW diagonal 41 42 real(dp),public,parameter :: GW_TOLQ =0.0001_dp 43 ! Tolerance below which two BZ points are considered equal within a RL vector: 44 ! for each red. direct. the abs value of the difference btw the two coord must be smaller that tolq. 45 46 real(dp),public,parameter :: GW_TOLQ0=0.001_dp 47 ! Tolerance below which a q-point is treated as zero (long wavelength limit) 48 49 real(dp),public,parameter :: GW_TOL_DOCC=0.01_dp 50 ! Tolerance on the difference between two occupation numbers. 51 ! below this value, the contribution of the transition is neglected in the evaluation of chi0 52 53 real(dp),public,parameter :: GW_TOL_W0=0.001_dp 54 ! Tolerance on the real part of the frequency appearing in the denominator of the 55 ! non-interacting Green function G0. Above this value, a small purely imaginary 56 ! complex shift is added to the denominator during the evaluation of chi0. 57 58 real(gwp),public,parameter :: one_gw = 1._gwp 59 real(gwp),public,parameter :: zero_gw = 0._gwp 60 61 complex(gwpc),public,parameter :: czero_gw = (0._gwp,0._gwp) 62 complex(gwpc),public,parameter :: cone_gw = (1._gwp,0._gwp) 63 complex(gwpc),public,parameter :: j_gw = (0._gwp,1._gwp) 64 65 !arrays 66 real(dp),public,parameter :: GW_Q0_DEFAULT(3) = [0.00001_dp, 0.00002_dp, 0.00003_dp] 67 68 ! Weights and nodes for Gauss-Kronrod integration rules 69 ! Gauss 7 Kronrod 15 70 real(dp),public,parameter :: Kron15N(15) = (/-0.991455371120812639206854697526328516642040_dp, & 71 & -0.949107912342758524526189684047851262400770_dp,-0.864864423359769072789712788640926201210972_dp, & 72 & -0.741531185599394439863864773280788407074150_dp,-0.586087235467691130294144838258729598436780_dp, & 73 & -0.405845151377397166906606412076961463347382_dp,-0.207784955007898467600689403773244913479780_dp, & 74 & 0.000000000000000000000000000000000000000000_dp, 0.207784955007898467600689403773244913479780_dp, & 75 & 0.405845151377397166906606412076961463347382_dp, 0.586087235467691130294144838258729598436780_dp, & 76 & 0.741531185599394439863864773280788407074150_dp, 0.864864423359769072789712788640926201210972_dp, & 77 & 0.949107912342758524526189684047851262400770_dp, 0.991455371120812639206854697526328516642040_dp/) 78 79 real(dp),public,parameter :: Kron15W(15) = (/0.022935322010529224963732008058969591993561_dp, & 80 & 0.063092092629978553290700663189204286665070_dp,0.104790010322250183839876322541518017443757_dp, & 81 & 0.140653259715525918745189590510237920399890_dp,0.169004726639267902826583426598550284106245_dp, & 82 & 0.190350578064785409913256402421013682826078_dp,0.204432940075298892414161999234649084716518_dp, & 83 & 0.209482141084727828012999174891714263697760_dp,0.204432940075298892414161999234649084716518_dp, & 84 & 0.190350578064785409913256402421013682826078_dp,0.169004726639267902826583426598550284106245_dp, & 85 & 0.140653259715525918745189590510237920399890_dp,0.104790010322250183839876322541518017443757_dp, & 86 & 0.063092092629978553290700663189204286665070_dp,0.022935322010529224963732008058969591993561_dp/) 87 88 real(dp),public,parameter :: Gau7W(7) = (/0.12948496616886969327061143267908201832859_dp, & 89 & 0.279705391489276667901467771423779582486925_dp,0.38183005050511894495036977548897513387837_dp, & 90 & 0.417959183673469387755102040816326530612245_dp,0.38183005050511894495036977548897513387837_dp, & 91 & 0.279705391489276667901467771423779582486925_dp,0.12948496616886969327061143267908201832859_dp/) 92 93 ! Gauss 11 Kronrod 23 94 real(dp),public,parameter :: Kron23N(23) = (/-0.996369613889542634360164573335160773367030_dp, & 95 & -0.978228658146056992803938001122857390771420_dp,-0.941677108578067946455396730352134962188750_dp, & 96 & -0.887062599768095299075157769303927266631680_dp,-0.816057456656220942392261355192625879277860_dp, & 97 & -0.730152005574049324093416252031153458049643_dp,-0.630599520161965092168263312405390318822425_dp, & 98 & -0.519096129206811815925725669458609554480227_dp,-0.397944140952377573675073943298232277259112_dp, & 99 & -0.269543155952344972331531985400861524679620_dp,-0.136113000799361815798364355994952934344660_dp, & 100 & 0.000000000000000000000000000000000000000000_dp, 0.136113000799361815798364355994952934344660_dp, & 101 & 0.269543155952344972331531985400861524679620_dp, 0.397944140952377573675073943298232277259112_dp, & 102 & 0.519096129206811815925725669458609554480227_dp, 0.630599520161965092168263312405390318822425_dp, & 103 & 0.730152005574049324093416252031153458049643_dp, 0.816057456656220942392261355192625879277860_dp, & 104 & 0.887062599768095299075157769303927266631680_dp, 0.941677108578067946455396730352134962188750_dp, & 105 & 0.978228658146056992803938001122857390771420_dp, 0.996369613889542634360164573335160773367030_dp/) 106 107 real(dp),public,parameter :: Kron23W(23) = (/0.00976544104596075802247917260964352216369_dp, & 108 & 0.027156554682104262051721401617851679412810_dp,0.04582937856442641598526161637154832107050_dp, & 109 & 0.063097424750374906584540530495371467781318_dp,0.07866457193222732928421712412387330212537_dp, & 110 & 0.092953098596900827769293667912429161939839_dp,0.10587207448138939648189189823198112055496_dp, & 111 & 0.116739502461047270810811060893282832324909_dp,0.125158799100319505060067189609770147044437_dp, & 112 & 0.131280684229805644255977537339957161670610_dp,0.135193572799884533184261853141533217156135_dp, & 113 & 0.136577794711118301018953895305516133510830_dp,0.135193572799884533184261853141533217156135_dp, & 114 & 0.131280684229805644255977537339957161670610_dp,0.125158799100319505060067189609770147044437_dp, & 115 & 0.116739502461047270810811060893282832324909_dp,0.10587207448138939648189189823198112055496_dp, & 116 & 0.092953098596900827769293667912429161939839_dp,0.07866457193222732928421712412387330212537_dp, & 117 & 0.063097424750374906584540530495371467781318_dp,0.04582937856442641598526161637154832107050_dp, & 118 & 0.027156554682104262051721401617851679412810_dp,0.00976544104596075802247917260964352216369_dp/) 119 120 real(dp),public,parameter :: Gau11W(11) = (/0.0556685671161736664827537204425485787285_dp, & 121 & 0.125580369464904624634694299223940100197616_dp,0.186290210927734251426097641431655891691285_dp, & 122 & 0.233193764591990479918523704843175139431800_dp,0.262804544510246662180688869890509195372765_dp, & 123 & 0.272925086777900630714483528336342189156042_dp,0.262804544510246662180688869890509195372765_dp, & 124 & 0.233193764591990479918523704843175139431800_dp,0.186290210927734251426097641431655891691285_dp, & 125 & 0.125580369464904624634694299223940100197616_dp,0.055668567116173666482753720442548578728500_dp/) 126 127 ! Gauss 15 Kronrod 31 128 real(dp),public,parameter :: Kron31N(31) = (/-0.998002298693397060285172840152271209073410_dp, & 129 & -0.987992518020485428489565718586612581146970_dp,-0.967739075679139134257347978784337225283360_dp, & 130 & -0.937273392400705904307758947710209471244000_dp,-0.897264532344081900882509656454495882831780_dp, & 131 & -0.848206583410427216200648320774216851366260_dp,-0.790418501442465932967649294817947346862140_dp, & 132 & -0.724417731360170047416186054613938009630900_dp,-0.650996741297416970533735895313274692546948_dp, & 133 & -0.570972172608538847537226737253910641238390_dp,-0.485081863640239680693655740232350612866339_dp, & 134 & -0.394151347077563369897207370981045468362750_dp,-0.299180007153168812166780024266388962661603_dp, & 135 & -0.201194093997434522300628303394596207812836_dp,-0.101142066918717499027074231447392338787451_dp, & 136 & 0.000000000000000000000000000000000000000000_dp, 0.101142066918717499027074231447392338787451_dp, & 137 & 0.201194093997434522300628303394596207812836_dp, 0.299180007153168812166780024266388962661603_dp, & 138 & 0.394151347077563369897207370981045468362750_dp, 0.485081863640239680693655740232350612866339_dp, & 139 & 0.570972172608538847537226737253910641238390_dp, 0.650996741297416970533735895313274692546948_dp, & 140 & 0.724417731360170047416186054613938009630900_dp, 0.790418501442465932967649294817947346862140_dp, & 141 & 0.848206583410427216200648320774216851366260_dp, 0.897264532344081900882509656454495882831780_dp, & 142 & 0.937273392400705904307758947710209471244000_dp, 0.967739075679139134257347978784337225283360_dp, & 143 & 0.987992518020485428489565718586612581146970_dp, 0.998002298693397060285172840152271209073410_dp/) 144 145 real(dp),public,parameter :: Kron31W(31) = (/0.0053774798729233489877920514301276498183100_dp, & 146 & 0.0150079473293161225383747630758072680946390_dp, 0.0254608473267153201868740010196533593972700_dp, & 147 & 0.0353463607913758462220379484783600481226300_dp, 0.0445897513247648766082272993732796902232570_dp, & 148 & 0.0534815246909280872653431472394302967715500_dp, 0.0620095678006706402851392309608029321904000_dp, & 149 & 0.0698541213187282587095200770991474757860450_dp, 0.0768496807577203788944327774826590067221100_dp, & 150 & 0.0830805028231330210382892472861037896015540_dp, 0.0885644430562117706472754436937743032122700_dp, & 151 & 0.0931265981708253212254868727473457185619300_dp, 0.0966427269836236785051799076275893351366570_dp, & 152 & 0.0991735987217919593323931734846031310595673_dp, 0.1007698455238755950449466626175697219163500_dp, & 153 & 0.1013300070147915490173747927674925467709270_dp, 0.1007698455238755950449466626175697219163500_dp, & 154 & 0.0991735987217919593323931734846031310595673_dp, 0.0966427269836236785051799076275893351366570_dp, & 155 & 0.0931265981708253212254868727473457185619300_dp, 0.0885644430562117706472754436937743032122700_dp, & 156 & 0.0830805028231330210382892472861037896015540_dp, 0.0768496807577203788944327774826590067221100_dp, & 157 & 0.0698541213187282587095200770991474757860450_dp, 0.0620095678006706402851392309608029321904000_dp, & 158 & 0.0534815246909280872653431472394302967715500_dp, 0.0445897513247648766082272993732796902232570_dp, & 159 & 0.0353463607913758462220379484783600481226300_dp, 0.0254608473267153201868740010196533593972700_dp, & 160 & 0.0150079473293161225383747630758072680946390_dp, 0.0053774798729233489877920514301276498183100_dp/) 161 162 real(dp),public,parameter :: Gau15W(15) = (/0.030753241996117268354628393577204417721700_dp, & 163 0.070366047488108124709267416450667338466710_dp, 0.107159220467171935011869546685869303415544_dp, & 164 0.139570677926154314447804794511028322520850_dp, 0.166269205816993933553200860481208811130900_dp, & 165 0.186161000015562211026800561866422824506226_dp, 0.198431485327111576456118326443839324818693_dp, & 166 0.202578241925561272880620199967519314838662_dp, 0.198431485327111576456118326443839324818693_dp, & 167 0.186161000015562211026800561866422824506226_dp, 0.166269205816993933553200860481208811130900_dp, & 168 0.139570677926154314447804794511028322520850_dp, 0.107159220467171935011869546685869303415544_dp, & 169 0.070366047488108124709267416450667338466710_dp, 0.030753241996117268354628393577204417721700_dp/) 170 171 172 ! Flags for self-consistent GW calculations used in gw_driver and for parsing the input file. 173 integer,public,parameter :: GWSC_one_shot =1 174 integer,public,parameter :: GWSC_only_W =2 175 integer,public,parameter :: GWSC_only_G =3 176 integer,public,parameter :: GWSC_both_G_and_W =4 177 178 ! Flags defining the approximation used for the self-energy (used in csigme). 179 integer,public,parameter :: SIG_GW_PPM =0 ! standard GW with PPM 180 integer,public,parameter :: SIG_GW_AC =1 ! standard GW without PPM (analytical continuation) 181 integer,public,parameter :: SIG_GW_CD =2 ! standard GW without PPM (contour deformation) 182 integer,public,parameter :: SIG_HF =5 ! Hartree-Fock calculation 183 integer,public,parameter :: SIG_SEX =6 ! Screened Exchange calculation 184 integer,public,parameter :: SIG_COHSEX =7 ! COHSEX calculation 185 integer,public,parameter :: SIG_QPGW_PPM =8 ! model GW with PPM 186 integer,public,parameter :: SIG_QPGW_CD =9 ! model GW without PPM 187 188 public :: sigma_type_from_key 189 public :: sigma_is_herm 190 public :: sigma_needs_w 191 public :: sigma_needs_ppm 192 !public :: sigma_sc_on_wfs 193 !public :: sigma_sc_on_ene 194 public :: g0g0w 195 196 ! Private variables 197 integer,private,parameter :: STR_LEN=500
m_gwdefs/em1params_free [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
em1params_free
FUNCTION
Free dynamic memory allocated in the structure.
INPUTS
OUTPUT
SOURCE
410 subroutine em1params_free(Ep) 411 412 !Arguments ------------------------------------ 413 class(em1params_t),intent(inout) :: Ep 414 415 ! ************************************************************************* 416 417 !@em1params_t 418 419 !real 420 ABI_SFREE(Ep%qcalc) 421 ABI_SFREE(Ep%qibz) 422 ABI_SFREE(Ep%qlwl) 423 ABI_SFREE(Ep%omegasf) 424 !complex 425 ABI_SFREE(Ep%omega) 426 427 end subroutine em1params_free
m_gwdefs/em1params_t [ Types ]
[ Top ] [ m_gwdefs ] [ Types ]
NAME
em1params_t
FUNCTION
For the GW part of ABINIT, the em1params_t structured datatype gather different parameters used to calculate the inverse dielectric matrices.
SOURCE
212 type,public :: em1params_t 213 214 !scalars 215 integer :: awtr ! If 1 the Adler-Wiser expression for Chi_0 is evaluated 216 ! taking advantage of time-reversal symmetry 217 integer :: gwcalctyp ! Calculation type (see input variable) 218 integer :: gwcomp ! 1 if extrapolar technique is used. 0 otherwise. 219 integer :: inclvkb ! Integer flag related to the evaluation of the commutator for q-->0 220 integer :: spmeth ! Method used to approximate the delta function in the expression for Im Chi_0 221 integer :: nI ! Number of components (rows) in the chi0 matrix. 222 integer :: nJ ! Number of components (columns) in the chi0 matrix. 223 integer :: npwvec ! Max between npwe and npwwfn, used to pass the dimension of arrays e.g gvec 224 integer :: npwwfn ! Number of planewaves for wavefunctions 225 integer :: npwe ! Number of planewaves for $\tilde \epsilon$ 226 integer :: npwepG0 ! Number of planewaves in the enlarged sphere G-G0, to account for umklapp G0 vectors 227 integer :: nbnds ! Number of bands used to evaluate $\tilde \epsilon$ 228 integer :: nkibz ! Number of k-points in the IBZ 229 integer :: nsppol ! 1 for spin unpolarized, 2 for collinear spin polarized 230 integer :: nqcalc ! Number of q-points that are calculated (subset of qibz) 231 integer :: nqibz ! Number of q-points in the IBZ 232 integer :: nqlwl ! Number of directions to analyze the non analytical behavior for q-->0 233 integer :: nomega ! Number of frequencies where evaluate $\tilde \epsilon (\omega)$ 234 integer :: nomegaer,nomegaei ! Number of real and imaginary frequencies, respectively 235 integer :: nomegaec ! Number of frequencies on a grid in the complex plane nomegaec = nomegaei*(nomegaer-1) 236 integer :: nomegasf ! Number of frequencies used for the spectral function 237 integer :: symchi ! 0 ==> do not use symmetries to reduce the k-points summed over in chi0 238 ! 1 ==> take advantage of point group symmetries as well as time-reversal 239 240 real(dp) :: gwencomp ! Extrapolar energy used if gwcomp==1. 241 real(dp) :: omegaermin ! Minimum real frequency used in the contour deformation method 242 real(dp) :: omegaermax ! Maximum real frequency used in the contour deformation method 243 real(dp) :: mbpt_sciss ! Scissor energy used in chi0 244 real(dp) :: spsmear ! Smearing of the delta in case of spmeth==2 245 real(dp) :: zcut ! Small imaginary shift to avoid poles in chi0 246 247 logical :: analytic_continuation ! if true calculate chi0 only along the imaginary axis 248 logical :: contour_deformation ! if true calculated chi0 both along the real and the imaginary axis 249 logical :: plasmon_pole_model ! if true a plasmonpole model is used (only 1 or 2 frequencies are calculated) 250 251 !arrays 252 integer :: mG0(3) 253 ! For each reduced direction gives the max G0 component to account for umklapp processes 254 255 real(dp),allocatable :: qcalc(:,:) 256 ! (3,nqcalc) 257 ! q-points that are explicitely calculated (subset of qibz). 258 259 real(dp),allocatable :: qibz(:,:) 260 ! (3,nqibz) 261 ! q-points in the IBZ. 262 263 real(dp),allocatable :: qlwl(:,:) 264 ! (3,nqlwl) 265 ! q-points used for the long-wavelength limit. 266 267 real(dp),allocatable :: omegasf(:) 268 ! (nomegasf) 269 ! real frequencies used to calculate the imaginary part of chi0. 270 271 complex(dpc),allocatable :: omega(:) 272 ! (nomegasf) 273 ! real and imaginary frequencies in chi0,epsilon and epsilonm1. 274 275 contains 276 procedure :: free => em1params_free 277 end type em1params_t
m_gwdefs/g0g0w [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
g0g0w
FUNCTION
Calculates the frequency-dependent part of the RPA polarizability G0G0.
INPUTS
OUTPUT
SOURCE
673 function g0g0w(omega, numerator, delta_ene, zcut, TOL_W0, opt_poles) 674 675 !Arguments ------------------------------------ 676 !scalars 677 integer,intent(in):: opt_poles 678 real(dp),intent(in) :: TOL_W0,delta_ene,numerator,zcut 679 complex(dpc) :: g0g0w 680 complex(dpc),intent(in) :: omega 681 682 !Local variables ------------------------------ 683 !scalars 684 real(dp) :: sgn 685 character(len=500) :: msg 686 687 !************************************************************************ 688 689 if (delta_ene**2 > tol14) then 690 sgn = SIGN(1.0_dp,delta_ene) 691 692 if (opt_poles == 2) then 693 ! Resonant and anti-resonant contributions. 694 if (DABS(REAL(omega)) > TOL_W0) then 695 ! omega on the real axis 696 g0g0w = numerator / (omega + delta_ene - j_dpc*sgn*zcut) & 697 -numerator / (omega - delta_ene + j_dpc*sgn*zcut) 698 else 699 ! omega on the imag axis (g0g0w is purely real) 700 g0g0w = numerator / (omega + delta_ene) & 701 -numerator / (omega - delta_ene) 702 end if 703 704 else if (opt_poles == 1) then 705 ! Only resonant contribution is included. 706 if (DABS(REAL(omega)) > TOL_W0) then 707 g0g0w = numerator / (omega + delta_ene - j_dpc*sgn*zcut) 708 else 709 g0g0w = numerator / (omega + delta_ene) 710 end if 711 712 else 713 write(msg,'(a,i0)')" Wrong value for opt_poles: ",opt_poles 714 ABI_ERROR(msg) 715 end if ! opt_poles 716 717 else 718 ! delta_ene**2 < tol14 719 g0g0w = czero 720 end if 721 722 end function g0g0w
m_gwdefs/sigijtab_free [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigijtab_free
FUNCTION
deallocate all memory in a sigijtab_t datatype.
INPUTS
OUTPUT
SOURCE
445 subroutine sigijtab_free(Sigijtab) 446 447 !Arguments ------------------------------------ 448 !scalars 449 type(sigijtab_t),intent(inout) :: Sigijtab(:,:) 450 451 !Local variables 452 integer :: ii,jj,kk,ilow,iup 453 ! ************************************************************************* 454 455 !@sigijtab_t 456 do jj=1,SIZE(Sigijtab,DIM=2) 457 do ii=1,SIZE(Sigijtab,DIM=1) 458 459 ilow=LBOUND(Sigijtab(ii,jj)%col,DIM=1) 460 iup =UBOUND(Sigijtab(ii,jj)%col,DIM=1) 461 do kk=ilow,iup 462 ABI_FREE(Sigijtab(ii,jj)%col(kk)%bidx) 463 end do 464 ABI_FREE(Sigijtab(ii,jj)%col) 465 466 end do 467 end do 468 469 end subroutine sigijtab_free
m_gwdefs/sigma_is_herm [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigma_is_herm
FUNCTION
Return .TRUE. if the approximated self-energy is hermitian.
INPUTS
Sigp<sigparams_t>=datatype gathering data and info on the self-energy run.
SOURCE
579 pure logical function sigma_is_herm(Sigp) 580 581 !Arguments ------------------------------ 582 class(sigparams_t),intent(in) :: Sigp 583 584 !Local variables ------------------------------ 585 integer :: mod10 586 587 !************************************************************************ 588 589 mod10=MOD(Sigp%gwcalctyp,10) 590 sigma_is_herm = ANY(mod10 == [SIG_HF, SIG_SEX, SIG_COHSEX]) 591 592 end function sigma_is_herm
m_gwdefs/sigma_needs_ppm [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigma_needs_ppm
FUNCTION
Return .TRUE. if the self-energy run requires a plasmon-pole model.
INPUTS
Sigp<sigparams_t>=datatype gathering data and info on the self-energy run.
SOURCE
640 pure logical function sigma_needs_ppm(Sigp) 641 642 !Arguments ------------------------------ 643 class(sigparams_t),intent(in) :: Sigp 644 645 !Local variables ------------------------------ 646 integer :: mod10 647 648 !************************************************************************ 649 650 mod10=MOD(Sigp%gwcalctyp,10) 651 sigma_needs_ppm = (ANY(mod10 == (/SIG_GW_PPM, SIG_QPGW_PPM/)) .or. & 652 Sigp%gwcomp==1 & 653 ) 654 655 end function sigma_needs_ppm
m_gwdefs/sigma_needs_w [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigma_needs_w
FUNCTION
Return .TRUE. if self-energy requires the screened interaction W. For example HF does not need the SCR file.
INPUTS
Sigp<sigparams_t>=datatype gathering data and info on the self-energy run.
SOURCE
610 pure logical function sigma_needs_w(Sigp) 611 612 !Arguments ------------------------------ 613 class(sigparams_t),intent(in) :: Sigp 614 615 !Local variables ------------------------------ 616 integer :: mod10 617 618 !************************************************************************ 619 620 mod10=MOD(Sigp%gwcalctyp,10) 621 sigma_needs_w = (mod10/=SIG_HF) 622 623 end function sigma_needs_w
m_gwdefs/sigma_type_from_key [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigma_type_from_key
FUNCTION
Return a string definining the particular approximation used for the self-energy. Stops if the key is not in the list of allowed possibilities.
INPUTS
key=Integer
OUTPUT
SOURCE
536 function sigma_type_from_key(key) result(sigma_type) 537 538 integer,intent(in) :: key 539 character(len=STR_LEN) :: sigma_type 540 541 !Local variables ------------------------------ 542 !scalars 543 character(len=500) :: msg 544 545 !************************************************************************ 546 547 sigma_type = "None" 548 if (key==SIG_GW_PPM ) sigma_type = ' standard GW with PPM' 549 if (key==SIG_GW_AC ) sigma_type = ' standard GW without PPM (analytical continuation)' 550 if (key==SIG_GW_CD ) sigma_type = ' standard GW without PPM (contour deformation)' 551 if (key==SIG_HF ) sigma_type = ' Hartree-Fock calculation' 552 if (key==SIG_SEX ) sigma_type = ' Screened Exchange calculation' 553 if (key==SIG_COHSEX ) sigma_type = ' COHSEX calculation' 554 if (key==SIG_QPGW_PPM) sigma_type = ' model GW with PPM' 555 if (key==SIG_QPGW_CD ) sigma_type = ' model GW without PPM' 556 557 if (sigma_type == "None") then 558 write(msg,'(a,i0)')" Unknown value for key= ",key 559 ABI_ERROR(msg) 560 end if 561 562 end function sigma_type_from_key
m_gwdefs/sigparams_free [ Functions ]
[ Top ] [ m_gwdefs ] [ Functions ]
NAME
sigparams_free
FUNCTION
Free dynamic memory allocated in the structure.
INPUTS
OUTPUT
SOURCE
487 subroutine sigparams_free(Sigp) 488 489 !Arguments ------------------------------------ 490 class(sigparams_t),intent(inout) :: Sigp 491 492 ! ************************************************************************* 493 494 !@sigparams_t 495 496 !integer 497 ABI_SFREE(Sigp%kptgw2bz) 498 ABI_SFREE(Sigp%minbnd) 499 ABI_SFREE(Sigp%maxbnd) 500 !real 501 ABI_FREE(Sigp%kptgw) 502 !complex 503 ABI_SFREE(Sigp%omegasi) 504 ABI_SFREE(Sigp%omega_r) 505 506 if (allocated(Sigp%Sigcij_tab)) then 507 call sigijtab_free(Sigp%Sigcij_tab) 508 ABI_FREE(Sigp%Sigcij_tab) 509 end if 510 511 if (allocated(Sigp%Sigxij_tab)) then 512 call sigijtab_free(Sigp%Sigxij_tab) 513 ABI_FREE(Sigp%Sigxij_tab) 514 end if 515 516 end subroutine sigparams_free
m_gwdefs/sigparams_t [ Types ]
[ Top ] [ m_gwdefs ] [ Types ]
NAME
sigparams_t
FUNCTION
For the GW part of ABINIT, the sigparams_t structured datatype gather different parameters that characterize the calculation of the matrix elements of the self-energy operator.
SOURCE
304 type,public :: sigparams_t 305 306 integer :: gwcalctyp ! Calculation type 307 308 integer :: gwgamma ! If 1 include vertex correction (GWGamma) 309 integer :: gwcomp ! 1 if the extrapolar technique is used. 310 311 integer :: minbdgw,maxbdgw ! Minimum and maximum band index (considering the spin) defining 312 ! The set of bands where GW corrections are evaluated 313 314 integer :: mG0(3) ! For each reduced direction gives the max G0 component 315 ! to account for umklapp processes 316 317 integer :: npwvec ! Max betwenn npwe and npwwfn, used to pass the dimension of arrays e.g gvec 318 integer :: npwwfn ! No. of planewaves for wavefunctions 319 integer :: npwx ! No. of planewaves for $\Sigma_x$ 320 integer :: npwc ! No. of planewaves for $\Sigma_c$ and W 321 integer :: nbnds ! No. of bands summed over. 322 integer :: nomegasr ! No. of frequencies on the real axis to evaluate the spectral function 323 integer :: nomegasrd ! No. of frequencies on the real axis to evaluate $\Sigma(E)$ 324 integer :: nomegasi ! No. of frequencies along the imaginary axis for Sigma in case of AC 325 integer :: nsig_ab ! No. of components in the self-energy operator (1 if nspinor==1, 4 if nspinor==2) 326 integer :: nspinor ! No. of spinorial components. 327 integer :: nsppol ! 1 for unpolarized, 2 for spin-polarized calculation 328 integer :: nkptgw ! No. of k-points where GW corrections have been calculated 329 integer :: ppmodel ! Integer defining the plasmon pole model used, 0 for None. 330 integer :: symsigma ! 0 ==> do not use symmetries to reduce the k-points summed over in sigma 331 ! 1 ==> take advantage of space group symmetries as well as time-reversal 332 integer :: use_sigxcore ! 1 if core contribution to sigma is estimated by using Hartree-Fock 333 334 real(dp) :: deltae ! Energy step used to evaluate numerically the derivative of the self energy 335 ! $\frac{\partial \Re \Sigma(E)}{\partial E_o}$ 336 real(dp) :: ecutwfn ! cutoff energy for the wavefunctions. 337 real(dp) :: ecutsigx ! cutoff energy for the the exchange parth of Sigma. 338 real(dp) :: ecuteps ! cutoff energy for W 339 340 real(dp) :: gwencomp ! Extrapolar energy used if gwcomp==1. 341 342 real(dp) :: mbpt_sciss ! Scissor energy used in G0 343 344 real(dp) :: minomega_r ! Minimum real frequency for the evaluation of the spectral function 345 real(dp) :: maxomega_r ! Maximum real frequency for the evaluation of the spectral function 346 real(dp) :: maxomega4sd ! Maximum displacement around the KS energy where evaluate the diagonal 347 ! Elements of $ \Sigma(E)$ 348 real(dp) :: omegasimax ! Max omega for Sigma along the imag axis in case of analytic continuation 349 real(dp) :: omegasimin ! min omega for Sigma along the imag axis in case of analytic continuation 350 351 real(dp) :: sigma_mixing ! Global factor that multiplies Sigma to give the final matrix element. 352 ! Usually one, except for the hybrid functionals. 353 354 real(dp) :: zcut ! Value of $\delta$ used to avoid the divergences (see related input variable) 355 356 integer,allocatable :: kptgw2bz(:) 357 ! (nkptgw) 358 ! For each k-point where GW corrections are calculated, the corresponding index in the BZ. 359 360 integer,allocatable :: minbnd(:,:), maxbnd(:,:) 361 ! (nkptgw, nsppol) 362 ! For each k-point at which GW corrections are calculated, the min and Max band index considered 363 ! (see also input variable dtset%bdgw). 364 365 real(dp),allocatable :: kptgw(:,:) 366 ! (3, nkptgw) 367 ! k-points for the GW corrections in reduced coordinates. 368 369 !TODO should be removed, everything should be in Sr% 370 371 complex(dpc),allocatable :: omegasi(:) 372 ! (nomegasi) 373 ! Frequencies along the imaginary axis used for the analytical continuation. 374 375 complex(dpc),allocatable :: omega_r(:) 376 ! (nomegasr) 377 ! Frequencies used to evaluate the spectral function. 378 379 type(sigijtab_t),allocatable :: Sigcij_tab(:,:) 380 ! (nkptgw, nsppol)%col(kb)%bidx(ii) gives the index of the left wavefunction. 381 ! in the <i,kgw,s|\Sigma_c|j,kgw,s> matrix elements that has to be calculated in cisgme. 382 ! in the case of self-consistent GW on wavefunctions. 383 384 type(sigijtab_t),allocatable :: Sigxij_tab(:,:) 385 ! Save as Sigcij_tab but for the Hermitian \Sigma_x where only the upper triangle is needed. 386 387 contains 388 procedure :: free => sigparams_free 389 end type sigparams_t