constrained_sample_nI Subroutine

    subroutine constrained_sample_nI(this, contain, renormalization, pos, tgt, prob)
        class(AliasSampler_t), intent(in) :: this
        integer, intent(in) :: contain(:)
        real(dp), intent(in) :: renormalization
        integer, intent(out) :: pos, tgt
        real(dp), intent(out) :: prob
        character(*), parameter :: this_routine = 'constrained_sample_nI'

        ASSERT(is_set(contain))
        ASSERT(1 <= contain(1) .and. contain(size(contain)) <= size(this%probs%ptr))
        ASSERT(renormalization .isclose. (sum(this%get_prob(contain))))

        if (near_zero(renormalization)) then
            tgt = 0
            prob = 1.0
            return
        end if

        if (renormalization < redrawing_cutoff) then
        block
            type(CDF_Sampler_t) :: cdf_sampler
            cdf_sampler = CDF_Sampler_t(this%probs%ptr(contain), renormalization)
            call cdf_sampler%sample(pos, prob)
            tgt = contain(pos)
            return
        end block
        end if

        tgt = this%table%sample()
        pos = int(binary_search_int(contain, tgt))
        do while (pos == -1)
            tgt = this%table%sample()
            pos = int(binary_search_int(contain, tgt))
        end do
        prob = this%probs%ptr(tgt) / renormalization
        ASSERT(prob .isclose. (this%probs%ptr(tgt) / sum(this%probs%ptr(contain))))
        ASSERT(contain(pos) == tgt)
    end subroutine constrained_sample_nI