calc_rho_ii_and_sum_n Subroutine

public subroutine calc_rho_ii_and_sum_n(one_rdms, norm_1rdm, SumN_Rho_ii)


Type IntentOptional Attributes Name
type(one_rdm_t), intent(inout) :: one_rdms(:)
real(kind=dp), intent(in) :: norm_1rdm(:)
real(kind=dp), intent(out) :: SumN_Rho_ii(:)


Source Code

Source Code

    subroutine calc_rho_ii_and_sum_n(one_rdms, norm_1rdm, SumN_Rho_ii)

        ! Calculate the rho_ii arrays, which hold the diagonals of the 1-RDMs,
        ! appropriately normalised and sorted in terms of the energy ordering
        ! of the orbitals. Then, calculate SumN_Rho_ii, which holds the sum of
        ! this diagonal for orbitals occupied in the HF determinant.

        use FciMCData, only: HFDet_True
        use LoggingData, only: tDiagRDM
        use rdm_data, only: tOpenShell
        use RotateOrbsData, only: SymLabelListInv_rot, NoOrbs
        use SystemData, only: BRR, nel
        use UMatCache, only: gtID

        type(one_rdm_t), intent(inout) :: one_rdms(:)
        real(dp), intent(in) :: norm_1rdm(:)
        real(dp), intent(out) :: SumN_Rho_ii(:)

        integer :: irdm, i, HFDet_ID, BRR_ID

        ! Want to sum the diagonal elements of the 1-RDM for the HF orbitals.
        ! Given the HF orbitals, SymLabelListInv_rot tells us their position
        ! in the 1-RDM.
        SumN_Rho_ii = 0.0_dp

        do irdm = 1, size(one_rdms)
            do i = 1, NoOrbs
                ! rho_ii is the diagonal elements of the 1-RDM. We want this
                ! ordered according to the energy of the orbitals. Brr has the
                ! orbital numbers in order of energy... i.e Brr(2) = the orbital
                ! index with the second lowest energy. BRR is always in spin
                ! orbitals. i gives the energy level, BRR gives the orbital,
                ! SymLabelListInv_rot gives the position of  this orbital in
                ! one_rdm.

                associate(ind => SymLabelListInv_rot)
                    if (tDiagRDM) then
                        if (tOpenShell) then
                            one_rdms(irdm)%rho_ii(i) = one_rdms(irdm)%matrix(ind(BRR(i)), ind(BRR(i))) * norm_1rdm(irdm)
                            BRR_ID = gtID(BRR(2 * i))
                            one_rdms(irdm)%rho_ii(i) = one_rdms(irdm)%matrix(ind(BRR_ID), ind(BRR_ID)) * norm_1rdm(irdm)
                        end if
                    end if

                    if (i <= nel) then
                        if (tOpenShell) then
                            SumN_Rho_ii(irdm) = SumN_Rho_ii(irdm) + &
                                                (one_rdms(irdm)%matrix(ind(HFDet_True(i)), ind(HFDet_True(i))) * norm_1rdm(irdm))
                            HFDet_ID = gtID(HFDet_True(i))
                            SumN_Rho_ii(irdm) = SumN_Rho_ii(irdm) + &
                                                (one_rdms(irdm)%matrix(ind(HFDet_ID), ind(HFDet_ID)) * norm_1rdm(irdm)) / 2.0_dp
                        end if
                    end if
                end associate
            end do
        end do

    end subroutine calc_rho_ii_and_sum_n