IffGit has a new shared runner for building Docker images in GitLab CI. Visit https://iffgit.fz-juelich.de/examples/ci-docker-in-docker for more details.

atom_input.f 38.5 KB
Newer Older
1
2
3
4
5
6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

7
8
9
      MODULE m_atominput
      use m_juDFT

10
      INTEGER, PARAMETER  :: dbgfh=6, errfh=6, warnfh=6
11
12
13
14
15
16
17
18
      REAL, PARAMETER     :: eps=0.00000001

      CONTAINS
!***********************************************************************
!     reads in the parameters associated with atom types from input
!     file. Part of inp-generator
!***********************************************************************
      SUBROUTINE atom_input(
19
     >                      infh,xl_buffer,bfh,buffer,
20
     >                      input,idlist,xmlCoreRefOccs,
21
     X                      nline,xmlElectronStates,
22
     X                      xmlPrintCoreStates,xmlCoreOccs,
23
     >                      atomTypeSpecies,numSpecies,
24
25
26
27
     <                      nel,atoms,enpara )

      USE m_types
      USE m_juDFT_init
28
29
      USE m_readrecord
      USE m_setatomcore, ONLY : setatom_bystr, setcore_bystr
30
      USE m_constants
31
    
32
33
      IMPLICIT NONE

34
      TYPE(t_input),INTENT(INOUT)    :: input
35
      TYPE(t_enpara),INTENT(OUT)     :: enpara
36
37
      TYPE(t_atoms),INTENT(INOUT)    :: atoms

38
39
! ... Arguments ...
      INTEGER, INTENT (IN)    :: infh  ! file number of input-file
40
      INTEGER, INTENT (IN)    :: bfh
41
      INTEGER, INTENT (IN)    :: numSpecies
42
43
44
      INTEGER, INTENT (INOUT) :: nline ! current line in this file
      INTEGER, INTENT (INOUT) :: nel   ! number of valence electrons

45
      INTEGER, INTENT (IN)     :: xl_buffer
46
      INTEGER, INTENT (IN)     :: atomTypeSpecies(atoms%ntype)
47
      REAL   , INTENT (IN)     :: idlist(atoms%ntype)
48
      REAL   , INTENT (IN)     :: xmlCoreRefOccs(29)
49
      REAL, INTENT (INOUT)     :: xmlCoreOccs(2,29,atoms%ntype)
50
      INTEGER, INTENT (INOUT)  :: xmlElectronStates(29,atoms%ntype)
51
      LOGICAL, INTENT (INOUT)  :: xmlPrintCoreStates(29,atoms%ntype)
52
53
54
55
56
57
58
59
60
61
62
63
      CHARACTER(len=xl_buffer) :: buffer

!===> data
      INTEGER, PARAMETER ::  l_buffer=512   ! maximum length of e-config string
      INTEGER, PARAMETER ::  nwdd=2         ! maximum number of windows
      INTEGER, PARAMETER ::  nstd=31        ! maximum number of core states

!===> Local Variables
      INTEGER :: nbuffer,ios,n,i,j,l,d1,d10,aoff,up,dn
      INTEGER :: lmax0_def,lnonsph0_def,jri0_def,ncst0_def
      INTEGER :: lmax0,lnonsph0,jri0,ncst0,nlod0,llod
      INTEGER :: natomst,ncorest,nvalst,z,nlo0
64
      INTEGER :: xmlCoreStateNumber, lmaxdTemp
65
      REAL    :: rmt0_def,dx0_def,bmu0_def, upReal, dnReal
66
      REAL    :: rmt0,dx0,bmu0,zat0,id,electronsOnAtom, electronsLeft
67
      LOGICAL :: fatalerror, h_atom, h_allatoms
68
69
70
71
72
      LOGICAL :: idone(atoms%ntype) 
      INTEGER :: lonqn(atoms%nlod,atoms%ntype),z_int(atoms%ntype)
      INTEGER :: coreqn(2,nstd,atoms%ntype),lval(nstd,atoms%ntype)
      INTEGER :: llo0(atoms%nlod)
      REAL    :: nelec(0:nwdd),coreocc(nstd,atoms%ntype)
73
74
75


      CHARACTER(len= l_buffer) :: econfig0_def,econfig0
76
77
      CHARACTER(len= l_buffer) :: econfig(atoms%ntype) ! verbose electronic config
      CHARACTER(len=80) :: lo(atoms%ntype), lo0
78
      CHARACTER(len=13) :: fname
79
      CHARACTER(LEN=20) :: speciesName0
80
81

      CHARACTER(len=1) :: lotype(0:3)
82

83
84
85
86
87
88
89
      DATA lotype /'s','p','d','f'/

!---> initialize some variables

      fatalerror = .false.
      h_atom=.false.;h_allatoms=.false.

90
91
92
93
      idone(1:atoms%ntype) = .false.
      z_int(1:atoms%ntype) = NINT(atoms%zatom(1:atoms%ntype))
      lo = ' ' ; nlod0 = 0 ; atoms%nlo = 0 ; llod = 0 ; lonqn = 0
      atoms%ncst = 0 ; econfig(1:atoms%ntype) = ' '
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
!
      lmax0_def    = -9999  
      lnonsph0_def = -9999
      rmt0_def     = -9999.9
      dx0_def      = -9999.9
      jri0_def     = -9999
      ncst0_def    = -9999
      econfig0_def = 'NONE'
      bmu0_def     = -9999.9

      WRITE(6,*)
      WRITE(6,'(a50)') '==============================================='
      WRITE(6,'(a50)') '===  modifying atomic input for &(all)atom  ==='
      WRITE(6,'(a50)') '==============================================='
      WRITE(6,*)

!===> continue reading input
      
      nbuffer = len_trim(buffer)

      IF ((buffer(1:9)=='&allatoms') .OR. 
     &    (buffer(1:5)=='&atom') .OR.
!    resetting nbuffer for &qss or &soc interferes with the lapw
!    namelist, therefore, its contributions are also checked for.
!    might interfere with other namelists, too.
!    Klueppelberg Jul 2012
     &    (buffer(1:5)=='&comp') .OR.
     &    (buffer(1:5)=='&exco') .OR.
     &    (buffer(1:5)=='&film') .OR.
     &    (buffer(1:4)=='&kpt') ) THEN
      ELSE
        nbuffer = 0 ! reset, to read in after &qss or &soc 
      ENDIF
      
      loop: DO

      IF (nbuffer == 0) THEN
        DO
132
          CALL read_record(infh,xl_buffer,bfh,nline,nbuffer,buffer,ios)
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
          IF (ios==1) GOTO 999
          IF (ios == 2)  CALL juDFT_error
     +         ("end of file while reading a record",
     +         calledby ="atom_input")
          IF (buffer(1:1)=='&') EXIT
          CALL err(0)
          fatalerror = .true.
        ENDDO
      ENDIF

!===> allatoms

      IF (buffer(1:9)=='&allatoms') THEN
        IF (h_allatoms) CALL err(1)
        h_allatoms = .true.
        IF (h_atom) then
          WRITE (errfh,*)
          WRITE (errfh,*) 'atom_input: ERROR',
     &     'namelist &allatoms must appear before namelist(s) &atom.'
          WRITE (errfh,*)
          fatalerror = .true.
        ELSE
!--->     read defaults for atom defaults
          CALL read_allatoms(
157
     >                       bfh,l_buffer,
158
159
160
161
162
163
     <                       rmt0_def,dx0_def,jri0_def,lmax0_def,
     <                       lnonsph0_def,ncst0_def,econfig0_def,
     <                       bmu0_def,ios)

          IF (ios.NE.0) GOTO 912
          IF (rmt0_def > -9999.8) THEN
164
            atoms%rmt     = rmt0_def
165
166
167
            WRITE (6,'(a25,f12.6)') 'globally changed rmt to',rmt0_def
          ENDIF
          IF (dx0_def  > -9999.8)   THEN
168
            atoms%dx      = dx0_def
169
170
171
            WRITE (6,'(a25,f12.6)') 'globally changed dx  to',dx0_def
          ENDIF
          IF (jri0_def > -9998  )   THEN
172
            atoms%jri     = jri0_def
173
174
175
            WRITE (6,'(a25,i12)') 'globally changed jri to',jri0_def
          ENDIF
          IF (lmax0_def > -9998 )   THEN
176
            atoms%lmax    = lmax0_def
177
178
179
180
            WRITE (6,'(a26,i12)') 'globally changed lmax to',
     &                                                       lmax0_def
          ENDIF
          IF (lnonsph0_def > -9998) THEN
181
            atoms%lnonsph = lnonsph0_def
182
183
184
185
            WRITE (6,'(a28,i12)') 'globally changed lnonsph to ',
     &                                                    lnonsph0_def
          ENDIF
          IF (ncst0_def > -9998 )   THEN
186
            atoms%ncst    = ncst0_def
187
188
189
190
191
192
193
194
195
            WRITE (6,'(a26,i12)') 'globally changed ncst to',
     &                                                       ncst0_def
          ENDIF
          IF (econfig0_def.NE.'NONE') THEN
            econfig = econfig0_def
            WRITE (6,'(a26,a80)') 'globally set econfig to ',
     &                                                    econfig0_def
          ENDIF
          IF (bmu0_def > -9999.8)   THEN
196
            atoms%bmu     = bmu0_def
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
            WRITE (6,'(a25,f12.6)') 'globally changed bmu to',bmu0_def
          ENDIF
        ENDIF

!===> atom
      ELSEIF (buffer(1:5)=='&atom') THEN
        h_atom=.true.

!--->   set atom defaults
        lmax0    = -9999  
        lnonsph0 = -9999
        rmt0     = -9999.9
        dx0      = -9999.9
        jri0     = -9999
        ncst0    = -9999
        econfig0 = 'NONE'
        bmu0     = -9999.9
        lo0      = ' '
215
        speciesName0 = ''
216
217
218

!--->   read namelist
        CALL read_atom(
219
     >                 bfh,l_buffer,lotype,
220
     <                 id,zat0,rmt0,jri0,dx0,lmax0,lnonsph0,
221
222
     <                 ncst0,econfig0,speciesName0,bmu0,lo0,nlod0,llod,
     <                 ios)
223
224
225
226
        IF (ios.ne.0) THEN
          CALL err(3)
        ELSE
!--->     put the data into the correct place
227
          DO n = 1, atoms%ntype
228
229
230
231
232
233
            IF (abs( id - idlist(n) ) > 0.001) CYCLE
            IF (idone(n)) then
              WRITE (errfh,*) 'atom_input: ERROR. did that one already'
              fatalerror=.true.
              EXIT
            ELSE
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
              IF (speciesName0.NE.'') THEN
                 atoms%speciesName(atomTypeSpecies(n)) = 
     &              TRIM(ADJUSTL(speciesName0))
                 DO i = 1, numSpecies
                    IF (i.NE.atomTypeSpecies(n)) THEN
                       IF((TRIM(ADJUSTL(speciesName0))).EQ.
     &                    (TRIM(ADJUSTL(atoms%speciesName(i))))) THEN
                          WRITE(*,*) ''
                          WRITE(*,*) 'Error for species name'
                          WRITE(*,*) TRIM(ADJUSTL(speciesName0))
                          WRITE(*,*) ''
                          CALL juDFT_error
     +                       ("Same name for different species",
     +                        calledby ="atom_input")
                       END IF
                    END IF
                 END DO
              END IF
252
              IF (rmt0 > -9999.8) THEN
253
                atoms%rmt(n)  = rmt0
254
255
256
257
                WRITE (6,'(a9,i4,2a2,a16,f12.6)') 'for atom ',n,
     &                ' (',namat_const(z_int(n)),') changed rmt to',rmt0
              ENDIF
              IF (dx0 > -9999.8) THEN
258
                atoms%dx(n)  = dx0
259
260
261
262
                WRITE (6,'(a9,i4,2a2,a16,f12.6)') 'for atom ',n,
     &                ' (',namat_const(z_int(n)),') changed dx  to', dx0
              ENDIF
              IF (jri0 > -9998  ) THEN
263
                atoms%jri(n)  = jri0
264
265
266
267
                WRITE (6,'(a9,i4,2a2,a16,i12)') 'for atom ',n,
     &                ' (',namat_const(z_int(n)),') changed jri to',jri0
              ENDIF
              IF (lmax0 > -9998  ) THEN
268
                atoms%lmax(n)  = lmax0
269
270
271
272
                WRITE (6,'(a9,i4,2a2,a17,i12)') 'for atom ',n,
     &              ' (',namat_const(z_int(n)),') changed lmax to',lmax0
              ENDIF
              IF (lnonsph0 > -9998  ) THEN
273
                atoms%lnonsph(n)  = lnonsph0
274
275
276
277
                WRITE (6,'(a9,i4,2a2,a20,i12)') 'for atom ',n,
     &        ' (',namat_const(z_int(n)),') changed lnonsph to',lnonsph0
              ENDIF
              IF (bmu0 > -9999.8  ) THEN
278
                atoms%bmu(n)  = bmu0
279
280
281
282
                WRITE (6,'(a9,i4,2a2,a16,f12.6)') 'for atom ',n,
     &              ' (',namat_const(z_int(n)),  ') changed bmu to',bmu0
              ENDIF
              IF (ncst0 > -9998  ) THEN
283
                atoms%ncst(n)  = ncst0
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
                WRITE (6,'(a9,i4,2a2,a17,i12)') 'for atom ',n,
     &             ' (',namat_const(z_int(n)), ') changed ncst to',ncst0
              ENDIF
! ===> electronic configuration
              IF (econfig0.NE.'NONE') THEN
                 econfig(n) = econfig0
                 WRITE (6,'(a9,i4,2a2,a17,a80)') 'for atom ',n,
     &           ' (',namat_const(z_int(n)),') set econfig to ',econfig0
                 CALL setatom_bystr(
     >                              l_buffer,nwdd,econfig(n),
     <                          natomst,ncorest,nvalst,nelec)
                 WRITE (6,'("   corestates =",i3," with",f6.1,
     &                                 " electrons")')  ncorest,nelec(0)
                 WRITE (6,'("   valence st.=",i3," with",f6.1,
     &                                 " electrons")')   nvalst,nelec(1)
                 IF (nelec(2) /= 0) THEN
                 WRITE (6,'("second window found!")')
                 WRITE (6,'("   valence st.=",i3," with",f6.1,
     &                                 " electrons")')   nvalst,nelec(2)
                 ENDIF
304
305
                 IF (nelec(0)+nelec(1)+nelec(2)-
     &               atoms%zatom(n)>0.01) THEN
306
307
308
309
310
311
312
                    CALL juDFT_error
     +                   ("econfig does not fit to this atom type!"
     +                   ,calledby ="atom_input")
                 ENDIF
                 IF (ncst0 > -9998  ) THEN
                   IF (ncorest /= ncst0) THEN
                     WRITE (6,'("  ==> core-states (ncst):",i3,
313
314
     &                              " =/= (econfig):",i3)') 
     &                        atoms%ncst,ncorest
315
316
317
318
319
                     CALL juDFT_error
     +                    ("econfig does not fit to the specified ncst"
     +                    ,calledby ="atom_input")
                   ENDIF
                 ELSE
320
                   atoms%ncst(n) = ncorest
321
322
323
324
                 ENDIF
              ENDIF
! ===> local orbitals
              IF (lo0 /= ' ') THEN
325
                WRITE (6,'(a6,i3,a7,i3,a3,a80)')
326
327
     &                     "nlod =",nlod0," llod =",llod," : ",lo0
                lo(n)      = lo0
328
329
                IF (nlod0 > atoms%nlod)  
     &             CALL juDFT_error("atom_input: too "
330
331
     &                              //"many lo",calledby="atom_input")

332
333
                atoms%nlo(n) = len_trim(lo(n))/2
                DO i = 1, atoms%nlo(n)
334
335
336
                  j = 2*i
                  DO l = 0, 3
                    IF (lo(n)(j:j) == lotype(l)) THEN
337
                      atoms%llo(i,n) = l
338
339
340
341
342
343
                    ENDIF
                  ENDDO
                  j = j - 1
                  READ (lo(n)(j:j),*) lonqn(i,n)
                ENDDO
                WRITE (6,'("   nlo(",i3,") = ",i2," llo = ",8i2)') n,
344
345
346
     &                    atoms%nlo(n),(atoms%llo(i,n),i=1,atoms%nlo(n))
                WRITE (6,'("   lonqn = ",8i2)') 
     &                    (lonqn(i,n),i=1,atoms%nlo(n))
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
              ENDIF
              idone(n)   = .true.
            ENDIF
          ENDDO
        ENDIF

!===> not an atom related namelist, we are done
      ELSE
        exit loop
      ENDIF

      nbuffer = 0
      ENDDO loop

 999  CONTINUE
      IF (fatalerror) 
     &   CALL juDFT_error("ERROR(S) reading input. Check output for "
     &                  //"details.",calledby="atom_input")

!----------- adjust the core-levels, lo's and the energy parameters ----

368
369
      coreqn(1:2,1:nstd,1:atoms%ntype) = 0
      coreocc(1:nstd,1:atoms%ntype) = -1.0
370

371
      nel = 0
372
373
374

      IF ( ANY(atoms%bmu(:) > 0.0) ) input%jspins=2 

375
376
      lmaxdTemp = atoms%lmaxd
      atoms%lmaxd = 3
377
      call enpara%init(atoms,input%jspins)
378
      DO n = 1, atoms%ntype
379
380

        CALL setcore_bystr(
381
     >                      n,nstd,atoms%ntype,l_buffer,
382
383
384
385
386
387
388
389
390
391
392
393
     X                      econfig,natomst,ncorest,
     <                      coreqn,coreocc)

        IF ( coreqn(1,1,n) /= 0 ) THEN
          DO i = 1, natomst
            IF (coreqn(2,i,n) < 0) THEN
               lval(i,n) = - coreqn(2,i,n) - 1
            ELSE
               lval(i,n) = coreqn(2,i,n)
            ENDIF
          ENDDO 

394
395
           d1  = mod(nint(atoms%zatom(n)),10)
           d10 = int( (nint(atoms%zatom(n)) + 0.5)/10 )
396
          aoff = iachar('1')-1
397
398
399
400
401
          IF(.NOT.input%l_inpXML) THEN
             fname = 'corelevels.'//achar(d10+aoff)//achar(d1+aoff)
             OPEN (27,file=fname,form='formatted')
             write(27,'(i3)') natomst
          END IF
402
403

          WRITE (6,*) '----------'
404
          electronsOnAtom = 0
405
406
407
408
409
410
411
          DO i = 1, ncorest
            WRITE(6,'("     core :",2i3,f6.1)') 
     &             coreqn(1,i,n),coreqn(2,i,n),coreocc(i,n)
            j = INT(coreocc(i,n) / 2)
            IF (coreocc(i,n) > 2*j) THEN
              j = - coreocc(i,n)
            ENDIF
412
413
414
            IF(.NOT.input%l_inpXML) THEN
               write(27,'(4i3)') coreqn(1,i,n),coreqn(2,i,n),j,j
            END IF
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
            xmlCoreStateNumber = 0
            SELECT CASE(coreqn(1,i,n))
               CASE (1)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 1   !(1s1/2)
               CASE (2)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 2   !(2s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 3    !(2p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 4   !(2p3/2)
               CASE (3)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 5   !(3s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 6    !(3p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 7   !(3p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 9    !(3d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 10  !(3d5/2)
               CASE (4)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 8   !(4s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 11   !(4p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 12  !(4p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 14   !(4d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 15  !(4d5/2)
                  IF(coreqn(2,i,n).EQ.3) xmlCoreStateNumber = 19   !(4f5/2)
                  IF(coreqn(2,i,n).EQ.-4) xmlCoreStateNumber = 20  !(4f7/2)
               CASE (5)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 13  !(5s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 16   !(5p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 17  !(5p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 21   !(5d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 22  !(5d5/2)
                  IF(coreqn(2,i,n).EQ.3) xmlCoreStateNumber = 26   !(5f5/2)
                  IF(coreqn(2,i,n).EQ.-4) xmlCoreStateNumber = 27  !(5f7/2)
               CASE (6)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 18  !(6s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 23   !(6p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 24  !(6p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 28   !(6d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 29  !(6d5/2)
               CASE (7)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 25  !(7s1/2)
            END SELECT
            IF(xmlCoreStateNumber.EQ.0) STOP 'Invalid core state!'
455
            xmlElectronStates(xmlCoreStateNumber,n) = coreState_const
456
457
458
459
460
461
462
463
464
465
466
467
468
            xmlPrintCoreStates(xmlCoreStateNumber,n) = 
     +         coreocc(i,n).NE.xmlCoreRefOccs(xmlCoreStateNumber)
            SELECT CASE(xmlCoreStateNumber)
               CASE (9:10,14:15,19:22,26:29)
                  up = MIN((xmlCoreRefOccs(xmlCoreStateNumber)/2),
     +                          coreocc(i,n))
                  dn = MAX(0.0,coreocc(i,n)-up)
               CASE DEFAULT
                  up = CEILING(coreocc(i,n)/2)
                  dn = FLOOR(coreocc(i,n)/2)
            END SELECT
            xmlCoreOccs(1,xmlCoreStateNumber,n) = up
            xmlCoreOccs(2,xmlCoreStateNumber,n) = dn
469
            electronsOnAtom = electronsOnAtom + up + dn
470
471
472
473
474
          ENDDO
          DO i = ncorest+1, natomst
            WRITE(6,'("  valence :",2i3,f6.1,i4,a1)') 
     &             coreqn(1,i,n),coreqn(2,i,n),coreocc(i,n),
     &                      coreqn(1,i,n),lotype(lval(i,n))
475
            nel = nel + coreocc(i,n) * atoms%neq(n)
476
            electronsOnAtom = electronsOnAtom + coreocc(i,n)
477

478
479
480
481
c           In d and f shells a magnetic alignment of the spins
c           is preferred in the valence bands
c           Hence the up and down occupation is chosen such that
c           the total spin is maximized
482

483
484
485
486
487
488
489
490
            IF ( abs(coreqn(2,i,n)+0.5) > 2.499 )
     +      THEN
              IF ( coreocc(i,n) > abs(coreqn(2,i,n)) ) THEN
                up = abs(coreqn(2,i,n))
                dn = coreocc(i,n) - abs(coreqn(2,i,n))
              ELSE
                up = coreocc(i,n)
                dn = 0
491
492
493
             END IF
             upreal=up
             dnreal=dn
494

495
c           in s and p states equal occupation of up and down states
496

497
498
499
500
501
502
503
            ELSE
              j = INT(coreocc(i,n) / 2)
              IF (coreocc(i,n) > 2*j) THEN
                j = - coreocc(i,n)
              ENDIF
              up = j
              dn = j
504
505
              upReal = coreocc(i,n) / 2.0
              dnReal = coreocc(i,n) / 2.0
506
            END IF
507
508
509
510
            IF(.NOT.input%l_inpXML) THEN
               WRITE(27,'(4i3,i4,a1)') coreqn(1,i,n),coreqn(2,i,n),
     &                             up,dn,coreqn(1,i,n),lotype(lval(i,n))
            END IF
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
            xmlCoreStateNumber = 0
            SELECT CASE(coreqn(1,i,n))
               CASE (1)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 1   !(1s1/2)
               CASE (2)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 2   !(2s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 3    !(2p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 4   !(2p3/2)
               CASE (3)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 5   !(3s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 6    !(3p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 7   !(3p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 9    !(3d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 10  !(3d5/2)
               CASE (4)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 8   !(4s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 11   !(4p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 12  !(4p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 14   !(4d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 15  !(4d5/2)
                  IF(coreqn(2,i,n).EQ.3) xmlCoreStateNumber = 19   !(4f5/2)
                  IF(coreqn(2,i,n).EQ.-4) xmlCoreStateNumber = 20  !(4f7/2)
               CASE (5)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 13  !(5s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 16   !(5p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 17  !(5p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 21   !(5d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 22  !(5d5/2)
                  IF(coreqn(2,i,n).EQ.3) xmlCoreStateNumber = 26   !(5f5/2)
                  IF(coreqn(2,i,n).EQ.-4) xmlCoreStateNumber = 27  !(5f7/2)
               CASE (6)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 18  !(6s1/2)
                  IF(coreqn(2,i,n).EQ.1) xmlCoreStateNumber = 23   !(6p1/2)
                  IF(coreqn(2,i,n).EQ.-2) xmlCoreStateNumber = 24  !(6p3/2)
                  IF(coreqn(2,i,n).EQ.2) xmlCoreStateNumber = 28   !(6d3/2)
                  IF(coreqn(2,i,n).EQ.-3) xmlCoreStateNumber = 29  !(6d5/2)
               CASE (7)
                  IF(coreqn(2,i,n).EQ.-1) xmlCoreStateNumber = 25  !(7s1/2)
            END SELECT
            IF(xmlCoreStateNumber.EQ.0) STOP 'Invalid valence state!'
            xmlElectronStates(xmlCoreStateNumber,n) = valenceState_const
            xmlPrintCoreStates(xmlCoreStateNumber,n) = 
     +         coreocc(i,n).NE.xmlCoreRefOccs(xmlCoreStateNumber)
!            SELECT CASE(xmlCoreStateNumber)
!               CASE (9:10,14:15,19:22,26:29)
!                  up = MIN((xmlCoreRefOccs(xmlCoreStateNumber)/2),
!     +                          coreocc(i,n))
!                  dn = MAX(0.0,coreocc(i,n)-up)
!               CASE DEFAULT
!                  up = CEILING(coreocc(i,n)/2)
!                  dn = FLOOR(coreocc(i,n)/2)
!            END SELECT
563
564
            xmlCoreOccs(1,xmlCoreStateNumber,n) = upReal
            xmlCoreOccs(2,xmlCoreStateNumber,n) = dnReal
565
566
          ENDDO
          WRITE (6,*) '----------'
567
568

5392  FORMAT (' atom type: ',i5,' protons: ',f0.8,' electrons: ',f0.8)
569
          IF (ABS(electronsOnAtom-atoms%zatom(n)).GT.1e-13) THEN
570
             WRITE(*,*) 'Note: atom is charged. Is this Intended?'
571
             WRITE(*,5392) n, atoms%zatom(n), electronsOnAtom
572
             WRITE(6,*) 'Note: atom is charged. Is this Intended?'
573
             WRITE(6,5392) n, atoms%zatom(n), electronsOnAtom
574
575
          END IF

576
          IF(.NOT.input%l_inpXML) CLOSE(27)
577
578

          DO i = natomst,1,-1                    ! determine valence states
579
580
            IF (enpara%el0(lval(i,n),n,1) < -9999.8) THEN ! not processed already
              enpara%el0(lval(i,n),n,:) = REAL(coreqn(1,i,n))
581
              IF (i <= ncorest) THEN
582
                enpara%el0(lval(i,n),n,:) = coreqn(1,i,n) + 1.0 ! was already in the core
583
584
585
586
              ENDIF
            ENDIF
          ENDDO
          DO j = 0,3
587
588
            IF (enpara%el0(j,n,1) < -9999.8) THEN
              enpara%el0(j,n,:) = REAL(j+1)
589
590
591
592
593
            ENDIF
          ENDDO

        ELSE  ! determine defauts  as usual

594
595
596
          z = NINT(atoms%zatom(n))
          nlo0 = atoms%nlo(n)
          llo0 = atoms%llo(:,n)
597
          CALL atom_defaults(
598
599
     >                       n,atoms%ntype,atoms%nlod,z,atoms%neq,
     X                       ncst0,nel,atoms%nlo,atoms%llo)
600

601
          IF (atoms%ncst(n) == 0) atoms%ncst(n) = ncst0
602
          IF (lonqn(1,n) /= 0) THEN ! already set before
603
604
605
606
            DO i = 1,atoms%nlo(n)                       ! subtract lo-charge
              nel = nel - 2*(2*atoms%llo(i,n)+1)*atoms%neq(n)
              IF (atoms%llo(i,n) == 0) atoms%ncst(n) = atoms%ncst(n) + 1
              IF (atoms%llo(i,n) >  0) atoms%ncst(n) = atoms%ncst(n) + 2
607
            ENDDO
608
609
610
611
612
613
            atoms%nlo(n) = nlo0                         ! set old values
            atoms%llo(:,n) = llo0 
            DO i = 1,atoms%nlo(n)                       ! add old lo-charge
              nel = nel + 2*(2*atoms%llo(i,n)+1)*atoms%neq(n)   
              IF (atoms%llo(i,n) == 0) atoms%ncst(n) = atoms%ncst(n) - 1
              IF (atoms%llo(i,n) >  0) atoms%ncst(n) = atoms%ncst(n) - 2
614
615
            ENDDO
          ELSE
616
617
             lonqn(1:atoms%nlo(n),n) = 0 !LO check below should not be needed
                                         !for default setting of enparas
618
619
          ENDIF

620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
          IF(juDFT_was_argument("-electronConfig")) THEN
             electronsLeft = NINT(atoms%zatom(n))
             DO i = 1, 29
                electronsLeft = electronsLeft - xmlCoreRefOccs(i)
                IF(electronsLeft.GT.-xmlCoreRefOccs(i)+eps) THEN
                   natomst = i
                   IF(electronsLeft.LT.-eps) THEN
                      xmlPrintCoreStates(i,n) = .TRUE.
                      SELECT CASE(i)
                         CASE (9:10,14:15,19:22,26:29)
                            up = MIN((xmlCoreRefOccs(i)/2),
     +                               -electronsLeft)
                            dn = MAX(0.0,(-electronsLeft)-up)
                         CASE DEFAULT
                            up = CEILING((-electronsLeft)/2)
                            dn = FLOOR((-electronsLeft)/2)
                      END SELECT
                      xmlCoreOccs(1,i,n) = up
                      xmlCoreOccs(2,i,n) = dn
                   END IF
                END IF
             END DO
             xmlElectronStates(1:atoms%ncst(n),n) = coreState_const
             xmlElectronStates(atoms%ncst(n)+1:natomst,n) = 
     +          valenceState_const
          END IF
646
647
648

        ENDIF

649
650
651
        IF (atoms%nlo(n) /= 0) THEN                    ! check for local orbitals
          DO i = 1, atoms%nlo(n)
            enpara%ello0(i,n,:) = REAL(lonqn(i,n))
652
653
654
            IF (lonqn(i,n) == enpara%qn_el(atoms%llo(i,n),n,1)) THEN  ! increase qn
              enpara%qn_el(atoms%llo(i,n),n,:) = 
     &           enpara%qn_el(atoms%llo(i,n),n,1) + 1          ! in LAPW's by 1
655
656
657
            ENDIF
          ENDDO
        ENDIF
658
659
        enpara%skiplo(n,:) = 0
        DO i = 1, atoms%nlo(n)
660
          enpara%qn_ello(i,n,:) = enpara%qn_el(atoms%llo(i,n),n,:) - 1
661
          enpara%skiplo(n,:) = enpara%skiplo(n,1) + (2*atoms%llo(i,n)+1)
662
663
        ENDDO

664
665
666
667
668
669
670
671
672
      ENDDO

      DO n = 1, atoms%ntype
! correct valence charge
         DO i = 1,atoms%nlo(n)
            IF (atoms%llo(i,n).GT.3) THEN
               nel = nel - 2*(2*atoms%llo(i,n)+1)*atoms%neq(n)   
               IF (atoms%llo(i,n) == 0) atoms%ncst(n) = atoms%ncst(n)+1
               IF (atoms%llo(i,n) >  0) atoms%ncst(n) = atoms%ncst(n)+2
673
674
            ELSE IF (enpara%qn_ello(i,n,1).GE.
     &               enpara%qn_el(atoms%llo(i,n),n,1)) THEN
675
676
677
               nel = nel - 2*(2*atoms%llo(i,n)+1)*atoms%neq(n)   
               IF (atoms%llo(i,n) == 0) atoms%ncst(n) = atoms%ncst(n)+1
               IF (atoms%llo(i,n) >  0) atoms%ncst(n) = atoms%ncst(n)+2
678
679
680
681
            END IF
         ENDDO
      ENDDO

682
683
684
      WRITE (6,'("Valence Electrons =",i5)') nel


685
686
687
688
      enpara%llochg = .FALSE.
      enpara%lchange = .FALSE.
      enpara%enmix = 1.0
      enpara%lchg_v = .TRUE.
689
      IF(juDFT_was_argument("-genEnpara")) THEN
690
         CALL enpara%write(atoms,input%jspins,input%film)
691
      END IF
692
      atoms%lmaxd = lmaxdTemp
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
      RETURN

!===> error handling

 911  CONTINUE
      WRITE (errfh,*) 'atom_input: ERROR reading input. ios  =',ios,
     &               ', line =',nline
      CALL juDFT_error("atom_input: ERROR reading input",
     &                calledby="atom_input")

 912  CONTINUE
      WRITE (errfh,*) 'atom_input: ERROR reading namelist.',
     &               ' ios =',ios,
     &               ' line =',nline
      WRITE (errfh,*) buffer(1:nbuffer)
      WRITE (errfh,*) 'The cause of this error may be ...'
      WRITE (errfh,*) '        a variable not defined in this namelist,'
      WRITE (errfh,*) '        wrong type of data for a variable.'
      CALL juDFT_error("atom_input: ERROR reading input",
     &                calledby="atom_input")

 913  CONTINUE
      WRITE (errfh,*) 'atom_input: ERROR reading record.',
     &               ' ios =',ios,
     &               ' line =',nline
      WRITE (errfh,*) buffer(1:nbuffer)
      CALL juDFT_error("atom_input: ERROR reading input",
     &                calledby="atom_input")

!----------------------------------------------------------------
      CONTAINS   ! INTERNAL subroutines
!----------------------------------------------------------------
      SUBROUTINE err( n )

      INTEGER, INTENT (IN) :: n

      WRITE(errfh,*)
      IF (n==1) THEN
        WRITE (errfh,*) 'atom_input: ERROR multiple namelists.',
     &               ' line =',nline
      ELSEIF (n==2) THEN
        WRITE (errfh,*) 'atom_input: ERROR unknown namelist.',
     &               ' line =',nline
      ELSEIF (n==3) THEN
        WRITE (errfh,*) 'atom_input: ERROR reading namelist.',
     &               ' line =',nline
      ELSE
        WRITE (errfh,*) 'atom_input: ERROR reading input.',
     &               ' line =',nline
      ENDIF
      WRITE (errfh,*) buffer(1:nbuffer)
      WRITE (errfh,*)
      fatalerror = .true.
      RETURN
      END SUBROUTINE err
!----------------------------------------------------------------
      END SUBROUTINE atom_input
!----------------------------------------------------------------
!================================================================
      SUBROUTINE read_allatoms(
753
     >                         bfh,l_buffer,
754
755
756
757
758
759
760
761
     X                         rmt,dx,jri,lmax,lnonsph,ncst,econfig,
     <                         bmu,ios)
!****************************************************************
!     reads in defaults for muffin-tin radius, mesh, etc.
!****************************************************************

      IMPLICIT NONE

762
      INTEGER, INTENT (IN)    :: bfh,l_buffer
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
      INTEGER, INTENT (INOUT) :: jri     ! mt radial mesh points
      INTEGER, INTENT (INOUT) :: lmax    ! max. l to include for density, overlap etc.
      INTEGER, INTENT (INOUT) :: lnonsph ! max. l for nonspherical MT-contributions
      INTEGER, INTENT (INOUT) :: ncst    ! # of core levels
      INTEGER, INTENT (OUT)   :: ios
      REAL, INTENT (INOUT)    :: rmt, dx ! muffin-tin radius and log. spacing
      REAL, INTENT (INOUT)    :: bmu     ! magnetic moment
      CHARACTER(len=l_buffer) :: econfig ! verbose electronic config

      NAMELIST /allatoms/ rmt,dx,jri,lmax,lnonsph,ncst,econfig,
     &                    bmu

      READ (bfh,allatoms,err=911,end=911,iostat=ios)

 911  CONTINUE
      END SUBROUTINE read_allatoms
!================================================================
      SUBROUTINE read_atom(
781
     >                     bfh,l_buffer,lotype,
782
     X                     id,z,rmt,jri,dx,lmax,lnonsph,ncst,econfig,
783
     <                     speciesName,bmu,lo,nlod,llod,ios )
784
785
786
787
788
789
790
791
!***********************************************************************
!     reads in muffin-tin radius, mesh, etc.
!***********************************************************************

      USE m_element, ONLY : z_namat
      IMPLICIT NONE

! ... arguments ...
792
793
794
795
796
797
798
799
      INTEGER, INTENT (IN)           :: bfh,l_buffer
      REAL, INTENT (OUT)             :: id,z,rmt,dx,bmu
      INTEGER                        :: lmax,lnonsph,ncst,jri,nlod,llod
      CHARACTER(len=l_buffer)        :: econfig
      CHARACTER(len=80)              :: lo
      CHARACTER(len=20), INTENT(OUT) :: speciesName
      INTEGER, INTENT (OUT)          :: ios
      CHARACTER(len=1), INTENT (IN)  :: lotype(0:3)
800
801
802
803
804

! ... internal variables ...
      INTEGER                  :: i,j,k,l,n
      REAL                     :: zz
      CHARACTER(len=2)         :: element
805
      CHARACTER(len=20)        :: name
806
807
808
809
810
811
      CHARACTER(len=80)        :: lo1

      CHARACTER(len=2) :: lotype2(0:3)
      DATA lotype2 /'sS','pP','dD','fF'/

      NAMELIST /atom/ id,z,rmt,dx,jri,lmax,lnonsph,ncst,
812
     &                econfig,bmu,lo,element,name
813
814
815
816

      id = -9999.9
      z  = -9999.9
      element = ' '
817
818
      speciesName = ''
      name = ''
819
820
821

      READ (bfh,atom,err=911,end=911,iostat=ios)

822
823
      speciesName = TRIM(ADJUSTL(name))

824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
! -> determine which atom we are concerned with ...

      IF ((z < -9999.8).AND.(element.EQ.' ')) THEN
        WRITE (errfh,*)
     &       'ERROR! No element specified  in namelist atom...'
        WRITE (errfh,*) 'use z=.. or element=.. to define it!'
        ios = 3001
        RETURN
      ENDIF
      IF (id < -9999.8) THEN     ! if no id specified
        zz = REAL(z_namat(element))
        IF (z < 0.00) THEN
          IF (zz > eps) THEN
            id = zz              ! use element name
          ELSE
            id = z               ! or use "z" for id
          ENDIF
        ELSE
          IF (zz > 0 .AND. abs(zz-z)>0.5) THEN
            WRITE (warnfh,*)
            WRITE (warnfh,*) 'atom_input: WARNING! ',
     &       'z and z of specified element differ by more than 0.5. '
            WRITE (warnfh,*) '  z = ', z, 'element =',element
            WRITE (warnfh,*)
          ENDIF
        ENDIF
      ENDIF

!---> order local orbitals determine nlod, llod
      lo1 = adjustl(lo)
      lo = ' '
      i = 0
      n = 0
      DO l = 0, 3
        DO
          j = SCAN(lo1,lotype2(l))  ! search for 's' or 'S', 'p' or 'P' etc.
          IF (j > 0) THEN
            lo1(j:j) = ' '
            n = n + 1
            i = i + 1
            IF (j > 1) THEN
              k = SCAN(lo1(j-1:j-1),'123456') ! determine principal quantum number
              IF (k > 0) THEN
                lo(i:i) = lo1(j-1:j-1)
                lo1(j-1:j-1) = ' '
              ELSE
                lo(i:i) = '0'
              ENDIf
            ELSE
              lo(i:i) = '0'
            ENDIF
            i = i + 1
            lo(i:i) = lotype(l)
            nlod = max( nlod, n )
            llod = max( llod, l )
          ELSE
            EXIT
          ENDIF
        ENDDO
      ENDDO
      IF (len_trim(lo1) > 0) then
        WRITE (errfh,*) 'ERROR reading local orbital input...',lo1
        ios = 3002
      ENDIF

 911  CONTINUE
      END SUBROUTINE read_atom
!================================================================

      SUBROUTINE atom_defaults(
     >                         n,ntype,nlod,z,neq,
     X                         ncst2,nel,nlo,llo)
896
      USE m_juDFT
897
898
899
900
901
902
903
904
      IMPLICIT NONE

      INTEGER, INTENT (IN)    :: n,ntype,nlod,z
      INTEGER, INTENT (IN)    :: neq(ntype)
      INTEGER, INTENT (INOUT) :: nel,ncst2
      INTEGER, INTENT (INOUT) :: nlo(ntype),llo(nlod,ntype)
      

905
906
907
908
909
910
911
      INTEGER locore,lo
      INTEGER ncst1(0:103),nce(0:24)

!
! electrons associated with a given number of core-levels
!
      nce=-1;
912
913
      nce(0) = 0  ; nce(1) = 2  ; nce(2)= 4   ; nce(4) = 10  
      nce(5) = 12 ; nce(7) =18 ;nce(8)=20; 
914
915
916
917
      nce(9) = 28 ; nce(12) = 36; nce(14) = 46; nce(17) = 54
      nce(19) = 68; nce(21) = 78; nce(24) = 86


918
919
!
! number of core levels for each element; the INT(ncst1/100) number
920
921
922
923
924
925
926
! provides information about possible local orbitals: 
!    0  no LO
!    1  s-LO
!    2  p-LO
!    4  d-LO
!    8  f-LO
! Sums are allowed, i.e. 3 (s,p)-LO
927
!
928
929
930

!Defaults
      ncst1 =(/0,0,                                            0,        ! Va,H,He
931
932
933
     + 01, 01,                                         1, 1, 1, 1, 1, 1,        ! Li - Ne
     + 304,304,                                        4, 4, 4, 4, 4, 4,       ! Na - Ar
     + 307,307,307,307,307,307,307,307,207,207, 7,  7,409,409,409,409,
934
     +                                                          409, 9,  ! K - Kr
935
     + 312,312,312,312,312,312,312,212,212,212,312,212,414,414,414,414,
936
937
     +                                                         414,414,  ! Rb - Xe
     + 317,317,217,217,217,217,217,217,217,217,217,17, 17,17,17,17,17,    ! Cs - Lu
938
     +    1119,1119,319,319,219,219,219,219,219,421,421,421,421,421,421,  ! Hf - Rn
939
940
     + 324,324,224,224,224,24,24,24,24,24,24,24, 24,24,24,24,24/)   ! Fr - Lw

941
      if (judft_was_argument("-fast_defaults")) 
942
     + ncst1 =(/0,0,                                                0,  ! Va,H,He
943
944
     +     01, 01,                                  1, 1, 1, 1, 1, 1,  ! Li - Ne
     +     04, 04,                                  4, 4, 4, 4, 4, 4,  ! Na - Ar
945
946
947
948
949
950
     +    307,307,207,207, 7, 7, 7, 7, 7, 7, 7, 7,409, 9, 9, 9, 9, 9,  ! K - Kr
     +    312,312,212,212,12,12,12,12,12,12,12,12,414,14,14,14,14,14,  ! Rb - Xe
     +    317,317,217,217,17,17,17,17,17,17,17,17, 17,17,17,17,17,     ! Cs - Lu
     +                219,19,19,19,19,19,19,19,19,421,21,21,21,21,21,  ! Hf - Rn
     +    324,324,224,224,224,24,24,24,24,24,24,24, 24,24,24,24,24/)   ! Fr - Lw

951
952
953
954

!
!--> determine core levels
!
955
956
957
958
959
960
961
962
963
964
      ncst2 = mod(ncst1( z ),100)
      lo=int(ncst1(z)/100)
      
      nel=nel+(z - nce(ncst2))*neq(n)
      nlo(n) = 0 ; locore = 0

      IF (btest(lo,0)) THEN !s-lo
         locore=locore+2
         ncst2=ncst2-1
         nlo(n)=nlo(n)+1
965
         llo(nlo(n),n)=0
966
      ENDIF
967
968
969
970
      IF (btest(lo,1)) THEN !p-lo
         locore=locore+6
         ncst2=ncst2-2
         nlo(n)=nlo(n)+1
971
         llo(nlo(n),n)=1
972
973
974
975
976
      ENDIF
      IF (btest(lo,2)) THEN !d-lo
         locore=locore+10
         ncst2=ncst2-2
         nlo(n)=nlo(n)+1
977
         llo(nlo(n),n)=2
978
979
980
981
982
      ENDIF
      IF (btest(lo,3)) THEN !f-lo
         locore=locore+14
         ncst2=ncst2-2
         nlo(n)=nlo(n)+1
983
         llo(nlo(n),n)=3
984
985
986
987
      ENDIF
      
      nel = nel +  locore  * neq(n)
      
988
989
990
      END SUBROUTINE atom_defaults

      END MODULE m_atominput