Commit 3032c4e5 authored by Gregor Michalicek's avatar Gregor Michalicek

And more and more fixes to integrate the Wannier code

parent 5be524db
This diff is collapsed.
This diff is collapsed.
......@@ -22,14 +22,18 @@ c*****************************************c
> n3d,nv_b,nv_b2,nbnd,neigd,
> nslibd_b,nslibd_b2,nbasfcn,
> addnoco,addnoco2,
> k1_b ,k2_b ,k3_b, gb,
> k1_b ,k2_b ,k3_b, gb,
> k1_b2,k2_b2,k3_b2,gb2,
> bkpt,bbmat,vpw,z_b,z_b2,
> bkpt,bbmat,vpw,zMat_b,zMat_b2,
> rgphs,ustep,ig,l_kin,sign,uHu)
#include "cpp_double.h"
USE m_types
implicit none
TYPE(t_zmat), INTENT(IN) :: zMat_b, zMat_b2
c ..arguments..
logical, intent(in) :: l_kin
integer, intent(in) :: gb(3),gb2(3)
......@@ -44,36 +48,35 @@ c ..arguments..
real, intent(in) :: bkpt(3),bbmat(3,3)
complex, intent(in) :: ustep(n3d),vpw(n3d)
complex, intent(in) :: chi
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
complex, intent(in) :: z_b(nbasfcn,neigd),z_b2(nbasfcn,neigd)
#else
real, intent(in) :: z_b(nbasfcn,neigd),z_b2(nbasfcn,neigd)
#endif
complex, intent(inout) :: uHu(nbnd,nbnd)
c ..local variables..
#if (!defined(CPP_INVERSION)||defined(CPP_SOC))
complex,allocatable :: vstep(:,:)
complex,allocatable :: mat(:,:)
complex :: th
#else
real,allocatable :: vstep(:,:)
real,allocatable :: mat(:,:)
complex,allocatable :: vstep_c(:,:)
complex,allocatable :: mat_c(:,:)
complex :: th_c
real,allocatable :: vstep_r(:,:)
real,allocatable :: mat_r(:,:)
real,allocatable :: uHu_tmp(:,:)
real :: th
#endif
real :: th_r
real,allocatable :: rk_b(:),rk_b2(:)
real :: phase,phase2,ekk,s(3)
integer :: i,j,j1,j2,j3,i1,i2,i3,in,ind
allocate( vstep(nv_b2,nv_b) )
allocate( mat(nv_b,nslibd_b2) )
allocate( rk_b (nv_b) )
allocate( rk_b2(nv_b2) )
#if (!defined(CPP_INVERSION)||defined(CPP_SOC))
#else
allocate( uHu_tmp(nslibd_b,nslibd_b2) )
#endif
IF(zMat_b%l_real) THEN
allocate( uHu_tmp(nslibd_b,nslibd_b2) )
allocate( vstep_r(nv_b2,nv_b) )
allocate( mat_r(nv_b,nslibd_b2) )
vstep_r(:,:) = 0.0
ELSE
allocate( vstep_c(nv_b2,nv_b) )
allocate( mat_c(nv_b,nslibd_b2) )
vstep_c(:,:) = CMPLX(0.0,0.0)
END IF
! set up |k+G-G(k+b1)|^2
do i=1,nv_b
......@@ -95,7 +98,6 @@ c ..local variables..
! construct vstep(g,g') ~ V(g-g')
! + Theta(g-g')*[rk_b+rk_b2]
vstep(:,:) = 0.0
do i=1,nv_b
j1 = -k1_b(i) + gb(1) - gb2(1)
j2 = -k2_b(i) + gb(2) - gb2(2)
......@@ -113,42 +115,47 @@ c ..local variables..
endif
ekk = rk_b(i) + rk_b2(j)
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
th = phase*conjg(vpw(in))
if(l_kin) th = th + phase*0.25*ekk*conjg(ustep(in))
#else
th = phase*real(vpw(in))
if(l_kin) th = th + phase*0.25*ekk*real(ustep(in))
#endif
vstep(j,i) = th
IF(zMat_b%l_real) THEN
th_r = phase*real(vpw(in))
if(l_kin) th_r = th_r + phase*0.25*ekk*real(ustep(in))
vstep_r(j,i) = th_r
ELSE
th_c = phase*conjg(vpw(in))
if(l_kin) th_c = th_c + phase*0.25*ekk*conjg(ustep(in))
vstep_c(j,i) = th_c
END IF
enddo
enddo
! complex conjugate of (z(k+b1,g))^* vstep(g,g') z(k+b2,g')
#if (!defined(CPP_INVERSION) || defined(CPP_SOC))
call CPP_BLAS_cgemm('T','N',nv_b,nslibd_b2,nv_b2,cmplx(1.0),
> vstep,nv_b2,z_b2(1+addnoco2,1),nbasfcn,
> cmplx(0.0),mat,nv_b)
mat = conjg(mat)
call CPP_BLAS_cgemm('T','N',nslibd_b,nslibd_b2,nv_b,
> chi,z_b(1+addnoco,1),nbasfcn,
> mat,nv_b,cmplx(1.0),uHu,nbnd)
#else
call CPP_BLAS_sgemm('T','N',nv_b,nslibd_b2,nv_b2,real(1.0),
> vstep,nv_b2,z_b2(1+addnoco2,1),nbasfcn,
> real(0.0),mat,nv_b)
call CPP_BLAS_sgemm('T','N',nslibd_b,nslibd_b2,nv_b,
> real(1.0),z_b(1+addnoco,1),nbasfcn,
> mat,nv_b,real(0.0),uHu_tmp,nslibd_b)
uHu(1:nslibd_b,1:nslibd_b2) = uHu(1:nslibd_b,1:nslibd_b2)
IF(zMat_b%l_real) THEN
call CPP_BLAS_sgemm('T','N',nv_b,nslibd_b2,nv_b2,real(1.0),
> vstep_r,nv_b2,zMat_b2%z_r(1+addnoco2,1),nbasfcn,
> real(0.0),mat_r,nv_b)
call CPP_BLAS_sgemm('T','N',nslibd_b,nslibd_b2,nv_b,
> real(1.0),zMat_b%z_r(1+addnoco,1),nbasfcn,
> mat_r,nv_b,real(0.0),uHu_tmp,nslibd_b)
uHu(1:nslibd_b,1:nslibd_b2) = uHu(1:nslibd_b,1:nslibd_b2)
> + uHu_tmp(1:nslibd_b,1:nslibd_b2)*chi
#endif
ELSE
call CPP_BLAS_cgemm('T','N',nv_b,nslibd_b2,nv_b2,cmplx(1.0),
> vstep_c,nv_b2,zMat_b2%z_c(1+addnoco2,1),nbasfcn,
> cmplx(0.0),mat_c,nv_b)
mat_c = conjg(mat_c)
call CPP_BLAS_cgemm('T','N',nslibd_b,nslibd_b2,nv_b,
> chi,zMat_b%z_c(1+addnoco,1),nbasfcn,
> mat_c,nv_b,cmplx(1.0),uHu,nbnd)
END IF
deallocate( vstep, mat )
deallocate( rk_b, rk_b2 )
#if (!defined(CPP_INVERSION) || defined(CPP_SOC))
#else
deallocate( uHu_tmp )
#endif
IF(zMat_b%l_real) THEN
deallocate( uHu_tmp )
deallocate( vstep_r, mat_r )
ELSE
deallocate( vstep_c, mat_c )
END IF
end subroutine
end module m_wann_uHu_int
......@@ -17,7 +17,7 @@ c***************************************c
> nslibd,nslibd_b,jspin,jspin_b,ico,
> k1,k2,k3,k1_b,k2_b,k3_b,
> jspd,nvd,area,nbasfcn,neigd,
> z,z_b,nv,nv_b,sk2,phi2,omtil,gb,gb2,qss,sign2,
> zMat,zMat_b,nv,nv_b,sk2,phi2,omtil,gb,gb2,qss,sign2,
< uHu)
use m_constants
......@@ -35,6 +35,8 @@ c***************************************c
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_zmat), INTENT(IN) :: zMat, zMat_b
c .. scalar Arguments..
logical, intent (in) :: l_noco,l_soc
......@@ -60,11 +62,6 @@ c ..array arguments..
integer, intent (in) :: k1_b(:,:),k2_b(:,:),k3_b(:,:) !k1_b(nvd,jspd),k2_b(nvd,jspd),k3_b(nvd,jspd)
complex, intent (inout) :: uHu(:,:) !uHu(nbnd,nbnd)
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
complex, intent (in):: z(:,:),z_b(:,:) !z(nbasfcn,neigd),z_b(nbasfcn,neigd)
#else
real, intent (in):: z(:,:),z_b(:,:) !z(nbasfcn,neigd),z_b(nbasfcn,neigd)
#endif
complex, intent (in):: vxy(nmzxyd,odi%n2d-1,2)
......@@ -228,15 +225,26 @@ c...for the k-point
bvac = exp(-cmplx(0.0,m*arg))*(ic**m)*
* cmplx(-duz(l,m)*bess(m) +
- uz(l,m)*sk2(irec2)*dbss(m),0.0)/
/ ((wronk1)*sqrt(omtil))
do n = 1,nslibd
acof(l,m,n) = acof(l,m,n) +
+ z(k+addnoco,n)*avac
c + conjg(z(k,n))*avac
bcof(l,m,n) = bcof(l,m,n) +
+ z(k+addnoco,n)*bvac
c + conjg(z(k,n))*bvac
enddo
/ ((wronk1)*sqrt(omtil))
IF(zMat%l_real) THEN
do n = 1,nslibd
acof(l,m,n) = acof(l,m,n) +
+ zMat%z_r(k+addnoco,n)*avac
c + conjg(zMat%z_r(k,n))*avac
bcof(l,m,n) = bcof(l,m,n) +
+ zMat%z_r(k+addnoco,n)*bvac
c + conjg(zMat%z_r(k,n))*bvac
enddo
ELSE
do n = 1,nslibd
acof(l,m,n) = acof(l,m,n) +
+ zMat%z_c(k+addnoco,n)*avac
c + conjg(zMat%z_c(k,n))*avac
bcof(l,m,n) = bcof(l,m,n) +
+ zMat%z_c(k+addnoco,n)*bvac
c + conjg(zMat%z_c(k,n))*bvac
enddo
END IF
enddo ! -mb:mb
endif
enddo
......@@ -273,15 +281,27 @@ c...for the b-point
bvac = exp(-cmplx(0.0,m*arg))*(ic**m)*
& cmplx(-duz_b(l,m)*bess(m) +
- uz_b(l,m)*sk2(irec2)*dbss(m),0.0)/
/ ((wronk1)*sqrt(omtil))
do n = 1,nslibd_b
acof_b(l,m,n) = acof_b(l,m,n) +
+ z_b(k+addnoco2,n)*avac
c + conjg(z(k,n))*avac
bcof_b(l,m,n) = bcof_b(l,m,n) +
+ z_b(k+addnoco2,n)*bvac
c + conjg(z(k,n))*bvac
enddo
/ ((wronk1)*sqrt(omtil))
IF(zMat_b%l_real) THEN
do n = 1,nslibd_b
acof_b(l,m,n) = acof_b(l,m,n) +
+ zMat_b%z_r(k+addnoco2,n)*avac
c + conjg(zMat_b%z_r(k,n))*avac
bcof_b(l,m,n) = bcof_b(l,m,n) +
+ zMat_b%z_r(k+addnoco2,n)*bvac
c + conjg(zMat_b%z_r(k,n))*bvac
enddo
ELSE
do n = 1,nslibd_b
acof_b(l,m,n) = acof_b(l,m,n) +
+ zMat_b%z_c(k+addnoco2,n)*avac
c + conjg(zMat_b%z_c(k,n))*avac
bcof_b(l,m,n) = bcof_b(l,m,n) +
+ zMat_b%z_c(k+addnoco2,n)*bvac
c + conjg(zMat_b%z_c(k,n))*bvac
enddo
END IF
enddo ! -mb:mb
endif
enddo ! k = 1,nv
......
......@@ -16,13 +16,17 @@ c***************************************c
> nmz,delz,ig2,nq2,kv2,area,bmat,bbmat,evac,evac_b,
> bkpt,bkpt_b,vzxy,vz,nslibd,nslibd_b,jspin,jspin_b,
> ico,k1,k2,k3,k1_b,k2_b,k3_b,jspd,nvd,nbasfcn,neigd,
> z,z_b,nv,nv_b,omtil,gb,gb2,sign2,uHu)
> zMat,zMat_b,nv,nv_b,omtil,gb,gb2,sign2,uHu)
#include "cpp_double.h"
USE m_types
use m_constants, only : pimach
use m_intgr, only : intgz0
use m_d2fdz2cmplx
implicit none
TYPE(t_zmat), INTENT(IN) :: zMat, zMat_b
c .. scalar Arguments..
logical, intent (in) :: l_noco,l_soc,zrfs
integer, intent (in) :: nlotot,jspin_b,n2d,jspins,ico
......@@ -49,12 +53,6 @@ c ..array arguments..
complex, intent (in) :: vzxy(nmzxyd,n2d-1,2)
complex, intent (inout) :: uHu(nbnd,nbnd)
#if ( !defined(CPP_INVERSION) || defined(CPP_SOC) )
complex, intent (in):: z(nbasfcn,neigd),z_b(nbasfcn,neigd)
#else
real, intent (in):: z(nbasfcn,neigd),z_b(nbasfcn,neigd)
#endif
c ..local arrays..
complex, allocatable :: ac(:,:),bc(:,:),ac_b(:,:),bc_b(:,:)
complex, allocatable :: vv(:),fac(:),fac_b(:),fac1(:)
......@@ -275,10 +273,18 @@ c if (nvac==1) symvacvac=2
av = -c_1 * cmplx( dte(l),zks*te(l) )
bv = c_1 * cmplx( dt(l),zks* t(l) )
c-----> loop over basis functions
do n = 1,nslibd
ac(l,n) = ac(l,n) + z(k+addnoco,n)*av
bc(l,n) = bc(l,n) + z(k+addnoco,n)*bv
enddo
IF(zMat%l_real) THEN
do n = 1,nslibd
ac(l,n) = ac(l,n) + zMat%z_r(k+addnoco,n)*av
bc(l,n) = bc(l,n) + zMat%z_r(k+addnoco,n)*bv
enddo
ELSE
do n = 1,nslibd
ac(l,n) = ac(l,n) + zMat%z_c(k+addnoco,n)*av
bc(l,n) = bc(l,n) + zMat%z_c(k+addnoco,n)*bv
enddo
END IF
enddo
......@@ -314,10 +320,18 @@ c-----> construct a and b coefficients for the k+b point
av = -c_1 * cmplx( dte_b(l),zks*te_b(l) )
bv = c_1 * cmplx( dt_b(l),zks*t_b(l) )
c-----> loop over basis functions
do n = 1,nslibd_b
ac_b(l,n) = ac_b(l,n) + z_b(k+addnoco2,n)*av
bc_b(l,n) = bc_b(l,n) + z_b(k+addnoco2,n)*bv
enddo
IF(zMat_b%l_real) THEN
do n = 1,nslibd_b
ac_b(l,n) = ac_b(l,n) + zMat_b%z_r(k+addnoco2,n)*av
bc_b(l,n) = bc_b(l,n) + zMat_b%z_r(k+addnoco2,n)*bv
enddo
ELSE
do n = 1,nslibd_b
ac_b(l,n) = ac_b(l,n) + zMat_b%z_c(k+addnoco2,n)*av
bc_b(l,n) = bc_b(l,n) + zMat_b%z_c(k+addnoco2,n)*bv
enddo
END IF
enddo
......
......@@ -4,7 +4,7 @@ c composite matrices as c
c input for wannier90 c
c*************************c
module m_wann_gwf_commat2
USE m_fleurenv
USE m_juDFT
implicit none
contains
......@@ -101,7 +101,7 @@ c*************************c
inquire(file='proj',exist=l_proj)
if(.not.l_proj) write(*,*)'missing: proj'
if(l_miss.or.(.not.l_exist).or.(.not.l_proj))
> call fleur_err("missing file(s) for HDWFs")
> call juDFT_error("missing file(s) for HDWFs")
! get number of bands and wfs from proj
open(405,file='proj',status='old')
......@@ -238,7 +238,7 @@ c*****************************c
> gb_q(1:3,1:nntot_q,q))
mmn(:,:,nn,ikqpt)=mmnq(:,:,nq,k)
else ! otherwise problem
call fleur_err("problem in overlap mmn gwf",
call juDFT_error("problem in overlap mmn gwf",
> calledby="wann_gwf_commat")
endif
enddo!nn
......@@ -318,7 +318,7 @@ c*****************************c
enddo
enddo
else ! otherwise problem
call fleur_err("problem in overlap mmn gwf",
call juDFT_error("problem in overlap mmn gwf",
> calledby="wann_gwf_commat")
endif
enddo!nn
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment