Something went wrong on our end
-
Philipp Rüssmann authoredPhilipp Rüssmann authored
readinput12.f90 35.81 KiB
SUBROUTINE READINPUT(BRAVAIS,LCARTESIAN,RBASIS,ABASIS,BBASIS,CBASIS, &
& DX,DY,DZ, &
& ALATC,BLATC,CLATC, &
& IRNS,NAEZ,NVAC,NEMB,KAOEZ,IRM,ZAT,SITEAT, &
& INS,KSHAPE, &
& LMAX,LMMAX,LPOT, &
& NATYP,NSPIN, &
& NMIN,NRAD,NSMALL, &
& TOLHS,TOLVDIST,TOLAREA, &
& KXC,TXC, &
& I13, &
& NLBASIS,NRBASIS,NLEFT,NRIGHT,ZPERLEFT,ZPERIGHT, &
& TLEFT,TRIGHT,LINTERFACE,RCUTZ,RCUTXY,RMTCORE, &
& LMTREF,RMTREF,SIZEFAC,NFACELIM, EFSET, AOUT_ALL, NPOI)
use mod_version_info, only: serialnr
!#@# KKRtags: VORONOI input-output
implicit none
include 'inc.geometry'
! .. Local Arrays ..
CHARACTER*4 TSPIN(2)
CHARACTER*8 TKWS(3)
CHARACTER*43 TINS(0:3),TKCOR(0:3),TVREL(0:2)
! ..
! .. External Subroutines ..
EXTERNAL RCSTOP
! ..
! .. Intrinsic Functions ..
INTRINSIC MIN
! ..
! .. Array Arguments ..
INTEGER IRNS(*),KFG(4,NATYPD),LMXC(NATYPD)
INTEGER NTCELL(NATYPD),CLS(NATYPD)
INTEGER REFPOT(NATYPD+NEMBD)
INTEGER KAOEZ(*)
INTEGER NBSHIFT ! How many basis atoms have been shifted from polyhedron center
INTEGER SITEAT(NATYPD),NOQ(NAEZD)
REAL*8 ZAT(*),MTFAC(NATYPD),BRAVAIS(3,3),RBASIS(3,*),RMTCORE(*),RMTREF(*)
REAL*8 DX(*),DY(*),DZ(*),RBUNSHIFT(3,NATYPD),RAUX(3) ! Unshifted positions
REAL*8 TRIGHT(3,*),TLEFT(3,*),ZPERLEFT(3),ZPERIGHT(3)
REAL*8 MTWGHT(0:NTOTD)
REAL*8 SIZEFAC(-NLEMBD*NEMBD:NTOTD)
REAL*8 AOUT_ALL(NATYPD)
CHARACTER*124 TXC(3)
CHARACTER*256 UIO
CHARACTER*40 I12,I13,I19,I25,I40
character*80 text
! ..
! .. Scalar Arguments ..
REAL*8 ALAT,E1,E2,ESHIFT,FCM,HFIELD,MIXING,QBOUND,TK, &
& VCONST,ABASIS,BBASIS,CBASIS,RCUTZ,RCUTXY,RMTREFDEF
REAL*8 TOLHS,TOLVDIST,TOLAREA
REAL*8 EFSET ! set Fermi level to this value
INTEGER ICC,ICST,IFILE,IGF,IMIX,INS, &
& IPE,IPF,IPFE,IPOTOU,IPRCOR,ISITE, &
& IRM,IRNUMX,ISHIFT, &
& KPRE,KSCOEF,KSHAPE, &
& KVREL,KWS,KXC,LMAX,LMMAX,LMPOT,LPOT,MD, &
& NATYP,NPNT1,NPNT2,NPNT3,NPOL,NSPIN,INDX,IAT, NPOI
INTEGER NMIN,NSMALL,NRAD,NFACELIM,NBR
INTEGER NSTEPS,KMT,NAEZ,NVAC,NEMB
INTEGER NINEQ,NEMBZ,NZ,CENTEROFINV(3)
REAL*8 ALATC,BLATC,CLATC
INTEGER MMIN,MMAX,SINN,SOUT,RIN,ROUT
INTEGER INTERVX,INTERVY,INTERVZ,NREF,NCLS
INTEGER NLBASIS,NRBASIS,NLEFT,NRIGHT
LOGICAL COMPLX,LINTERFACE,LMTREF,LNEW,LCARTESIAN
! ..
! .. Local Scalars ..
REAL*8 BRYMIX,STRMIX,E3,TX,TY,TZ
INTEGER I,IL,IP,J,IER,I1,IC,II,M2,IQ CHARACTER*43 TSHAPE
!
CHARACTER*8 TESTC(16),OPTC(8)
COMMON /TESTC/TESTC
COMMON /OPTC/OPTC
! ..
! .. Data statements ..
DATA RMTREFDEF/2.0/ ! default rmt-ref
DATA TSPIN/'non-',' '/
DATA TSHAPE/' exact cell treatment (shape correction) '/
DATA TINS/' spherical averaged input potential ', &
& ' non spherical input potential for cluster ', &
& ' non spherical input potential for cluster ', &
& ' non spherical input potential '/
DATA TKWS/' full mt',' ws ',' full ws'/
! ..
!
!------------ array set up and definition of input parameter -----------
!
TXC(1) = ' Morruzi,Janak,Williams #serial: ' // serialnr
TXC(2) = ' von Barth,Hedin #serial: ' // serialnr
TXC(3) = ' Vosko,Wilk,Nusair #serial: ' // serialnr
OPEN(111,FILE='inputcard_generated.txt') ! Write out found or assumed values
IL=1
I=1
KXC = 1
! =============================================================================
! Begin Structure
WRITE(*,*) 'Begin Structure'
! set interface/bulk (2D/3D) mode
LINTERFACE = .FALSE.
CALL IoInput('INTERFACE ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) LINTERFACE
WRITE(111,*) 'INTERFACE= ',LINTERFACE
ENDIF
! Read in the bravais vectors (normalized to alatc)
! Notation: BRAVAIS(J,I) J=x,y,z I=1,2,3
! If the third bravais vector is zero, then surface (2-dimentional) geometry
! is implied.
!
CALL IoInput('BRAVAIS ',UIO,I,7,IER)
IF (IER.NE.0) THEN
WRITE(*,*) 'readinput: error: BRAVAIS not found, stopping.'
STOP 'readinput: error: BRAVAIS not found'
ENDIF
NBR = 3
IF (LINTERFACE) NBR = 2
BRAVAIS(:,:) = 0.D0
DO I=1,NBR
CALL IoInput('BRAVAIS ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) (BRAVAIS(J,I), J=1,3)
ENDDO
IF (BRAVAIS(1,3).EQ.0.D0.AND.BRAVAIS(2,3).EQ.0.D0.AND.BRAVAIS(3,3).EQ.0.D0) THEN
LINTERFACE=.TRUE.
WRITE(*,*) 'Surface geometry'
ENDIF
WRITE(111,FMT='(A10)') 'BRAVAIS '
DO I=1,NBR
WRITE(111,FMT='(3E16.8)') BRAVAIS(1:3,I)
ENDDO
LCARTESIAN= .FALSE. ! defalt is false, then bravais lattice coordinates are used
CALL IoInput('CARTESIAN ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) LCARTESIAN
WRITE(*,*) 'readinput: CARTESIAN=',LCARTESIAN
WRITE(111,*) 'CARTESIAN= ',LCARTESIAN
CALL IoInput('NVAC ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NVAC
WRITE(*,*) 'NVAC=', NVAC
ELSE
NVAC = 0
ENDIF
CALL IoInput('NAEZ ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NAEZ
WRITE(*,*) 'NAEZ=', NAEZ
ELSE
WRITE(*,*) 'readinput12: NAEZ not found, stopping.'
STOP 'readinput12: NAEZ not found, stopping.'
ENDIF
IF(NAEZ.GT.NAEZD) THEN
WRITE(6,*) 'Please, increase the parameter naezd (',naezd,') in inc.p to',naez
STOP 'readinput12: ERROR in NAEZD.'
ENDIF
! Basis atoms
WRITE(111,FMT='(A16)') '<RBASIS> '
DO I=1,NAEZ
CALL IoInput('<RBASIS> ',UIO,I,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) (RBASIS(J,I), J=1,3)
WRITE(111,FMT='(3E24.12)') (RBASIS(J,I), J=1,3)
ELSE
IER=0
CALL IoInput('RBASIS ',UIO,I,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) (RBASIS(J,I), J=1,3)
WRITE(111,FMT='(3E24.12)') (RBASIS(J,I), J=1,3)
ELSE
WRITE(*,*) 'RINPUT13: Keyword <RBASIS> or RBASIS not found. Stopping.'
STOP 'RINPUT13: RBASIS'
ENDIF
ENDIF
ENDDO ! I=1,NAEZ
NFACELIM = NFACED ! Upper limit of allowed number of faces can be smaller than dimension for speedup
CALL IoInput('<NFACELIM> ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NFACELIM
WRITE(111,*) '<NFACELIM>=',NFACELIM
ENDIF
IF (NFACELIM.GT.NFACED) THEN
WRITE(*,*) 'readinput: NFACELIM.GT.NFACED',NFACELIM,NFACED
STOP 'readinput: NFACELIM.GT.NFACED'
ENDIF
RBUNSHIFT(1:3,1:NAEZ) = RBASIS(1:3,1:NAEZ)
DX(1:NAEZ) = 0.D0
DY(1:NAEZ) = 0.D0
DZ(1:NAEZ) = 0.D0
CALL IoInput('NBSHIFT ',UIO,0,7,IER)
IF (IER.EQ.0) THEN
WRITE(*,*) 'Reading unshifted positions...'
READ (UNIT=UIO,FMT=*) NBSHIFT
IF (NBSHIFT.GT.NAEZ) STOP 'NBSHIFT > NAEZ'
WRITE(*,*) '... for',NBSHIFT,' basis atoms.'
DO I = 1,NBSHIFT
CALL IoInput('RBUNSHIFT ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) INDX,(RBUNSHIFT(J,INDX), J=1,3)
ENDDO
! Now put the unshifted positions in array rbasis,
! so that the voronoi cells are centered in these.
! The extra shift is saved in DX,DY,DZ.
DO I = 1,NAEZ
DX(I) = RBASIS(1,I) - RBUNSHIFT(1,I)
DY(I) = RBASIS(2,I) - RBUNSHIFT(2,I)
DZ(I) = RBASIS(3,I) - RBUNSHIFT(3,I)
RBASIS(1:3,I) = RBUNSHIFT(1:3,I)
ENDDO
ENDIF
DO I=1,NAEZ
KAOEZ(I) = I
END DO
ABASIS = 1.D0
BBASIS = 1.D0
CBASIS = 1.D0
CALL IoInput('BASISCALE ',UIO,IL,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) ABASIS,BBASIS,CBASIS
CALL IoInput('ALATBASIS ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) ALATC
ELSE
WRITE(*,*) 'readinput12: ALATBASIS not found, stopping.'
STOP 'readinput12: ALATBASIS not found, stopping.'
ENDIF
BLATC = 1.D0
CLATC = 1.D0
! =============================================================================
! Start Interface mode
!
! --- > Read Left and Right host and set up the embending positions
!
!
NLEFT = 0 ; NRIGHT = 0 ; NEMB = 0 ! in 3d geometry
IF (LINTERFACE) THEN
WRITE(6,9410)
IER = 0
! Check if the keywords exist for old/new treatment of left and right host
CALL IoInput('LEFTBASIS ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
LNEW = .FALSE.
ELSE
LNEW = .TRUE. ! New type input, KAOEZ not needed, keywords <RBLEFT>,<RBRIGHT> used. IER = 0
CALL IoInput('<RBLEFT> ',UIO,1,7,IER)
ENDIF
IF (IER.NE.0) THEN
WRITE(*,*) 'readinput12: LEFTBASIS or <RBLEFT> not found in inputcard'
STOP 'rinput13: LEFTBASIS or <RBLEFT> not found in inputcard'
ENDIF
IER = 0
CALL IoInput('RIGHBASIS ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
LNEW = .FALSE.
ELSE
LNEW = .TRUE.
IER = 0
CALL IoInput('<RBRIGHT> ',UIO,1,7,IER)
ENDIF
IF (IER.NE.0) THEN
WRITE(*,*) 'readinput12: RIGHBASIS or <RBRIGHT> not found in inputcard'
STOP 'readinput12: RIGHBASIS or <RBRIGHT> not found in inputcard'
ENDIF
NRIGHT = 10
CALL IoInput('NRIGHTHO ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NRIGHT
WRITE(111,*) 'NRIGHTHO=',NRIGHT
ELSE
WRITE(111,*) 'Default NRIGHTHO=',NRIGHT
ENDIF
NLEFT = 10
CALL IoInput('NLEFTHOS ',UIO,1,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NLEFT
WRITE(111,*) 'NLEFTHOS=',NLEFT
ELSE
WRITE(111,*) 'Default NLEFTHOS=',NLEFT
ENDIF
IF (NLEFT+NRIGHT.GT.NLEMBD) THEN
WRITE(*,*) 'readinput12: NLEFTHOS+NRIGHTHO=',NLEFT+NRIGHT,' exceeds dimension NLEMBD=',NLEMBD
STOP 'readinput12: NLEFTHOS+NRIGHTHO exceeds dimension NLEMBD'
ENDIF
CALL IoInput('<NLBASIS> ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NLBASIS
ELSE
IER=0
CALL IoInput('NLBASIS ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NLBASIS
ELSE
WRITE(*,*) 'readinput12:<NLBASIS> or NLBASIS not found in inputcard'
STOP 'readinput12: NLBASIS not found in inputcard'
ENDIF
ENDIF
CALL IoInput('<NRBASIS> ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NRBASIS
ELSE
IER=0
CALL IoInput('NRBASIS ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) NRBASIS
ELSE
WRITE(*,*) 'readinput12:<NRBASIS> or NRBASIS not found in inputcard' STOP 'readinput12: NLBASIS not found in inputcard'
ENDIF
ENDIF
! Information is enought to define NEMB
NEMB = NLBASIS + NRBASIS
IF(NEMB.GT.NEMBD) THEN
write(6,*) 'Please, increase the parameter nembd (',nembd,') in inc.p to',nemb
STOP 'ERROR in NEMBD.'
ENDIF
IF (LNEW) THEN
DO I=1,NLBASIS
CALL IoInput('<RBLEFT> ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) (TLEFT(I1,I),I1=1,3)
KAOEZ(NAEZ+I) = I ! Default
CALL IoInput('<KAOEZL> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) KAOEZ(NAEZ+I)
CALL IoInput('<RMTREFL> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTREF(NAEZ+I)
WRITE (111,FMT='(3E20.12,3X,F9.6,3X,I5)') (TLEFT(I1,I),I1=1,3),RMTREF(NAEZ+I),KAOEZ(NAEZ+I)
ENDDO
WRITE(111,FMT='(A82)') '<RBRIGHT> <RMTREFL> <KAOEZL>'
DO I=1,NRBASIS
CALL IoInput('<RBRIGHT> ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) (TRIGHT(I1,I),I1=1,3)
KAOEZ(NAEZ+NLBASIS+I) = I ! Default
CALL IoInput('<KAOEZR> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) KAOEZ(NAEZ+NLBASIS+I)
CALL IoInput('<RMTREFR> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTREF(NAEZ+NLBASIS+I)
WRITE (111,FMT='(3E20.12,3X,F9.6,3X,I5)') (TRIGHT(I1,I),I1=1,3),RMTREF(NAEZ+NLBASIS+I),KAOEZ(NAEZ+NLBASIS+I)
ENDDO
ELSE
DO I=1,NLBASIS
TLEFT(:,I) = 0.D0
CALL IoInput('LEFTBASIS ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) (TLEFT(I1,I),I1=1,3),II
KAOEZ(NAEZ+I) = II ! changed 1.11.99
write(6,*) 'this must be 1',KAOEZ(NAEZ+I),NAEZ+I
RMTREF(NAEZ+I) = RMTREFDEF
CALL IoInput('<RMTREFL> ',UIO,I,7,IER) ! referen-pot rmt
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTREF(NAEZ+I)
END DO
DO I=1,NRBASIS
TRIGHT(:,I) = 0.D0
CALL IoInput('RIGHBASIS ',UIO,I,7,IER)
READ (UNIT=UIO,FMT=*) (TRIGHT(I1,I),I1=1,3),II
KAOEZ(NAEZ+NLBASIS+I) = II ! changed 1.11.99
RMTREF(NAEZ+NLBASIS+I) = RMTREFDEF
CALL IoInput('<RMTREFR> ',UIO,I,7,IER) ! referen-pot rmt
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTREF(NAEZ+NLBASIS+I)
END DO
ENDIF
! Put The additional atoms in the "embending" positions
DO I=1,NLBASIS DO I1=1,3
RBASIS(I1,NAEZ+I) = TLEFT(I1,I)
END DO
END DO
DO I=1,NRBASIS
DO I1=1,3
RBASIS(I1,NAEZ+NLBASIS+I) = TRIGHT(I1,I)
END DO
END DO
!ccccccccccccccccccccccccccccccccccccccccccccccccccccc
! In RBASIS we have first the basis atoms or the interface
! atoms then the left host then the right host the host
! goes in the NEMB positions
!
!ccccccccccccccccccccccccccccccccccccccccccccccccccccc
CALL IoInput('ZPERIODL ',UIO,IL,7,IER)
READ (UNIT=UIO,FMT=*) (ZPERLEFT(i1),I1=1,3)
CALL IoInput('ZPERIODR ',UIO,IL,7,IER)
READ (UNIT=UIO,FMT=*) (ZPERIGHT(i1),I1=1,3)
WRITE(6,9430) NLEFT,NLBASIS
WRITE(6,9440) NRIGHT,NRBASIS
WRITE(6,9450) (ZPERLEFT(i1),I1=1,3)
WRITE(6,9460) (ZPERIGHT(i1),I1=1,3)
WRITE(6,9465)
WRITE(6,9470)
DO I=NLEFT,1,-1
DO I1=NLBASIS,1,-1
tx = TLEFT(1,i1) + (I-1)*ZPERLEFT(1)
ty = TLEFT(2,i1) + (I-1)*ZPERLEFT(2)
tz = TLEFT(3,i1) + (I-1)*ZPERLEFT(3)
WRITE(6,9420) (I-1)*NLBASIS+i1, tx,ty,tz
END DO
END DO
WRITE(6,9475)
DO I=1,NAEZ
WRITE(6,9420) I, (RBASIS(I1,I),I1=1,3)
END DO
WRITE(6,9480)
DO I=1,NRIGHT
DO I1=1,NRBASIS
tx = TRIGHT(1,i1) + (I-1)*ZPERIGHT(1)
ty = TRIGHT(2,i1) + (I-1)*ZPERIGHT(2)
tz = TRIGHT(3,i1) + (I-1)*ZPERIGHT(3)
WRITE(6,9420) (I-1)*NRBASIS+i1,tx,ty,tz
END DO
END DO
END IF ! LINTERFACE
! End Interface mode
! =============================================================================
! =============================================================================
! Begin define volume weights
WRITE(*,*) 'Begin define volume weights'
! returned to main program in array sizefac that is also read in from atominfo below
! sizefac runs over all atoms, incl. left/right embedded atoms in all layers
! up to nleft/nright, and incl. all imp. atoms (set in subroutine that sets the
! impurity details).
! However, mtwght goes not over all left/right embedded atoms in all layers
! but only over all left/right embedded basis atoms.
SIZEFAC(:) = 1.D0 ! Default
! Wished muffin-tin (units of alat)
DO ISITE=1,NAEZ
MTWGHT(ISITE) = DSQRT(SIZEFAC(ISITE)) ! Initialize
CALL IoInput('<MTWAL> ',UIO,ISITE,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) MTWGHT(ISITE) IF (IER.EQ.0) WRITE(*,*) 'Read in wished MT radius in ALAT:',MTWGHT(ISITE)
ENDDO
! Wished muffin-tin (atomic units, overrides units of alat)
DO ISITE=1,NAEZ
CALL IoInput('<MTWAU> ',UIO,ISITE,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) MTWGHT(ISITE)
MTWGHT(ISITE) = MTWGHT(ISITE)/ALATC
WRITE(*,*) 'Read in wished MT radius in AU:',MTWGHT(ISITE)
ENDIF
ENDDO
! Same for left/right embedded atoms
IF (LINTERFACE) THEN
IL = 0
DO ISITE = NAEZ + 1,NAEZ + NLBASIS
MTWGHT(ISITE) = 1.D0 ! Initialize
IL = IL + 1
CALL IoInput('<LFMTWAL> ',UIO,IL,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) MTWGHT(ISITE)
IF (IER.EQ.0) WRITE(*,*) 'Read in wished MT radius in ALAT:',MTWGHT(ISITE)
ENDDO
IL = 0
DO ISITE = NAEZ + NLBASIS + 1,NAEZ + NLBASIS + NRBASIS
MTWGHT(ISITE) = 1.D0 ! Initialize
IL = IL + 1
CALL IoInput('<RTMTWAL> ',UIO,IL,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) MTWGHT(ISITE)
IF (IER.EQ.0) WRITE(*,*) 'Read in wished MT radius in ALAT:',MTWGHT(ISITE)
ENDDO
! Wished muffin-tin (atomic units, overrides units of alat)
IL = 0
DO ISITE = NAEZ + 1,NAEZ + NLBASIS
IL = IL + 1
CALL IoInput('<LFMTWAU> ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) MTWGHT(ISITE)
MTWGHT(ISITE) = MTWGHT(ISITE)/ALATC
WRITE(*,*) 'Read in wished MT radius in AU:',MTWGHT(ISITE)
ENDIF
ENDDO
IL = 0
DO ISITE = NAEZ + NLBASIS + 1,NAEZ + NLBASIS + NRBASIS
IL = IL + 1
CALL IoInput('<RTMTWAU> ',UIO,IL,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) MTWGHT(ISITE)
MTWGHT(ISITE) = MTWGHT(ISITE)/ALATC
WRITE(*,*) 'Read in wished MT radius in AU:',MTWGHT(ISITE)
ENDIF
ENDDO
IL = 1
ENDIF
! Re-define size-factor as muffin-tin squared (units of alat)
DO ISITE=1,NAEZ
SIZEFAC(ISITE) = MTWGHT(ISITE)**2
ENDDO
! Define weights for left/right embedded atoms
IF (LINTERFACE) THEN
ISITE = 0
DO II = 1,NLEFT
DO I1 = NAEZ + 1,NAEZ + NLBASIS
ISITE = ISITE + 1
SIZEFAC(-ISITE) = MTWGHT(I1)**2
ENDDO ENDDO
DO II = 1,NRIGHT
DO I1 = NAEZ + NLBASIS + 1,NAEZ + NLBASIS + NRBASIS
ISITE = ISITE + 1
SIZEFAC(-ISITE) = MTWGHT(I1)**2
ENDDO
ENDDO
ENDIF
WRITE(*,*) 'End Define volume weights'
! End Define volume weights
! =============================================================================
WRITE(*,*) 'End Structure'
! End Structure
! =============================================================================
! =============================================================================
! Start clusters setup
CALL IoInput('RCLUSTZ ',UIO,IL,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RCUTZ
IF (IER.NE.0) THEN
WRITE(*,*) 'readinput12: RCLUSTZ not found, stopping.'
STOP 'readinput12: RCLUSTZ not found, stopping.'
END IF
RCUTXY = RCUTZ
CALL IoInput('RCLUSTXY ',UIO,IL,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RCUTXY
! RMT of ref-pot in TB-KKR to determine TB-KKR clusters
LMTREF = .FALSE.
DO ISITE=1,NAEZ
RMTREF(ISITE) = RMTREFDEF ! Initialize
CALL IoInput('<RMTREF> ',UIO,ISITE,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTREF(ISITE)
IF (IER.EQ.0) LMTREF = .TRUE.
RMTREFDEF = RMTREF(1) ! Re-define default in case it is only included for 1st atom.
ENDDO
WRITE(6,*) 'Parameters Used for the cluster calculation'
if (abs(rcutz-rcutxy).lt.1.d-4) then
write(6,*) 'Clusters inside spheres with radius R = ',rcutz
else
write(6,*) 'Clusters inside cylindels with '
write(6,*) 'Rz = ',rcutz,' Rxy = ',rcutxy
end if
write(6,2018) ! rbasis
write(6,2025) ((rbasis(j,i),j=1,3),i,i=1,naez)
if (nemb.gt.0) write(6,*)
write(6,2031) ((rbasis(j,i),j=1,3),i,refpot(kaoez(i)),i=naez+1,naez+nemb)
! End clusters setup
! =============================================================================
! =============================================================================
! Start chemistry
IL=1
NSPIN = 1
CALL IoInput('NSPIN ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NSPIN
IF (IER.NE.0) WRITE(*,*) 'NSPIN not found, setting NSPIN=1.'
IF (NSPIN.NE.1.AND.NSPIN.NE.2) THEN
WRITE(*,*) 'readinput: NSPIN not 1 or 2',NSPIN
STOP 'readinput: incomprehensible value for NSPIN'
ENDIF
NATYP = NAEZ
CALL IoInput('NATYP ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NATYP
IF (IER.NE.0) WRITE(*,*) 'NATYP not found, setting NATYP=NAEZ.'
IF (NATYP.LT.NAEZ) THEN
WRITE(*,*) 'readinput: NATYP.LT.NAEZ',NATYP,NAEZ
STOP 'readinput: NATYP.LT.NAEZ'
ENDIF
WRITE(6,8000) NAEZ,NATYP,NSPIN
8000 FORMAT('NAEZ=',I5,' ; NATYP=',I5' ; NSPIN=',I2)
RMTCORE(1:NAEZ) = 1.D10 ! A large number to be reduced to touching MT later
ZAT(1:NAEZ) = 29.D0 ! Default is copper.
CALL IoInput('ATOMINFO ',UIO,I+1,7,IER)
IF (IER.EQ.0) THEN
DO I=1,NAEZ
CALL IoInput('ATOMINFO ',UIO,I+1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) ZAT(I), &
& LMXC(I), &
& (KFG(J,I),J=1,4), &
& CLS(I), &
& REFPOT(I), &
& NTCELL(I), &
& MTFAC(I), &
& IRNS(I),RMTCORE(I),SIZEFAC(I)
END DO
ENDIF
! Read ZAT (overrides ATOMINFO read-in)
DO I = 1,NATYP
CALL IoInput('<ZATOM> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) ZAT(I)
ENDDO
! Read RMTCORE (overrides ATOMINFO read-in)
DO I = 1,NAEZ
CALL IoInput('<RMTCORE> ',UIO,I,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) RMTCORE(I)
ENDDO
DO IAT = 1,NATYPD
SITEAT(IAT) = IAT ! For mapping in jellstart if there is no cpa mode
ENDDO
! CPA mode:
IF (NATYP.GT.NAEZ) THEN ! CPA mode, read site for each atom
SITEAT(1:NATYPD) = 0 ! reset the site index per atom, it will be read in in case of cpa.
IER = 0
CALL IoInput('<SITE> ',UIO,1,7,IER)
IF (IER.NE.0) THEN
WRITE(*,*) 'readinput: <SITE> in CPA mode not found, stopping.'
STOP 'readinput: <SITE> in CPA mode not found.'
ENDIF
WRITE(111,FMT='(A18)') '<SITE> '
DO IAT = 1,NATYP
CALL IoInput('<SITE> ',UIO,IAT,7,IER)
READ (UNIT=UIO,FMT=*) SITEAT(IAT)
WRITE(111,FMT='(I5)') SITEAT(IAT)
IF (SITEAT(IAT).GT.NAEZ.OR.SITEAT(IAT).LT.1) THEN
WRITE(*,*) 'readinput: CPA SITE ill-defined',IAT,SITEAT(IAT)
STOP 'readinput: CPA SITE ill-defined'
ENDIF
ENDDO
NOQ(1:NAEZD) = 0
DO IAT = 1,NATYP
ISITE = SITEAT(IAT)
NOQ(ISITE) = NOQ(ISITE) + 1
ENDDO
DO ISITE=1,NAEZ
IF (NOQ(ISITE).LT.1) THEN
WRITE(6,*) 'readinput: CPA: SITE',ISITE,'HAS NO ASSIGNED ATOM'
STOP 'readinput: CPA, atoms not properly assigned.'
ENDIF
END DO
ENDIF
! End CPA mode
! read in whished value of Fermi level (core state energies of starting potential is shifted accordingly)
EFSET = -1.0d0 ! default value -1 signals no shift (EF=0.4...)
CALL IoInput('EFSET ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) EFSET
WRITE(*,*) 'readinput: EFSET=',EFSET
WRITE(111,*) 'EFSET= ',EFSET
WRITE(6,2028) NATYP
WRITE(6,2104)
WRITE(6,1029) ( &
& ZAT(I), &
& LMXC(I), &
& (KFG(J,I),J=1,4), &
& CLS(I), &
& REFPOT(I), &
& NTCELL(I), &
& MTFAC(I), &
& IRNS(I),RMTCORE(I),SIZEFAC(I),I=1,NATYP)
WRITE(6,2108)
WRITE(6,2104)
! Read in `a` factor that defines radial mesh
! r(i) = b*(exp(a*(i-1))-1) where b is determined from a and rmt automatically
! setting to a negative number uses default values (chosen in maindriver)
WRITE(111,FMT='(A16)') '<AFAC_RAD> '
AOUT_ALL(1:NATYPD) = -1.0d0
DO I=1,NATYP
CALL IoInput('<AFAC_RAD> ',UIO,I,7,IER)
IF (IER.EQ.0) THEN
READ (UNIT=UIO,FMT=*) AOUT_ALL(I)
ENDIF
WRITE(111,FMT='(1E24.12)') AOUT_ALL(I)
ENDDO ! I=1,NATYPD
! End chemistry
! =============================================================================
! =============================================================================
! Start control
LMAX = 3
CALL IoInput('LMAX ',UIO,0,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) LMAX
IF (IER.NE.0) WRITE(*,*) 'LMAX not found, setting LMAX=3.'
KSHAPE = 2 ! Default
CALL IoInput('KSHAPE ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) KSHAPE
IRM = 484 ! Default
CALL IoInput('IRM ',UIO,1,7,IER) IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) IRM
INS = 1 ! Default
CALL IoInput('INS ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) INS
CALL IoInput('NMIN ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NMIN
CALL IoInput('NRAD ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NRAD
NRAD = MAX(NMIN,NRAD)
CALL IoInput('NSMALL ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NSMALL
NSMALL = MAX(NMIN,NSMALL)
! read in the total number of shape function points
NPOI = 125 ! Default, see also maindriver.f
CALL IoInput('NPOI ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) NPOI
! Tolerance for voronoi construction, defaults in maindriver data
CALL IoInput('<TOLHS> ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) TOLHS
CALL IoInput('<TOLVD> ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) TOLVDIST
CALL IoInput('<TOLAREA> ',UIO,1,7,IER)
IF (IER.EQ.0) READ (UNIT=UIO,FMT=*) TOLAREA
IF (INS.GT.0) THEN
WRITE (6,FMT=9130)
WRITE (6,FMT=9140)
DO I = 1,NATYP
WRITE (6,FMT=9150) I,IRNS(I),IRNSD
IF (IRNS(I).GT.IRNSD) CALL RCSTOP('19 ')
ENDDO
IF (LMAX.GT.LMAXD) THEN
WRITE (6,FMT=9120)
CALL RCSTOP('20 ')
END IF
END IF
LMMAX = (LMAX+1)**2
LPOT = MIN(2*LMAX,LPOTD)
LMPOT = (LPOT+1)* (LPOT+1)
! End control
! =============================================================================
write(6,2101)
write(6,2012) (kaoez(i),i=1,naez)
!
! ------------------------------------------------------------------------
do i=1,naez
if (kaoez(i).lt.1) STOP 'Error in KAOEZ'
enddo
NCLS = 0
NREF = 0
!
DO I=1,NATYP
NCLS = MAX(NCLS,CLS(I)) ENDDO
DO I=1,NATYP
NREF = MAX(NREF,REFPOT(I))
ENDDO
WRITE(6,2016) NCLS,NREF
! ------------------------------------------------------------------------
! =============================================================================
! Start read run and test options and files
TESTC(1:16) = ' '
CALL IoInput('TESTOPT ',UIO,1,7,IER)
IF(IER.EQ.0) READ(UNIT=UIO,FMT=980)(TESTC(i),i=1,8)
CALL IoInput('TESTOPT ',UIO,2,7,IER)
IF(IER.EQ.0) READ(UNIT=UIO,FMT=980)(TESTC(8+i),i=1,8)
WRITE(6,52) (TESTC(I),I=1,16)
52 FORMAT(79('-')/' TEST OPTIONS:'/2(1X,A8,7('//',A8)/)/79('-'))
CALL IoInput('RUNOPT ',UIO,1,7,IER)
READ (UNIT=UIO,FMT=980)(OPTC(i),i=1,8)
WRITE(6,62) (OPTC(i),i=1,8)
62 FORMAT(79('-')/' EXECUTION OPTIONS:'/1X,A8,7('//',A8)/79('-'))
980 FORMAT(8A8)
IL=1
I12=' '
CALL IoInput('FILES ',UIO,IL,7,IER)
IF(IER.EQ.0) READ (UNIT=UIO,FMT='(A40)') I12
I13=' '
CALL IoInput('FILES ',UIO,IL+1,7,IER)
IF(IER.EQ.0) READ (UNIT=UIO,FMT='(A40)') I13
I40=' '
CALL IoInput('FILES ',UIO,IL+2,7,IER)
IF(IER.EQ.0) READ (UNIT=UIO,FMT='(A40)') I40
I19=' '
CALL IoInput('FILES ',UIO,IL+3,7,IER)
IF(IER.EQ.0) READ (UNIT=UIO,FMT='(A40)') I19
write(6,*) 'I12="',I12,'"'
write(6,*) 'I13="',I13,'"'
write(6,*) 'I19="',I19,'"'
! End read run and test options and files
! =============================================================================
CLOSE(111)
WRITE(6,2110)
WRITE(6,*) ' >>>>>>>>> RINPUT99 EXITS NOW <<<<<<<<<< '
RETURN
! *********************************************Input-End ********
1000 FORMAT(10I4)
1001 FORMAT(3F12.9)
1002 FORMAT(80A1)
1003 FORMAT(10L4)
1004 FORMAT(I4)
1005 FORMAT(4I4)
1006 FORMAT(3F12.9,I4)
1029 FORMAT((F4.0,I4,4x,4I1,3I4,F8.4,I4,E12.4,F8.4))
! ------------------------------------------------------------------------
2001 FORMAT(/,' ATOM',I3/)
2002 FORMAT(/,1X,80(A1)///' INPUT DATA : (NBG=1: LDA; NBG=2: LSDA)')
2004 FORMAT( 79(1H=)/ &
& 'I',77X,'I'/ &
& 'I',4X,' Screened Korringa-Kohn-Rostoker ', &
& 'Electronic Structure Code',10X,'I'/ & & 'I',4X,'for Bulk and Interfaces',20X,'I'/ &
& 'I',77X,'I',/,'I',4X,' ',A8,57X,'I'/ &
& 'I',77X,'I',/,'I',4X,'Version ',A8,57X,'I'/ &
& 'I',77X,'I'/ &
& 'I',77X,'I',/,79(1H=),/,/,1X,80A1 )
2010 FORMAT(' NSPIN '/I4)
2011 FORMAT(' NSTEPS'/I4)
2012 FORMAT(' KAOEZ '/,(10I4))
2013 FORMAT(' M2 MMIN MMAX SINN', &
& ' SOUT RIN ROUT'/7I8)
2014 FORMAT(' ALATC BLATC CLATC'/3F15.8)
2015 FORMAT(' INTERVX INTERVY INTERVZ'/3I10)
2016 FORMAT(' NCLS NREF '/,3I8)
2017 FORMAT(' COMPLX'/(3X,L1))
2018 FORMAT(' RBASIS')
2019 FORMAT(' ABASIS BBASIS CBASIS'/3F15.8)
2020 FORMAT(' NPRINCD NLAYER'/,2I8)
2021 FORMAT(' INIPOL'/,(10I4))
2022 FORMAT(' IXIPOL'/,(10I4))
2023 FORMAT(' NAEZ NEMB NEMBZ'/,3I8)
2024 FORMAT(' NZ '/,I4)
2025 FORMAT((3F15.8,I6))
2026 FORMAT(' COFINV'/,3F15.8)
2027 FORMAT(' LATT '/,I4)
2028 FORMAT &
&(' NATYP '/,I4/,' Z lmx KFG cls pot ntc MTFAC irns MT Size')
2030 FORMAT(' EQINV '/,(10I4))
2031 FORMAT((3F15.8,2I6))
2040 FORMAT(' IFULLD ISPARD ISLABD'/,4I8)
2050 FORMAT(' Dimension and Input Data CHECK')
! ------------------------------------------------------------------------
2100 FORMAT(79(1H-))
2101 format( 3(1H-),1H+ , 3(14(1H-),1H+), 30(1H-))
2102 format( 3(9(1H-),1H+) ,49(1H-))
2103 FORMAT(10(3(1H-),1H+) ,39(1H-))
2104 format( 3(1H-),1H+ ,75(1H-))
2105 format( 3(3(1H-),1H+) ,67(1H-))
2106 format( 5(3(1H-),1H+) ,59(1H-))
2107 format( 3(14(1H-),1H+),34(1H-))
2108 format( 2(3(1H-),1H+), 7(1H-),1H+, 3(3(1H-),1H+), &
& 7(1H-),1H+, 3(1H-),1H+, 39(1H-))
2109 format( 5(7(1H-),1H+) ,39(1H-))
2110 format( 3(7(1H-),1H+) ,55(1H-))
2111 format( 7(7(1H-),1H+) ,23(1H-))
2112 format( 2(7(1H-),1H+) ,63(1H-))
9000 FORMAT (I2,3X,4I5)
9010 FORMAT (1X,I1,1X,4I1)
9020 FORMAT (/,33x,'check of dimension-data consistency',/,33x, &
& 35 ('-'),/,40x,'lmax : (',i6,',',i6,')',/,40x, &
& 'natyp : (',i6,',',i6,')',/,40x,'irm : (',i6,',',i6, &
& ')',/,40x,'nspin : (',i6,',',i6,')',/)
9030 FORMAT (1x,10 ('*'),' external magnetic field applied hfield=', f8.5)
9040 FORMAT (3f12.7)
9050 FORMAT (20x,a4,'spin polarized calculation')
9060 FORMAT (8i4)
9070 FORMAT (1x,20x,' calculation with',a8,'-potential')
9080 FORMAT (1x,79 ('*'))
9090 FORMAT (' mixing factor used :',f15.6,/, &
& ' convergence quality required :',1p,d15.2)
9100 FORMAT (1x,20x,a24,'exchange-correlation potential')
9110 FORMAT (/,20x,'broyden"s method # :',i3, &
& ' is used up to iteration- ',/,20x,'depth :',i3, &
& ' then jacobian is fixed and potential ',/,20x, &
& 'is updated using that jacobian')
9120 FORMAT (13x,' in case of calculating non - spherical wavefcts ', &
& 'the parameter lmaxd has to be less-equal lmax !')
9130 FORMAT (/)
9140 FORMAT (20x,'full potential calculation ', &
& '- cut off of non spherical potential',/,' >',/)
9150 FORMAT (31x,'representive atom no.',i3,' irns :',i5,' irnsd :',i5) 9160 FORMAT (21x,a43)
9170 FORMAT (21x,a43)
9180 FORMAT (2i5)
9190 FORMAT (3f12.7,/,4i4)
9200 FORMAT (3i4,1f12.7)
9210 FORMAT (' lmax'/,i4)
9220 FORMAT (' E1 E2 TK'/,3f12.6)
9230 FORMAT (' NPOL NPNT1 NPNT2 NPNT3'/,4i7)
9240 FORMAT (' IRNUMX ITCCOR IPRCOR'/,3i7)
9250 FORMAT (' IFILE IPE ISHIFT ESHIFT'/,3i7,f12.6)
9260 FORMAT (' KSHAPE IRM INS '/,3i7)
9280 FORMAT (' external magnetic hfield :',f15.4/, &
& ' VCONST :',f15.6)
9290 FORMAT (' IMIX IPOTOU IGF ICC'/,4i7)
9300 FORMAT (' ITDBRY'/,i7)
9310 FORMAT (' STRMIX FCM QBOUND'/,3f12.6)
9320 FORMAT (' BRYMIX'/,f12.6)
9330 FORMAT (' KTE KPRE KEFG KVMAD KSCOEF'/,5i7)
9301 format( 3(1H-),1H+ ,75(1H-))
9302 format( 3(11(1H-),1H+),43(1H-))
9303 format(3(6(1H-),1H+) ,58(1H-))
9304 format(4(6(1H-),1H+) ,51(1H-))
9305 format(3(6(1H-),1H+),11(1H-),1H+ ,46(1H-))
9306 format(6(6(1H-),1H+) ,37(1H-))
9307 format(6(1H-),1H+,72(1H-))
9308 format(11(1H-),1H+,67(1H-))
9309 format(5(6(1H-),1H+) ,44(1H-))
9410 format('*** SLAB - INTERFACE CALCULATION ***'/)
9420 format(I5,3F12.6)
9430 format('Number of LEFT Host Layers : ',I5,' with ',I5,' basis')
9440 format('Number of RIGHT Host Layers : ',I5,' with ',I5,' basis')
9450 format('Left side periodicity : ',3F10.5)
9460 format('Right side periodicity : ',3F10.5)
9465 format(' Geommetry used : '/, &
& ' ATOM TX TY TZ ')
9470 format('--------------- Left Host -------------- ')
9475 format('--------------- S L A B -------------- ')
9480 format('--------------- Right Host -------------- ')
END SUBROUTINE READINPUT