setup_exchange_matrix Subroutine

    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