Calculate the Wannier Function spread !
Centre constraint contribution. Zero if slwf_constrain=false Contribution from constrains on centres wann_spread%om_c = wann_spread%om_iod + wann_spread%om_d + wann_spread%om_nu
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| complex(kind=dp), | intent(in) | :: | csheet(:,:,:) | |||
| real(kind=dp), | intent(in) | :: | sheet(:,:,:) | |||
| real(kind=dp), | intent(out) | :: | rave(:,:) | |||
| real(kind=dp), | intent(out) | :: | r2ave(:) | |||
| real(kind=dp), | intent(out) | :: | rave2(:) | |||
| type(localisation_vars), | intent(out) | :: | wann_spread |
subroutine wann_omega(csheet, sheet, rave, r2ave, rave2, wann_spread)
!==================================================================!
! !
!! Calculate the Wannier Function spread !
! !
! Modified by Valerio Vitale for the SLWF+C method (PRB 90, 165125)!
! Jun 2018, based on previous work by Charles T. Johnson and !
! Radu Miron at Implerial College London
!===================================================================
use w90_parameters, only: num_wann, m_matrix, nntot, wb, bk, num_kpts, &
omega_invariant, timing_level, &
selective_loc, slwf_constrain, slwf_num, &
ccentres_cart
use w90_io, only: io_stopwatch
implicit none
complex(kind=dp), intent(in) :: csheet(:, :, :)
real(kind=dp), intent(in) :: sheet(:, :, :)
real(kind=dp), intent(out) :: rave(:, :)
real(kind=dp), intent(out) :: r2ave(:)
real(kind=dp), intent(out) :: rave2(:)
type(localisation_vars), intent(out) :: wann_spread
!local variables
real(kind=dp) :: summ, mnn2
real(kind=dp) :: brn
integer :: ind, nkp, nn, m, n, iw, nkp_loc
if (timing_level > 1 .and. on_root) call io_stopwatch('wann: omega', 1)
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
do n = 1, num_wann
! Note that this ln_tmp is defined differently wrt the one in wann_domega
ln_tmp_loc(n, nn, nkp_loc) = (aimag(log(csheet(n, nn, nkp) &
*m_matrix_loc(n, n, nn, nkp_loc))) - sheet(n, nn, nkp))
end do
end do
end do
rave = 0.0_dp
do iw = 1, num_wann
do ind = 1, 3
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
rave(ind, iw) = rave(ind, iw) + wb(nn)*bk(ind, nn, nkp) &
*ln_tmp_loc(iw, nn, nkp_loc)
enddo
enddo
enddo
enddo
call comms_allreduce(rave(1, 1), num_wann*3, 'SUM')
rave = -rave/real(num_kpts, dp)
rave2 = 0.0_dp
do iw = 1, num_wann
rave2(iw) = sum(rave(:, iw)*rave(:, iw))
enddo
! aam: is this useful?
!~ rtot=0.0_dp
!~ do ind = 1, 3
!~ do loop_wann = 1, num_wann
!~ rtot (ind) = rtot (ind) + rave (ind, loop_wann)
!~ enddo
!~ enddo
r2ave = 0.0_dp
do iw = 1, num_wann
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
mnn2 = real(m_matrix_loc(iw, iw, nn, nkp_loc)*conjg(m_matrix_loc(iw, iw, nn, nkp_loc)), kind=dp)
r2ave(iw) = r2ave(iw) + wb(nn)*(1.0_dp - mnn2 + ln_tmp_loc(iw, nn, nkp_loc)**2)
enddo
enddo
enddo
call comms_allreduce(r2ave(1), num_wann, 'SUM')
r2ave = r2ave/real(num_kpts, dp)
!~ wann_spread%om_1 = 0.0_dp
!~ do nkp = 1, num_kpts
!~ do nn = 1, nntot
!~ do loop_wann = 1, num_wann
!~ wann_spread%om_1 = wann_spread%om_1 + wb(nn) * &
!~ ( 1.0_dp - m_matrix(loop_wann,loop_wann,nn,nkp) * &
!~ conjg(m_matrix(loop_wann,loop_wann,nn,nkp)) )
!~ enddo
!~ enddo
!~ enddo
!~ wann_spread%om_1 = wann_spread%om_1 / real(num_kpts,dp)
!~
!~ wann_spread%om_2 = 0.0_dp
!~ do loop_wann = 1, num_wann
!~ sqim = 0.0_dp
!~ do nkp = 1, num_kpts
!~ do nn = 1, nntot
!~ sqim = sqim + wb(nn) * &
!~ ( (aimag(log(csheet(loop_wann,nn,nkp) * &
!~ m_matrix(loop_wann,loop_wann,nn,nkp))) - &
!~ sheet(loop_wann,nn,nkp))**2 )
!~ enddo
!~ enddo
!~ sqim = sqim / real(num_kpts,dp)
!~ wann_spread%om_2 = wann_spread%om_2 + sqim
!~ enddo
!~
!~ wann_spread%om_3 = 0.0_dp
!~ do loop_wann = 1, num_wann
!~ bim = 0.0_dp
!~ do ind = 1, 3
!~ do nkp = 1, num_kpts
!~ do nn = 1, nntot
!~ bim(ind) = bim(ind) &
!~ + wb(nn) * bk(ind,nn,nkp) &
!~ * ( aimag(log(csheet(loop_wann,nn,nkp) &
!~ * m_matrix(loop_wann,loop_wann,nn,nkp))) &
!~ - sheet(loop_wann,nn,nkp) )
!~ enddo
!~ enddo
!~ enddo
!~ bim = bim/real(num_kpts,dp)
!~ bim2 = 0.0_dp
!~ do ind = 1, 3
!~ bim2 = bim2 + bim (ind) * bim (ind)
!~ enddo
!~ wann_spread%om_3 = wann_spread%om_3 - bim2
!~ enddo
!jry: Either the above (om1,2,3) or the following is redundant
! keep it in the code base for testing
if (selective_loc) then
wann_spread%om_iod = 0.0_dp
do nkp_loc = 1, counts(my_node_id)
do nn = 1, nntot
summ = 0.0_dp
do n = 1, slwf_num
summ = summ &
+ real(m_matrix_loc(n, n, nn, nkp_loc) &
*conjg(m_matrix_loc(n, n, nn, nkp_loc)), kind=dp)
if (slwf_constrain) then
!! Centre constraint contribution. Zero if slwf_constrain=false
summ = summ - lambda_loc*ln_tmp_loc(n, nn, nkp_loc)**2
end if
enddo
wann_spread%om_iod = wann_spread%om_iod &
+ wb(nn)*(real(slwf_num, dp) - summ)
enddo
enddo
call comms_allreduce(wann_spread%om_iod, 1, 'SUM')
wann_spread%om_iod = wann_spread%om_iod/real(num_kpts, dp)
wann_spread%om_d = 0.0_dp
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
do n = 1, slwf_num
brn = sum(bk(:, nn, nkp)*rave(:, n))
wann_spread%om_d = wann_spread%om_d + (1.0_dp - lambda_loc)*wb(nn) &
*(ln_tmp_loc(n, nn, nkp_loc) + brn)**2
enddo
enddo
enddo
call comms_allreduce(wann_spread%om_d, 1, 'SUM')
wann_spread%om_d = wann_spread%om_d/real(num_kpts, dp)
wann_spread%om_nu = 0.0_dp
!! Contribution from constrains on centres
if (slwf_constrain) then
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
do n = 1, slwf_num
wann_spread%om_nu = wann_spread%om_nu + 2.0_dp*wb(nn)* &
ln_tmp_loc(n, nn, nkp_loc)*lambda_loc* &
sum(bk(:, nn, nkp)*ccentres_cart(n, :))
enddo
enddo
enddo
call comms_allreduce(wann_spread%om_nu, 1, 'SUM')
wann_spread%om_nu = wann_spread%om_nu/real(num_kpts, dp)
do n = 1, slwf_num
wann_spread%om_nu = wann_spread%om_nu + lambda_loc*sum(ccentres_cart(n, :)**2)
end do
end if
wann_spread%om_tot = wann_spread%om_iod + wann_spread%om_d + wann_spread%om_nu
!! wann_spread%om_c = wann_spread%om_iod + wann_spread%om_d + wann_spread%om_nu
else
if (first_pass) then
wann_spread%om_i = 0.0_dp
nkp = nkp_loc + displs(my_node_id)
do nkp_loc = 1, counts(my_node_id)
do nn = 1, nntot
summ = 0.0_dp
do m = 1, num_wann
do n = 1, num_wann
summ = summ &
+ real(m_matrix_loc(n, m, nn, nkp_loc)*conjg(m_matrix_loc(n, m, nn, nkp_loc)), kind=dp)
enddo
enddo
wann_spread%om_i = wann_spread%om_i &
+ wb(nn)*(real(num_wann, dp) - summ)
enddo
enddo
call comms_allreduce(wann_spread%om_i, 1, 'SUM')
wann_spread%om_i = wann_spread%om_i/real(num_kpts, dp)
first_pass = .false.
else
wann_spread%om_i = omega_invariant
endif
wann_spread%om_od = 0.0_dp
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
do m = 1, num_wann
do n = 1, num_wann
if (m .ne. n) wann_spread%om_od = wann_spread%om_od &
+ wb(nn)*real(m_matrix_loc(n, m, nn, nkp_loc) &
*conjg(m_matrix_loc(n, m, nn, nkp_loc)), kind=dp)
enddo
enddo
enddo
enddo
call comms_allreduce(wann_spread%om_od, 1, 'SUM')
wann_spread%om_od = wann_spread%om_od/real(num_kpts, dp)
wann_spread%om_d = 0.0_dp
do nkp_loc = 1, counts(my_node_id)
nkp = nkp_loc + displs(my_node_id)
do nn = 1, nntot
do n = 1, num_wann
brn = sum(bk(:, nn, nkp)*rave(:, n))
wann_spread%om_d = wann_spread%om_d + wb(nn) &
*(ln_tmp_loc(n, nn, nkp_loc) + brn)**2
enddo
enddo
enddo
call comms_allreduce(wann_spread%om_d, 1, 'SUM')
wann_spread%om_d = wann_spread%om_d/real(num_kpts, dp)
wann_spread%om_tot = wann_spread%om_i + wann_spread%om_d + wann_spread%om_od
end if
if (timing_level > 1 .and. on_root) call io_stopwatch('wann: omega', 2)
return
end subroutine wann_omega