get_offdiag_helement_rs_hub Function

    function get_offdiag_helement_rs_hub(nI, ex, tpar) result(hel)
        integer, intent(in) :: nI(nel), ex(2)
        logical, intent(in) :: tpar
        HElement_t(dp) :: hel

        integer(n_int) :: ilut(0:NIfTot), ilutJ(0:NIfTot)
        real(dp) :: n_i, n_j
        type(ExcitationInformation_t) :: excitInfo

        if (tGUGA) then
            call EncodeBitDet(nI, ilut)
            ilutJ = make_ilutJ(ilut, ex, 1)

            call calc_guga_matrix_element(ilut, CSF_Info_t(ilut), ilutJ, CSF_Info_t(ilutJ), excitInfo, hel, .true.)

            if (tpar) hel = -hel
            return
        end if

        ! in case we need it, the off-diagonal, except parity is just
        ! -t if the hop is possible
        hel = GetTMatEl(ex(1), ex(2))

        ! like niklas, choose the alpha spins to be the correlated ones
        if (t_spin_dependent_transcorr .and. is_alpha(ex(1))) then
            hel = hel + get_1_body_contrib_spin_transcorr(nI, ex)
        end if

        if (t_trans_corr_hop) then
            hel = hel + get_2_body_contrib_transcorr_hop(nI, ex)
        end if

        if (tpar) hel = -hel

        ! put the transcorrelated stuff here for now, alhtough it would be
        ! better to do it as a procedure pointer..
        if (t_trans_corr) then
            call EncodeBitDet(nI, ilut)
            if (t_trans_corr_new) then
                n_j = get_opp_spin(ilut, ex(1))
                n_i = get_opp_spin(ilut, ex(2))

                ! try to go one step further before the transformation to
                ! the momentum space.. and check if the results are still
                ! correct..

                hel = hel * (1.0_dp + (exp(trans_corr_param) - 1.0_dp) * n_j + &
                             (exp(-trans_corr_param) - 1.0_dp) * n_i - &
                             2.0_dp * (cosh(trans_corr_param) - 1.0_dp) * n_i * n_j)

            else

                hel = hel * exp(trans_corr_param * &
                                (get_opp_spin(ilut, ex(1)) - get_opp_spin(ilut, ex(2))))
            end if

        end if

        ! it is not really a 2-body trans, but still use the flag and
        ! parameter
        if (t_trans_corr_2body) then
            call EncodeBitDet(nI, ilut)
            hel = hel * exp(trans_corr_param_2body * &
                            (get_spin_opp_neighbors(ilut, ex(1)) - get_spin_opp_neighbors(ilut, ex(2))))
        end if
    end function get_offdiag_helement_rs_hub