1! ---
2! Copyright (C) 1996-2016	The SIESTA group
3!  This file is distributed under the terms of the
4!  GNU General Public License: see COPYING in the top directory
5!  or http://www.gnu.org/copyleft/gpl.txt .
6! See Docs/Contributors.txt for a list of contributors.
7! ---
8
9      SUBROUTINE READSTS(VOLUME,EMIN, NE, DELTAE, ETA, XSTS, BFUNC,
10     .                   IUNITCD,  ARMUNI)
11
12
13C **********************************************************************
14C Reads from the input file the data to calculate STS spectra
15C
16C Coded by P. Ordejon, July 2004
17C **********************************************************************
18
19      USE FDF
20
21      IMPLICIT NONE
22
23      DOUBLE PRECISION, INTENT(IN) ::
24     .  VOLUME
25
26      INTEGER, INTENT(OUT) ::
27     .  NE, IUNITCD, BFUNC
28
29      DOUBLE PRECISION, INTENT(OUT) ::
30     .  EMIN, DELTAE, ETA, XSTS(3), ARMUNI
31
32
33C ****** INPUT  ********************************************************
34C REAL*8 VOLUME          : Volume of the unit cell
35C ****** OUTPUT ********************************************************
36C REAL*8  EMIN           : Minimum energy in STS spectrum
37C INTEGER NE             : Number of energy points
38C REAL*8  DELTAE         : Step in energy plot
39C REAL*8  ETA            : Lorenzian broadening of levels
40C REAL*8  XSTS(3)        : Position where LDOS is computed
41C INTEGER BFUNC          : Broadening function for STS spectra
42C                          1 = Lorentzian, 0 = Gaussian
43C INTEGER IUNITCD        : Units for the electron density
44C                          IUNITCD = 1 => Ele/(bohr)**3
45C                          IUNITCD = 2 => Ele/(Ang)**3
46C                          IUNITCD = 3 => Ele/(unitcell)
47C REAL*8  ARMUNI         : Conversion factors for the charge density
48C **********************************************************************
49
50C Internal variables ---------------------------------------------------
51C INTEGER ISCALE         : Units for the atomic positions
52C                          (ISCALE = 1 => Bohrs, ISCALE = 2 => Ang)
53      CHARACTER
54     .  CPF*22, CPF_DEFECT*22,
55     .  UCD*22, UCD_DEFECT*22,
56     .  BROAD*22, BROAD_DEFECT*22
57
58      INTEGER
59     .  IUNIT, IX, JX, ISCALE
60
61      DOUBLE PRECISION EMAX
62
63      LOGICAL
64     .  LEQI
65
66      EXTERNAL
67     .  LEQI
68C ----
69
70      CPF_DEFECT = 'Bohr'
71
72      CPF = FDF_STRING('Denchar.CoorUnits',CPF_DEFECT)
73
74      IF (LEQI(CPF,'bohr')) then
75        ISCALE = 1
76      ELSEIF (LEQI(CPF,'ang')) then
77        ISCALE = 2
78      ELSE
79        WRITE(6,*)'readsts: ERROR Denchar.CoorUnits must be Ang or Bohr'
80        STOP
81      ENDIF
82
83      UCD_DEFECT = 'Ele/bohr**3'
84      UCD = FDF_STRING('Denchar.DensityUnits',UCD_DEFECT)
85      IF (LEQI(UCD,'ele/bohr**3')) then
86        IUNITCD = 1
87      ELSEIF (LEQI(UCD,'ele/ang**3')) then
88        IUNITCD = 2
89      ELSEIF (LEQI(UCD,'ele/unitcell')) then
90        IUNITCD = 3
91      ELSE
92       WRITE(6,'(A)')' readsts: ERROR   Wrong Option in Units of      '
93       WRITE(6,'(A)')' readsts:  Charge Density                       '
94       WRITE(6,'(A)')' readsts:  You must choose one of the following:'
95       WRITE(6,'(A)')' readsts:                                       '
96       WRITE(6,'(A)')' readsts:      - Ele/bohr**3                    '
97       WRITE(6,'(A)')' readsts:      - Ele/ang**3                     '
98       WRITE(6,'(A)')' readsts:      - Ele/unitcell                   '
99       STOP
100      ENDIF
101
102
103      BROAD_DEFECT = 'Lorentzian'
104
105      BROAD = FDF_STRING('Denchar.STSBroadeningFunc',BROAD_DEFECT)
106
107      IF (LEQI(BROAD,'lorentzian')) then
108        BFUNC  = 1
109      ELSEIF (LEQI(BROAD,'gaussian')) then
110        BFUNC  = 2
111      ELSE
112        WRITE(6,*)'readsts: ERROR Denchar.STSBroadeningFunc must be'
113        WRITE(6,*)'readsts:       Lorentzian or Gaussian'
114        STOP
115      ENDIF
116
117      NE = FDF_INTEGER('Denchar.STSEnergyPoints',100)
118
119      IF ( FDF_BLOCK('Denchar.STSPosition',IUNIT) ) THEN
120        READ(IUNIT,*)(XSTS(IX),IX=1,3)
121      ELSE
122        WRITE(6,*)'To calculate STS spectra, you must'
123        WRITE(6,*)'specify a position (Denchar.STSPosition)'
124        WRITE(6,*)'in the input file'
125        STOP
126      ENDIF
127
128      EMIN = FDF_PHYSICAL('Denchar.STSEmin',-1.0D0,'eV')
129      EMAX = FDF_PHYSICAL('Denchar.STSEmax',1.0D0,'eV')
130      ETA  = FDF_PHYSICAL('Denchar.STSEta',0.1D0,'eV')
131
132      DELTAE = (EMAX-EMIN)/DFLOAT(NE-1)
133
134C Scale points coordinates
135C   Iscale = 1 => Do nothing
136C   Iscale = 2 => Multiply by 1./0.529177 (Ang --> Bohr)
137
138      IF (ISCALE .EQ. 2) THEN
139        DO IX = 1,3
140          XSTS(IX) = 1.D0 / 0.529177D0 * XSTS(IX)
141        ENDDO
142      ENDIF
143
144C Units of Charge Density
145C   Iunitcd = 1 => Do nothing
146C   Iunitcd = 2 => Multiply by (1.d0 / 0.529177d0) **3 (bohr**3 --> Ang**3)
147C   Iunitcd = 3 => Multiply by volume unit cell (in bohrs**3)
148
149      IF (IUNITCD .EQ. 1) THEN
150        ARMUNI = 1.D0
151      ELSEIF( IUNITCD .EQ. 2 ) THEN
152        ARMUNI = (1.D0 / 0.529177D0)**3
153      ELSEIF( IUNITCD .EQ. 3 ) THEN
154        ARMUNI = VOLUME
155      ENDIF
156
157      END
158
159