Statistics
| Revision:

root / src / ReadInput_siesta.f90 @ 5

History | View | Annotate | Download (11.4 kB)

1
 SUBROUTINE ReadInput_siesta
2

    
3
! This routine reads an input template for Siesta
4

    
5
  use VarTypes
6
  use Path_module
7
  use Io_module
8

    
9
  IMPLICIT NONE
10

    
11
  INTERFACE
12
     function valid(string) result (isValid)
13
       CHARACTER(*), intent(in) :: string
14
       logical                  :: isValid
15
     END function VALID
16

    
17
    FUNCTION SearchInput(Input,String,Line,Clean)  Result (Found)
18

    
19
      Use Vartypes
20
      use io_module
21
! Input
22
      TYPE (Input_line), POINTER, INTENT(IN) :: Input
23
      CHARACTER(*), INTENT(IN) :: String
24
      CHARACTER(*), OPTIONAL, INTENT(IN) :: Clean
25

    
26
! Output
27
      TYPE (Input_line), POINTER, INTENT(OUT) :: Line
28

    
29
      LOGICAL :: Found
30

    
31
    END FUNCTION SearchInput
32

    
33
    FUNCTION InString(Line,String,Case,Clean,Back)  Result(Pos)
34

    
35
      Use VarTypes
36

    
37
      implicit none
38
! Input
39
      CHARACTER(*), INTENT(IN) :: Line
40
      CHARACTER(*), INTENT(IN) :: String
41
      LOGICAL, OPTIONAL, INTENT(IN) :: Case
42
      LOGICAL, OPTIONAL, INTENT(IN) :: Back
43
      CHARACTER(*), OPTIONAL, INTENT(IN) :: Clean
44

    
45
! Output
46
! the position of String in Line (the first one) unless Back is present
47
      INTEGER(KINT) :: Pos
48
    END FUNCTION InString
49

    
50
    SUBROUTINE die(routine, msg, file, line, unit)
51

    
52
      Use VarTypes
53
      Use io_module
54

    
55
      implicit none
56
!--------------------------------------------------------------- Input Variables
57
      character(len=*), intent(in)           :: routine, msg
58
      character(len=*), intent(in), optional :: file
59
      integer(KINT), intent(in), optional      :: line, unit
60

    
61
    END SUBROUTINE die
62

    
63
    SUBROUTINE Warning(routine, msg, file, line, unit)
64

    
65
      Use VarTypes
66
      Use io_module
67

    
68
      implicit none
69
!--------------------------------------------------------------- Input Variables
70
      character(len=*), intent(in)           :: routine, msg
71
      character(len=*), intent(in), optional :: file
72
      integer(KINT), intent(in), optional      :: line, unit
73

    
74
    END SUBROUTINE Warning
75

    
76

    
77
 SUBROUTINE WriteList(Input,Unit)
78

    
79
! This routine reads an input template for Siesta
80

    
81
  use VarTypes
82
  use Io_module
83

    
84
  IMPLICIT NONE
85

    
86
! Input
87
      TYPE (Input_line), POINTER, INTENT(IN) :: Input
88
      INTEGER(KINT), OPTIONAL, INTENT(IN) :: Unit
89

    
90
    END SUBROUTINE WriteList
91
  END INTERFACE
92

    
93

    
94
  CHARACTER(132) ::  Line,Line2
95
  INTEGER(KINT) ::  Idx
96
  INTEGER(KINT) ::  IAt,I
97
  INTEGER(KINT) :: IoRead, ITmp
98
  REAL(KREAL) :: Xtmp, Ytmp, Ztmp
99
  
100
  LOGICAL :: Debug
101
  LOGICAL :: FSpecies, FCoord
102

    
103
  TYPE(Input_Line), POINTER :: Search,Bla
104

    
105
  Debug=Valid("readinput").OR.Valid("readinput_siesta")
106

    
107
 if (debug) Call Header("Entering ReadInput_Siesta")
108

    
109
  ! We read the Siesta input file
110

    
111
     ALLOCATE(Siesta_Input)
112
     NULLIFY(Siesta_Input%next)
113
     NULLIFY(Siesta_Input%prev)
114
     Current => Siesta_Input
115

    
116
     FSpecies=.FALSE.
117
     FCoord=.FALSE.
118

    
119
     READ(IOIN,'(A)',iostat=IoRead) Line
120

    
121
     DO WHILE (IoRead==0)
122
        Line=AdjustL(Line)
123
        if (debug) WRITE(*,*) 'Line:', Line
124
        current%Line=TRIM(Line)
125
        ALLOCATE(current%next)
126
        Current%next%prev => Current
127
        Current => Current%next
128
        current%line="toto"
129
        Nullify(Current%next)
130

    
131
        READ(IOIN,'(A)',iostat=IoRead) Line
132
     END DO
133
! With this procedure, the last current is not valid and should be deleted.
134
       Bla => Current%prev
135
       nullify(Bla%next)
136
       deallocate(Current)
137
       current => Bla
138

    
139
     if (debug) THEN
140
        WRITE(*,*) 'Input just read'
141
        Call WriteList(Siesta_input,Unit=IOOUT)
142
     END IF
143

    
144
! We analyse the input
145

    
146
! We look for the NumberofAtoms
147
     If (SearchInput(Siesta_Input,"NUMBEROFATOMS",Search,Clean=".-_")) THEN
148
        Line=AdjustL(Search%line)
149
        Idx=Index(Line," ")
150
        Line2=Trim(AdjustL(Line(Idx+1:)))
151
        Read(Line2,*) IAt
152
        if (Iat/=Nat) THEN
153
           Call Die('ReadInput_siesta','Nat in FDF sample different from  Nat Path input', Unit=IOOUT)
154
        END IF
155
     ELSE
156
! There is no atom number defined !!!
157
           Call Warning('ReadInput_siesta','No NumberOfAtoms in FDF sample', Unit=IOOUT)
158
           WRITE(current%Line,'(1X,"NumberOfAtoms ",I5)') Nat
159
           ALLOCATE(current%next)
160
           Current%next%prev => Current
161
           Current => Current%next
162
           Nullify(Current%next)
163
         
164
        END IF
165

    
166
! We look for the NumberofSpecies
167
        IF (SearchInput(Siesta_Input,"NUMBEROFSPECIES",Search,Clean=".-_")) THEN
168
           Line=AdjustL(Search%line)
169
           Idx=Index(Line," ")
170
           Line2=Trim(AdjustL(Line(Idx+1:)))
171
           Read(Line2,*) Siesta_NbSpecies
172
           ALLOCATE(ListSpecies(Siesta_NbSpecies))
173
           ALLOCATE(Siesta_SpeciesName(Siesta_NbSpecies))
174
        END IF
175

    
176
! We look for SystemLabel
177
        If (SearchInput(Siesta_Input,"SYSTEMLABEL",Search,Clean=".-_")) THEN
178
           Line=AdjustL(Search%line)
179
           Idx=Index(Line," ")
180
           Siesta_Label=Trim(adjustl(Line(Idx+1:)))
181
        ELSE
182
           Siesta_label='siesta'
183
        END IF
184

    
185
! We look for the ChemicalSpeciesLabel block
186
        IF (SearchInput(Siesta_Input,"CHEMICALSPECIESLABEL",Search,Clean=".-_")) THEN
187
           if (.NOT.ALLOCATED(IdxSpecies)) ALLOCATE(IdxSpecies(NAt))
188
           I=0
189
           Search => Search%next
190
           DO WHILE (InString(Search%line,'ENDBLOCK',Case=.FALSE.,Clean=".-_")==0)
191
              I=I+1
192
              if (I>Siesta_NbSpecies) THEN
193
                 Call Die('ReadInput_siesta:Reading ChemicalSpeciesLabel','Found more line in this block than NbSpecies !')
194
              END IF
195
              Line=AdjustL(Search%Line)
196
              Read(Line,*) ITmp
197
              Idx=Index(Line,' ')
198
              Line=AdjustL(Line(Idx+1:))
199
              Read(Line,*) Ztmp
200
              ListSpecies(ITmp)=ZTmp
201
              Idx=Index(Line,' ')
202
              Line=AdjustL(Line(Idx+1:))
203
              Idx=Index(Line,' ')
204
              Siesta_SpeciesName(ITmp)=Line(1:Idx-1)
205
              Search => Search%next
206
           END DO
207
           IF (I/=Siesta_NbSpecies) Call Die('ReadInput_siesta:Reading ChemicalSpeciesLabel', &
208
                'Number of lines in this block different from NbSpecies')
209
        ELSE
210
           Call Die('ReadInput_siesta:Reading ChemicalSpeciesLabel', &
211
                'Block ChemicalSpeciesLabel is mandatory !')
212
        END IF
213

    
214

    
215
! We look for the AtomicCoordinatesAndAtomicSpecies  block
216
        IF (SearchInput(Siesta_Input,"ATOMICCOORDINATESANDATOMICSPECIES",Search,Clean=".-_")) THEN
217
           ALLOCATE(Siesta_Paste(Nat))
218
           Current=>Search
219
           DO I=1,NAt
220
              Search => Search%next
221
              Read(Search%line,*) XTmp,YTmp,ZTmp,IdxSpecies(I)
222
! We save everything but the x,y,z, and species description              
223
              Line=AdjustL(search%line)
224
! We skip x
225
              Idx=Index(Line,' ')
226
              Line=AdjustL(Line(Idx+1:))
227
! we skip y
228
              Idx=Index(Line,' ')
229
              Line=AdjustL(Line(Idx+1:))
230
! we skip z
231
              Idx=Index(Line,' ')
232
              Line=AdjustL(Line(Idx+1:))
233
! we skip species
234
              Idx=Index(Line,' ')
235
              Line=AdjustL(Line(Idx+1:))
236
              Siesta_Paste(I)=TRIM(Line)
237
           END DO
238
           if (debug) THEN
239
              WRITE(*,*) 'Input before deletion of coord block'
240
              Call WriteList(Siesta_input)
241
           END IF
242
! We will now delete this block from our input sample as it will then be
243
! written directly by Opt'n Path
244
! Search%next point on %endblock
245
           search => search%next
246
           IF (ASSOCIATED(Search%next)) THEN
247
! we are not at the end of the input file
248
              if (ASSOCIATED(Current%Prev)) THEN
249
                 Search%next%prev => current%prev
250
                 current%prev%next => search%next
251
              ELSE
252
! the coordinate block is at the begining of the input
253
                 siesta_input => Search%next
254
                 Nullify(Siesta_Input%Prev)
255
              END IF
256
           ELSE
257
! the coordinate block is the last part of the input
258
              nullify(current%prev%next)
259
           END IF
260
           if (debug) THEN
261
              WRITE(*,*) 'Input after deletion of coord block 1'
262
              Call WriteList(Siesta_input)
263
           END IF
264

    
265
           DO I=1,Nat+2
266
              bla=>current
267
              current=> current%next
268
              deallocate(bla)
269
           END DO
270
           if (debug) THEN
271
              WRITE(*,*) 'Input after deletion of coord block 2'
272
              Call WriteList(Siesta_input)
273
           END IF
274
        ELSE
275
           IF (SearchInput(Siesta_Input,"ZMATRIX",Search,Clean=".-_")) THEN
276
              call Die('ReadInput_Siesta','For now, I need the full block'//&
277
                   'AtomicCoordinatesAndAtomicSpecies, ZMatrix not yet handled.')
278
           ELSE
279
              call Die('ReadInput_Siesta','For now, I need the full block' //&
280
              'AtomicCoordinatesAndAtomicSpecies.')
281
           END IF
282
        END IF
283

    
284
        IF (SearchInput(Siesta_Input,"ATOMICCOORDINATESFORMAT",Search,Clean=".-_")) THEN
285
           Line=Adjustl(Search%Line)
286
           Idx=Index(Line,' ')
287
           Line=AdjustL(Line(Idx+1:))
288
           Call UpCase(Line)
289
           SELECT CASE (Line) 
290
              CASE ('ANG','NOTSCALEDCARTESIANANG')
291
                 Siesta_Unit_Read=1.d0
292
              CASE ('FRACTIONAL','SCALEDBYLATTICEVECTORS')
293
              CASE ('BOHR','NOTSCALEDCARTESIANBOHR')
294
                 Siesta_Unit_Read=a0
295
                 IF (INPUT=='SIESTA') THEN
296
! We have read the coordinates, but not the unit. This is corrected here
297
                    XyZGeomI=XyzGeomI*a0
298
                 END IF
299
            END SELECT
300
         END IF
301

    
302
        IF (SearchInput(Siesta_Input,"ATOMICCOORFORMATOUT",Search,Clean=".-_")) THEN
303
           Line=Adjustl(Search%Line)
304
           Idx=Index(Line,' ')
305
           Line=AdjustL(Line(Idx+1:))
306
           Call UpCase(Line)
307
           SELECT CASE (Line) 
308
              CASE ('ANG','NOTSCALEDCARTESIANANG')
309
                 Siesta_Unit_Write=1.d0
310
              CASE ('FRACTIONAL','SCALEDBYLATTICEVECTORS')
311
              CASE ('BOHR','NOTSCALEDCARTESIANBOHR')
312
                 Siesta_Unit_Write=Unit
313
            END SELECT
314
         END IF
315

    
316
         IF (SearchInput(Siesta_Input,"LATTICECONSTANT",Search,Clean=".-_")) THEN
317
           Line=Adjustl(Search%Line)
318
! We discar the label
319
           Idx=Index(Line,' ')
320
           Line=AdjustL(Line(Idx+1:))
321
! we read the value
322
           Read(Line,*) Siesta_LatticeConstant
323
! We discard the value
324
           Idx=Index(Line,' ')
325
           Line=AdjustL(Line(Idx+1:))
326
! We read the unit
327
           Call UpCase(Line)
328
           SELECT CASE (Line) 
329
              CASE ('ANG')
330
                 Siesta_Lat_Unit=1.d0
331
              CASE ('BOHR')
332
                 Siesta_Lat_Unit=Unit
333
            END SELECT
334
! for now!
335
            Call Die('ReadInput_siesta:LatticeConstant',"For now, periodic calculations are NOT possible in Opt'n Path")
336

    
337
         END IF
338

    
339
         IF (SearchInput(Siesta_Input,"LATTICEVECTORS",Search,Clean=".-_")) THEN
340

    
341
            Call Die('ReadInput_siesta:LatticeVectors',"For now, periodic calculations are NOT possible in Opt'n Path")
342

    
343
         END IF
344

    
345
         IF (SearchInput(Siesta_Input,"LATTICEPARAMETERS",Search,Clean=".-_")) THEN
346

    
347
            Call Die('ReadInput_siesta:LatticeParameters',"For now, periodic calculations are NOT possible in Opt'n Path")
348

    
349
         END IF
350

    
351
         IF (SearchInput(Siesta_Input,"SUPERCELL",Search,Clean=".-_")) THEN
352

    
353
            Call Die('ReadInput_siesta:SuperCell',"SuperCell  NOT possible in Opt'n Path")
354

    
355
         END IF
356

    
357
!         Call Die('Readinput_siesta:fin',' Lecture finie')
358

    
359
     if (debug) Call Header("Exiting ReadInput_Siesta")
360

    
361
   END SUBROUTINE READINPUT_SIESTA