w_inpXML.f90 27.8 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 10 11 12 13 14 15 16 17 18 19 20 21 22
MODULE m_winpXML

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!
!!!   XML input file generator
!!!
!!!   This subroutine is supposed to write out a file inp.xml
!!!   containing all required input data.
!!!                                         GM'16
!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
CONTAINS
SUBROUTINE w_inpXML(&
&                   atoms,obsolete,vacuum,input,stars,sliceplot,banddos,&
&                   cell,sym,xcpot,noco,jij,oneD,hybrid,kpts,div,l_gamma,&
&                   noel,namex,relcor,a1,a2,a3,scale,dtild_opt,name_opt,&
23
&                   xmlElectronStates,xmlPrintCoreStates,xmlCoreOccs,&
24
&                   atomTypeSpecies,speciesRepAtomType,l_outFile,numSpecies,&
25
&                   enpara)
26 27 28

   USE m_types
   USE m_juDFT_init
29
   USE m_constants
30
   USE m_xmlOutput
31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49

   IMPLICIT NONE

! arguments

   TYPE(t_input),INTENT(IN)   :: input
   TYPE(t_sym),INTENT(IN)     :: sym
   TYPE(t_stars),INTENT(IN)   :: stars 
   TYPE(t_atoms),INTENT(IN)   :: atoms
   TYPE(t_vacuum),INTENT(IN)   :: vacuum
   TYPE(t_obsolete),INTENT(IN) :: obsolete
   TYPE(t_kpts),INTENT(IN)     :: kpts
   TYPE(t_oneD),INTENT(IN)     :: oneD
   TYPE(t_hybrid),INTENT(IN)   :: hybrid
   TYPE(t_Jij),INTENT(IN)      :: Jij
   TYPE(t_cell),INTENT(IN)     :: cell
   TYPE(t_banddos),INTENT(IN)  :: banddos
   TYPE(t_sliceplot),INTENT(IN):: sliceplot
   TYPE(t_xcpot),INTENT(IN)    :: xcpot
50
   TYPE(t_noco),INTENT(IN)     :: noco
51
   TYPE(t_enpara),INTENT(IN)   :: enpara
52
   INTEGER, INTENT (IN)        :: numSpecies
53 54
   INTEGER, INTENT (IN)        :: div(3)
   INTEGER, INTENT (IN)        :: atomTypeSpecies(atoms%ntype)
55 56
   INTEGER, INTENT (IN)        :: speciesRepAtomType(numSpecies)
   LOGICAL, INTENT (IN)        :: l_gamma, l_outFile
57 58
   REAL,    INTENT (IN)        :: a1(3),a2(3),a3(3),scale
   REAL, INTENT (IN)     :: xmlCoreOccs(2,29,atoms%ntype)
59
   INTEGER, INTENT (IN)  :: xmlElectronStates(29,atoms%ntype)
60 61
   LOGICAL, INTENT (IN)  :: xmlPrintCoreStates(29,atoms%ntype)
   CHARACTER(len=3),INTENT(IN) :: noel(atoms%ntypd)
62
   CHARACTER(len=4),INTENT(IN) :: namex
63 64 65 66 67
   CHARACTER(len=12),INTENT(IN):: relcor
   REAL,INTENT(IN),OPTIONAL    :: dtild_opt
   CHARACTER(len=8),INTENT(IN),OPTIONAL:: name_opt(10)


68
   INTEGER :: iSpecies, fileNum
69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110
   CHARACTER(len=8) :: name(10)

!+lda+u
   REAL    u,j
   INTEGER l
   LOGICAL l_amf
   CHARACTER(len=3) ch_test
   NAMELIST /ldaU/ l,u,j,l_amf
!-lda+u
!+odim
   INTEGER MM,vM,m_cyl
   LOGICAL invs1,zrfs1
   INTEGER chi,rot
   LOGICAL d1,band
   NAMELIST /odim/ d1,MM,vM,m_cyl,chi,rot,invs1,zrfs1
!-odim
! ..
! ..  Local Variables
   REAL     ::dtild ,scpos, zc, sumWeight
   INTEGER  ::nw,idsprs, n1, n2
   INTEGER ieq,i,k,na,n,ilo
   REAL s3,ah,a,hs2,rest
   LOGICAL l_hyb,l_sym,ldum
   INTEGER :: ierr
! ..
!...  Local Arrays
   CHARACTER :: helpchar(atoms%ntypd)
   CHARACTER(len=  4) :: chntype
   CHARACTER(len= 41) :: chform
   CHARACTER(len=100) :: line

!     added for HF and hybrid functionals
   REAL                  ::  aMix,omega
   INTEGER               :: idum
   CHARACTER (len=1)     ::  check

   CHARACTER(len=20) :: tempNumberString, speciesName
   CHARACTER(len=150) :: format
   CHARACTER(len=20) :: mixingScheme
   CHARACTER(len=10) :: loType
   CHARACTER(len=10) :: bzIntMode
   CHARACTER(len=200) :: symFilename
111
   LOGICAL :: kptGamma, l_relcor, l_explicit
112 113
   INTEGER :: iAtomType, startCoreStates, endCoreStates
   CHARACTER(len=100) :: xPosString, yPosString, zPosString
114
   CHARACTER(len=200) :: coreStatesString, valenceStatesString
115 116 117 118 119 120 121
   REAL :: tempTaual(3,atoms%nat)
   REAL :: a1Temp(3),a2Temp(3),a3Temp(3)
   REAL :: amatTemp(3,3), bmatTemp(3,3)
   CHARACTER(len=7) :: coreStateList(29) !'(1s1/2)'
   CHARACTER(len=4) :: nobleGasConfigList(6) !'[He]'

   DATA coreStateList / '(1s1/2)','(2s1/2)','(2p1/2)','(2p3/2)','(3s1/2)',&
122
&                       '(3p1/2)','(3p3/2)','(3d3/2)','(3d5/2)','(4s1/2)',&
123 124 125 126 127 128 129 130 131 132 133 134 135 136
&                       '(4p1/2)','(4p3/2)','(5s1/2)','(4d3/2)','(4d5/2)',&
&                       '(5p1/2)','(5p3/2)','(6s1/2)','(4f5/2)','(4f7/2)',&
&                       '(5d3/2)','(5d5/2)','(6p1/2)','(6p3/2)','(7s1/2)',&
&                       '(5f5/2)','(5f7/2)','(6d3/2)','(6d5/2)' /

   DATA nobleGasConfigList / '[He]','[Ne]','[Ar]','[Kr]','[Xe]','[Rn]' /

   IF (PRESENT(dtild_opt)) dtild=dtild_opt
   IF (PRESENT(name_opt)) name=name_opt

   symFilename = 'sym.out'
   kptGamma = l_gamma
   band = .false.
   nw=1
137 138
   IF (TRIM(ADJUSTL(namex)).EQ.'hf'.OR.TRIM(ADJUSTL(namex)).EQ.'exx'.OR.&
       TRIM(ADJUSTL(namex)).EQ.'hse'.OR.TRIM(ADJUSTL(namex)).EQ.'vhse') l_hyb = .true.
139 140 141 142 143 144 145 146 147 148 149 150 151
   l_relcor=.true.
   IF(relcor.EQ.'relativi') THEN
      l_relcor=.true.
   ELSE 
      l_relcor=.false.
   END IF

   DO i = 1, 3
      a1Temp(i) = a1(i)
      a2Temp(i) = a2(i)
      a3Temp(i) = a3(i)
   END DO

152 153 154
   fileNum = -1
   IF(l_outFile) THEN
      fileNum = getXMLOutputUnitNumber()
155
      CALL openXMLElementNoAttributes('inputData')
156 157 158 159
   ELSE
      fileNum = 5
      OPEN (fileNum,file='inp.xml',form='formatted',status='unknown')
      REWIND (fileNum)
160

161 162 163
      WRITE (fileNum,'(a)') '<?xml version="1.0" encoding="UTF-8" standalone="no"?>'
      WRITE (fileNum,'(a)') '<fleurInput fleurInputVersion="0.27">'
   END IF
164 165

   IF(PRESENT(name_opt)) THEN
166 167 168
      WRITE (fileNum,'(a)') '   <comment>'
      WRITE (fileNum,'(a6,10a8)') '      ',name
      WRITE (fileNum,'(a)') '   </comment>'
169 170
   END IF

171
   WRITE (fileNum,'(a)') '   <calculationSetup>'
172 173 174

!      <cutoffs Kmax="3.60000" Gmax="11.000000" GmaxXC="9.200000" numbands="0"/>
   110 FORMAT('      <cutoffs Kmax="',f0.8,'" Gmax="',f0.8,'" GmaxXC="',f0.8,'" numbands="',i0,'"/>')
175
   WRITE (fileNum,110) input%rkmax,stars%gmax,xcpot%gmaxxc,input%gw_neigd
176 177 178 179 180 181 182 183 184 185 186 187 188 189 190

!      <scfLoop itmax="9" maxIterBroyd="99" imix="Anderson" alpha="0.05" spinf="2.00"/>
   120 FORMAT('      <scfLoop itmax="',i0,'" maxIterBroyd="',i0,'" imix="',a,'" alpha="',f0.8,'" spinf="',f0.8,'"/>')
   SELECT CASE (input%imix)
      CASE (1) 
         mixingScheme='straight'
      CASE (3) 
         mixingScheme='Broyden1'
      CASE (5) 
         mixingScheme='Broyden2'
      CASE (7) 
         mixingScheme='Anderson'
      CASE DEFAULT 
         mixingScheme='errorUnknownMixing'
   END SELECT
191
   WRITE (fileNum,120) input%itmax,input%maxiter,TRIM(mixingScheme),input%alpha,input%spinf
192 193 194

!      <coreElectrons ctail="T" frcor="F" kcrel="0"/>
   130 FORMAT('      <coreElectrons ctail="',l1,'" frcor="',l1,'" kcrel="',i0,'"/>')
195
   WRITE (fileNum,130) input%ctail,input%frcor,input%kcrel
196 197 198

!      <magnetism jspins="1" l_noco="F" l_J="F" swsp="F" lflip="F"/>
   140 FORMAT('      <magnetism jspins="',i0,'" l_noco="',l1,'" l_J="',l1,'" swsp="',l1,'" lflip="',l1,'"/>')
199
   WRITE (fileNum,140) input%jspins,noco%l_noco,jij%l_J,input%swsp,input%lflip
200 201 202

!      <soc theta="0.00000" phi="0.00000" l_soc="F" spav="F" off="F" soc66="F"/>
   150 FORMAT('      <soc theta="',f0.8,'" phi="',f0.8,'" l_soc="',l1,'" spav="',l1,'" off="',l1,'" soc66="',l1,'"/>')
203
   WRITE (fileNum,150) noco%theta,noco%phi,noco%l_soc,noco%soc_opt(atoms%ntype+2),noco%soc_opt(atoms%ntype+1),obsolete%eig66(2)
204 205 206 207 208 209 210 211

   IF (noco%l_noco) THEN
      160 FORMAT('      <nocoParams l_ss="',l1,'" l_mperp="',l1,'" l_constr="',l1,'" l_disp="',l1,'" sso_opt="',a3,'" mix_b="',f0.8,'" thetaJ="',f0.8,'" nsh="',i0,'"/>')
      STOP 'Output of Noco input not yet implemented!'
   END IF

   IF (oneD%odd%d1) THEN
      170 FORMAT('      <oneDParams d1="',l1,'" MM="',i0,'" vM="',i0,'" m_cyl="',i0,'" chi="',i0,'" rot="',i0,'" invs1="',l1,'" zrfs1="',l1,'"/>')
212
      WRITE (fileNum,170) oneD%odd%d1,oneD%odd%M,oneD%odd%mb,oneD%odd%m_cyl,oneD%odd%chi,oneD%odd%rot,oneD%odd%invs,oneD%odd%zrfs
213 214 215 216
   END IF

!      <expertModes gw="0" pot8="F" eig66="F" lpr="0" isec1="99" secvar="F" />
   180 FORMAT('      <expertModes gw="',i0,'" pot8="',l1,'" eig66="',l1,'" lpr="',i0,'" isec1="',i0,'" secvar="',l1,'"/>')
217
   WRITE (fileNum,180) input%gw,obsolete%pot8,obsolete%eig66(1),obsolete%lpr,input%isec1,input%secvar
218 219 220

!      <geometryOptimization l_f="F" xa="2.00000" thetad="330.00000" epsdisp="0.00001" epsforce="0.00001"/>
   190 FORMAT('      <geometryOptimization l_f="',l1,'" xa="',f0.8,'" thetad="',f0.8,'" epsdisp="',f0.8,'" epsforce="',f0.8,'"/>')
221
   WRITE (fileNum,190) input%l_f,input%xa,input%thetad,input%epsdisp,input%epsforce
222 223 224 225 226 227 228 229 230 231 232 233 234

   IF(input%gauss.AND.input%tria) THEN
      STOP 'Error: bz integration modes gauss AND tria selected!'
   END IF

   bzIntMode = 'hist'
   IF(input%gauss) THEN
      bzIntMode = 'gauss'
   ELSE IF(input%tria) THEN
      bzIntMode = 'tria'
   END IF
!      <bzIntegration valenceElectrons="8.00000" mode="hist" fermiSmearingEnergy="0.00100">
   200 FORMAT('      <bzIntegration valenceElectrons="',f0.8,'" mode="',a,'" fermiSmearingEnergy="',f0.8,'">')
235
   WRITE (fileNum,200) input%zelec,TRIM(ADJUSTL(bzIntMode)),input%tkb
236

237 238
   l_explicit = juDFT_was_argument("-explicit").OR.l_outFile
   IF(l_explicit) THEN
239
      sumWeight = 0.0
240
      DO i = 1, kpts%nkpt
241 242 243
         sumWeight = sumWeight + kpts%weight(i)
      END DO
      205 FORMAT('         <kPointList posScale="',f0.8,'" weightScale="',f0.8,'" count="',i0,'">')
244
      WRITE (fileNum,205) kpts%posScale, sumWeight, kpts%nkpt
245 246
      DO i = 1, kpts%nkpt
         206 FORMAT('            <kPoint weight="',f12.6,'">',f12.6,' ',f12.6,' ',f12.6,'</kPoint>')
247
         WRITE (fileNum,206) kpts%weight(i), kpts%bk(1,i), kpts%bk(2,i), kpts%bk(3,i)
248
      END DO
249
      WRITE (fileNum,'(a)')('         </kPointList>')
250 251 252
   ELSE IF( (div(1) == 0).OR.(div(2) == 0) ) THEN
!            <kPointCount count="100" gamma="F"/>
      208 FORMAT('         <kPointCount count="',i0,'" gamma="',l1,'"/>')
253
      WRITE (fileNum,208) kpts%nkpt,kptGamma
254 255 256
   ELSE
!            <kPointMesh nx="10" ny="10" nz="10" gamma="F"/>
      210 FORMAT('         <kPointMesh nx="',i0,'" ny="',i0,'" nz="',i0,'" gamma="',l1,'"/>')
257
      WRITE (fileNum,210) div(1),div(2),div(3),kptGamma
258
   END IF
259
   WRITE (fileNum,'(a)') '      </bzIntegration>'
260 261 262

!      <energyParameterLimits ellow="-2.00000" elup="2.00000"/>
   220 FORMAT('      <energyParameterLimits ellow="',f0.8,'" elup="',f0.8,'"/>')
263
   WRITE (fileNum,220) input%ellow,input%elup
264

265 266
   WRITE (fileNum,'(a)') '   </calculationSetup>'
   WRITE (fileNum,'(a)') '   <cell>'
267

268 269
   IF(l_explicit) THEN
      WRITE(fileNum,'(a)') '      <symmetryOperations>'
270
      DO i = 1, sym%nop
271
      WRITE(fileNum,'(a)') '         <symOp>'
272
      224 FORMAT('            <row-1>',i0,' ',i0,' ',i0,' ',f0.15,'</row-1>')
273
      WRITE(fileNum,224) sym%mrot(1,1,i), sym%mrot(1,2,i), sym%mrot(1,3,i), sym%tau(1,i)
274
      225 FORMAT('            <row-2>',i0,' ',i0,' ',i0,' ',f0.15,'</row-2>')
275
      WRITE(fileNum,225) sym%mrot(2,1,i), sym%mrot(2,2,i), sym%mrot(2,3,i), sym%tau(2,i)
276
      226 FORMAT('            <row-3>',i0,' ',i0,' ',i0,' ',f0.15,'</row-3>')
277 278
      WRITE(fileNum,226) sym%mrot(3,1,i), sym%mrot(3,2,i), sym%mrot(3,3,i), sym%tau(3,i)
      WRITE(fileNum,'(a)') '         </symOp>'
279
      END DO
280
      WRITE(fileNum,'(a)') '      </symmetryOperations>'
281 282
   ELSE IF(TRIM(ADJUSTL(sym%namgrp)).EQ.'any') THEN
      228 FORMAT('      <symmetryFile filename="',a,'"/>')
283
      WRITE(fileNum,228) TRIM(ADJUSTL(symFilename))
284 285 286
   ELSE
!      <symmetry spgrp="any" invs="T" zrfs="F"/>
      230 FORMAT('      <symmetry spgrp="',a,'" invs="',l1,'" zrfs="',l1,'"/>')
287
      WRITE (fileNum,230) TRIM(ADJUSTL(sym%namgrp)),sym%invs,sym%zrfs
288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309

   END IF
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!! Note: Different options for the cell definition!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

   IF (cell%latnam.EQ.'c-b') THEN
      a1Temp(1) = sqrt(2.)* a1Temp(1)
   END IF
   IF (cell%latnam.EQ.'hex') THEN
      s3 = sqrt(3.)
      a1Temp(1) = 2*a1Temp(1)/sqrt(3.)
   END IF
   IF (cell%latnam.EQ.'hx3') THEN
      a1Temp(1) = 2*a1Temp(1)
   END IF

   IF (input%film) THEN
!      <xsd:attribute name="dVac" type="xsd:double" use="required"/>
!      <xsd:attribute name="dTilda" type="xsd:double" use="required"/>
!      <filmLattice ...>
      241 FORMAT('      <filmLattice scale="',f0.8,'" latnam="',a,'" dVac="',f0.8,'" dTilda="',f0.8,'">')
310
      WRITE(fileNum,241) scale, TRIM(ADJUSTL(cell%latnam)), vacuum%dvac, dtild
311
      IF (cell%latnam.EQ.'any') THEN
312
         WRITE (fileNum,'(a)') '         <bravaisMatrix>'
313
         255 FORMAT('            <row-1>',f0.12,' ',f0.12,' ',f0.12,'</row-1>')
314
         WRITE (fileNum,255) a1Temp(1),a1Temp(2),a1Temp(3)
315
         265 FORMAT('            <row-2>',f0.12,' ',f0.12,' ',f0.12,'</row-2>')
316
         WRITE (fileNum,265) a2Temp(1),a2Temp(2),a2Temp(3)
317
         275 FORMAT('            <row-3>',f0.12,' ',f0.12,' ',f0.12,'</row-3>')
318 319
         WRITE (fileNum,275) a3Temp(1),a3Temp(2),a3Temp(3)
         WRITE (fileNum,'(a)') '         </bravaisMatrix>'
320 321 322 323 324
      ELSE
         IF ((cell%latnam.EQ.'squ').OR.(cell%latnam.EQ.'hex').OR.&
     &       (cell%latnam.EQ.'c-b').OR.(cell%latnam.EQ.'hx3').OR.&
     &       (cell%latnam.EQ.'c-r').OR.(cell%latnam.EQ.'p-r')) THEN
            256 FORMAT('         <a1>',f0.12,'</a1>')
325
            WRITE (fileNum,256) a1Temp(1)
326 327 328
         END IF
         IF ((cell%latnam.EQ.'c-r').OR.(cell%latnam.EQ.'p-r')) THEN
            266 FORMAT('         <a2>',f0.12,'</a2>')
329
            WRITE (fileNum,266) a2Temp(2)
330 331 332 333
         END IF

         IF (cell%latnam.EQ.'obl') THEN
            257 FORMAT('         <row-1>',f0.12,' ',f0.12,'</row-1>')
334
            WRITE (fileNum,257) a1Temp(1), a1Temp(2)
335
            267 FORMAT('         <row-2>',f0.12,' ',f0.12,'</row-2>')
336
            WRITE (fileNum,267) a2Temp(1), a2Temp(2)
337 338 339
         END IF
      END IF

340
      268 FORMAT('         <vacuumEnergyParameters vacuum="',i0,'" spinUp="',f0.8,'" spinDown="',f0.8,'"/>')
341 342 343 344
      DO i = 1, vacuum%nvac
         WRITE(fileNum,268) i, enpara%evac0(i,1), enpara%evac0(i,input%jspins)
      END DO

345
      WRITE (fileNum,'(a)') '      </filmLattice>'
346 347 348
   ELSE

      242 FORMAT('      <bulkLattice scale="',f0.12,'" latnam="',a,'">')
349
      WRITE (fileNum,242) scale, TRIM(ADJUSTL(cell%latnam))
350 351 352 353

      IF (cell%latnam.EQ.'any') THEN

!         <bravaisMatrix scale="1.0000000">
354
         WRITE (fileNum,'(a)') '         <bravaisMatrix>'
355 356 357

!            <row-1>0.00000 5.13000 5.13000</row-1>
         250 FORMAT('            <row-1>',f0.12,' ',f0.12,' ',f0.12,'</row-1>')
358
         WRITE (fileNum,250) a1Temp(1),a1Temp(2),a1Temp(3)
359 360
!            <row-2>5.13000 0.00000 5.13000</row-2>
         260 FORMAT('            <row-2>',f0.12,' ',f0.12,' ',f0.12,'</row-2>')
361
         WRITE (fileNum,260) a2Temp(1),a2Temp(2),a2Temp(3)
362 363
!            <row-3>5.13000 5.13000 0.00000</row-3>
         270 FORMAT('            <row-3>',f0.12,' ',f0.12,' ',f0.12,'</row-3>')
364
         WRITE (fileNum,270) a3Temp(1),a3Temp(2),a3Temp(3)
365

366
         WRITE (fileNum,'(a)') '         </bravaisMatrix>'
367 368 369 370 371 372
      END IF

      IF ((cell%latnam.EQ.'squ').OR.(cell%latnam.EQ.'hex').OR.&
     &    (cell%latnam.EQ.'c-b').OR.(cell%latnam.EQ.'hx3').OR.&
     &    (cell%latnam.EQ.'c-r').OR.(cell%latnam.EQ.'p-r')) THEN
         252 FORMAT('         <a1>',f0.12,'</a1>')
373
         WRITE (fileNum,252) a1Temp(1)
374 375 376

         IF ((cell%latnam.EQ.'c-r').OR.(cell%latnam.EQ.'p-r')) THEN
            262 FORMAT('         <a2>',f0.12,'</a2>')
377
            WRITE (fileNum,262) a2Temp(2)
378 379 380
         END IF

         272 FORMAT('         <c>',f0.12,'</c>')
381
         WRITE (fileNum,272) dtild
382 383 384 385
      END IF

      IF (cell%latnam.EQ.'obl') THEN
         254 FORMAT('         <row-1>',f0.12,' ',f0.12,'</row-1>')
386
         WRITE (fileNum,254) a1Temp(1), a1Temp(2)
387 388

         264 FORMAT('         <row-2>',f0.12,' ',f0.12,'</row-2>')
389
         WRITE (fileNum,264) a2Temp(1), a2Temp(2)
390 391

         274 FORMAT('         <c>',f0.12,'</c>')
392
         WRITE (fileNum,274) dtild
393 394
      END IF

395
      WRITE (fileNum,'(a)') '      </bulkLattice>'
396
   END IF
397
   WRITE (fileNum,'(a)') '   </cell>'
398 399 400

!   <xcFunctional name="pbe" relativisticCorrections="F">
   280 FORMAT('   <xcFunctional name="',a,'" relativisticCorrections="',l1,'"/>')
401
   WRITE (fileNum,280) TRIM(namex), l_relcor
402 403 404 405

!      <xcParams igrd="1" lwb="F" ndvgrd="6" idsprs="0" chng="-0.100e-11"/>

!   290 FORMAT('      <xcParams igrd="',i0,'" lwb="',l1,'" ndvgrd="',i0,'" idsprs="',i0,'" chng="',e,'"/>')
406 407
!   WRITE (fileNum,290) xcpot%igrd,obsolete%lwb,obsolete%ndvgrd,0,obsolete%chng
!   WRITE (fileNum,'(a)') '   </xcFunctional>'
408

409 410
   WRITE (fileNum,'(a)') '   <atomSpecies>'
   DO iSpecies=1, numSpecies
411 412 413 414 415 416 417
      iAtomType = speciesRepAtomType(iSpecies)
      IF(iAtomType.EQ.-1) THEN
         EXIT
      END IF
!      <species name="Si-1" element="Si" atomicNumber="14" coreStates="4" magMom="0.0" flipSpin="F">
      300 FORMAT('      <species name="',a,'" element="',a,'" atomicNumber="',i0,'" coreStates="',i0,'" magMom="',f0.8,'" flipSpin="',l1,'">')
      tempNumberString = ''
418
      WRITE(tempNumberString,'(i0)') iSpecies
419
      speciesName = TRIM(ADJUSTL(noel(iAtomType))) // '-' // TRIM(ADJUSTL(tempNumberString))
420
      WRITE (fileNum,300) TRIM(ADJUSTL(speciesName)),TRIM(ADJUSTL(noel(iAtomType))),atoms%nz(iAtomType),atoms%ncst(iAtomType),atoms%bmu(iAtomType),atoms%nflip(iAtomType)
421 422 423

!         <mtSphere radius="2.160000" gridPoints="521" logIncrement="0.022000"/>
      310 FORMAT('         <mtSphere radius="',f0.8,'" gridPoints="',i0,'" logIncrement="',f0.8,'"/>')
424
      WRITE (fileNum,310) atoms%rmt(iAtomType),atoms%jri(iAtomType),atoms%dx(iAtomType)
425 426 427

!         <atomicCutoffs lmax="8" lnonsphr="6"/>
      320 FORMAT('         <atomicCutoffs lmax="',i0,'" lnonsphr="',i0,'"/>')
428
      WRITE (fileNum,320) atoms%lmax(iAtomType),atoms%lnonsph(iAtomType)
429

430
      IF (ALL((enpara%el0(0:3,iAtomType,1)-INT(enpara%el0(0:3,iAtomType,1))).LE.0.00000001)) THEN
431
!         <energyParameters s="3" p="3" d="3" f="4"/>
432
         321 FORMAT('         <energyParameters s="',i0,'" p="',i0,'" d="',i0,'" f="',i0,'"/>')
433 434
         WRITE (fileNum,321) INT(enpara%el0(0,iAtomType,1)),INT(enpara%el0(1,iAtomType,1)),&
                             INT(enpara%el0(2,iAtomType,1)),INT(enpara%el0(3,iAtomType,1))
435 436
      END IF

437
      IF(ANY(xmlElectronStates(:,iAtomType).NE.noState_const)) THEN
438 439 440
         endCoreStates = 1
         startCoreStates = 1
         coreStatesString = ''
441
         valenceStatesString = ''
442
         DO i = 1, 29
443
            IF (xmlElectronStates(i,iAtomType).EQ.coreState_const) endCoreStates = i
444 445
         END DO
         IF ((endCoreStates.GE.24).AND.&
446
&            (ALL(xmlPrintCoreStates(1:24,iAtomType).EQV..FALSE.)).AND.&
447
&            (ALL(xmlElectronStates(1:24,iAtomType).EQ.coreState_const)) ) THEN
448 449 450
            coreStatesString = nobleGasConfigList(6)
            startCoreStates = 25
         ELSE IF ((endCoreStates.GE.17).AND.&
451
&                 (ALL(xmlPrintCoreStates(1:17,iAtomType).EQV..FALSE.)).AND.&
452
&                 (ALL(xmlElectronStates(1:17,iAtomType).EQ.coreState_const))) THEN
453 454 455
            coreStatesString = nobleGasConfigList(5)
            startCoreStates = 18
         ELSE IF ((endCoreStates.GE.12).AND.&
456
&                 (ALL(xmlPrintCoreStates(1:12,iAtomType).EQV..FALSE.)).AND.&
457
&                 (ALL(xmlElectronStates(1:12,iAtomType).EQ.coreState_const))) THEN
458 459 460
            coreStatesString = nobleGasConfigList(4)
            startCoreStates = 13
         ELSE IF ((endCoreStates.GE.7).AND.&
461
&                 (ALL(xmlPrintCoreStates(1:7,iAtomType).EQV..FALSE.)).AND.&
462
&                 (ALL(xmlElectronStates(1:7,iAtomType).EQ.coreState_const))) THEN
463 464 465
            coreStatesString = nobleGasConfigList(3)
            startCoreStates = 8
         ELSE IF ((endCoreStates.GE.4).AND.&
466
&                 (ALL(xmlPrintCoreStates(1:4,iAtomType).EQV..FALSE.)).AND.&
467
&                 (ALL(xmlElectronStates(1:4,iAtomType).EQ.coreState_const))) THEN
468 469 470
            coreStatesString = nobleGasConfigList(2)
            startCoreStates = 5
         ELSE IF ((endCoreStates.GE.1).AND.&
471
&                 (ALL(xmlPrintCoreStates(1:1,iAtomType).EQV..FALSE.)).AND.&
472 473
&                 (ALL(xmlElectronStates(1:1,iAtomType).EQ.coreState_const))) THEN
            coreStatesString = nobleGasConfigList(1)
474 475 476
            startCoreStates = 2
         END IF
         DO i = startCoreStates, endCoreStates
477
            IF(xmlElectronStates(i,iAtomType).EQ.coreState_const) THEN
478 479 480
               coreStatesString = TRIM(ADJUSTL(coreStatesString)) // ' ' // coreStateList(i)
            END IF
         END DO
481 482 483 484 485
         DO i = 1, 29
            IF(xmlElectronStates(i,iAtomType).EQ.valenceState_const) THEN
               valenceStatesString = TRIM(ADJUSTL(valenceStatesString)) // ' ' // coreStateList(i)
            END IF
         END DO
486
         WRITE (fileNum,'(a)') '         <electronConfig>'
487
!         <coreConfig>[He] (2s1/2) (2p1/2) (2p3/2)</coreConfig>
488
         322 FORMAT('            <coreConfig>',a,'</coreConfig>')
489
         WRITE(fileNum,322) TRIM(ADJUSTL(coreStatesString))
490
         323 FORMAT('            <valenceConfig>',a,'</valenceConfig>')
491
         WRITE(fileNum,323) TRIM(ADJUSTL(valenceStatesString))
492 493
         DO i = startCoreStates, 29
            IF ((xmlElectronStates(i,iAtomType).NE.noState_const).AND.(xmlPrintCoreStates(i,iAtomType))) THEN
494
!         <coreStateOccupation state="(2s1/2)" spinUp="1.0" spinDown="1.0"/>
495
               325 FORMAT('            <stateOccupation state="',a,'" spinUp="',f0.8,'" spinDown="',f0.8,'"/>')
496
               WRITE(fileNum,325) coreStateList(i), xmlCoreOccs(1,i,iAtomType), xmlCoreOccs(2,i,iAtomType)
497 498
            END IF
         END DO
499
         WRITE (fileNum,'(a)') '         </electronConfig>'
500 501 502 503 504
      END IF

      DO ilo = 1, atoms%nlo(iAtomType)
!         <lo type="HELO" l="0" n="4"/>
         l = atoms%llo(ilo,iAtomType)
505 506 507 508
         n = INT(enpara%ello0(ilo,iAtomType,1))
         loType = 'SCLO'
         IF(n.LT.0) THEN
            loType = 'HELO'
509
         END IF
510 511 512
         n = ABS(n)
         324 FORMAT('         <lo type="',a,'" l="',i0,'" n="',i0,'" eDeriv="',i0,'"/>')
         WRITE (fileNum,324) TRIM(ADJUSTL(loType)), l, n, atoms%ulo_der(ilo,iAtomType)
513 514
      END DO

515
      WRITE (fileNum,'(a)') '      </species>'
516
   END DO
517 518
   WRITE (fileNum,'(a)') '   </atomSpecies>'
   WRITE (fileNum,'(a)') '   <atomGroups>'
519 520 521 522 523 524
   na = 0
   DO iAtomType=1, atoms%ntype
      iSpecies = atomTypeSpecies(iAtomType)
!      <atomGroup species="Si-1">
      330 FORMAT('      <atomGroup species="',a,'">')
      tempNumberString = ''
525
      WRITE(tempNumberString,'(i0)') iSpecies
526
      speciesName = TRIM(ADJUSTL(noel(iAtomType))) // '-' // TRIM(ADJUSTL(tempNumberString))
527
      WRITE (fileNum,330) TRIM(ADJUSTL(speciesName))
528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544 545 546 547

      DO ieq=1,atoms%neq(iAtomType)
         na = na + 1
         tempTaual(1,na) = atoms%taual(1,na)
         tempTaual(2,na) = atoms%taual(2,na)
         tempTaual(3,na) = atoms%taual(3,na)
         DO i = 2,9
            rest = ABS(i*tempTaual(1,na) - NINT(i*tempTaual(1,na)) ) + ABS(i*tempTaual(2,na) - NINT(i*tempTaual(2,na)))
            IF (.not.input%film) THEN
               rest = rest + ABS(i*tempTaual(3,na) - NINT(i*tempTaual(3,na)) )
            END IF
            IF (rest.LT.(i*0.000001)) EXIT
         END DO
         scpos = 1.0
         IF (i.LT.10) scpos = real(i)  ! common factor found (x,y)
         DO i = 1,2
            tempTaual(i,na) = tempTaual(i,na)*scpos
         ENDDO
         IF (.not.input%film) tempTaual(3,na) = tempTaual(3,na)*scpos
         IF (input%film) THEN
548
            tempTaual(3,na) = dtild*tempTaual(3,na)/scale
549 550 551 552 553 554 555 556 557 558 559 560 561 562 563 564 565 566 567 568 569 570 571
         END IF
!+odim in 1D case all the coordinates are given in cartesian YM
         IF (oneD%odd%d1) THEN
            tempTaual(1,na) = tempTaual(1,na)*a1(1)
            tempTaual(2,na) = tempTaual(2,na)*a2(2)
         END IF
!-odim
         IF (oneD%odd%d1) THEN
            STOP '1D position output not implemented!'
         ELSE IF (input%film) THEN
!         <filmPos> x/myConstant  y/myConstant  1/myConstant</filmPos>
            340 FORMAT('         <filmPos>',a,' ',a,' ',a,'</filmPos>')
            xPosString = ''
            yPosString = ''
            zPosString = ''
            IF((scpos.NE.1.0).AND.((tempTaual(1,na).NE.0.0).OR.(tempTaual(2,na).NE.0.0).OR.(tempTaual(3,na).NE.0.0))) THEN
               WRITE(xPosString,'(f0.12,a1,f0.12)') tempTaual(1,na), '/', scpos
               WRITE(yPosString,'(f0.12,a1,f0.12)') tempTaual(2,na), '/', scpos
            ELSE
               WRITE(xPosString,'(f0.12)') tempTaual(1,na)
               WRITE(yPosString,'(f0.12)') tempTaual(2,na)
            END IF
            WRITE(zPosString,'(f0.12)') tempTaual(3,na)
572
            WRITE (fileNum,340) TRIM(ADJUSTL(xPosString)),TRIM(ADJUSTL(yPosString)),TRIM(ADJUSTL(zPosString))
573 574 575 576 577 578 579 580 581 582 583 584 585 586 587
         ELSE
!         <relPos> x/myConstant  y/myConstant  z/myConstant</relPos>
            350 FORMAT('         <relPos>',a,' ',a,' ',a,'</relPos>')
            xPosString = ''
            yPosString = ''
            zPosString = ''
            IF((scpos.NE.1.0).AND.((tempTaual(1,na).NE.0.0).OR.(tempTaual(2,na).NE.0.0).OR.(tempTaual(3,na).NE.0.0))) THEN
               WRITE(xPosString,'(f0.12,a1,f0.12)') tempTaual(1,na), '/', scpos
               WRITE(yPosString,'(f0.12,a1,f0.12)') tempTaual(2,na), '/', scpos
               WRITE(zPosString,'(f0.12,a1,f0.12)') tempTaual(3,na), '/', scpos
            ELSE
               WRITE(xPosString,'(f0.12)') tempTaual(1,na)
               WRITE(yPosString,'(f0.12)') tempTaual(2,na)
               WRITE(zPosString,'(f0.12)') tempTaual(3,na)
            END IF
588
            WRITE (fileNum,350) TRIM(ADJUSTL(xPosString)),TRIM(ADJUSTL(yPosString)),TRIM(ADJUSTL(zPosString))
589 590 591 592
         END IF
      END DO
!         <force calculate="F" relaxX="T" relaxY="T" relaxZ="T"/>
      360 FORMAT('         <force calculate="',l1,'" relaxXYZ="',3l1,'"/>')
593
      WRITE (fileNum,360) atoms%l_geo(iAtomType),atoms%relax(1,iAtomType),atoms%relax(2,iAtomType),atoms%relax(3,iAtomType)
594

595
      WRITE (fileNum,'(a)') '      </atomGroup>'
596
   END DO
597
   WRITE (fileNum,'(a)') '   </atomGroups>'
598 599

   368 FORMAT('   <output dos="',l1,'" band="',l1,'" vacdos="',l1,'" slice="',l1,'">')
600
   WRITE (fileNum,368) banddos%dos,band,banddos%vacdos,sliceplot%slice
601 602 603

!      <checks vchk="F" cdinf="F" disp="F"/>
   370 FORMAT('      <checks vchk="',l1,'" cdinf="',l1,'" disp="',l1,'"/>')
604
   WRITE (fileNum,370) input%vchk,input%cdinf,obsolete%disp
605 606 607

!      <densityOfStates ndir="0" minEnergy="-0.50000" maxEnergy="0.50000" sigma="0.01500"/>  
   380 FORMAT('      <densityOfStates ndir="',i0,'" minEnergy="',f0.8,'" maxEnergy="',f0.8,'" sigma="',f0.8,'"/>')
608
   WRITE (fileNum,380) banddos%ndir,banddos%e2_dos,banddos%e1_dos,banddos%sig_dos
609 610 611

!      <vacuumDOS layers="0" integ="F" star="F" nstars="0" locx1="0.00" locy1="0.00" locx2="0.00" locy2="0.00" nstm="0" tworkf="0.000000"/>
   390 FORMAT('      <vacuumDOS layers="',i0,'" integ="',l1,'" star="',l1,'" nstars="',i0,'" locx1="',f0.8,'" locy1="',f0.8,'" locx2="',f0.8,'" locy2="',f0.8,'" nstm="',i0,'" tworkf="',f0.8,'"/>')
612
   WRITE (fileNum,390) vacuum%layers,input%integ,vacuum%starcoeff,vacuum%nstars,vacuum%locx(1),vacuum%locy(1),vacuum%locx(2),vacuum%locy(2),vacuum%nstm,vacuum%tworkf
613 614 615

!      <plotting iplot="F" score="F" plplot="F"/>
   400 FORMAT('      <plotting iplot="',l1,'" score="',l1,'" plplot="',l1,'"/>')
616
   WRITE (fileNum,400) sliceplot%iplot,input%score,sliceplot%plpot
617 618 619

!      <chargeDensitySlicing numkpt="0" minEigenval="0.000000" maxEigenval="0.000000" nnne="0" pallst="F"/>
   410 FORMAT('      <chargeDensitySlicing numkpt="',i0,'" minEigenval="',f0.8,'" maxEigenval="',f0.8,'" nnne="',i0,'" pallst="',l1,'"/>')
620
   WRITE (fileNum,410) sliceplot%kk,sliceplot%e1s,sliceplot%e2s,sliceplot%nnne,input%pallst
621 622 623

!      <specialOutput form66="F" eonly="F" bmt="F"/>
   420 FORMAT('      <specialOutput form66="',l1,'" eonly="',l1,'" bmt="',l1,'"/>')
624
   WRITE (fileNum,420) obsolete%form66,input%eonly,input%l_bmt
625

626
   WRITE (fileNum,'(a)') '   </output>'
627 628 629
   IF(l_outFile) THEN
      CALL closeXMLElement('inputData')
   ELSE
630 631 632
      WRITE (fileNum,'(a)') '</fleurInput>'
      CLOSE (fileNum)
   END IF
633 634 635

END SUBROUTINE w_inpXML
END MODULE m_winpXML