TABLE OF CONTENTS


ABINIT/dos_hdr_write [ Functions ]

[ Top ] [ Functions ]

NAME

 dos_hdr_write

FUNCTION

 Write the header of the DOS files, for both smearing and tetrahedron methods.

COPYRIGHT

 Copyright (C) 1998-2017 ABINIT group (XG, AF)
 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

 deltaene=increment of DOS energy arguments
 enemax=maximal value of the DOS energy argument
 enemin=minimal value of the DOS energy argument
 nene=number of DOS energy argument
 buffer=approximative buffer energy for the output of the DOS
  (beyond the max and min energy values).
 dosdeltae=DOS delta of Energy (if zero, take default values)
 eigen(mband*nkpt*nsppol)=eigenvalues (input or init to large number), hartree
 fermie=fermi energy useful for band alignment...
 mband=maximum number of bands
 nband(nkpt*nsppol)=number of bands at each k point
 nkpt=number of k points
 nsppol=1 for unpolarized, 2 for spin-polarized
 occopt=option for occupancies, or re-smearing scheme if dblsmr /= 0
 prtdos=1 for smearing technique, 2 or 3 for tetrahedron technique
 tphysel="physical" electronic temperature with FD occupations
 tsmear=smearing width (or temperature)
 unitdos=unit number of output of the DOS.

OUTPUT

   Only writing.

PARENTS

      getnel,m_epjdos

CHILDREN

      wrtout

SOURCE

 48 #if defined HAVE_CONFIG_H
 49 #include "config.h"
 50 #endif
 51 
 52 #include "abi_common.h"
 53 
 54 
 55 subroutine dos_hdr_write(buffer,deltaene,dosdeltae,&
 56 &  eigen,enemax,enemin,fermie,mband,nband,nene,&
 57 &  nkpt,nsppol,occopt,prtdos,tphysel,tsmear,unitdos)
 58 
 59  use defs_basis
 60  use m_profiling_abi
 61 
 62 !This section has been created automatically by the script Abilint (TD).
 63 !Do not modify the following lines by hand.
 64 #undef ABI_FUNC
 65 #define ABI_FUNC 'dos_hdr_write'
 66  use interfaces_14_hidewrite
 67 !End of the abilint section
 68 
 69  implicit none
 70 
 71 !Arguments ------------------------------------
 72 !scalars
 73  integer,intent(in) :: mband,nkpt,nsppol,occopt,prtdos,unitdos,nene
 74  real(dp),intent(in) :: buffer,dosdeltae,fermie,tphysel,tsmear
 75  real(dp),intent(in) :: deltaene,enemax,enemin
 76 !arrays
 77  integer,intent(in) :: nband(nkpt*nsppol)
 78  real(dp),intent(in) :: eigen(mband*nkpt*nsppol)
 79 
 80 !Local variables-------------------------------
 81 !scalars
 82  character(len=500) :: msg
 83 
 84 ! *************************************************************************
 85 
 86 !Write the DOS file
 87  write(msg, '(7a,i2,a,i5,a,i4)' ) "#",ch10, &
 88 & '# ABINIT package : DOS file  ',ch10,"#",ch10,&
 89 & '# nsppol =',nsppol,', nkpt =',nkpt,', nband(1)=',nband(1)
 90  call wrtout(unitdos,msg,'COLL')
 91 
 92  if (any(prtdos== [1,4])) then
 93    write(msg, '(a,i2,a,f6.3,a,f6.3,a)' )  &
 94 &   '# Smearing technique, occopt =',occopt,', tsmear=',tsmear,' Hartree, tphysel=',tphysel,' Hartree'
 95  else
 96    write(msg, '(a)' ) '# Tetrahedron method '
 97  end if
 98  call wrtout(unitdos,msg,'COLL')
 99 
100  if (mband*nkpt*nsppol>=3) then
101    write(msg, '(a,3f8.3,2a)' )'# For identification : eigen(1:3)=',eigen(1:3),ch10,"#"
102  else
103    write(msg, '(a,3f8.3)' ) '# For identification : eigen=',eigen
104    write(msg, '(3a)')trim(msg),ch10,"#"
105  end if
106  call wrtout(unitdos,msg,'COLL')
107 
108  write(msg, '(a,f16.8)' ) '# Fermi energy : ', fermie
109  call wrtout(unitdos,msg,'COLL')
110 
111  if (prtdos==1) then
112    write(msg, '(5a)' ) "#",ch10,&
113 &   '# The DOS (in electrons/Hartree/cell) and integrated DOS (in electrons/cell),',&
114 &   ch10,'# as well as the DOS with tsmear halved and doubled, are computed,'
115 
116  else if (prtdos==2)then
117    write(msg, '(3a)' ) "#",ch10,&
118 &   '# The DOS (in electrons/Hartree/cell) and integrated DOS (in electrons/cell) are computed,'
119 
120  else if (any(prtdos == [3, 4])) then
121    write(msg, '(5a)' ) "#",ch10,&
122 &   '# The local DOS (in electrons/Hartree for one atomic sphere)',ch10,&
123 &   '# and integrated local DOS (in electrons for one atomic sphere) are computed.'
124 
125  else if (prtdos==5)then
126    write(msg, '(9a)' ) "#",ch10,&
127 &   '# The spin component DOS (in electrons/Hartree/cell)',ch10,&
128 &   '# and integrated spin component DOS (in electrons/cell) are computed.',ch10,&
129 &   '# Remember that the wf are eigenstates of S_z and S^2, not S_x and S_y',ch10,&
130 &   '#   so the latter will not always sum to 0 for paired electronic states.'
131  end if
132  call wrtout(unitdos,msg,'COLL')
133 
134  write(msg, '(a,i5,a,a,a,f9.4,a,f9.4,a,f8.5,a,a,a)' )&
135 & '# at ',nene,' energies (in Hartree) covering the interval ',ch10,&
136 & '# between ',enemin,' and ',enemax,' Hartree by steps of ',deltaene,' Hartree.',ch10,"#"
137  call wrtout(unitdos,msg,'COLL')
138 
139  if (prtdos==1) then
140    write(msg, '(a,a)' )&
141 &   '#       energy        DOS       Integr. DOS   ','     DOS           DOS    '
142    call wrtout(unitdos,msg,'COLL')
143 
144    write(msg, '(a)' )&
145 &   '#                                              (tsmear/2)    (tsmear*2) '
146    call wrtout(unitdos,msg,'COLL')
147  else
148    write(msg, '(a)' ) '#       energy        DOS '
149  end if
150 
151 end subroutine dos_hdr_write