wann_get_mp.f 2.01 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 56 57
      module m_wann_get_mp
      use m_juDFT
      contains
      subroutine wann_get_mp(
     >               nkpts,kpoints,
     <               num)
c**************************************
c     Determine the structure of the
c     Monkhorst-Pack mesh.
c     
c     Frank Freimuth
c**************************************
      implicit none
      integer,intent(in)  :: nkpts
      real,intent(in)     :: kpoints(:,:)
      integer,intent(out) :: num(:)
      
      integer :: dim,iter
      real    :: maxi,mini,increm,compare

      IF(SIZE(kpoints,1)/=3)      CALL juDFT_error("wann_get_mp: 1"
     +     ,calledby ="wann_get_mp")
      IF(SIZE(kpoints,2)/=nkpts)  CALL juDFT_error("wann_get_mp: 2"
     +     ,calledby ="wann_get_mp")
      IF(SIZE(num,1)/=3)          CALL juDFT_error("wann_get_mp: 3"
     +     ,calledby ="wann_get_mp")

      do dim=1,3
         maxi=maxval(kpoints(dim,:))
         mini=minval(kpoints(dim,:))
         if(mini==maxi)then
            num(dim)=1
         else   
            increm=maxi-mini
            do iter=1,nkpts
               compare=maxi-kpoints(dim,iter)
               if(abs(compare).lt.1e-6)cycle
               if(compare.lt.increm) then
                  increm=compare
               endif   
            enddo
            num(dim)=(maxi-mini)/increm+1.01
         endif   
      enddo
      write(6,*)"wann_get_mp: determination of mp-grid parameters:"
      write(6,*)"mp_1=",num(1),"mp_2=",num(2),"mp_3=",num(3)
      IF(num(1)*num(2)*num(3)/=nkpts)  CALL juDFT_error
     +     ("mysterious kpoints",calledby ="wann_get_mp")

      end subroutine wann_get_mp
      end module m_wann_get_mp