spin_get_moment Subroutine

public subroutine spin_get_moment()

Uses

  • proc~~spin_get_moment~~UsesGraph proc~spin_get_moment spin_get_moment module~w90_postw90_common w90_postw90_common proc~spin_get_moment->module~w90_postw90_common module~w90_constants w90_constants proc~spin_get_moment->module~w90_constants module~w90_get_oper w90_get_oper proc~spin_get_moment->module~w90_get_oper module~w90_io w90_io proc~spin_get_moment->module~w90_io module~w90_comms w90_comms proc~spin_get_moment->module~w90_comms module~w90_parameters w90_parameters proc~spin_get_moment->module~w90_parameters module~w90_postw90_common->module~w90_constants module~w90_postw90_common->module~w90_comms module~w90_get_oper->module~w90_constants module~w90_io->module~w90_constants module~w90_comms->module~w90_constants module~w90_comms->module~w90_io module~w90_parameters->module~w90_constants module~w90_parameters->module~w90_io

Computes the spin magnetic moment by Wannier interpolation

Arguments

None

Calls

proc~~spin_get_moment~~CallsGraph proc~spin_get_moment spin_get_moment proc~get_ss_r get_SS_R proc~spin_get_moment->proc~get_ss_r proc~get_hh_r get_HH_R proc~spin_get_moment->proc~get_hh_r proc~io_file_unit io_file_unit proc~get_ss_r->proc~io_file_unit proc~get_win_min get_win_min proc~get_hh_r->proc~get_win_min proc~fourier_q_to_r fourier_q_to_R proc~get_hh_r->proc~fourier_q_to_r proc~io_error io_error proc~get_hh_r->proc~io_error proc~get_hh_r->proc~io_file_unit

Contents

Source Code


Source Code

  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