inpeig.f90 5.54 KB
Newer Older
1 2 3
      MODULE m_inpeig
      CONTAINS
      SUBROUTINE inpeig(&
4
     &                  atoms,cell,input,l_is_oneD,kpts,enpara,kptsFilename,latnam)
5 6 7 8 9 10 11 12 13 14 15 16 17 18
!*********************************************************************
!     inputs the necessary quantities for the eigenvalue part (energy
!     parameters, k-points, wavefunction cutoffs, etc.).
!                  m. weinert   jan. 1987
!     modification dec. 1990:
!     dummyline before reading l-dependent energies to make reading of
!     input easier (e.g. insert name of atom etc.)
!     modification dec. 93:
!     for step-forward diagonalization a la wu in case of more
!     than 1 window we read now
!     number of occupied states for EACH window
!*********************************************************************
      USE m_gkptwgt
      USE m_constants
19 20 21 22 23
      USE m_types_atoms
      USE m_types_cell
      USE m_types_input
      USE m_types_kpts
      USE m_types_enpara
24 25
      USE m_juDFT
     
26 27 28 29 30 31 32
      IMPLICIT NONE
!     ..
      TYPE(t_atoms),INTENT(IN)     :: atoms
      TYPE(t_cell),INTENT(IN)      :: cell
      TYPE(t_input),INTENT(IN)     :: input
      LOGICAL,INTENT(IN)           :: l_is_oneD
      TYPE(t_kpts),INTENT(INOUT)   :: kpts
33
      TYPE(t_enpara),OPTIONAL,INTENT(INOUT) :: enpara
34
      
35
      CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: kptsFilename
36
      CHARACTER(len=*),INTENT(IN)    :: latnam
37 38 39 40

!     ..
!     .. Local Scalars ..
      REAL      :: wt,scale
41
      INTEGER   :: i,j,nk,jsp,n
42
      LOGICAL   :: xyu,l_enpara,l_clf,l_k
43
      CHARACTER(LEN=255) :: fname
44 45 46 47 48 49 50
!     ..
!
     
!---> input energy parameters for each atom.
!---> the energy parameters for l.ge.3 have the same value
!---> read from file 40='enpara'  shz Jan.96
!
51

52 53 54
      IF(PRESENT(enpara)) THEN
         IF (.NOT.input%l_inpXML) THEN
            !read enpara file if present!
Daniel Wortmann's avatar
Daniel Wortmann committed
55
            CALL enpara%init_enpara(atoms,input%jspins,input%film)
56
         END IF
57 58 59 60
      END IF
!
!---> read k-points from file 41='kpts'
!
61 62 63 64 65 66 67

      IF(PRESENT(kptsFilename)) THEN
         fname = TRIM(ADJUSTL(kptsFilename))
      ELSE
         fname = 'kpts'
      END IF

68 69 70
      INQUIRE(file=TRIM(ADJUSTL(fname)),exist=l_k)
      if (.not.l_k) return 

71
      OPEN (41,file=TRIM(ADJUSTL(fname)),form='formatted',status='old')
72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87
!
!---> k-mesh: given in units of the reciprocal lattice basis vectors
!---> scale is a factor to make input easier (default=1.0). k-pt
!---> weights can be relative weights since they are renormalized.
!---> input: for bulk - k1,k2,k3,wtkpt
!--->        for film - k1,k2,wtkpt
!--->           weights are calculated for films, if wtkpt=0
!     for film calculation k1,k2 may also be read in xy - units : xyu=T
!     1 = boundery of BZ on kx/ky axis
!                                                  shz Feb.96
         READ (41,FMT=8110,ERR=911,END=911) kpts%nkpt,scale,xyu
         GOTO 912
  911    CONTINUE
         xyu = .false.
  912    CONTINUE
         
Daniel Wortmann's avatar
Daniel Wortmann committed
88
         IF (kpts%nkpt.GT.kpts%nkpt)  THEN
89 90 91 92 93 94 95 96 97 98 99 100 101
           CALL juDFT_error('nkptd too small',calledby='inpeig')
         ENDIF
 8100    FORMAT (i5,f20.10)
 8110    FORMAT (i5,f20.10,3x,l1)
         IF (scale.EQ.0.0) scale = 1.0
         DO nk = 1,kpts%nkpt
            READ (41,FMT=8040) (kpts%bk(i,nk),i=1,3),kpts%wtkpt(nk)
 8040       FORMAT (4f10.5)
            IF (input%film .AND. .NOT.l_is_oneD) THEN
               kpts%wtkpt(nk) = kpts%bk(3,nk)
               kpts%bk(3,nk) = 0.0
               IF (xyu) THEN
!           transform to cartesian coordinates
102
                  IF (latnam.EQ.'hex') THEN
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 132 133 134 135 136 137 138 139 140 141
                     kpts%bk(1,nk) = kpts%bk(1,nk)*tpi_const/cell%amat(2,2)
                     kpts%bk(2,nk) = kpts%bk(2,nk)*pi_const/cell%amat(1,1)
                  ELSE
                     kpts%bk(1,nk) = kpts%bk(1,nk)*pi_const/cell%amat(1,1)
                     kpts%bk(2,nk) = kpts%bk(2,nk)*pi_const/cell%amat(2,2)
                  END IF
!           transform to internal coordinates
                  kpts%bk(1:2,nk)=matmul(kpts%bk(1:2,nk),cell%amat(1:2,1:2))/tpi_const
               END IF
            ELSEIF (.NOT.input%film) THEN
              IF (xyu) THEN
                call juDFT_warn("The xyu feature is not tested",calledby="inpeig")
                kpts%bk(:,nk) = kpts%bk(:,nk)!*2.0/sc !TODO what is this scaling?
                kpts%bk(:,nk) = matmul( cell%amat, kpts%bk(:,nk) )
              ENDIF
            END IF
            DO  i = 1,3
               kpts%bk(i,nk) = kpts%bk(i,nk)/scale
            ENDDO
!-odim
            IF (l_is_oneD) THEN
!--> trapezoidal
               IF (kpts%bk(3,nk).EQ.0. .OR. kpts%bk(3,nk).EQ.0.5) THEN
                  kpts%wtkpt(nk) = 1.
               ELSE
                  kpts%wtkpt(nk) = 2.
               END IF
            END IF
!-odim
         ENDDO
         wt = sum(kpts%wtkpt(:kpts%nkpt))

         IF (wt.EQ.0.0) THEN
            IF (input%film) THEN
!
!---> generate k-point weights for 2d BZ: squ, rec, cen, hex
!     determine new wt
!
               CALL gkptwgt(&
142
     &                      kpts,cell,latnam)
Daniel Wortmann's avatar
Daniel Wortmann committed
143
               wt=sum(kpts%wtkpt)
144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161
             ELSE
               CALL juDFT_error("wtkpts",calledby ="inpeig",hint&
     &              ="The sum of weights in the kpts file is zero")
            END IF
         END IF
         IF (l_is_oneD)  kpts%bk(1:2,:) = 0.0    
         kpts%wtkpt(:) = kpts%wtkpt(:)/wt
     
         WRITE (6,FMT=8120)  kpts%nkpt
         DO  nk = 1,kpts%nkpt
            WRITE (6,FMT=8040)  kpts%bk(:,nk),kpts%wtkpt(nk)
         ENDDO
 8120    FORMAT (1x,/,' number of k-points for this window =',i5,/,t12,&
     &          'coordinates',t34,'weights')
      CLOSE (41)

      END SUBROUTINE inpeig
      END MODULE m_inpeig