tran_reduce_hr Subroutine

  subroutine tran_reduce_hr(real_space_ham, ham_r, hr_one_dim, real_lattice, irvec, mp_grid, &
                            irvec_max, nrpts, nrpts_one_dim, num_wann, one_dim_vec, timing_level, &
                            stdout, timer, error, comm)
    !================================================!
    !
    ! reduce ham_r from 3-d to 1-d
    !
    !================================================!

    use w90_constants, only: dp, eps8
    use w90_io, only: io_stopwatch_start, io_stopwatch_stop
    use w90_wannier90_types, only: real_space_ham_type
    use w90_error, only: w90_error_type, set_error_alloc, set_error_fatal
    use w90_types, only: timer_list_type

    implicit none

    ! passed vars
    type(real_space_ham_type), intent(in) :: real_space_ham
    type(timer_list_type), intent(inout) :: timer
    type(w90_error_type), allocatable, intent(out) :: error
    type(w90_comm_type), intent(in) :: comm

    integer, intent(in) :: irvec(:, :)
    integer, intent(inout) :: irvec_max ! limits of hr_one_dim final dim
    integer, intent(in) :: mp_grid(3)
    integer, intent(in) :: nrpts
    integer, intent(in) :: num_wann
    integer, intent(inout) :: nrpts_one_dim
    integer, intent(inout) :: one_dim_vec
    integer, intent(in) :: timing_level
    integer, intent(in) :: stdout

    real(kind=dp), allocatable, intent(inout) :: hr_one_dim(:, :, :)
    real(kind=dp), intent(in) :: real_lattice(3, 3)

    complex(kind=dp), intent(in) :: ham_r(:, :, :)

    ! local variables
    integer :: ierr
    integer :: irvec_tmp(3), two_dim_vec(2)
    integer :: i, j
    integer :: i1, i2, i3, n1, nrpts_tmp, loop_rpt

    if (timing_level > 1) call io_stopwatch_start('tran: reduce_hr', timer)

    ! Find one_dim_vec which is parallel to one_dim_dir
    ! two_dim_vec - the other two lattice vectors
    j = 0
    do i = 1, 3
      if (abs(abs(real_lattice(real_space_ham%one_dim_dir, i)) &
              - sqrt(dot_product(real_lattice(:, i), real_lattice(:, i)))) .lt. eps8) then
        one_dim_vec = i
        j = j + 1
      end if
    end do
    if (j .ne. 1) then
      write (stdout, '(i3,a)') j, ' : 1-D LATTICE VECTOR NOT DEFINED'
      call set_error_fatal(error, 'Error: 1-d lattice vector not defined in tran_reduce_hr', comm)
      return
    end if

    j = 0
    do i = 1, 3
      if (i .ne. one_dim_vec) then
        j = j + 1
        two_dim_vec(j) = i
      end if
    end do

    ! starting H matrix should include all W-S supercell where
    ! the center of the cell spans the full space of the home cell
    ! adding one more buffer layer when mp_grid(one_dim_vec) is an odd number

    !irvec_max = (mp_grid(one_dim_vec)+1)/2
    irvec_tmp = maxval(irvec, DIM=2) + 1
    irvec_max = irvec_tmp(one_dim_vec)
    nrpts_one_dim = 2*irvec_max + 1

    allocate (hr_one_dim(num_wann, num_wann, -irvec_max:irvec_max), stat=ierr)
    if (ierr /= 0) then
      call set_error_alloc(error, 'Error in allocating hr_one_dim in tran_reduce_hr', comm)
      return
    end if
    hr_one_dim = 0.0_dp

    ! check imaginary part
    write (stdout, '(1x,a,F12.6)') 'Maximum imaginary part of the real-space Hamiltonian: ', &
      maxval(abs(aimag(ham_r)))

    ! select a subset of ham_r, where irvec is 0 along the two other lattice vectors

    nrpts_tmp = 0
    loop_n1: do n1 = -irvec_max, irvec_max
      do loop_rpt = 1, nrpts
        i1 = mod(n1 - irvec(one_dim_vec, loop_rpt), mp_grid(one_dim_vec))
        i2 = irvec(two_dim_vec(1), loop_rpt)
        i3 = irvec(two_dim_vec(2), loop_rpt)
        if (i1 .eq. 0 .and. i2 .eq. 0 .and. i3 .eq. 0) then
          nrpts_tmp = nrpts_tmp + 1
          hr_one_dim(:, :, n1) = real(ham_r(:, :, loop_rpt), dp)
          cycle loop_n1
        end if
      end do
    end do loop_n1

    if (nrpts_tmp .ne. nrpts_one_dim) then
      write (stdout, '(a)') 'FAILED TO EXTRACT 1-D HAMILTONIAN'
      call set_error_fatal(error, 'Error: cannot extract 1d hamiltonian in tran_reduce_hr', comm)
      return
    end if

    if (timing_level > 1) call io_stopwatch_stop('tran: reduce_hr', timer)

    return

  end subroutine tran_reduce_hr