subroutine check_and_sort_similar_centres(signatures, num_G)
!=======================================================!
! Here, we consider the possiblity of wannier functions !
! with similar centres, such as a set of d-orbitals !
! centred on an atom. We use tran_group_threshold to !
! identify them, if they exist, then use the signatures !
! to dishinguish and sort then consistently from unit !
! cell to unit cell. !
! !
! MS: For 2two-shot and beyond, some parameters, !
! eg, first_group_element will need changing to consider!
! the geometry of the new systems. !
!=======================================================!
use w90_constants, only: dp
use w90_io, only: stdout, io_stopwatch, io_error
use w90_parameters, only: tran_num_ll, num_wann, tran_num_cell_ll, iprint, timing_level, &
tran_group_threshold, write_xyz
use w90_hamiltonian, only: wannier_centres_translated
implicit none
integer, intent(in) :: num_G
real(dp), intent(in), dimension(:, :) :: signatures
integer :: i, j, k, l, ierr, group_iterator, coord_iterator, num_wf_iterator, &
num_wann_cell_ll, iterator, max_position(1), p, num_wf_cell_iter
integer, allocatable, dimension(:) :: idx_similar_wf, group_verifier, sorted_idx, centre_id
real(dp), allocatable, dimension(:) :: dot_p
integer, allocatable, dimension(:, :) :: tmp_wf_verifier, wf_verifier, first_group_element, &
ref_similar_centres, unsorted_similar_centres
integer, allocatable, dimension(:, :, :) :: wf_similar_centres
logical, allocatable, dimension(:) :: has_similar_centres
if (timing_level > 2) call io_stopwatch('tran: lcr_2c2_sort: similar_centres', 1)
num_wann_cell_ll = tran_num_ll/tran_num_cell_ll
allocate (wf_similar_centres(tran_num_cell_ll*4, num_wann_cell_ll, num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating wf_similar_centre in check_and_sort_similar_centres')
allocate (idx_similar_wf(num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating idx_similar_wf in check_and_sort_similar_centres')
allocate (has_similar_centres(num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating has_similar_centres in check_and_sort_similar_centres')
allocate (tmp_wf_verifier(4*tran_num_cell_ll, num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating tmp_wf_verifier in check_and_sort_similar_centres')
allocate (group_verifier(4*tran_num_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating group_verifier in check_and_sort_similar_centres')
allocate (first_group_element(4*tran_num_cell_ll, num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating first_group_element in check_and_sort_similar_centres')
allocate (centre_id(num_wann_cell_ll), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating centre_id in check_and_sort_similar_centres')
!
! First find WFs with similar centres: store in wf_similar_centres(cell#,group#,WF#)
!
group_verifier = 0
tmp_wf_verifier = 0
first_group_element = 0
centre_id = 0
!
! Loop over unit cells in PL1,PL2,PL3 and PL4
!
do i = 1, 4*tran_num_cell_ll
group_iterator = 0
has_similar_centres = .false.
!
! Loops over wannier functions in present unit cell
!
num_wf_cell_iter = 0
do j = 1, num_wann_cell_ll
num_wf_iterator = 0
!
! 2nd Loop over wannier functions in the present unit cell
!
do k = 1, num_wann_cell_ll
if ((.not. has_similar_centres(k)) .and. (j .ne. k)) then
coord_iterator = 0
!
! Loop over x,y,z to find similar centres
!
do l = 1, 3
if (i .le. 2*tran_num_cell_ll) then
if (abs(wannier_centres_translated(coord(l), tran_sorted_idx(j + (i - 1)*num_wann_cell_ll)) - &
wannier_centres_translated(coord(l), tran_sorted_idx(k + (i - 1)*num_wann_cell_ll))) &
.le. tran_group_threshold) then
coord_iterator = coord_iterator + 1
else
exit
endif
else
if (abs(wannier_centres_translated(coord(l), tran_sorted_idx(num_wann - 2*tran_num_ll + &
j + (i - 2*tran_num_cell_ll - 1)*num_wann_cell_ll)) - &
wannier_centres_translated(coord(l), tran_sorted_idx(num_wann - 2*tran_num_ll + &
k + (i - 2*tran_num_cell_ll - 1)*num_wann_cell_ll))) &
.le. tran_group_threshold) then
coord_iterator = coord_iterator + 1
else
exit
endif
endif
enddo
if (coord_iterator .eq. 3) then
if (.not. has_similar_centres(j)) then
num_wf_iterator = num_wf_iterator + 1
if (i .le. 2*tran_num_cell_ll) then
idx_similar_wf(num_wf_iterator) = tran_sorted_idx(j + (i - 1)*num_wann_cell_ll)
else
idx_similar_wf(num_wf_iterator) = tran_sorted_idx(j + num_wann - 2*tran_num_ll + &
(i - 2*tran_num_cell_ll - 1)*num_wann_cell_ll)
endif
if (i .le. 2*tran_num_cell_ll) then
first_group_element(i, j) = j + (i - 1)*num_wann_cell_ll
else
first_group_element(i, j) = num_wann - 2*tran_num_ll + &
j + (i - 2*tran_num_cell_ll - 1)*num_wann_cell_ll
endif
num_wf_cell_iter = num_wf_cell_iter + 1
centre_id(num_wf_cell_iter) = j
endif
has_similar_centres(k) = .true.
has_similar_centres(j) = .true.
num_wf_iterator = num_wf_iterator + 1
if (i .le. 2*tran_num_cell_ll) then
idx_similar_wf(num_wf_iterator) = tran_sorted_idx(k + (i - 1)*num_wann_cell_ll)
else
idx_similar_wf(num_wf_iterator) = tran_sorted_idx(k + num_wann - 2*tran_num_ll + &
(i - 2*tran_num_cell_ll - 1)*num_wann_cell_ll)
endif
endif
endif
enddo ! loop over k
if (num_wf_iterator .gt. 0) then
group_iterator = group_iterator + 1
wf_similar_centres(i, group_iterator, :) = idx_similar_wf(:)
!
!Save number of WFs in each group
!
tmp_wf_verifier(i, group_iterator) = num_wf_iterator
endif
enddo
if ((count(has_similar_centres) .eq. 0) .and. (i .eq. 1)) then
write (stdout, '(a)') ' No wannier functions found with similar centres: sorting completed'
exit
elseif (i .eq. 1) then
write (stdout, *) ' Wannier functions found with similar centres: '
write (stdout, *) ' -> using signatures to complete sorting '
endif
!
!Save number of group of WFs in each unit cell and compare to previous unit cell
!
group_verifier(i) = group_iterator
if (iprint .ge. 4) write (stdout, '(a11,i4,a13,i4)') ' Unit cell:', i, ' Num groups:', group_verifier(i)
if (i .ne. 1) then
if (group_verifier(i) .ne. group_verifier(i - 1)) then
if (write_xyz) call tran_write_xyz()
call io_error('Inconsitent number of groups of similar centred wannier functions between unit cells')
elseif (i .eq. 4*tran_num_cell_ll) then
write (stdout, *) ' Consistent groups of similar centred wannier functions between '
write (stdout, *) ' unit cells found'
write (stdout, *) ' '
endif
endif
enddo !Loop over all unit cells in PL1,PL2,PL3,PL4
!
! Perform check to ensure consistent number of WFs between equivalent groups in different unit cells
!
if (any(has_similar_centres)) then
!
!
allocate (wf_verifier(4*tran_num_cell_ll, group_verifier(1)), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating wf_verifier in check_and_sort_similar_centres')
!
!
if (iprint .ge. 4) write (stdout, *) 'Unit cell Group number Num WFs'
wf_verifier = 0
wf_verifier = tmp_wf_verifier(:, 1:group_verifier(1))
do i = 1, 4*tran_num_cell_ll
do j = 1, group_verifier(1)
if (iprint .ge. 4) write (stdout, '(a3,i4,a9,i4,a7,i4)') ' ', i, ' ', j, ' ', wf_verifier(i, j)
if (i .ne. 1) then
if (wf_verifier(i, j) .ne. wf_verifier(i - 1, j)) &
call io_error('Inconsitent number of wannier &
&functions between equivalent groups of similar &
¢red wannier functions')
endif
enddo
enddo
write (stdout, *) ' Consistent number of wannier functions between equivalent groups of similar'
write (stdout, *) ' centred wannier functions'
write (stdout, *) ' '
!
write (stdout, *) ' Fixing order of similar centred wannier functions using parity signatures'
!
do i = 2, 4*tran_num_cell_ll
do j = 1, group_verifier(1)
!
! Make array of WF numbers which act as a reference to sort against
! and an array which need sorting
!
allocate (ref_similar_centres(group_verifier(1), wf_verifier(1, j)), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating ref_similar_centres in check_and_sort_similar_centres')
allocate (unsorted_similar_centres(group_verifier(1), wf_verifier(1, j)), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating unsorted_similar_centres in check_and_sort_similar_centres')
allocate (sorted_idx(wf_verifier(1, j)), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating sorted_idx in check_and_sort_similar_centres')
allocate (dot_p(wf_verifier(1, j)), stat=ierr)
if (ierr /= 0) call io_error('Error in allocating dot_p in check_and_sort_similar_centres')
!
do k = 1, wf_verifier(1, j)
ref_similar_centres(j, k) = wf_similar_centres(1, j, k)
unsorted_similar_centres(j, k) = wf_similar_centres(i, j, k)
enddo
!
sorted_idx = 0
do k = 1, wf_verifier(1, j)
dot_p = 0.0_dp
!
! building the array of positive dot products of signatures between unsorted_similar_centres(j,k)
! and all the ref_similar_centres(j,:)
!
do l = 1, wf_verifier(1, j)
do p = 1, num_G
dot_p(l) = dot_p(l) + abs(signatures(p, unsorted_similar_centres(j, k)))* &
abs(signatures(p, ref_similar_centres(j, l)))
enddo
enddo
!
max_position = maxloc(dot_p)
!
sorted_idx(max_position(1)) = unsorted_similar_centres(j, k)
enddo
!
! we have the properly ordered indexes for group j in unit cell i, now we need
! to overwrite the tran_sorted_idx array at the proper position
!
tran_sorted_idx(first_group_element(i, centre_id(j)):first_group_element(i, centre_id(j)) +&
&wf_verifier(i, j) - 1) = sorted_idx(:)
!
deallocate (dot_p, stat=ierr)
if (ierr /= 0) call io_error('Error in deallocating dot_p in check_and_sort_similar_centres')
deallocate (sorted_idx, stat=ierr)
if (ierr /= 0) call io_error('Error in deallocating sorted_idx in check_and_sort_similar_centres')
deallocate (unsorted_similar_centres, stat=ierr)
if (ierr /= 0) call io_error('Error in deallocating unsorted_similar_centres in check_and_sort_similar_centres')
deallocate (ref_similar_centres, stat=ierr)
if (ierr /= 0) call io_error('Error in deallocating ref_similar_centres in check_and_sort_similar_centres')
enddo
enddo
!
! checking that all the indices of WFs in the new tran_sorted_idx are distinct
! Remark: physically, no two WFs with similar centres can have the same type so we should expect
! this check to always pass unless the signatures/wf are very weird !!
!
do k = 1, num_wann
iterator = 0
do l = 1, num_wann
if (tran_sorted_idx(l) .eq. k) then
iterator = iterator + 1
endif
enddo
!
if ((iterator .ge. 2) .or. (iterator .eq. 0)) call io_error( &
'A Wannier Function appears either zero times or twice after sorting, this may be due to a &
&poor wannierisation and/or disentanglement')
!write(stdout,*) ' WF : ',k,' appears ',iterator,' time(s)'
enddo
deallocate (wf_verifier, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating wf_verifier in check_and_sort_similar_centres')
endif
deallocate (centre_id, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating centre_id in check_and_sort_similar_centres')
deallocate (first_group_element, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating first_group_element in check_and_sort_similar_centres')
deallocate (group_verifier, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating group_verifier in check_and_sort_similar_centres')
deallocate (tmp_wf_verifier, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating tmp_wf_verifier in check_and_sort_similar_centres')
deallocate (has_similar_centres, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating has_similar_centres in check_and_sort_similar_centres')
deallocate (idx_similar_wf, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating idx_similar_wf in check_and_sort_similar_centres')
deallocate (wf_similar_centres, stat=ierr)
if (ierr /= 0) call io_error('Error deallocating wf_similar_centre in check_and_sort_similar_centres')
if (timing_level > 2) call io_stopwatch('tran: lcr_2c2_sort: similar_centres', 2)
return
end subroutine check_and_sort_similar_centres