wann_read_umatrix.F 9.94 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 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55
      MODULE m_wann_read_umatrix
      use m_juDFT
c************************************************************
c     Read in the transformation matrix between Bloch
c     and Wannier functions which is calculated by 
c     wannier90.
c     Frank Freimuth, October 2006
c************************************************************
      CONTAINS
      SUBROUTINE wann_read_umatrix2(
     >               fullnkpts,num_wann,num_bands,
     >               um_format,jspin,wan90version,
     <               have_disentangled,
     <               lwindow,ndimwin,
     <               u_matrix_opt,u_matrix_tmp,m_matrix)
      implicit none
      integer,intent(in)  :: fullnkpts,jspin
      integer,intent(in)  :: num_wann,num_bands
      logical,intent(in)  :: um_format
      integer,intent(in)  :: wan90version
      logical,intent(out) :: have_disentangled
      logical,intent(out) :: lwindow(num_bands,fullnkpts)
      integer,intent(out) :: ndimwin(fullnkpts)
      complex,intent(out) :: u_matrix_opt(num_bands,num_wann,fullnkpts)
      complex,intent(out) :: u_matrix_tmp(num_wann,num_wann,fullnkpts)
      complex,intent(out),optional :: m_matrix(:,:,:,:)

      integer,allocatable :: exclude_bands(:)
      integer             :: num_kpts,chk_unit,mp_grid(3),l
      integer             :: i,j,k,nkp,ntmp,num_nnmax,nntot
      logical             :: l_chk,l_umdat
      character(len=3)    :: spin12(2)
      character(len=33)   :: header
      real                :: tmp_latt(3,3)
      real,allocatable    :: tmp_kpt_latt(:,:)
      character(len=20)   :: checkpoint
      real                :: tmp_omi
      real                :: omega_invariant
      character(len=2)    :: spinspin12(2)
      data spin12/'WF1' , 'WF2'/
      data spinspin12/'.1','.2'/

c**************************************************************
c     read in chk
c**************************************************************

      write(6,*)"read in chk"
      l_chk=.false.
      inquire (file=spin12(jspin)//'.chk',exist=l_chk)
Daniel Wortmann's avatar
Daniel Wortmann committed
56 57 58
      IF(.NOT.l_chk)  then
         write(*,*)spin12(jspin)//'.chk'
         CALL juDFT_error("file chk not found",calledby
59
     +     ="wann_read_umatrix")
Daniel Wortmann's avatar
Daniel Wortmann committed
60
      ENDIF
61 62 63 64 65

      chk_unit=152
      open(chk_unit,file=spin12(jspin)//'.chk',status='old',
     &     form='unformatted')

Daniel Wortmann's avatar
Daniel Wortmann committed
66
      !write(*,*)'wan90version',wan90version
67 68 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 111 112 113
      if(wan90version.eq.1)then !wannier90 version 1.1
        ! Read comment line
        read(chk_unit) header
        write(6,*)header
        ! Real and reciprocal lattice (units: Angstroem)
        read(chk_unit) ((tmp_latt(i,j),i=1,3),j=1,3)  ! Real lattice
        write(6,*)tmp_latt
        read(chk_unit) ((tmp_latt(i,j),i=1,3),j=1,3)  ! Reciprocal lattice
        write(6,*)tmp_latt
 
        read(chk_unit) num_kpts  ! K-points
        write(6,*)"num_kpts=",num_kpts
        if (num_kpts.ne.fullnkpts) CALL 
     &     juDFT_error("num_kpts.ne.fullnkpts",calledby
     +     ="wann_read_umatrix")
      elseif((wan90version.eq.2).or.(wan90version.eq.3)) then 
        !wannier90 version 1.2 or wannier90 version 2.0

        ! Read comment line
        read(chk_unit) header
        write(6,*)header

        read(chk_unit) ntmp 

        read(chk_unit) ntmp 
        allocate(exclude_bands(ntmp))
        read(chk_unit) (exclude_bands(i),i=1,ntmp) ! Excluded bands

        ! Real and reciprocal lattice (units: Angstroem)
        read(chk_unit) ((tmp_latt(i,j),i=1,3),j=1,3)  ! Real lattice
        write(6,*)tmp_latt
        read(chk_unit) ((tmp_latt(i,j),i=1,3),j=1,3)  ! Reciprocal lattice
        write(6,*)tmp_latt
 
        read(chk_unit) num_kpts                ! K-points
        write(6,*)"num_kpts=",num_kpts
        if (num_kpts.ne.fullnkpts) stop "num_kpts.ne.fullnkpts"
        read(chk_unit) (mp_grid(i),i=1,3)         ! M-P grid
      else
          CALL juDFT_error("unknown wan90version",calledby
     +     ="wann_read_umatrix")
      endif

      allocate(tmp_kpt_latt(3,fullnkpts))  
      read(chk_unit) ((tmp_kpt_latt(i,nkp),i=1,3),nkp=1,fullnkpts)
      deallocate(tmp_kpt_latt)

Daniel Wortmann's avatar
Daniel Wortmann committed
114
      read(chk_unit) nntot                ! nntot
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 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167
      read(chk_unit) ntmp                ! num_wann
      IF (ntmp/=num_wann)  CALL juDFT_error("mismatch in num_wann"
     +     ,calledby ="wann_read_umatrix")

      read(chk_unit) checkpoint             ! checkpoint
      checkpoint=adjustl(trim(checkpoint))
      write(6,*)checkpoint





      read(chk_unit) have_disentangled      
! whether a disentanglement has been performed

      if (have_disentangled) then
         write(6,*)"You used disentangling"
         write(6,*)"Reading in disentangling information"
         read(chk_unit) omega_invariant     ! omega invariant
         write(6,*)omega_invariant


       ! U matrix opt
         read(chk_unit) ((lwindow(i,nkp),i=1,num_bands),nkp=1,num_kpts)
         read(chk_unit) (ndimwin(nkp),nkp=1,num_kpts)
         read(chk_unit)(((u_matrix_opt(i,j,nkp),i=1,num_bands)
     &  ,j=1,num_wann),nkp=1,num_kpts)

      endif



      if(.not.wan90version.eq.0)then

         read(chk_unit) (((u_matrix_tmp(i,j,k),i=1,num_wann),
     &                       j=1,num_wann),k=1,num_kpts)

      endif

      if(present(m_matrix))then
         write(*,*)"nntot=",nntot
         write(*,*)"num_kpts=",num_kpts
         write(*,*)"num_wann=",num_wann
               read(chk_unit) ((((m_matrix(i,j,k,l),i=1,num_wann),
     &                   j=1,num_wann),k=1,nntot),l=1,num_kpts)
      endif   

      close(chk_unit)

      if(wan90version.eq.0)then
c****************************************************************
c                     read in _um.dat (old version of wannier90)
c****************************************************************
Daniel Wortmann's avatar
Daniel Wortmann committed
168 169 170 171 172 173 174 175 176 177 178 179
        l_umdat=.false.
        inquire (file=spin12(jspin)//'_um.dat',exist=l_umdat)
      IF(.NOT.l_umdat)  CALL juDFT_error("where is your um_dat?",
     +     calledby ="wann_read_umatrix")

        open(111,file=spin12(jspin)//'_um.dat',form='unformatted')
        write(6,*)"read in um_data"
        read(111)header
        read(111)tmp_omi
        read(111) ntmp,num_kpts,num_nnmax
      IF(ntmp/=num_wann)  CALL judft_error("mismatch in num_wann",
     +     calledby ="wann_read_umatrix")
180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248

        read(111)(((u_matrix_tmp(i,j,k),i=1,num_wann),j=1,num_wann),
     &    k=1,num_kpts)
        close(111)
      endif !wannier90 version

c$$$      if (um_format) then
c$$$          write(6,*)"write um_data to formatted file WF1.umn/WF2.umn"
c$$$          open(222,
c$$$     &       file=spin12(jspin)//'.umn',form='formatted')
c$$$c          write(222,*)header
c$$$c          write(222,*)tmp_omi
c$$$c          write(222,*)num_wann,num_kpts,num_nnmax
c$$$          write(222,*)"transformation between Bloch and Wannier"
c$$$          write(222,*)num
c$$$          do k=1,num_kpts
c$$$             do i=1,num_wann
c$$$                do j=1,num_wann
c$$$                  write(222,'(i3,3x,i3,3x,i3,3x,f20.16,f20.16)')
c$$$     &                                 i,j,k,u_matrix_tmp(i,j,k)
c$$$                enddo
c$$$             enddo
c$$$          enddo
c$$$          close(222)
c$$$      endif !um_format       
      END SUBROUTINE wann_read_umatrix2

      SUBROUTINE wann_read_umatrix(
     >               fullnkpts,num_wann,num_bands,
     >               um_format,jspin,wan90version,
     <               have_disentangled,
     <               lwindow,ndimwin,u_matrix)
      implicit none
      integer,intent(in)  :: fullnkpts,jspin
      integer,intent(in)  :: num_wann,num_bands
      logical,intent(in)  :: um_format
      integer,intent(in)  :: wan90version
      logical,intent(out) :: have_disentangled
      logical,intent(out) :: lwindow(num_bands,fullnkpts)
      integer,intent(out) :: ndimwin(fullnkpts)
      complex,intent(out) :: u_matrix(num_bands,num_wann,fullnkpts)

      complex             :: u_matrix_opt(num_bands,num_wann,fullnkpts)
      complex             :: u_matrix_tmp(num_wann,num_wann,fullnkpts)
      integer             :: i,j,k
      character(len=3)    :: spin12(2)
      data spin12/'WF1' , 'WF2'/

      call wann_read_umatrix2(
     >         fullnkpts,num_wann,num_bands,
     >         um_format,jspin,wan90version,
     <         have_disentangled,
     <         lwindow,ndimwin,u_matrix_opt,
     <         u_matrix_tmp)
      
      u_matrix(:,:,:)=0.0
      if (have_disentangled) then 
c**************************************************************
c     calculate u_matrix for the entangled case
c*************************************************************
         do k=1,fullnkpts
            do i=1,num_wann
               do j=1,num_wann
           u_matrix(:,j,k)= u_matrix(:,j,k)+
     &   u_matrix_opt(:,i,k)*u_matrix_tmp(i,j,k)
               enddo
            enddo
         enddo
      else   
249
         IF(num_bands/=num_wann) CALL juDFT_error("num_bands/=num_wann"
250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276
     +        ,calledby ="wann_read_umatrix")
        u_matrix(:,:,:)=u_matrix_tmp(:,:,:)
      endif   

      if (um_format) then
          write(6,*)"write um_data to formatted file WF1.umn/WF2.umn"
          open(222,
     &       file=spin12(jspin)//'.umn',form='formatted')
c          write(222,*)header
c          write(222,*)tmp_omi
c          write(222,*)num_wann,num_kpts,num_nnmax
          write(222,*)"transformation between Bloch and Wannier"
          write(222,*)num_bands,fullnkpts,num_wann
          do k=1,fullnkpts
             do i=1,num_wann
                do j=1,num_bands
                  write(222,'(i5,1x,i5,1x,i5,1x,f18.12,1x,f18.12)')
     &                                 j,i,k,u_matrix(j,i,k)
                enddo
             enddo
          enddo
          close(222)
      endif !um_format       


      END SUBROUTINE wann_read_umatrix
      END MODULE m_wann_read_umatrix