mixedFullStart Subroutine

private subroutine mixedFullStart(ilut, csf_i, excitInfo, plusWeight, minusWeight, zeroWeight, tempExcits, nExcits)

Arguments

Type IntentOptional Attributes Name
integer(kind=n_int), intent(in) :: ilut(0:nifguga)
type(CSF_Info_t), intent(in) :: csf_i
type(ExcitationInformation_t), intent(in) :: excitInfo
real(kind=dp), intent(in) :: plusWeight
real(kind=dp), intent(in) :: minusWeight
real(kind=dp), intent(in) :: zeroWeight
integer(kind=n_int), intent(out), allocatable :: tempExcits(:,:)
integer, intent(out) :: nExcits

Contents

Source Code


Source Code

    subroutine mixedFullStart(ilut, csf_i, excitInfo, plusWeight, minusWeight, &
                              zeroWeight, tempExcits, nExcits)
        ! remember full-start matrix element are stored in the same row
        ! as deltaB = -1 mixed ones... so access matrix element below with
        ! deltaB = -1 !!
        integer(n_int), intent(in) :: ilut(0:nifguga)
        type(CSF_Info_t), intent(in) :: csf_i
        type(ExcitationInformation_t), intent(in) :: excitInfo
        real(dp), intent(in) :: plusWeight, minusWeight, zeroWeight
        integer(n_int), intent(out), allocatable :: tempExcits(:, :)
        integer, intent(out) :: nExcits
        character(*), parameter :: this_routine = "mixedFullStart"

        integer :: st, nmax, ierr
        integer(n_int) :: t(0:nifguga)
        real(dp) :: tempWeight, tempWeight_1, bVal

        ASSERT(.not. isZero(ilut, excitInfo%fullStart))
        ASSERT(isProperCSF_ilut(ilut))

        ! depending on the stepvector create certain starting branches

        st = excitInfo%fullStart
        bVal = csf_i%B_real(st)

        ! determine worst case amount of excitations:
        nMax = 2 + 2**count_open_orbs_ij(csf_i, st, excitInfo%fullEnd)
        allocate(tempExcits(0:nifguga, nMax), stat=ierr)

        ! assert that at least one of the weights is non-zero
        ! otherwise checkCompatibility() has done smth wrong

        t = ilut

        select case (csf_i%stepvector(st))
        case (3)
            ! only deltaB 0 branch possible, and even no change in stepvector
            call encode_matrix_element(t, -Root2, 1)
            call encode_matrix_element(t, 0.0_dp, 2)

            call setDeltaB(0, t)
            tempExcits(:, 1) = t
            nExcits = 1

            ! do the mixed fullstart new with all weights contributed for
        case (1)
            ASSERT(zeroWeight + plusWeight > 0.0_dp)
            if (zeroWeight > 0.0_dp .and. plusWeight > 0.0_dp) then

                ! depending on weights, and b value maybe 2 excitations possible
                ! zero weight always > 0, so only check other
                ! both excitations possible
                ! first do 1->1
                call setDeltaB(0, t)

                call getDoubleMatrixElement(1, 1, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, tempWeight, tempWeight_1)

                call encode_matrix_element(t, tempWeight, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                ! the  change it to the 1->2 start

                clr_orb(t, 2 * st - 1)
                set_orb(t, 2 * st)

                call setDeltaB(2, t)

                call getDoubleMatrixElement(2, 1, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, x1_element=tempWeight_1)

                call encode_matrix_element(t, 0.0_dp, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 2) = t

                nExcits = 2

            else if (near_zero(plusWeight)) then
                ! only 0 branch possible
                call setDeltaB(0, t)

                call getDoubleMatrixElement(1, 1, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, tempWeight, tempWeight_1)

                call encode_matrix_element(t, tempWeight, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                nExcits = 1

            else if (near_zero(zeroWeight)) then
                ! only the switch to the +2 branch valid

                clr_orb(t, 2 * st - 1)
                set_orb(t, 2 * st)

                call setDeltaB(2, t)

                call getDoubleMatrixElement(2, 1, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, x1_element=tempWeight_1)

                call encode_matrix_element(t, 0.0_dp, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                nExcits = 1

            else
                ! something went wrong probably... 0 branch always
                ! possible remember, so only 0 branch not possible..
                call stop_all(this_routine, "something went wrong. should not be here!")
            end if

        case (2)
            ! here b value is never a problem so i just have to check for
            ! minusWeight
            ASSERT(minusWeight + zeroWeight > 0.0_dp)

            if (near_zero(minusWeight)) then
                ! only 0 branch possible
                call setDeltaB(0, t)

                call getDoubleMatrixElement(2, 2, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, tempWeight, tempWeight_1)

                call encode_matrix_element(t, tempWeight, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                nExcits = 1

            else if (minusWeight > 0.0_dp .and. zeroWeight > 0.0_dp) then
                ! both branches possible, ! do 0 first
                call setDeltaB(0, t)

                call getDoubleMatrixElement(2, 2, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, tempWeight, tempWeight_1)

                call encode_matrix_element(t, tempWeight, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                ! then change that to 2->1 start
                set_orb(t, 2 * st - 1)
                clr_orb(t, 2 * st)

                call setDeltaB(-2, t)

                call getDoubleMatrixElement(1, 2, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, x1_element=tempWeight_1)

                call encode_matrix_element(t, 0.0_dp, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 2) = t

                nExcits = 2

            else if (near_zero(zeroWeight)) then
                ! only -2 start possible

                ! then change that to 2->1 start
                set_orb(t, 2 * st - 1)
                clr_orb(t, 2 * st)

                call setDeltaB(-2, t)

                call getDoubleMatrixElement(1, 2, -1, gen_type%L, gen_type%R, &
                                            bVal, 1.0_dp, x1_element=tempWeight_1)

                call encode_matrix_element(t, 0.0_dp, 1)
                call encode_matrix_element(t, tempWeight_1, 2)

                tempExcits(:, 1) = t

                nExcits = 1

            else
                ! smth went wrong
                call stop_all(this_routine, "both starting branch weights are 0")

            end if
        end select

    end subroutine mixedFullStart