subroutine setup_exchange_matrix(in_lat)
! by convention encode the exchange matrix element by the two
! involved electrons with opposite spins! so we can encode this
! as a 2D matrix
class(lattice), intent(in), pointer :: in_lat
character(*), parameter :: this_routine = "setup_exchange_matrix"
integer :: i, ind
ASSERT(associated(in_lat))
! create the exchange matrix from the given lattice
! connections
if (allocated(exchange_matrix)) deallocate(exchange_matrix)
allocate(exchange_matrix(nbasis, nbasis))
exchange_matrix = 0.0_dp
ASSERT(in_lat%get_nsites() == nbasis / 2)
do i = 1, in_lat%get_nsites()
ind = in_lat%get_site_index(i)
! print *, "ind, next:", ind, in_lat%get_neighbors(ind)
associate(next => in_lat%get_neighbors(ind))
exchange_matrix(2 * ind - 1, 2 * next) = exchange_j / 2.0_dp
exchange_matrix(2 * ind, 2 * next - 1) = exchange_j / 2.0_dp
ASSERT(all(next > 0))
ASSERT(all(next <= nbasis / 2))
end associate
ASSERT(ind > 0)
ASSERT(ind <= nbasis / 2)
end do
end subroutine setup_exchange_matrix