TABLE OF CONTENTS


ABINIT/m_gwdefs [ Modules ]

[ Top ] [ 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-2018 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 .

PARENTS

SOURCE

 19 #if defined HAVE_CONFIG_H
 20 #include "config.h"
 21 #endif
 22 
 23 #include "abi_common.h"
 24 
 25 MODULE m_gwdefs
 26 
 27  use defs_basis
 28  use m_abicore
 29  use m_errors
 30 
 31  implicit none
 32 
 33  private
 34 
 35 ! Unit number for formatted files produced by GW calculations.
 36 ! These files are not supposed to be read by abinit therefore
 37 ! their names and unit numbers are not defined in the Dtfil% structure.
 38  integer,public,parameter :: unt_gw  = 21  ! GW corrections
 39  integer,public,parameter :: unt_sig = 22  ! Self-energy as a function of frequency
 40  integer,public,parameter :: unt_sgr = 23  ! Derivative wrt omega of the Self-energy
 41  integer,public,parameter :: unt_sgm = 20  ! Sigma on the Matsubara axis
 42  integer,public,parameter :: unt_gwdiag  = 40 ! GW diagonal
 43 
 44  real(dp),public,parameter :: GW_TOLQ =0.0001_dp
 45  ! Tolerance below which two BZ points are considered equal within a RL vector:
 46  ! for each red. direct. the abs value of the difference btw the two coord must be smaller that tolq.
 47 
 48  real(dp),public,parameter :: GW_TOLQ0=0.001_dp
 49  ! Tolerance below which a q-point is treated as zero (long wavelength limit)
 50 
 51  real(dp),public,parameter :: GW_TOL_DOCC=0.01_dp
 52  ! Tolerance on the difference between two occupation numbers.
 53  ! below this value, the contribution of the transition is neglected in the evaluation of chi0
 54 
 55  real(dp),public,parameter :: GW_TOL_W0=0.001_dp
 56  ! Tolerance on the real part of the frequency appearing in the denominator of the
 57  ! non-interacting Green function G0. Above this value, a small purely imaginary
 58  ! complex shift is added to the denominator during the evaluation of chi0.
 59 
 60  real(gwp),public,parameter :: one_gw =1._gwp
 61  real(gwp),public,parameter :: zero_gw=0._gwp
 62 
 63  complex(gwpc),public,parameter :: czero_gw=(0._gwp,0._gwp)
 64  complex(gwpc),public,parameter :: cone_gw =(1._gwp,0._gwp)
 65  complex(gwpc),public,parameter :: j_gw    =(0._gwp,1._gwp)
 66 
 67 !arrays
 68  real(dp),public,parameter :: GW_Q0_DEFAULT(3) = (/0.00001_dp, 0.00002_dp, 0.00003_dp/)
 69 
 70 ! Weights and nodes for Gauss-Kronrod integration rules
 71 ! Gauss 7 Kronrod 15
 72  real(dp),public,parameter :: Kron15N(15) = (/-0.991455371120812639206854697526328516642040_dp, &
 73 & -0.949107912342758524526189684047851262400770_dp,-0.864864423359769072789712788640926201210972_dp, &
 74 & -0.741531185599394439863864773280788407074150_dp,-0.586087235467691130294144838258729598436780_dp, &
 75 & -0.405845151377397166906606412076961463347382_dp,-0.207784955007898467600689403773244913479780_dp, &
 76 &  0.000000000000000000000000000000000000000000_dp, 0.207784955007898467600689403773244913479780_dp, &
 77 &  0.405845151377397166906606412076961463347382_dp, 0.586087235467691130294144838258729598436780_dp, &
 78 &  0.741531185599394439863864773280788407074150_dp, 0.864864423359769072789712788640926201210972_dp, &
 79 &  0.949107912342758524526189684047851262400770_dp, 0.991455371120812639206854697526328516642040_dp/)
 80 
 81  real(dp),public,parameter :: Kron15W(15) = (/0.022935322010529224963732008058969591993561_dp, &
 82 & 0.063092092629978553290700663189204286665070_dp,0.104790010322250183839876322541518017443757_dp, &
 83 & 0.140653259715525918745189590510237920399890_dp,0.169004726639267902826583426598550284106245_dp, &
 84 & 0.190350578064785409913256402421013682826078_dp,0.204432940075298892414161999234649084716518_dp, &
 85 & 0.209482141084727828012999174891714263697760_dp,0.204432940075298892414161999234649084716518_dp, &
 86 & 0.190350578064785409913256402421013682826078_dp,0.169004726639267902826583426598550284106245_dp, &
 87 & 0.140653259715525918745189590510237920399890_dp,0.104790010322250183839876322541518017443757_dp, &
 88 & 0.063092092629978553290700663189204286665070_dp,0.022935322010529224963732008058969591993561_dp/)
 89 
 90  real(dp),public,parameter :: Gau7W(7) = (/0.12948496616886969327061143267908201832859_dp, &
 91 & 0.279705391489276667901467771423779582486925_dp,0.38183005050511894495036977548897513387837_dp, &
 92 & 0.417959183673469387755102040816326530612245_dp,0.38183005050511894495036977548897513387837_dp, &
 93 & 0.279705391489276667901467771423779582486925_dp,0.12948496616886969327061143267908201832859_dp/)
 94 
 95 ! Gauss 11 Kronrod 23
 96  real(dp),public,parameter :: Kron23N(23) = (/-0.996369613889542634360164573335160773367030_dp, &
 97 & -0.978228658146056992803938001122857390771420_dp,-0.941677108578067946455396730352134962188750_dp, &
 98 & -0.887062599768095299075157769303927266631680_dp,-0.816057456656220942392261355192625879277860_dp, &
 99 & -0.730152005574049324093416252031153458049643_dp,-0.630599520161965092168263312405390318822425_dp, &
100 & -0.519096129206811815925725669458609554480227_dp,-0.397944140952377573675073943298232277259112_dp, &
101 & -0.269543155952344972331531985400861524679620_dp,-0.136113000799361815798364355994952934344660_dp, &
102 &  0.000000000000000000000000000000000000000000_dp, 0.136113000799361815798364355994952934344660_dp, &
103 &  0.269543155952344972331531985400861524679620_dp, 0.397944140952377573675073943298232277259112_dp, &
104 &  0.519096129206811815925725669458609554480227_dp, 0.630599520161965092168263312405390318822425_dp, &
105 &  0.730152005574049324093416252031153458049643_dp, 0.816057456656220942392261355192625879277860_dp, &
106 &  0.887062599768095299075157769303927266631680_dp, 0.941677108578067946455396730352134962188750_dp, &
107 &  0.978228658146056992803938001122857390771420_dp, 0.996369613889542634360164573335160773367030_dp/)
108 
109  real(dp),public,parameter :: Kron23W(23) = (/0.00976544104596075802247917260964352216369_dp, &
110 & 0.027156554682104262051721401617851679412810_dp,0.04582937856442641598526161637154832107050_dp, &
111 & 0.063097424750374906584540530495371467781318_dp,0.07866457193222732928421712412387330212537_dp, &
112 & 0.092953098596900827769293667912429161939839_dp,0.10587207448138939648189189823198112055496_dp, &
113 & 0.116739502461047270810811060893282832324909_dp,0.125158799100319505060067189609770147044437_dp, &
114 & 0.131280684229805644255977537339957161670610_dp,0.135193572799884533184261853141533217156135_dp, &
115 & 0.136577794711118301018953895305516133510830_dp,0.135193572799884533184261853141533217156135_dp, &
116 & 0.131280684229805644255977537339957161670610_dp,0.125158799100319505060067189609770147044437_dp, &
117 & 0.116739502461047270810811060893282832324909_dp,0.10587207448138939648189189823198112055496_dp, &
118 & 0.092953098596900827769293667912429161939839_dp,0.07866457193222732928421712412387330212537_dp, &
119 & 0.063097424750374906584540530495371467781318_dp,0.04582937856442641598526161637154832107050_dp, &
120 & 0.027156554682104262051721401617851679412810_dp,0.00976544104596075802247917260964352216369_dp/)
121 
122  real(dp),public,parameter :: Gau11W(11) = (/0.0556685671161736664827537204425485787285_dp, &
123 & 0.125580369464904624634694299223940100197616_dp,0.186290210927734251426097641431655891691285_dp, &
124 & 0.233193764591990479918523704843175139431800_dp,0.262804544510246662180688869890509195372765_dp, &
125 & 0.272925086777900630714483528336342189156042_dp,0.262804544510246662180688869890509195372765_dp, &
126 & 0.233193764591990479918523704843175139431800_dp,0.186290210927734251426097641431655891691285_dp, &
127 & 0.125580369464904624634694299223940100197616_dp,0.055668567116173666482753720442548578728500_dp/)
128 
129 ! Gauss 15 Kronrod 31
130  real(dp),public,parameter :: Kron31N(31) = (/-0.998002298693397060285172840152271209073410_dp, &
131 & -0.987992518020485428489565718586612581146970_dp,-0.967739075679139134257347978784337225283360_dp, &
132 & -0.937273392400705904307758947710209471244000_dp,-0.897264532344081900882509656454495882831780_dp, &
133 & -0.848206583410427216200648320774216851366260_dp,-0.790418501442465932967649294817947346862140_dp, &
134 & -0.724417731360170047416186054613938009630900_dp,-0.650996741297416970533735895313274692546948_dp, &
135 & -0.570972172608538847537226737253910641238390_dp,-0.485081863640239680693655740232350612866339_dp, &
136 & -0.394151347077563369897207370981045468362750_dp,-0.299180007153168812166780024266388962661603_dp, &
137 & -0.201194093997434522300628303394596207812836_dp,-0.101142066918717499027074231447392338787451_dp, &
138 &  0.000000000000000000000000000000000000000000_dp, 0.101142066918717499027074231447392338787451_dp, &
139 &  0.201194093997434522300628303394596207812836_dp, 0.299180007153168812166780024266388962661603_dp, &
140 &  0.394151347077563369897207370981045468362750_dp, 0.485081863640239680693655740232350612866339_dp, &
141 &  0.570972172608538847537226737253910641238390_dp, 0.650996741297416970533735895313274692546948_dp, &
142 &  0.724417731360170047416186054613938009630900_dp, 0.790418501442465932967649294817947346862140_dp, &
143 &  0.848206583410427216200648320774216851366260_dp, 0.897264532344081900882509656454495882831780_dp, &
144 &  0.937273392400705904307758947710209471244000_dp, 0.967739075679139134257347978784337225283360_dp, &
145 &  0.987992518020485428489565718586612581146970_dp, 0.998002298693397060285172840152271209073410_dp/)
146 
147  real(dp),public,parameter :: Kron31W(31) = (/0.0053774798729233489877920514301276498183100_dp, &
148 & 0.0150079473293161225383747630758072680946390_dp, 0.0254608473267153201868740010196533593972700_dp, &
149 & 0.0353463607913758462220379484783600481226300_dp, 0.0445897513247648766082272993732796902232570_dp, &
150 & 0.0534815246909280872653431472394302967715500_dp, 0.0620095678006706402851392309608029321904000_dp, &
151 & 0.0698541213187282587095200770991474757860450_dp, 0.0768496807577203788944327774826590067221100_dp, &
152 & 0.0830805028231330210382892472861037896015540_dp, 0.0885644430562117706472754436937743032122700_dp, &
153 & 0.0931265981708253212254868727473457185619300_dp, 0.0966427269836236785051799076275893351366570_dp, &
154 & 0.0991735987217919593323931734846031310595673_dp, 0.1007698455238755950449466626175697219163500_dp, &
155 & 0.1013300070147915490173747927674925467709270_dp, 0.1007698455238755950449466626175697219163500_dp, &
156 & 0.0991735987217919593323931734846031310595673_dp, 0.0966427269836236785051799076275893351366570_dp, &
157 & 0.0931265981708253212254868727473457185619300_dp, 0.0885644430562117706472754436937743032122700_dp, &
158 & 0.0830805028231330210382892472861037896015540_dp, 0.0768496807577203788944327774826590067221100_dp, &
159 & 0.0698541213187282587095200770991474757860450_dp, 0.0620095678006706402851392309608029321904000_dp, &
160 & 0.0534815246909280872653431472394302967715500_dp, 0.0445897513247648766082272993732796902232570_dp, &
161 & 0.0353463607913758462220379484783600481226300_dp, 0.0254608473267153201868740010196533593972700_dp, &
162 & 0.0150079473293161225383747630758072680946390_dp, 0.0053774798729233489877920514301276498183100_dp/)
163 
164  real(dp),public,parameter :: Gau15W(15) = (/0.030753241996117268354628393577204417721700_dp, &
165 0.070366047488108124709267416450667338466710_dp, 0.107159220467171935011869546685869303415544_dp, &
166 0.139570677926154314447804794511028322520850_dp, 0.166269205816993933553200860481208811130900_dp, &
167 0.186161000015562211026800561866422824506226_dp, 0.198431485327111576456118326443839324818693_dp, &
168 0.202578241925561272880620199967519314838662_dp, 0.198431485327111576456118326443839324818693_dp, &
169 0.186161000015562211026800561866422824506226_dp, 0.166269205816993933553200860481208811130900_dp, &
170 0.139570677926154314447804794511028322520850_dp, 0.107159220467171935011869546685869303415544_dp, &
171 0.070366047488108124709267416450667338466710_dp, 0.030753241996117268354628393577204417721700_dp/)
172 
173 
174 ! Flags for self-consistent GW calculations used in gw_driver and for parsing the input file.
175  integer,public,parameter :: GWSC_one_shot      =1
176  integer,public,parameter :: GWSC_only_W        =2
177  integer,public,parameter :: GWSC_only_G        =3
178  integer,public,parameter :: GWSC_both_G_and_W  =4
179 
180 ! Flags defining the approximation used for the self-energy (used in csigme).
181  integer,public,parameter :: SIG_GW_PPM      =0  ! standard GW with PPM
182  integer,public,parameter :: SIG_GW_AC       =1  ! standard GW without PPM (analytical continuation)
183  integer,public,parameter :: SIG_GW_CD       =2  ! standard GW without PPM (contour deformation)
184  integer,public,parameter :: SIG_HF          =5  ! Hartree-Fock calculation
185  integer,public,parameter :: SIG_SEX         =6  ! Screened Exchange calculation
186  integer,public,parameter :: SIG_COHSEX      =7  ! COHSEX calculation
187  integer,public,parameter :: SIG_QPGW_PPM    =8  ! model GW with PPM
188  integer,public,parameter :: SIG_QPGW_CD     =9  ! model GW without PPM
189 
190  public :: sigma_type_from_key
191  public :: sigma_is_herm
192  public :: sigma_needs_w
193  public :: sigma_needs_ppm
194  !public :: sigma_sc_on_wfs
195  !public :: sigma_sc_on_ene
196  public :: g0g0w
197 
198 ! Private variables
199  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

PARENTS

      screening

CHILDREN

      sigijtab_free

SOURCE

418 subroutine em1params_free(Ep)
419 
420 
421 !This section has been created automatically by the script Abilint (TD).
422 !Do not modify the following lines by hand.
423 #undef ABI_FUNC
424 #define ABI_FUNC 'em1params_free'
425 !End of the abilint section
426 
427  implicit none
428 
429 !Arguments ------------------------------------
430 !scalars
431  type(em1params_t),intent(inout) :: Ep
432 
433 ! *************************************************************************
434 
435  !@em1params_t
436 
437 !real
438  if (allocated(Ep%qcalc)) then
439    ABI_FREE(Ep%qcalc)
440  end if
441  if (allocated(Ep%qibz)) then
442    ABI_FREE(Ep%qibz)
443  end if
444  if (allocated(Ep%qlwl)) then
445    ABI_FREE(Ep%qlwl)
446  end if
447  if (allocated(Ep%omegasf)) then
448    ABI_FREE(Ep%omegasf)
449  end if
450 
451 !complex
452  if (allocated(Ep%omega)) then
453    ABI_FREE(Ep%omega)
454  end if
455 
456 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

214  type,public :: em1params_t
215 
216 !scalars
217   integer :: awtr                   ! If 1 the Adler-Wiser expression for Chi_0 is evaluated
218                                     !  taking advantage of time-reversal symmetry
219   integer :: gwcalctyp              ! Calculation type (see input variable)
220   integer :: gwcomp                 ! 1 if extrapolar technique is used. 0 otherwise.
221   integer :: inclvkb                ! Integer flag related to the evaluation of the commutator for q-->0
222   integer :: spmeth                 ! Method used to approximate the delta function in the expression for Im Chi_0
223   integer :: nI                     ! Number of components (rows) in the chi0 matrix.
224   integer :: nJ                     ! Number of components (columns) in the chi0 matrix.
225   integer :: npwvec                 ! Max between npwe and npwwfn, used to pass the dimension of arrays e.g gvec
226   integer :: npwwfn                 ! Number of planewaves for wavefunctions
227   integer :: npwe                   ! Number of planewaves for $\tilde \epsilon$
228   integer :: npwepG0                ! Number of planewaves in the enlarged sphere G-G0, to account for umklapp G0 vectors
229   integer :: nbnds                  ! Number of bands used to evaluate $\tilde \epsilon$
230   integer :: nkibz                  ! Number of k-points in the IBZ
231   integer :: nsppol                 ! 1 for spin unpolarized, 2 for collinear spin polarized
232   integer :: nqcalc                 ! Number of q-points that are calculated (subset of qibz)
233   integer :: nqibz                  ! Number of q-points in the IBZ
234   integer :: nqlwl                  ! Number of directions to analyze the non analytical behavior for q-->0
235   integer :: nomega                 ! Number of frequencies where evaluate $\tilde \epsilon (\omega)$
236   integer :: nomegaer,nomegaei      ! Number of real and imaginary frequencies, respectively
237   integer :: nomegaec               ! Number of frequencies on a grid in the complex plane nomegaec = nomegaei*(nomegaer-1)
238   integer :: nomegasf               ! Number of frequencies used for the spectral function
239   integer :: symchi                 ! 0 ==> do not use symmetries to reduce the k-points summed over in chi0
240                                     ! 1 ==> take advantage of point group symmetries as well as time-reversal
241 
242   real(dp) :: gwencomp              ! Extrapolar energy used if gwcomp==1.
243   real(dp) :: omegaermin            ! Minimum real frequency used in the contour deformation method
244   real(dp) :: omegaermax            ! Maximum real frequency used in the contour deformation method
245   real(dp) :: mbpt_sciss              ! Scissor energy used in chi0
246   real(dp) :: spsmear               ! Smearing of the delta in case of spmeth==2
247   real(dp) :: zcut                  ! Small imaginary shift to avoid poles in chi0
248 
249   logical :: analytic_continuation  ! if true calculate chi0 only along the imaginary axis
250   logical :: contour_deformation    ! if true calculated chi0 both along the real and the imaginary axis
251   logical :: plasmon_pole_model     ! if true a plasmonpole model is used (only 1 or 2 frequencies are calculated)
252 
253 !arrays
254   integer :: mG0(3)
255   ! For each reduced direction gives the max G0 component to account for umklapp processes
256 
257   real(dp),allocatable :: qcalc(:,:)
258   ! qcalc(3,nqcalc)
259   ! q-points that are explicitely calculated (subset of qibz).
260 
261   real(dp),allocatable :: qibz(:,:)
262   ! qibz(3,nqibz)
263   ! q-points in the IBZ.
264 
265   real(dp),allocatable :: qlwl(:,:)
266   ! qlwl(3,nqlwl)
267   ! q-points used for the long-wavelength limit.
268 
269   real(dp),allocatable :: omegasf(:)
270   ! omegasf(nomegasf)
271   ! real frequencies used to calculate the imaginary part of chi0.
272 
273   complex(dpc),allocatable :: omega(:)
274   ! omega(nomegasf)
275   ! real and imaginary frequencies in chi0,epsilon and epsilonm1.
276 
277  end type em1params_t
278 
279  public :: em1params_free

m_gwdefs/g0g0w [ Functions ]

[ Top ] [ m_gwdefs ] [ Functions ]

NAME

 g0g0w

FUNCTION

  Calculates the frequency-dependent part of the RPA polarizability G0G0.

INPUTS

OUTPUT

SOURCE

805 function g0g0w(omega,numerator,delta_ene,zcut,TOL_W0,opt_poles)
806 
807 
808 !This section has been created automatically by the script Abilint (TD).
809 !Do not modify the following lines by hand.
810 #undef ABI_FUNC
811 #define ABI_FUNC 'g0g0w'
812 !End of the abilint section
813 
814  implicit none
815 
816 !Arguments ------------------------------------
817 !scalars
818  integer,intent(in):: opt_poles
819  real(dp),intent(in) :: TOL_W0,delta_ene,numerator,zcut
820  complex(dpc) :: g0g0w
821  complex(dpc),intent(in) :: omega
822 
823 !Local variables ------------------------------
824 !scalars
825  real(dp) :: sgn
826  character(len=500) :: msg
827 
828 !************************************************************************
829 
830  if (delta_ene**2>tol14) then
831    sgn=delta_ene/ABS(delta_ene)
832    !
833    if (opt_poles == 2) then ! Resonant and anti-resonant contributions.
834      if (DABS(REAL(omega))>TOL_W0) then
835        g0g0w =  numerator / (omega + delta_ene - j_dpc*sgn*zcut)&
836 &              -numerator / (omega - delta_ene + j_dpc*sgn*zcut)
837      else
838        g0g0w =  numerator / (omega + delta_ene)&
839 &              -numerator / (omega - delta_ene)
840      end if
841 
842    else if (opt_poles == 1) then ! Only resonant contribution is included.
843      if (DABS(REAL(omega))>TOL_W0) then
844        g0g0w =  numerator / (omega + delta_ene - j_dpc*sgn*zcut)
845      else
846        g0g0w =  numerator / (omega + delta_ene)
847      end if
848 
849    else
850      write(msg,'(a,i0)')" Wrong value for opt_poles: ",opt_poles
851      MSG_ERROR(msg)
852    end if ! opt_poles
853 
854  else ! delta_ene**2<tol14
855    g0g0w = czero
856  end if
857 
858 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

PARENTS

      m_gwdefs,setup_sigma

CHILDREN

      sigijtab_free

SOURCE

480 subroutine sigijtab_free(Sigijtab)
481 
482 
483 !This section has been created automatically by the script Abilint (TD).
484 !Do not modify the following lines by hand.
485 #undef ABI_FUNC
486 #define ABI_FUNC 'sigijtab_free'
487 !End of the abilint section
488 
489  implicit none
490 
491 !Arguments ------------------------------------
492 !scalars
493  type(sigijtab_t),intent(inout) :: Sigijtab(:,:)
494 
495 !Local variables
496  integer :: ii,jj,kk,ilow,iup
497 ! *************************************************************************
498 
499  !@sigijtab_t
500   do jj=1,SIZE(Sigijtab,DIM=2)
501     do ii=1,SIZE(Sigijtab,DIM=1)
502 
503       ilow=LBOUND(Sigijtab(ii,jj)%col,DIM=1)
504       iup =UBOUND(Sigijtab(ii,jj)%col,DIM=1)
505       do kk=ilow,iup
506         ABI_FREE(Sigijtab(ii,jj)%col(kk)%bidx)
507       end do
508       ABI_DT_FREE(Sigijtab(ii,jj)%col)
509 
510     end do
511   end do
512 
513 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.

PARENTS

CHILDREN

SOURCE

673 pure logical function sigma_is_herm(Sigp)
674 
675 
676 !This section has been created automatically by the script Abilint (TD).
677 !Do not modify the following lines by hand.
678 #undef ABI_FUNC
679 #define ABI_FUNC 'sigma_is_herm'
680 !End of the abilint section
681 
682  implicit none
683 
684 !Arguments ------------------------------
685 !scalars
686  type(sigparams_t),intent(in) :: Sigp
687 
688 !Local variables ------------------------------
689  integer :: mod10
690 
691 !************************************************************************
692 
693  mod10=MOD(Sigp%gwcalctyp,10)
694  sigma_is_herm = ANY(mod10 == [SIG_HF, SIG_SEX, SIG_COHSEX])
695 
696 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.

PARENTS

CHILDREN

SOURCE

762 pure logical function sigma_needs_ppm(Sigp)
763 
764 
765 !This section has been created automatically by the script Abilint (TD).
766 !Do not modify the following lines by hand.
767 #undef ABI_FUNC
768 #define ABI_FUNC 'sigma_needs_ppm'
769 !End of the abilint section
770 
771  implicit none
772 
773 !Arguments ------------------------------
774 !scalars
775  type(sigparams_t),intent(in) :: Sigp
776 
777 !Local variables ------------------------------
778  integer :: mod10
779 
780 !************************************************************************
781 
782  mod10=MOD(Sigp%gwcalctyp,10)
783  sigma_needs_ppm = ( ANY(mod10 == (/SIG_GW_PPM, SIG_QPGW_PPM/)) .or. &
784 &                   Sigp%gwcomp==1                                   &
785 &                  )
786 
787 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.

PARENTS

CHILDREN

SOURCE

718 pure logical function sigma_needs_w(Sigp)
719 
720 
721 !This section has been created automatically by the script Abilint (TD).
722 !Do not modify the following lines by hand.
723 #undef ABI_FUNC
724 #define ABI_FUNC 'sigma_needs_w'
725 !End of the abilint section
726 
727  implicit none
728 
729 !Arguments ------------------------------
730 !scalars
731  type(sigparams_t),intent(in) :: Sigp
732 
733 !Local variables ------------------------------
734  integer :: mod10
735 
736 !************************************************************************
737 
738  mod10=MOD(Sigp%gwcalctyp,10)
739  sigma_needs_w = (mod10/=SIG_HF)
740 
741 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

PARENTS

CHILDREN

SOURCE

617 function sigma_type_from_key(key) result(sigma_type)
618 
619 
620 !This section has been created automatically by the script Abilint (TD).
621 !Do not modify the following lines by hand.
622 #undef ABI_FUNC
623 #define ABI_FUNC 'sigma_type_from_key'
624 !End of the abilint section
625 
626  implicit none
627 
628  integer,intent(in) :: key
629  character(len=STR_LEN) :: sigma_type
630 
631 !Local variables ------------------------------
632 !scalars
633  character(len=500) :: msg
634 
635 !************************************************************************
636 
637  sigma_type = "None"
638  if (key==SIG_GW_PPM  )  sigma_type = ' standard GW with PPM'
639  if (key==SIG_GW_AC   )  sigma_type = ' standard GW without PPM (analytical continuation)'
640  if (key==SIG_GW_CD   )  sigma_type = ' standard GW without PPM (contour deformation)'
641  if (key==SIG_HF      )  sigma_type = ' Hartree-Fock calculation'
642  if (key==SIG_SEX     )  sigma_type = ' Screened Exchange calculation'
643  if (key==SIG_COHSEX  )  sigma_type = ' COHSEX calculation'
644  if (key==SIG_QPGW_PPM)  sigma_type = ' model GW with PPM'
645  if (key==SIG_QPGW_CD )  sigma_type = ' model GW without PPM'
646 
647  if (sigma_type == "None") then
648    write(msg,'(a,i0)')" Unknown value for key= ",key
649    MSG_ERROR(msg)
650  end if
651 
652 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

PARENTS

      sigma

CHILDREN

      sigijtab_free

SOURCE

537 subroutine sigparams_free(Sigp)
538 
539 
540 !This section has been created automatically by the script Abilint (TD).
541 !Do not modify the following lines by hand.
542 #undef ABI_FUNC
543 #define ABI_FUNC 'sigparams_free'
544 !End of the abilint section
545 
546  implicit none
547 
548 !Arguments ------------------------------------
549 !scalars
550  type(sigparams_t),intent(inout) :: Sigp
551 
552 ! *************************************************************************
553 
554  !@sigparams_t
555 
556 !integer
557  if (allocated(Sigp%kptgw2bz)) then
558    ABI_FREE(Sigp%kptgw2bz)
559  end if
560  if (allocated(Sigp%minbnd)) then
561    ABI_FREE(Sigp%minbnd)
562  end if
563  if (allocated(Sigp%maxbnd)) then
564    ABI_FREE(Sigp%maxbnd)
565  end if
566 
567 !real
568  if (allocated(Sigp%kptgw)) then
569    ABI_FREE(Sigp%kptgw)
570  end if
571 
572 !complex
573  if (allocated(Sigp%omegasi)) then
574    ABI_FREE(Sigp%omegasi)
575  end if
576  if (allocated(Sigp%omega_r)) then
577    ABI_FREE(Sigp%omega_r)
578  end if
579 
580 !logical
581 
582 !types
583  if (allocated(Sigp%Sigcij_tab)) then
584    call sigijtab_free(Sigp%Sigcij_tab)
585    ABI_DT_FREE(Sigp%Sigcij_tab)
586  end if
587 
588  if (allocated(Sigp%Sigxij_tab)) then
589    call sigijtab_free(Sigp%Sigxij_tab)
590    ABI_DT_FREE(Sigp%Sigxij_tab)
591  end if
592 
593 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

306  type,public :: sigparams_t
307 
308   integer :: gwcalctyp                   ! Calculation type
309 
310   integer :: gwgamma                     ! If 1 include vertex correction (GWGamma)
311   integer :: gwcomp                      ! 1 if the extrapolar technique is used.
312 
313   integer :: minbdgw,maxbdgw             ! Minimum and maximum band index (considering the spin) defining
314                                          ! The set of bands where GW corrections are evaluated
315 
316   integer :: mG0(3)                      ! For each reduced direction gives the max G0 component
317                                          ! to account for umklapp processes
318 
319   integer :: npwvec                      ! Max betwenn npwe and npwwfn, used to pass the dimension of arrays e.g gvec
320   integer :: npwwfn                      ! No. of planewaves for wavefunctions
321   integer :: npwx                        ! No. of planewaves for $\Sigma_x$
322   integer :: npwc                        ! No. of planewaves for $\Sigma_c$ and W
323   integer :: nbnds                       ! No. of bands summed over.
324   integer :: nomegasr                    ! No. of frequencies on the real axis to evaluate the spectral function
325   integer :: nomegasrd                   ! No. of frequencies on the real axis to evaluate $\Sigma(E)$
326   integer :: nomegasi                    ! No. of frequencies along the imaginary axis for Sigma in case of AC
327   integer :: nsig_ab                     ! No. of components in the self-energy operator (1 if nspinor==1, 4 if nspinor==2)
328   integer :: nspinor                     ! No. of spinorial components.
329   integer :: nsppol                      ! 1 for unpolarized, 2 for spin-polarized calculation
330   integer :: nkptgw                      ! No. of k-points where GW corrections have been calculated
331   integer :: ppmodel                     ! Integer defining the plasmon pole model used, 0 for None.
332   integer :: symsigma                    ! 0 ==> do not use symmetries to reduce the k-points summed over in sigma
333                                          ! 1 ==> take advantage of space group symmetries as well as time-reversal
334   integer :: use_sigxcore                ! 1 if core contribution to sigma is estimated by using Hartree-Fock
335 
336   real(dp) :: deltae                     ! Energy step used to evaluate numerically the derivative of the self energy
337                                          ! $\frac{\partial \Re \Sigma(E)}{\partial E_o}$
338   real(dp) :: ecutwfn                    ! cutoff energy for the wavefunctions.
339   real(dp) :: ecutsigx                   ! cutoff energy for the the exchange parth of Sigma.
340   real(dp) :: ecuteps                    ! cutoff energy for W
341 
342   real(dp) :: gwencomp                   ! Extrapolar energy used if gwcomp==1.
343 
344   real(dp) :: mbpt_sciss                 ! Scissor energy used in G0
345 
346   real(dp) :: minomega_r                 ! Minimum real frequency for the evaluation of the spectral function
347   real(dp) :: maxomega_r                 ! Maximum real frequency for the evaluation of the spectral function
348   real(dp) :: maxomega4sd                ! Maximum displacement around the KS energy where evaluate the diagonal
349                                          ! Elements of $ \Sigma(E)$
350   real(dp) :: omegasimax                 ! Max omega for Sigma along the imag axis in case of analytic continuation
351   real(dp) :: omegasimin                 ! min omega for Sigma along the imag axis in case of analytic continuation
352 
353   real(dp) :: sigma_mixing               ! Global factor that multiplies Sigma to give the final matrix element. 
354                                          ! Usually one, except for the hybrid functionals.
355 
356   real(dp) :: zcut                       ! Value of $\delta$ used to avoid the divergences (see related input variable)
357 
358   integer,allocatable :: kptgw2bz(:)
359   ! kptgw2bz(nkptgw)
360   ! For each k-point where GW corrections are calculated, the corresponding index in the BZ.
361 
362   integer,allocatable :: minbnd(:,:), maxbnd(:,:)
363   ! minbnd(nkptgw,nsppol), maxbnd(nkptgw,nsppol)
364   ! For each k-point at which GW corrections are calculated, the min and Max band index considered
365   ! (see also input variable dtset%bdgw).
366 
367   real(dp),allocatable :: kptgw(:,:)
368   ! kptgw(3,nkptgw)
369   ! k-points for the GW corrections in reduced coordinates.
370 
371   !TODO should be removed, everything should be in Sr%
372 
373   complex(dpc),allocatable :: omegasi(:)
374   ! omegasi(nomegasi)
375   ! Frequencies along the imaginary axis used for the analytical continuation.
376 
377   complex(dpc),allocatable :: omega_r(:)
378   ! omega_r(nomegasr)
379   ! Frequencies used to evaluate the spectral function.
380 
381   type(sigijtab_t),allocatable :: Sigcij_tab(:,:)
382   ! Sigcij_tab(nkptgw,nsppol)%col(kb)%bidx(ii)  gived the index of the left wavefunction.
383   ! in the <i,kgw,s|\Sigma_c|j,kgw,s> matrix elements that has to be calculated in cisgme.
384   ! in the case of self-consistent GW on wavefunctions.
385 
386   type(sigijtab_t),allocatable :: Sigxij_tab(:,:)
387   ! Save as Sigcij_tab but for the Hermitian \Sigma_x where only the upper triangle is needed.
388 
389  end type sigparams_t
390 
391  public :: sigparams_free