subroutine sitesym_symmetrize_u_matrix(ndim, umat, lwindow_in)
!==========================================================================!
! !
! calculate U(Rk)=d(R,k)*U(k)*D^{\dagger}(R,k) in the following two cases: !
! !
! 1. "disentanglement" phase (present(lwindow)) !
! ndim=num_bands !
! !
! 2. Minimization of Omega_{D+OD} (.not.present(lwindow)) !
! ndim=num_wann, d=d_matrix_band !
! !
!==========================================================================!
use w90_parameters, only: num_wann, num_bands, num_kpts
implicit none
integer, intent(in) :: ndim
complex(kind=dp), intent(inout) :: umat(ndim, num_wann, num_kpts)
logical, optional, intent(in) :: lwindow_in(num_bands, num_kpts)
! local
integer :: ik, ir, isym, irk, n
logical :: ldone(num_kpts)
complex(kind=dp) :: cmat(ndim, num_wann)
if (present(lwindow_in) .and. (ndim .ne. num_bands)) call io_error('ndim!=num_bands')
if (.not. present(lwindow_in)) then
if (ndim .ne. num_wann) call io_error('ndim!=num_wann')
endif
ldone = .false.
do ir = 1, nkptirr
ik = ir2ik(ir)
ldone(ik) = .true.
if (present(lwindow_in)) then
n = count(lwindow_in(:, ik))
else
n = ndim
endif
if (present(lwindow_in)) then
call symmetrize_ukirr(ir, ndim, umat(:, :, ik), n)
else
call symmetrize_ukirr(ir, ndim, umat(:, :, ik))
endif
do isym = 2, nsymmetry
irk = kptsym(isym, ir)
if (ldone(irk)) cycle
ldone(irk) = .true.
! cmat = d(R,k) * U(k)
call zgemm('N', 'N', n, num_wann, n, cmplx_1, &
d_matrix_band(:, :, isym, ir), ndim, &
umat(:, :, ik), ndim, cmplx_0, cmat, ndim)
! umat(Rk) = cmat*D^{+}(R,k) = d(R,k) * U(k) * D^{+}(R,k)
call zgemm('N', 'C', n, num_wann, num_wann, cmplx_1, cmat, ndim, &
d_matrix_wann(:, :, isym, ir), num_wann, cmplx_0, umat(:, :, irk), ndim)
enddo
enddo
if (any(.not. ldone)) call io_error('error in sitesym_symmetrize_u_matrix')
return
end subroutine sitesym_symmetrize_u_matrix