wann_mmkb_int.F 3.88 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
      module m_wann_mmkb_int
      contains
      subroutine wann_mmkb_int(
10
     >               interchi,addnoco,addnoco2,nvd,k1d,k2d,k3d,
11
     >               n3d,k1,k2,k3,
12
     >               nv,neigd,nbasfcn,zMat,nslibd,
13
     >               k1_b,k2_b,k3_b,
14
     >               nv_b,zMat_b,nslibd_b,
15 16
     >               nbnd,rgphs,ustep,ig,gb,
     <               mmnk)
17 18 19

      USE m_types

20 21
#include "cpp_double.h"
      implicit none
22

23
      TYPE(t_mat), INTENT(IN) :: zMat, zMat_b
24

25 26 27 28 29 30
      integer, intent(in) :: addnoco,addnoco2
      integer, intent(in) :: nvd,n3d,k1(nvd),k2(nvd),k3(nvd)
      integer, intent(in) :: nv,neigd,nbasfcn,nslibd,nslibd_b
      integer, intent(in) :: nv_b,k1_b(nvd),k2_b(nvd),k3_b(nvd)
      integer, intent(in) :: nbnd
      integer, intent(in) :: k1d,k2d,k3d
31
      complex, intent(in) :: rgphs(-k1d:k1d,-k2d:k2d,-k3d:k3d)
32
      complex, intent(in) :: ustep(n3d),interchi
33 34
      integer, intent(in) :: ig(-k1d:k1d,-k2d:k2d,-k3d:k3d)
      integer, intent(in) :: gb(3)
35
      complex, intent(inout) :: mmnk(nbnd,nbnd)
36

37 38 39 40
      complex,allocatable::stepf_c(:,:)
      complex,allocatable::phasusbmat_c(:,:)
      real,allocatable::stepf_r(:,:)
      real,allocatable::phasusbmat_r(:,:)
41 42 43 44
      real,allocatable::mmnk_tmp(:,:)
      integer i,j1,j2,j3,i1,i2,i3,j,in,m,n
      real phase

45 46 47 48 49 50 51 52 53 54 55
      IF(zMat%l_real) THEN
         allocate(mmnk_tmp(nslibd,nslibd_b))
         allocate(phasusbmat_r(nv,nslibd_b))
         allocate(stepf_r(nv_b,nv))
         stepf_r = 0.0
      ELSE
         allocate(phasusbmat_c(nv,nslibd_b))
         allocate(stepf_c(nv_b,nv))
         stepf_c = CMPLX(0.0,0.0)
      END IF

56 57 58 59 60 61 62 63 64 65 66 67
      do i =1,nv
       j1 =-k1(i)  - gb(1)
       j2 =-k2(i)  - gb(2)
       j3 =-k3(i)  - gb(3)
       do j = 1,nv_b
c-->     determine index and phase factor
         i1 =j1 + k1_b(j) 
         i2 =j2 + k2_b(j) 
         i3 =j3 + k3_b(j) 
         in = ig(i1,i2,i3)
         if (in.eq.0) cycle
         phase = rgphs(i1,i2,i3)
68 69 70 71 72
         IF (zMat%l_real) THEN
            stepf_r(j,i) = phase*real(ustep(in))
         ELSE
            stepf_c(j,i) = conjg(phase*ustep(in))
         END IF
73 74
       enddo
      enddo
75 76
      IF(zMat%l_real) THEN
         call CPP_BLAS_sgemm('T','N',nv,nslibd_b,nv_b,real(1.0),
77 78
     &               stepf_r,nv_b,zMat_b%data_r(1+addnoco2,1),nbasfcn,
c     &               stepf_r,nv_b,zMat_b%data_r,nbasfcn,
79 80
     &               real(0.0),phasusbmat_r,nv)
         call CPP_BLAS_sgemm('T','N',nslibd,nslibd_b,nv,real(1.0),
81 82
     &               zMat%data_r(1+addnoco,1),nbasfcn,phasusbmat_r,nv,
c     &               zMat%data_r,nbasfcn,phasusbmat_r,nv,
83 84 85 86 87
     &               real(0.0),mmnk_tmp,nbnd)
         mmnk(1:nslibd,1:nslibd_b)=mmnk(1:nslibd,1:nslibd_b)+
     &               mmnk_tmp(1:nslibd,1:nslibd_b)*interchi
      ELSE
         call CPP_BLAS_cgemm('T','N',nv,nslibd_b,nv_b,cmplx(1.0),
88 89
     &               stepf_c,nv_b,zMat_b%data_c(1+addnoco2,1),
c     &               stepf_c,nv_b,zMat_b%data_c,
90 91 92 93
     &               nbasfcn,cmplx(0.0),
     &               phasusbmat_c,nv)
         phasusbmat_c=conjg(phasusbmat_c)            
         call CPP_BLAS_cgemm('T','N',nslibd,nslibd_b,nv,interchi, 
94 95
     &               zMat%data_c(1+addnoco,1),nbasfcn,phasusbmat_c,nv,
c     &               zMat%data_c,nbasfcn,phasusbmat_c,nv,
96 97
     &               cmplx(1.0),mmnk,nbnd)
      END IF
98

99 100 101 102 103 104
      IF(zMat%l_real) THEN
         deallocate(mmnk_tmp)
         deallocate(phasusbmat_r,stepf_r)
      ELSE
         deallocate(phasusbmat_c,stepf_c)
      END IF
105 106 107

      end subroutine
      end module m_wann_mmkb_int