TABLE OF CONTENTS


ABINIT/smeared_delta [ Functions ]

[ Top ] [ Functions ]

NAME

 smeared_delta

FUNCTION

 This subroutine calculates the smeared delta that weights matrix elements: 
 \delta (\epsilon_{kn}-\epsilon_{k+Q,n'})

COPYRIGHT

 Copyright (C) 1999-2018 ABINIT group (PB, XG)
 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 .

INPUTS

 eigen0(mband*nsppol) : Eigenvalues at point K 
 eigenq(mband*nsppol)  : Eigenvalues at point K+Q
 mband : maximum number of bands
 smdelta : Variable controlling the smearinf scheme

OUTPUT

 smdfunc(mband,mband) : Smeared delta function weight corresponding to \delta(\epsilon_{n,k} - \epsilon_{n',k+Q})

PARENTS

      eig2stern,eig2tot

CHILDREN

SOURCE

 33 #if defined HAVE_CONFIG_H
 34 #include "config.h"
 35 #endif
 36 
 37 #include "abi_common.h"
 38 
 39 
 40 subroutine smeared_delta(eigen0,eigenq,esmear,mband,smdelta,smdfunc)
 41 
 42  use defs_basis
 43  use m_errors
 44  use m_profiling_abi
 45 
 46 !This section has been created automatically by the script Abilint (TD).
 47 !Do not modify the following lines by hand.
 48 #undef ABI_FUNC
 49 #define ABI_FUNC 'smeared_delta'
 50 !End of the abilint section
 51 
 52  implicit none
 53 
 54 !Arguments ------------------------------------
 55 !scalars
 56  integer,intent(in) :: mband,smdelta
 57 !arrays
 58  real(dp),intent(in) :: eigen0(mband),eigenq(mband),esmear
 59  real(dp),intent(out) :: smdfunc(mband,mband)
 60 
 61 !Local variables-------------------------------
 62 !tolerance for non degenerated levels
 63 !scalars
 64  integer :: ii,jj
 65  real(dp) :: aa,dsqrpi,gauss,xx
 66  character(len=500) :: message
 67 
 68 ! *********************************************************************
 69 
 70 
 71 !---------------------------------------------------------
 72 !Ordinary (unique) smearing function
 73 !---------------------------------------------------------
 74 
 75  if(smdelta==1)then
 76 
 77 !  Fermi-Dirac
 78    do ii=1,mband
 79      do jj= 1,mband
 80        xx= ( eigen0(ii) - eigenq(jj) )/esmear
 81        smdfunc(ii,jj)=0.25_dp/esmear/(cosh(xx/2.0_dp))**2
 82      end do
 83    end do
 84 
 85  else if(smdelta==2 .or. smdelta==3)then
 86 
 87 !  Cold smearing of Marzari, two values of the "a" parameter being possible
 88 !  first value gives minimization of the bump
 89    if(smdelta==2)aa=-.5634
 90 !  second value gives monotonic occupation function
 91    if(smdelta==3)aa=-.8165
 92 
 93    dsqrpi=1.0_dp/sqrt(pi)
 94    do ii=1,mband
 95      do jj=1,mband
 96        xx= ( eigen0(ii) - eigenq(jj) ) / esmear
 97        gauss=dsqrpi*exp(-xx**2)/esmear
 98        smdfunc(ii,jj)=gauss*(1.5_dp+xx*(-aa*1.5_dp+xx*(-1.0_dp+aa*xx)))
 99      end do
100    end do
101 
102  else if(smdelta==4)then
103 
104 !  First order Hermite-Gaussian of Paxton and Methfessel
105    dsqrpi=1.0_dp/sqrt(pi)
106    do ii=1,mband
107      do jj=1,mband
108        xx= ( eigen0(ii) - eigenq (jj) ) / esmear
109        smdfunc(ii,jj)=dsqrpi*(1.5_dp-xx**2)*exp(-xx**2)/esmear
110      end do
111    end do
112 
113  else if(smdelta==5)then
114 
115 !  Gaussian smearing
116    dsqrpi=1.0_dp/sqrt(pi)
117    do ii=1,mband
118      do jj=1,mband
119        xx= ( eigen0(ii) - eigenq (jj) ) / esmear
120        smdfunc(ii,jj)=dsqrpi*exp(-xx**2)/esmear
121      end do
122    end do
123 
124  else
125    write(message, '(a,i0,a)' )'  Smdelta= ',smdelta,' is not allowed in smdfunc'
126    MSG_BUG(message)
127  end if
128 
129 end subroutine smeared_delta