create_neel_state_chain Function

    function create_neel_state_chain() result(neel_state)
        integer :: neel_state(nel)

        integer :: i

        neel_state = [(i, i=1, 2 * nel - 1, 2)]
        neel_state(2:nel:2) = neel_state(2:nel:2) + 1

    end function create_neel_state_chain