getExcitation_guga Subroutine

    subroutine getExcitation_guga(nI, nJ, ex)
        ! routine to determine excitation in guga basis
        ! for now to it very naively and unellegant by converting to ilut
        ! format and calculating occupation vectors
        integer, intent(in) :: nI(nEl), nJ(nEl)
        integer, intent(out) :: ex(2, 2)

        integer(n_int) :: ilutI(0:niftot), ilutJ(0:niftot)
        integer :: first, last, cnt_e, cnt_h, occ_diff(nSpatOrbs), i

        call EncodeBitDet_guga(nI, ilutI)
        call EncodeBitDet_guga(nJ, ilutJ)

        occ_diff = calcOcc_vector_int(ilutI(0:nifd)) - calcOcc_vector_int(ilutJ(0:nifd))

        select case (sum(abs(occ_diff)))

        case (0)
            ! there is a different definition of RDMs in the guga case or?
            ! because for _RL_ -> ^RL^ excitations or nI = nJ i end up here
            ! but a lot of orbital index combinations can lead to
            ! this type of excitation.. where should i assign it to..?
            ! read the R.Shephard paper and think of a new calculation of
            ! the RDM calculation in the GUGA case..

            ! ... for now, but in the indices of the first and last switch..

            first = findFirstSwitch(ilutI, ilutJ, 1, nSpatOrbs)

            last = findLastSwitch(ilutI, ilutJ, first, nSpatOrbs)

            ex(1, 1) = 2 * first
            ex(1, 2) = 2 * last - 1

            ex(2, 1) = 2 * first - 1
            ex(2, 2) = 2 * last

        case (2)

            ! this is a "normal" double excitation
            ! find the electron in nI which gets excited

            ex(1, 2) = 0
            ex(2, 2) = 0

            do i = 1, nSpatOrbs
                if (occ_diff(i) == 1) ex(1, 1) = 2 * i
                if (occ_diff(i) == -1) ex(2, 1) = 2 * i
            end do

        case (4)

            ! keep count of the already found electrons and holes
            cnt_e = 1
            cnt_h = 1

            do i = 1, nSpatOrbs

                select case (occ_diff(i))

                    ! this choice of default spin-orbitals below
                    ! makes certain two_rdm samplings as default alos..
                    ! not sure if this choice alone is valid..
                    ! also have to ensure i get the "spins" right so the
                    ! rest of the NECI RDM routines can handle that..
                case (2)
                    ! two eletrons get excited from orb i:
                    ex(1, 1) = 2 * i - 1
                    ex(1, 2) = 2 * i

                case (1)
                    ! one electron gets excited from i
                    ! at the first encountered electron cnt_e = 1
                    ! -> so this below gives me an alpha electron!
                    ! -> at the second it will be an beta to ensure
                    ! i get "correct" spins..
                    ex(1, cnt_e) = 2 * i - cnt_e + 1

                    cnt_e = cnt_e + 1

                case (-1)
                    ! one hole found
                    ex(2, cnt_h) = 2 * i - cnt_h + 1

                    cnt_h = cnt_h + 1

                case (-2)
                    ! two electron get excited to orb i
                    ex(2, 1) = 2 * i - 1
                    ex(2, 2) = 2 * i

                end select

            end do

        end select

    end subroutine getExcitation_guga