Only used when interfaced to the CP code Not sure why this is done here and not in CP
subroutine overlap_rotate
!%%%%%%%%%%%%%%%%%%%%%
!! Only used when interfaced to the CP code
!! Not sure why this is done here and not in CP
use w90_parameters, only: num_bands, a_matrix, m_matrix_orig, nntot, timing_level
use w90_io, only: io_file_unit, io_error, io_stopwatch
implicit none
integer :: lam_unit, info, inn, i, j
real(kind=DP) :: lambda(num_bands, num_bands)
real(kind=DP) :: AP(num_bands*(num_bands + 1)/2)
real(kind=DP) :: eig(num_bands), work(3*num_bands)
if (timing_level > 1) call io_stopwatch('overlap: rotate', 1)
lam_unit = io_file_unit()
open (unit=lam_unit, file='lambda.dat', &
form='unformatted', status='old', action='read')
!~ write(stdout,*) ' Reading lambda.dat...'
read (lam_unit) lambda
!~ write(stdout,*) ' done'
close (lam_unit)
do j = 1, num_bands
do i = 1, j
AP(i + (j - 1)*j/2) = 0.5_dp*(lambda(i, j) + lambda(j, i))
end do
end do
CALL DSPEV('V', 'U', num_bands, AP, eig, lambda, num_bands, work, info)
if (info .ne. 0) &
call io_error('Diagonalization of lambda in overlap_rotate failed')
! For debugging
!~ write(stdout,*) 'EIGENVALUES - CHECK WITH CP OUTPUT'
!~ do i=1,num_bands
!~ write(stdout,*) 13.6058*eig(i)
!~ end do
! Rotate M_mn
do inn = 1, nntot
m_matrix_orig(:, :, inn, 1) = &
matmul(transpose(lambda), matmul(m_matrix_orig(:, :, inn, 1), lambda))
end do
! Rotate A_mn
a_matrix(:, :, 1) = matmul(transpose(lambda), a_matrix(:, :, 1))
! For debugging
!~ ! Write rotated A and M
!~ do i=1,num_bands
!~ do j=1,num_wann
!~ write(12,'(2i5,a,2f18.12)') i,j,' 1',a_matrix(i,j,1)
!~ enddo
!~ enddo
!~ do inn=1,nntot
!~ do i=1,num_bands
!~ do j=1,num_wann
!~ write(11,'(2i5,2a)') i,j,' 1',' 1'
!~ write(11,'(2f18.12)') m_matrix_orig(i,j,inn,1)
!~ enddo
!~ enddo
!~ enddo
!~ stop
if (timing_level > 1) call io_stopwatch('overlap: rotate', 2)
return
end subroutine overlap_rotate