check_and_sort_similar_centres Subroutine

private subroutine check_and_sort_similar_centres(signatures, num_G)

proc~~check_and_sort_similar_centres~~UsesGraph proc~check_and_sort_similar_centres check_and_sort_similar_centres module~w90_parameters w90_parameters module~w90_parameters->proc~check_and_sort_similar_centres module~w90_hamiltonian w90_hamiltonian module~w90_hamiltonian->proc~check_and_sort_similar_centres module~w90_constants w90_constants module~w90_constants->proc~check_and_sort_similar_centres module~w90_constants->module~w90_parameters module~w90_constants->module~w90_hamiltonian module~w90_io w90_io module~w90_constants->module~w90_io module~w90_io->proc~check_and_sort_similar_centres module~w90_io->module~w90_parameters
Help

Arguments

Type IntentOptional AttributesName
real(kind=dp), intent(in), dimension(:,:):: signatures
integer, intent(in) :: num_G

Calls

proc~~check_and_sort_similar_centres~~CallsGraph proc~check_and_sort_similar_centres check_and_sort_similar_centres proc~tran_write_xyz tran_write_xyz proc~check_and_sort_similar_centres->proc~tran_write_xyz proc~io_stopwatch io_stopwatch proc~check_and_sort_similar_centres->proc~io_stopwatch proc~io_error io_error proc~check_and_sort_similar_centres->proc~io_error proc~io_file_unit io_file_unit proc~tran_write_xyz->proc~io_file_unit proc~io_date io_date proc~tran_write_xyz->proc~io_date proc~io_stopwatch->proc~io_error
Help

Called By

proc~~check_and_sort_similar_centres~~CalledByGraph proc~check_and_sort_similar_centres check_and_sort_similar_centres proc~tran_lcr_2c2_sort tran_lcr_2c2_sort proc~tran_lcr_2c2_sort->proc~check_and_sort_similar_centres proc~tran_main tran_main proc~tran_main->proc~tran_lcr_2c2_sort proc~wannier_run wannier_run proc~wannier_run->proc~tran_main program~wannier wannier program~wannier->proc~tran_main
Help


Source Code

 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 &
                        &centred 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


berry_get_imf_klist berry_get_imfgh_klist berry_get_kubo_k berry_main boltzwann_main calcTDFandDOS check_and_sort_similar_centres clean_ws_translate comms_allreduce comms_allreduce_cmplx comms_allreduce_real comms_array_split comms_barrier comms_bcast comms_bcast_char comms_bcast_cmplx comms_bcast_int comms_bcast_logical comms_bcast_real comms_end comms_gatherv comms_gatherv_real comms_recv comms_recv_char comms_recv_cmplx comms_recv_int comms_recv_logical comms_recv_real comms_reduce comms_reduce_cmplx comms_reduce_int comms_reduce_real comms_scatterv comms_scatterv_int_1 comms_scatterv_int_2 comms_scatterv_int_3 comms_scatterv_real comms_send comms_send_char comms_send_cmplx comms_send_int comms_send_logical comms_send_real comms_setup conv_get_seedname conv_read_chkpt conv_read_chkpt_fmt conv_write_chkpt conv_write_chkpt_fmt dis_extract dis_extract_gamma dis_main dis_proj_froz dis_project dis_windows dos_get_k dos_get_levelspacing dos_main fourier_q_to_R gauss_freq geninterp_main get_AA_R get_BB_R get_CC_R get_FF_R get_HH_R get_module_kmesh get_smearing_index get_SS_R get_win_min group hamiltonian_dealloc hamiltonian_get_hr hamiltonian_setup hamiltonian_wigner_seitz hamiltonian_write_hr hamiltonian_write_rmn hamiltonian_write_tb internal_maxloc internal_set_kmesh internal_write_header io_date io_error io_file_unit io_get_seedname io_print_timings io_stopwatch io_time k_path k_slice kmesh_dealloc kmesh_get kmesh_get_bvectors kmesh_shell_automatic kmesh_shell_fixed kmesh_shell_from_file kmesh_spacing_mesh kmesh_spacing_singleinteger kmesh_supercell_sort kmesh_write master_sort_and_group MinusFermiDerivative my_ICOPY orthogonalize_u overlap_dealloc overlap_project overlap_project_gamma overlap_read overlap_rotate param_dealloc param_get_atoms param_get_block_length param_get_keyword param_get_keyword_block param_get_keyword_kpath param_get_keyword_vector param_get_projections param_get_range_vector param_get_smearing_type param_get_vector_length param_in_file param_lib_set_atoms param_memory_estimate param_postw90_write param_read param_read_chkpt param_uppercase param_write param_write_chkpt param_write_header plot_fermi_surface plot_interpolate_bands plot_main plot_u_matrices plot_wannier print_usage pw90common_fourier_R_to_k pw90common_fourier_R_to_k_new pw90common_fourier_R_to_k_vec pw90common_get_occ pw90common_kmesh_spacing pw90common_wanint_data_dist pw90common_wanint_get_kpoint_file pw90common_wanint_param_dist pw90common_wanint_setup qe_erf qe_erfc R_wz_sc R_wz_sc_equiv script_common script_fermi_lines sitesym_dealloc sitesym_dis_extract_symmetry sitesym_read sitesym_replace_d_matrix_band sitesym_slim_d_matrix_band sitesym_symmetrize_gradient sitesym_symmetrize_rotation sitesym_symmetrize_u_matrix sitesym_symmetrize_zmatrix sort spin_get_moment spin_get_moment_k spin_get_nk symmetrize_ukirr TDF_kpt tran_bulk tran_cut_hr_one_dim tran_dealloc tran_find_integral_signatures tran_get_ht tran_green tran_lcr tran_lcr_2c2_build_ham tran_lcr_2c2_sort tran_main tran_parity_enforce tran_read_htC tran_read_htX tran_read_htXY tran_reduce_hr tran_transfer tran_write_xyz utility_cart_to_frac utility_commutator_diag utility_compar utility_diagonalize utility_frac_to_cart utility_im_tr utility_inv2 utility_inv3 utility_lowercase utility_matmul_diag utility_metric utility_re_tr utility_recip_lattice utility_rotate utility_rotate_diag utility_string_to_coord utility_strip utility_translate_home utility_w0gauss utility_wgauss utility_zgemm wann_calc_projection wann_check_unitarity wann_domega wann_main wann_main_gamma wann_omega wann_omega_gamma wann_phases wann_spread_copy wann_svd_omega_i wann_write_r2mn wann_write_vdw_data wann_write_xyz wannier_run wannier_setup wham_get_D_h wham_get_D_h_a wham_get_deleig_a wham_get_eig_deleig wham_get_eig_UU_HH_JJlist wham_get_JJm_list wham_get_JJp_list wham_get_occ_mat_list wigner_seitz ws_translate_dist ws_write_vec