wann_svd_omega_i Subroutine

private subroutine wann_svd_omega_i()

Uses

  • proc~~wann_svd_omega_i~~UsesGraph proc~wann_svd_omega_i wann_svd_omega_i module~w90_constants w90_constants proc~wann_svd_omega_i->module~w90_constants module~w90_parameters w90_parameters proc~wann_svd_omega_i->module~w90_parameters module~w90_io w90_io proc~wann_svd_omega_i->module~w90_io module~w90_parameters->module~w90_constants module~w90_parameters->module~w90_io module~w90_io->module~w90_constants

Arguments

None

Calls

proc~~wann_svd_omega_i~~CallsGraph proc~wann_svd_omega_i wann_svd_omega_i zgesvd zgesvd proc~wann_svd_omega_i->zgesvd proc~io_error io_error proc~wann_svd_omega_i->proc~io_error

Contents

Source Code


Source Code

  subroutine wann_svd_omega_i()
    !========================================!

    use w90_constants, only: dp, cmplx_0
    use w90_io, only: io_stopwatch, io_error, stdout
    use w90_parameters, only: num_wann, num_kpts, nntot, wb, &
      m_matrix, lenconfac, length_unit, &
      timing_level

    implicit none

    complex(kind=dp), allocatable  :: cv1(:, :), cv2(:, :)
    complex(kind=dp), allocatable  :: cw1(:), cw2(:)
    complex(kind=dp), allocatable  :: cpad1(:)
    real(kind=dp), allocatable  :: singvd(:)

    integer :: ierr, info
    integer :: nkp, nn, nb, na, ind
    real(kind=dp) :: omt1, omt2, omt3

    if (timing_level > 1 .and. on_root) call io_stopwatch('wann: svd_omega_i', 1)

    allocate (cw1(10*num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating cw1 in wann_svd_omega_i')
    allocate (cw2(10*num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating cw2 in wann_svd_omega_i')
    allocate (cv1(num_wann, num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating cv1 in wann_svd_omega_i')
    allocate (cv2(num_wann, num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating cv2 in wann_svd_omega_i')
    allocate (singvd(num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating singvd in wann_svd_omega_i')
    allocate (cpad1(num_wann*num_wann), stat=ierr)
    if (ierr /= 0) call io_error('Error in allocating cpad1 in wann_svd_omega_i')

    cw1 = cmplx_0; cw2 = cmplx_0; cv1 = cmplx_0; cv2 = cmplx_0; cpad1 = cmplx_0
    singvd = 0.0_dp

    ! singular value decomposition
    omt1 = 0.0_dp; omt2 = 0.0_dp; omt3 = 0.0_dp
    do nkp = 1, num_kpts
      do nn = 1, nntot
        ind = 1
        do nb = 1, num_wann
          do na = 1, num_wann
            cpad1(ind) = m_matrix(na, nb, nn, nkp)
            ind = ind + 1
          enddo
        enddo
        call zgesvd('A', 'A', num_wann, num_wann, cpad1, num_wann, singvd, cv1, &
                    num_wann, cv2, num_wann, cw1, 10*num_wann, cw2, info)
        if (info .ne. 0) then
          call io_error('ERROR: Singular value decomp. zgesvd failed')
        endif

        do nb = 1, num_wann
          omt1 = omt1 + wb(nn)*(1.0_dp - singvd(nb)**2)
          omt2 = omt2 - wb(nn)*(2.0_dp*log(singvd(nb)))
          omt3 = omt3 + wb(nn)*(acos(singvd(nb))**2)
        enddo
      enddo
    enddo
    omt1 = omt1/real(num_kpts, dp)
    omt2 = omt2/real(num_kpts, dp)
    omt3 = omt3/real(num_kpts, dp)
    if (on_root) then
      write (stdout, *) ' '
      write (stdout, '(2x,a,f15.9,1x,a)') 'Omega Invariant:   1-s^2 = ', &
        omt1*lenconfac**2, '('//trim(length_unit)//'^2)'
      write (stdout, '(2x,a,f15.9,1x,a)') '                 -2log s = ', &
        omt2*lenconfac**2, '('//trim(length_unit)//'^2)'
      write (stdout, '(2x,a,f15.9,1x,a)') '                  acos^2 = ', &
        omt3*lenconfac**2, '('//trim(length_unit)//'^2)'
    endif

    deallocate (cpad1, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating cpad1 in wann_svd_omega_i')
    deallocate (singvd, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating singvd in wann_svd_omega_i')
    deallocate (cv2, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating cv2 in wann_svd_omega_i')
    deallocate (cv1, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating cv1 in wann_svd_omega_i')
    deallocate (cw2, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating cw2 in wann_svd_omega_i')
    deallocate (cw1, stat=ierr)
    if (ierr /= 0) call io_error('Error in deallocating cw1 in wann_svd_omega_i')

    if (timing_level > 1 .and. on_root) call io_stopwatch('wann: svd_omega_i', 2)

    return

  end subroutine wann_svd_omega_i