subroutine return_mp1_amp_and_mp2_energy(nI, ilut, ex, tParity, amp, energy_contrib)
! For a given determinant (input as nI), find the amplitude of it in the MP1 wavefunction.
! Also return the contribution from this determinant in the MP2 energy.
! To use this routine, generate an excitation from the Hartree-Fock determinant using the
! GenExcitations3 routine. This will return nI, ex and tParity which can be input into this
! routine.
use Determinants, only: GetH0Element3, GetH0Element4
use FciMCData, only: ilutHF, HFDet, Fii
use SystemData, only: tUEG
use SystemData, only: tGUGA
use guga_matrixElements, only: calcDiagMatEleGUGA_nI, calc_guga_matrix_element
use guga_data, only: ExcitationInformation_t
integer, intent(in) :: nI(nel)
integer(n_int), intent(in) :: ilut(0:NIfTot)
integer, intent(in) :: ex(2, maxExcit)
logical, intent(in) :: tParity
real(dp), intent(out) :: amp, energy_contrib
integer :: ic
HElement_t(dp) :: hel, H0tmp, denom
type(ExcitationInformation_t) :: excitInfo
amp = 0.0_dp
energy_contrib = 0.0_dp
if (ex(1, 2) == 0) then
ic = 1
else
ic = 2
end if
if (tHPHF) then
! Assume since we are using HPHF that the alpha and
! beta orbitals of the same spatial orbital have the same
! fock energies, so can consider either.
hel = hphf_off_diag_helement(HFDet, nI, iLutHF, ilut)
else if (tGUGA) then
call calc_guga_matrix_element(&
ilut, CSF_Info_t(ilut), ilutHF, CSF_Info_t(ilutHF), excitInfo, hel, .true.)
else
hel = get_helement(HFDet, nI, ic, ex, tParity)
end if
if (tUEG) then
! This will calculate the MP2 energies without having to use the fock eigenvalues.
! This is done via the diagonal determinant hamiltonian energies.
H0tmp = getH0Element4(nI, HFDet)
else if (tGUGA) then
! do i have a routine to calculate the diagonal and double
! contributions for GUGA csfs? yes!
H0tmp = calcDiagMatEleGUGA_nI(nI)
else
H0tmp = getH0Element3(nI)
end if
! If the relevant excitation from the Hartree-Fock takes electrons from orbitals
! (i,j) to (a,b), then denom will be equal to
! \epsilon_a + \epsilon_b - \epsilon_i - \epsilon_j
! as required in the denominator of the MP1 amplitude and MP2 energy.
denom = Fii - H0tmp
if (.not. (abs(denom) > 0.0_dp)) then
call warning_neci("return_mp1_amp_and_mp2_energy", &
"One of the determinants under consideration for the MP1 wave function is degenerate &
&with the Hartree-Fock determinant. Degenerate perturbation theory has not been &
&considered, but the amplitude of this determinant will be returned as huge(0.0_dp) &
&so that it should be included in the space.")
amp = huge(0.0_dp)
else
amp = hel / denom
energy_contrib = (hel**2) / denom
end if
end subroutine return_mp1_amp_and_mp2_energy