Chimie Théorique » Opt'n Path

 SUBROUTINE ReadInput_Vasp

! This routine reads an input template for Vasp

!----------------------------------------------------------------------

!  Copyright 2003-2014 Ecole Normale Supérieure de Lyon, 

!  Centre National de la Recherche Scientifique,

!  Université Claude Bernard Lyon 1. All rights reserved.

!

!  This work is registered with the Agency for the Protection of Programs 

!  as IDDN.FR.001.100009.000.S.P.2014.000.30625

!

!  Authors: P. Fleurat-Lessard, P. Dayal

!  Contact: optnpath@gmail.com

!

! This file is part of "Opt'n Path".

!

!  "Opt'n Path" is free software: you can redistribute it and/or modify

!  it under the terms of the GNU Affero General Public License as

!  published by the Free Software Foundation, either version 3 of the License,

!  or (at your option) any later version.

!

!  "Opt'n Path" is distributed in the hope that it will be useful,

!  but WITHOUT ANY WARRANTY; without even the implied warranty of

!

!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the

!  GNU Affero General Public License for more details.

!

!  You should have received a copy of the GNU Affero General Public License

!  along with "Opt'n Path". If not, see <http://www.gnu.org/licenses/>.

!

! Contact The Office of Technology Licensing, valorisation@ens-lyon.fr,

! for commercial licensing opportunities.

!----------------------------------------------------------------------

  use VarTypes

  use Path_module

  use Io_module

  IMPLICIT NONE

  INTERFACE

     function valid(string) result (isValid)

       CHARACTER(*), intent(in) :: string

       logical                  :: isValid

     END function VALID

  END INTERFACE

  CHARACTER(LCHARS) ::  Line

  INTEGER(KINT) ::  Iat

  INTEGER(KINT) :: I,J

  INTEGER(KINT), ALLOCATABLE :: NbAtType(:) !na

  INTEGER(KINT) :: NbType, NbTypeUser

  REAL(KREAL) :: B(3),xtmp,ytmp,ztmp

  LOGICAL :: TChk

  LOGICAL :: Debug

  Debug=Valid("readinput").OR.Valid("readinput_vasp")

  if (debug) Call Header("Entering ReadInput_Vasp")

  FPBC=.TRUE.

 if (Input/="VASP") THEN

     ! Input was not Vasp, so many parameters are missing like lattice 

     ! constants...

     ! we read them now !

     ALLOCATE(FFF(3,nat))

     ! First geometry is a bit special for VASP as we have to set

     ! many things

     IF (DEBUG) WRITE(*,*) "Reading Vasp Parameters"

     READ(IOIN,'(A)') Vasp_Title

     READ(IOIN,*) Vasp_param

     READ(IOIN,*) Lat_a

     READ(IOIN,*) Lat_b

     READ(IOIN,*) Lat_c

     Lat_a=Lat_a*Vasp_param

     Lat_b=Lat_b*Vasp_param

     Lat_c=Lat_c*Vasp_param

     ALLOCATE(NbAtType(nat))

     READ(IOIN,'(A)') Vasp_types

     ! First, how many different types ?

     NbAtType=0

     READ(Vasp_types,*,END=998) NbAtType

998  CONTINUE

     NbType=0

     DO WHILE (NbAtType(NbType+1).NE.0)

        NbType=NbType+1

     END DO

     ! Do we know the atom types ?

     IF (AtTypes(1).EQ.'  ') THEN

        ! user has not provided atom types... we have to find them ourselves

        ! by looking into the POTCAR file...

        INQUIRE(File="POTCAR", EXIST=Tchk)

        IF (.NOT.Tchk) THEN

           WRITE(*,*) "ERROR! No AtTypes provided, and POTCAR file not found"

           STOP

        END IF

        OPEN(IOTMP,File="POTCAR")

        DO I=1,NbType

           Line='Empty'

           DO WHILE (Line(1:2).NE.'US')

              READ(IOTMP,'(A)') Line

              Line=AdjustL(Line)

           END DO

           Line=adjustl(Line(3:))

           AtTypes(I)=Line(1:2)

        END DO

        if (debug) WRITE(*,'(A,100(1X,A2))') "ReadG:VASP AtTypes",AtTypes(1:NbType)

        CLOSE(IOTMP)

     ELSE  !AtTypes(1).EQ.'  '

        ! user has provided atom types

        NbTypeUser=0

        DO WHILE (AtTypes(NbTypeUser+1).NE.'  ')

           NbTypeUser=NbTypeUser+1

        END DO

        IF (NbType.NE.NbTypeUser) THEN

           WRITE(*,*) "ERROR Read_Geom : NbType in POSCAR do not match AtTypes"

           STOP

        END IF

     END IF

     IAt=1

     DO I=1,NbType

        DO J=1,NbAtType(I)

           AtName(Iat)=AtTypes(I)

           Iat=Iat+1

        END DO

     END DO

     DEALLOCATE(NbAtType)

     NbTypes=NbType

     READ(IOIN,'(A)') Vasp_comment

     READ(IOIN,'(A)') Vasp_direct

     V_direct=Adjustl(Vasp_direct)

     Call UpCase(V_direct)

! PFL 2011 Mar 8 ->

! We have to read the FFF flags :

     DO I=1,Nat

        READ(IOIN,*) Xtmp,ytmp,ztmp,FFF(1:3,I)

        DO J=1,3

           FFF(J,I)=AdjustL(FFF(J,I))

           CALL Upcase(FFF(J,I))

        END DO

     END DO

! <- PFL 2011 Mar 8 

  END IF

  ! In the case of VASP there is always the problem  of moving from one side

  ! of the box to the other...

     Renum=.TRUE.

     ! V_direct has been set in Read_geom

     IF (V_direct(1:6).EQ.'DIRECT') THEN

        Latr(1:3,1)=Lat_a

        Latr(1:3,2)=Lat_b

        Latr(1:3,3)=Lat_c

        B=1.

! TO DO: replace by LAPACK

        CALL Gaussj(Latr,3,3,B,1,1)

     ELSE

        Latr=0.

        Latr(1,1)=1.d0

        Latr(2,2)=1.d0

        Latr(3,3)=1.d0

     END IF

     ! Actualization of Frozen using the FFFF... 

     ! Frozen has the advantage ie if given, it imposes _ALL_ the FFF flags.

     IF (Frozen(1).NE.0) THEN

        WRITE(IOOUT,*) "Frozen is given. Flags of the given POSCAR are overriden"

        FFF='T'

        NFroz=0

        DO WHILE (Frozen(NFroz+1).NE.0)

           NFroz=NFroz+1

           FFF(1:3,Frozen(NFroz))='F'

        END DO

     ELSE

        WRITE(IOOUT,*) "Frozen not given : using  Flags of the given POSCAR"

        NFroz=0

        Frozen=0

        DO I=1,Nat

           IF ((FFF(1,I).EQ.'F').OR.(FFF(2,I).EQ.'F').OR.(FFF(3,I).EQ.'F')) THEN

              FFF(1:3,I)='F'

              NFroz=NFroz+1

              Frozen(NFroz)=I

           END IF

        END DO

        WRITE(IOOUT,*) "Frozen atoms:",Frozen(1:NFroz)

     END IF

     IF (Vmd) THEN

        if (debug) WRITE(*,*) "DBG MAIN, L803, VMD=T,NbTypes,AtTypes",NbTypes,AtTypes(1:NbTypes)

        Line=""

        DO I=1,NbTypes

           Line=TRIM(Line) // " " // TRIM(AdjustL(AtTypes(I))) 

        END DO

        Vasp_Title=Trim(Line) // " " // Trim(adjustL(Vasp_Title))

        if (debug) WRITE(*,*) "DBG MAIN, L809, VMD=T, Vasp_Title=",TRIM(Vasp_Title)

     END IF

 if (debug) Call Header("Exiting ReadInput_Vasp")

   END SUBROUTINE READINPUT_VASP