wann_postproc_setup4.F 6.36 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2017 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
      MODULE m_wann_postproc_setup4
8
#ifdef CPP_WANN4
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 58 59 60 61 62 63 64 65 66 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 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 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163
      CONTAINS
      SUBROUTINE wann_postproc_setup4(
     >           natd,nkpts,kpoints,amat,bmat,
     >           num,num_bands,ntype,neq,
     >           zatom,taual,namat,seed,bkpts_filename)
c******************************************************
c     Call the Wannier90 4d setup routine.
c     Purpose: get the file of nearest-neighbor kpoints
c     J. HANKE
c******************************************************
      IMPLICIT NONE

      integer, intent(in)         :: natd
      integer, intent(in)         :: nkpts
      real,    intent(in)         :: kpoints(4,nkpts)
      real,    intent(in)         :: amat(4,4)
      real,    intent(in)         :: bmat(4,4)
      integer, intent(in)         :: num(4)
      integer, intent(in)         :: num_bands
      integer, intent(in)         :: ntype
      integer, intent(in)         :: neq(ntype)
      real,    intent(in)         :: zatom(ntype)
      real,    intent(in)         :: taual(:,:)
      character(len=2),intent(in) :: namat(0:103)
      character(len=*),intent(in) :: seed    ! QPOINTS
      character(len=*),intent(in) :: bkpts_filename ! QPOINTS

      integer             :: i,j,at
      character(len=50)   :: seedname
      integer             :: num_atoms
      character(len=2)    :: atom_symbols(natd)
      logical             :: gamma_only
      logical             :: spinors
      integer,parameter   :: num_nnmax=14
      integer             :: nntot
      integer             :: nnlist(nkpts,num_nnmax)
      integer             :: nncell(4,nkpts,num_nnmax)
      integer             :: num_bands2
      integer             :: num_wann2
      real                :: proj_site(3,num_bands)
      integer             :: proj_l(num_bands)
      integer             :: proj_m(num_bands)
      integer             :: proj_radial(num_bands)
      real                :: proj_z(3,num_bands)
      real                :: proj_x(3,num_bands)
      real                :: proj_zona(num_bands)
      integer             :: exclude_bands(num_bands)
      integer             :: nn,ikpt
      real                :: pos(4,natd)


      ! Taken from wannier90-1.2/src/wannier_lib.F90
      interface
        subroutine wannier_setup4(seed__name, mp_grid_loc,
     +    num_kpts_loc,
     +    real_lattice_loc, recip_lattice_loc, kpt_latt_loc,
     +    num_bands_tot, num_atoms_loc, atom_symbols_loc,
     +    atoms_cart_loc, gamma_only_loc, spinors_loc, nntot_loc,
     +    nnlist_loc, nncell_loc, num_bands_loc, num_wann_loc,
     +    proj_site_loc, proj_l_loc, proj_m_loc, proj_radial_loc,
     +    proj_z_loc, proj_x_loc, proj_zona_loc, exclude_bands_loc)

         implicit none
         integer, parameter :: dp = selected_real_kind(15,300)
         integer, parameter :: num_nnmax=14
         character(len=*), intent(in) :: seed__name
         integer, dimension(4), intent(in) :: mp_grid_loc
         integer, intent(in) :: num_kpts_loc
         real(kind=dp), dimension(4,4), intent(in) :: real_lattice_loc
         real(kind=dp), dimension(4,4), intent(in) :: recip_lattice_loc
         real(kind=dp), dimension(4,num_kpts_loc), intent(in) ::
     +      kpt_latt_loc
         integer, intent(in) :: num_bands_tot
         integer, intent(in) :: num_atoms_loc
         character(len=*), dimension(num_atoms_loc), intent(in) ::
     +      atom_symbols_loc
         real(kind=dp), dimension(4,num_atoms_loc), intent(in) ::
     +      atoms_cart_loc
         logical, intent(in) :: gamma_only_loc
         logical, intent(in) :: spinors_loc

         integer, intent(out) :: nntot_loc
         integer, dimension(num_kpts_loc,num_nnmax), intent(out) ::
     +      nnlist_loc
         integer,dimension(4,num_kpts_loc,num_nnmax), intent(out) ::
     +      nncell_loc
         integer, intent(out) :: num_bands_loc
         integer, intent(out) :: num_wann_loc
         real(kind=dp), dimension(3,num_bands_tot), intent(out) ::
     +      proj_site_loc
         integer, dimension(num_bands_tot), intent(out) :: proj_l_loc
         integer, dimension(num_bands_tot), intent(out) :: proj_m_loc
         integer, dimension(num_bands_tot), intent(out) ::
     +      proj_radial_loc
         real(kind=dp), dimension(3,num_bands_tot), intent(out) ::
     +      proj_z_loc
         real(kind=dp), dimension(3,num_bands_tot), intent(out) ::
     +      proj_x_loc
         real(kind=dp), dimension(num_bands_tot), intent(out) ::
     +      proj_zona_loc
         integer, dimension(num_bands_tot), intent(out) ::
     +      exclude_bands_loc
        end subroutine wannier_setup4
      end interface

      do j=1,natd
         pos(:,j)=matmul(amat(:,:),taual(:,j))
      enddo

      !seedname='WF1'
      seedname=seed
      gamma_only=.false.
      spinors=.false.

      num_atoms=0
      do i=1,ntype
         at=nint(zatom(i))
         do j=1,neq(i)
            num_atoms=num_atoms+1
            atom_symbols(num_atoms)=namat(at)
         enddo !j
      enddo !i

c**********************************************************
c     Call Wannier90 routine for preparation.
c**********************************************************
      call wannier_setup4(
     >     seedname,num,
     >     nkpts,
     >     transpose(amat),bmat,
     >     kpoints,num_bands, 
     >     num_atoms,atom_symbols,pos, 
     >     gamma_only,spinors, 
     >     nntot,nnlist,nncell,num_bands2,
     >     num_wann2,
     >     proj_site,proj_l,proj_m,
     >     proj_radial,proj_z, 
     >     proj_x,proj_zona,exclude_bands)

c******************************************************
c           write bkpts
c******************************************************
      open(202,file=bkpts_filename,form='formatted')
      write (202,'(i4)') nntot
      do ikpt=1,nkpts
        do nn=1,nntot
          write (202,'(2i6,3x,4i4)') 
c     &     ikpt,bpt(nn,ikpt),(gb(i,nn,ikpt),i=1,3)
     &     ikpt,nnlist(ikpt,nn),(nncell(i,ikpt,nn),i=1,4)
        enddo
      enddo
      close(202)


      END SUBROUTINE wann_postproc_setup4
164
#endif
165
      END MODULE m_wann_postproc_setup4