BB_a(R)=<0n|H(r-R)|Rm> is the Fourier transform of BB_a(k) = i (a=x,y,z)
subroutine get_BB_R
!=====================================================
!
!! BB_a(R)=<0n|H(r-R)|Rm> is the Fourier transform of
!! BB_a(k) = i<u|H|del_a u> (a=x,y,z)
!
!=====================================================
use w90_constants, only: dp, cmplx_0, cmplx_i
use w90_parameters, only: num_kpts, nntot, nnlist, num_wann, num_bands, &
ndimwin, eigval, wb, bk, have_disentangled, &
timing_level, nncell, scissors_shift
use w90_postw90_common, only: nrpts, v_matrix
use w90_io, only: stdout, io_file_unit, io_error, io_stopwatch, &
seedname
use w90_comms, only: on_root, comms_bcast
integer :: idir, n, m, nn, i, ii, j, jj, &
ik, ik2, inn, nnl, nnm, nnn, &
winmin_q, winmin_qb, ncount, &
nb_tmp, nkp_tmp, nntot_tmp, mmn_in
complex(kind=dp), allocatable :: S_o(:, :)
complex(kind=dp), allocatable :: BB_q(:, :, :, :)
complex(kind=dp), allocatable :: H_q_qb(:, :)
integer, allocatable :: num_states(:)
real(kind=dp) :: m_real, m_imag
logical :: nn_found
character(len=60) :: header
if (timing_level > 1 .and. on_root) call io_stopwatch('get_oper: get_BB_R', 1)
if (.not. allocated(BB_R)) then
allocate (BB_R(num_wann, num_wann, nrpts, 3))
else
if (timing_level > 1 .and. on_root) call io_stopwatch('get_oper: get_BB_R', 2)
return
end if
if (on_root) then
if (abs(scissors_shift) > 1.0e-7_dp) &
call io_error('Error: scissors correction not yet implemented for BB_R')
allocate (BB_q(num_wann, num_wann, num_kpts, 3))
allocate (S_o(num_bands, num_bands))
allocate (H_q_qb(num_wann, num_wann))
allocate (num_states(num_kpts))
do ik = 1, num_kpts
if (have_disentangled) then
num_states(ik) = ndimwin(ik)
else
num_states(ik) = num_wann
endif
enddo
mmn_in = io_file_unit()
open (unit=mmn_in, file=trim(seedname)//'.mmn', &
form='formatted', status='old', action='read', err=103)
write (stdout, '(/a)', advance='no') &
' Reading overlaps from '//trim(seedname)//'.mmn in get_BB_R : '
! Read the comment line (header)
read (mmn_in, '(a)', err=104, end=104) header
write (stdout, '(a)') trim(header)
! Read the number of bands, k-points and nearest neighbours
read (mmn_in, *, err=104, end=104) nb_tmp, nkp_tmp, nntot_tmp
! Checks
if (nb_tmp .ne. num_bands) &
call io_error(trim(seedname)//'.mmn has wrong number of bands')
if (nkp_tmp .ne. num_kpts) &
call io_error(trim(seedname)//'.mmn has wrong number of k-points')
if (nntot_tmp .ne. nntot) &
call io_error &
(trim(seedname)//'.mmn has wrong number of nearest neighbours')
BB_q = cmplx_0
do ncount = 1, num_kpts*nntot
!
!Read from .mmn file the original overlap matrix
! S_o=<u_ik|u_ik2> between ab initio eigenstates
!
read (mmn_in, *, err=104, end=104) ik, ik2, nnl, nnm, nnn
do n = 1, num_bands
do m = 1, num_bands
read (mmn_in, *, err=104, end=104) m_real, m_imag
S_o(m, n) = cmplx(m_real, m_imag, kind=dp)
enddo
enddo
nn = 0
nn_found = .false.
do inn = 1, nntot
if ((ik2 .eq. nnlist(ik, inn)) .and. &
(nnl .eq. nncell(1, ik, inn)) .and. &
(nnm .eq. nncell(2, ik, inn)) .and. &
(nnn .eq. nncell(3, ik, inn))) then
if (.not. nn_found) then
nn_found = .true.
nn = inn
else
call io_error('Error reading '//trim(seedname)//'.mmn.&
& More than one matching nearest neighbour found')
endif
endif
end do
if (nn .eq. 0) then
write (stdout, '(/a,i8,2i5,i4,2x,3i3)') &
' Error reading '//trim(seedname)//'.mmn:', &
ncount, ik, ik2, nn, nnl, nnm, nnn
call io_error('Neighbour not found')
end if
call get_win_min(ik, winmin_q)
call get_win_min(nnlist(ik, nn), winmin_qb)
call get_gauge_overlap_matrix( &
ik, num_states(ik), &
nnlist(ik, nn), num_states(nnlist(ik, nn)), &
S_o, H=H_q_qb)
do idir = 1, 3
BB_q(:, :, ik, idir) = BB_q(:, :, ik, idir) &
+ cmplx_i*wb(nn)*bk(idir, nn, ik)*H_q_qb(:, :)
enddo
enddo !ncount
close (mmn_in)
call fourier_q_to_R(BB_q(:, :, :, 1), BB_R(:, :, :, 1))
call fourier_q_to_R(BB_q(:, :, :, 2), BB_R(:, :, :, 2))
call fourier_q_to_R(BB_q(:, :, :, 3), BB_R(:, :, :, 3))
endif !on_root
call comms_bcast(BB_R(1, 1, 1, 1), num_wann*num_wann*nrpts*3)
if (timing_level > 1 .and. on_root) call io_stopwatch('get_oper: get_BB_R', 2)
return
103 call io_error &
('Error: Problem opening input file '//trim(seedname)//'.mmn')
104 call io_error &
('Error: Problem reading input file '//trim(seedname)//'.mmn')
end subroutine get_BB_R