Computes the spin magnetic moment by Wannier interpolation
subroutine spin_get_moment
!============================================================!
! !
!! Computes the spin magnetic moment by Wannier interpolation
! !
!============================================================!
use w90_constants, only: dp, pi, cmplx_i
use w90_comms, only: on_root, my_node_id, num_nodes, comms_reduce
use w90_io, only: io_error, stdout
use w90_postw90_common, only: num_int_kpts_on_node, int_kpts, weight
use w90_parameters, only: spin_kmesh, wanint_kpoint_file, &
nfermi, fermi_energy_list
use w90_get_oper, only: get_HH_R, get_SS_R
integer :: loop_x, loop_y, loop_z, loop_tot
real(kind=dp) :: kweight, kpt(3), spn_k(3), spn_all(3), &
spn_mom(3), magnitude, theta, phi, conv
if (nfermi > 1) call io_error('Routine spin_get_moment requires nfermi=1')
call get_HH_R
call get_SS_R
if (on_root) then
write (stdout, '(/,/,1x,a)') '------------'
write (stdout, '(1x,a)') 'Calculating:'
write (stdout, '(1x,a)') '------------'
write (stdout, '(/,3x,a)') '* Spin magnetic moment'
end if
spn_all = 0.0_dp
if (wanint_kpoint_file) then
if (on_root) then
write (stdout, '(/,1x,a)') 'Sampling the irreducible BZ only'
write (stdout, '(5x,a)') &
'WARNING: - IBZ implementation is currently limited to simple cases:'
write (stdout, '(5x,a)') &
' Check results against a full BZ calculation!'
end if
!
! Loop over k-points on the irreducible wedge of the Brillouin zone,
! read from file 'kpoint.dat'
!
do loop_tot = 1, num_int_kpts_on_node(my_node_id)
kpt(:) = int_kpts(:, loop_tot)
kweight = weight(loop_tot)
call spin_get_moment_k(kpt, fermi_energy_list(1), spn_k)
spn_all = spn_all + spn_k*kweight
end do
else
if (on_root) &
write (stdout, '(/,1x,a)') 'Sampling the full BZ (not using symmetry)'
kweight = 1.0_dp/real(PRODUCT(spin_kmesh), kind=dp)
do loop_tot = my_node_id, PRODUCT(spin_kmesh) - 1, num_nodes
loop_x = loop_tot/(spin_kmesh(2)*spin_kmesh(3))
loop_y = (loop_tot - loop_x*(spin_kmesh(2)*spin_kmesh(3)))/spin_kmesh(3)
loop_z = loop_tot - loop_x*(spin_kmesh(2)*spin_kmesh(3)) &
- loop_y*spin_kmesh(3)
kpt(1) = (real(loop_x, dp)/real(spin_kmesh(1), dp))
kpt(2) = (real(loop_y, dp)/real(spin_kmesh(2), dp))
kpt(3) = (real(loop_z, dp)/real(spin_kmesh(3), dp))
call spin_get_moment_k(kpt, fermi_energy_list(1), spn_k)
spn_all = spn_all + spn_k*kweight
end do
end if
! Collect contributions from all nodes
!
call comms_reduce(spn_all(1), 3, 'SUM')
! No factor of g=2 because the spin variable spans [-1,1], not
! [-1/2,1/2] (i.e., it is really the Pauli matrix sigma, not S)
!
spn_mom(1:3) = -spn_all(1:3)
if (on_root) then
write (stdout, '(/,1x,a)') 'Spin magnetic moment (Bohr magn./cell)'
write (stdout, '(1x,a,/)') '===================='
write (stdout, '(1x,a18,f11.6)') 'x component:', spn_mom(1)
write (stdout, '(1x,a18,f11.6)') 'y component:', spn_mom(2)
write (stdout, '(1x,a18,f11.6)') 'z component:', spn_mom(3)
! Polar and azimuthal angles of the magnetization (defined as in pwscf)
!
conv = 180.0_dp/pi
magnitude = sqrt(spn_mom(1)**2 + spn_mom(2)**2 + spn_mom(3)**2)
theta = acos(spn_mom(3)/magnitude)*conv
phi = atan(spn_mom(2)/spn_mom(1))*conv
write (stdout, '(/,1x,a18,f11.6)') 'Polar theta (deg):', theta
write (stdout, '(1x,a18,f11.6)') 'Azim. phi (deg):', phi
end if
end subroutine spin_get_moment