Write out the matrix elements of r
subroutine hamiltonian_write_rmn()
!! Write out the matrix elements of r
!============================================!
use w90_parameters, only: m_matrix, wb, bk, num_wann, num_kpts, kpt_latt, &
nntot, write_bvec
use w90_constants, only: twopi, cmplx_i
use w90_io, only: io_error, io_file_unit, seedname, io_date
implicit none
complex(kind=dp) :: fac
real(kind=dp) :: rdotk
integer :: loop_rpt, m, n, nkp, ind, nn, file_unit
complex(kind=dp) :: position(3)
character(len=33) :: header
character(len=9) :: cdate, ctime
file_unit = io_file_unit()
open (file_unit, file=trim(seedname)//'_r.dat', form='formatted', status='unknown', err=101)
call io_date(cdate, ctime)
header = 'written on '//cdate//' at '//ctime
write (file_unit, *) header ! Date and time
write (file_unit, *) num_wann
write (file_unit, *) nrpts
do loop_rpt = 1, nrpts
do m = 1, num_wann
do n = 1, num_wann
position(:) = 0._dp
do nkp = 1, num_kpts
rdotk = twopi*dot_product(kpt_latt(:, nkp), real(irvec(:, loop_rpt), dp))
fac = exp(-cmplx_i*rdotk)/real(num_kpts, dp)
do ind = 1, 3
do nn = 1, nntot
if (m .eq. n) then
! For loop_rpt==rpt_origin, this reduces to
! Eq.(32) of Marzari and Vanderbilt PRB 56,
! 12847 (1997). Otherwise, is is Eq.(44)
! Wang, Yates, Souza and Vanderbilt PRB 74,
! 195118 (2006), modified according to
! Eqs.(27,29) of Marzari and Vanderbilt
position(ind) = position(ind) - &
wb(nn)*bk(ind, nn, nkp)*aimag(log(m_matrix(n, m, nn, nkp)))*fac
else
! Eq.(44) Wang, Yates, Souza and Vanderbilt PRB 74, 195118 (2006)
position(ind) = position(ind) + &
cmplx_i*wb(nn)*bk(ind, nn, nkp)*m_matrix(n, m, nn, nkp)*fac
endif
end do
end do
end do
write (file_unit, '(5I5,6F12.6)') irvec(:, loop_rpt), n, m, position(:)
end do
end do
end do
close (file_unit)
return
101 call io_error('Error: hamiltonian_write_rmn: problem opening file '//trim(seedname)//'_r')
end subroutine hamiltonian_write_rmn