create_hilbert_space Subroutine

    subroutine create_hilbert_space(nI, n_states, state_list_ni, state_list_ilut, &
                                    gen_all_excits_opt)
        ! a really basic routine, which creates the whole hilbert space based
        ! on a input determinant and other quantities, like symmetry sectors,
        ! already set outside the routine. for now this is specifically
        ! implemented for the k-space hubbard model, since i still need to
        ! test the transcorrelated approach there!
        integer, intent(in) :: nI(nel)
        integer, intent(out) :: n_states
        integer, intent(out), allocatable :: state_list_ni(:, :)
        integer(n_int), intent(out), allocatable :: state_list_ilut(:, :)
        procedure(generate_all_excits_t), optional :: gen_all_excits_opt
        character(*), parameter :: this_routine = "create_hilbert_space"

        procedure(generate_all_excits_t), pointer :: gen_all_excits
        integer(n_int), allocatable :: excit_list(:, :), temp_list_ilut(:, :)
        integer, allocatable :: temp_list_ni(:, :)
        integer :: n_excits, n_total, tmp_n_states, cnt, i, j, pos
        integer(n_int) :: ilutI(0:niftot)

        ! determine the type of system by the gen_all_excits routine

        if (present(gen_all_excits_opt)) then
            gen_all_excits => gen_all_excits_opt
        else
            gen_all_excits => gen_all_excits_default
        end if

        ! estimate the total number of excitations
        n_total = int(choose_i64(nBasis / 2, nOccAlpha) * choose_i64(nBasis / 2, nOccBeta))

        n_states = 1
        allocate(temp_list_ilut(0:niftot, n_total))
        allocate(temp_list_ni(nel, n_total))

        call EncodeBitDet(nI, ilutI)

        temp_list_ilut(:, 1) = ilutI
        temp_list_ni(:, 1) = nI

        tmp_n_states = 1
        ! thats a really inefficient way to do it:
        ! think of smth better at some point!
        do while (.true.)

            ! and i need a counter, which counts the number of added
            ! excitations to the whole list.. i guess if this is 0
            ! the whole hilbert space is reached..
            cnt = 0

            ! i need a temporary loop variable
            do i = 1, tmp_n_states
                call gen_all_excits(temp_list_ni(:, i), n_excits, excit_list)

                ! now i have to check if those states are already in the list
                do j = 1, n_excits

                    pos = binary_search_ilut(temp_list_ilut(:, 1:(tmp_n_states + cnt)), &
                                        excit_list(:, j), nifd + 1)

                    ! if not yet found:
                    if (pos < 0) then
                        ! insert it at the correct place
                        ! does - pos give me the correct place then?
                        pos = -pos
                        ! lets try.. and temp_list is always big enough i think..
                        ! first move
                        temp_list_ilut(:, (pos + 1):tmp_n_states + cnt + 1) = &
                            temp_list_ilut(:, pos:(tmp_n_states + cnt))

                        temp_list_ni(:, (pos + 1):(tmp_n_states + cnt + 1)) = &
                            temp_list_ni(:, pos:(tmp_n_states + cnt))
                        ! then insert
                        temp_list_ilut(:, pos) = excit_list(:, j)

                        call decode_bit_det(temp_list_ni(:, pos), excit_list(:, j))

                        ! and increase the number of state counter
                        cnt = cnt + 1
                    else
                        ! if already found i do not need to do anything i
                        ! guess..
                    end if
                end do

            end do
            tmp_n_states = tmp_n_states + cnt

            ! and somehow i need an exit criteria, if we found all the states..
            if (cnt == 0) exit
        end do

        n_states = tmp_n_states

        ! it should be already sorted or?? i think so..
        ! or does binary_search not indicate the position
        allocate(state_list_ni(nel, n_states), source=temp_list_ni(:, 1:n_states))
        allocate(state_list_ilut(0:niftot, n_states), source=temp_list_ilut(:, 1:n_states))

    end subroutine create_hilbert_space