subroutine CalcPotentials()
! Only temporarily like this, can tidy it up majorly
integer :: i, j, k, l, Starti, Finishi
real(dp) :: MaxTerm
l = 0
if (tERLocalization .and. (.not. tStoreSpinOrbs)) then
ERPotEnergy = 0.0_dp
if (tRotateVirtOnly) then
Starti = NoOcc + 1
Finishi = NoOrbs
else if (tRotateOccOnly) then
Starti = 1
Finishi = NoOcc
else
Starti = 1
Finishi = NoOrbs
end if
CoulPotEnergy = 0.0_dp
OffDiagPotEnergy = 0.0_dp
do i = Starti, Finishi
ERPotEnergy = ERPotEnergy + FourIndIntsER(i)
if (FourIndIntsER(i) < 0) then
call neci_flush(stdout)
call Stop_All('CalcPotentials', 'A <ii|ii> value is less than 0.')
end if
PotEnergy = PotEnergy + FourIndIntsER(i)
TwoEInts = TwoEInts + FourIndIntsER(i)
PEInts = PEInts + FourIndIntsER(i)
end do
else if (tERLocalization) then
ERPotEnergy = 0.0_dp
PotEnergy = 0.0_dp
if (tRotateVirtOnly) then
Starti = NoOcc + 1
Finishi = NoOrbs
else if (tRotateOccOnly) then
Starti = 1
Finishi = NoOcc
else
Starti = 1
Finishi = NoOrbs
end if
do i = Starti, Finishi
if (tStoreSpinOrbs) then
if (MOD(i, 2) == 0) then
j = i - 1
else
j = i + 1
end if
ERPotEnergy = ERPotEnergy + FourIndInts(i, j, i, j)
if ((FourIndInts(i, j, i, j) < 0) .or. (FourIndInts(j, i, j, i) < 0)) then
call neci_flush(stdout)
call Stop_All('CalcPotentials', 'A <ii|ii> value is less than 0.')
end if
PotEnergy = PotEnergy + FourIndInts(i, j, i, j)
TwoEInts = TwoEInts + FourIndInts(i, j, i, j)
PEInts = PEInts + FourIndInts(i, j, i, j)
else
ERPotEnergy = ERPotEnergy + FourIndInts(i, i, i, i)
if ((FourIndInts(i, i, i, i) < 0)) then
call neci_flush(stdout)
call Stop_All('CalcPotentials', 'A <ii|ii> value is less than 0.')
end if
PotEnergy = PotEnergy + FourIndInts(i, i, i, i)
TwoEInts = TwoEInts + FourIndInts(i, i, i, i)
PEInts = PEInts + FourIndInts(i, i, i, i)
end if
end do
end if
if (tOffDiagSqrdMin .or. tOffDiagSqrdMax .or. tOffDiagMin .or. tOffdiagMax) then
do l = 1, NoOrbs
do j = 1, l - 1
do k = 1, j - 1
do i = 1, k - 1
if (tOffDiagSqrdMin .or. tOffDiagSqrdMax) then
if (((i /= j) .and. (j /= l)) .and. ((i /= k) .or. (j /= l))) then
PotEnergy = PotEnergy + (FourIndInts(i, j, k, l)**2)
TwoEInts = TwoEInts + (FourIndInts(i, j, k, l)**2)
PEInts = PEInts + (FourIndInts(i, j, k, l)**2)
end if
end if
if (tOffDiagMin .or. tOffDiagMax) then
if (.not. ((k == i) .or. (j == l))) then
PotEnergy = PotEnergy + FourIndInts(i, j, k, l)
TwoEInts = TwoEInts + FourIndInts(i, j, k, l)
PEInts = PEInts + FourIndInts(i, j, k, l)
end if
end if
end do
end do
end do
end do
end if
if (tDoubExcMin) then
do i = 1, NoOrbs
do j = 1, NoOrbs
do k = 1, i - 1
if ((k == l) .and. (k == i)) cycle
do l = 1, j - 1
if ((j == k) .and. (j == l)) cycle
if ((j == k) .and. (j == i)) cycle
if ((j == l) .and. (j == i)) cycle
PotEnergy = PotEnergy + (FourIndInts(i, j, k, l)) - FourIndInts(i, j, l, k)
TwoEInts = TwoEInts + (FourIndInts(i, j, k, l)) - FourIndInts(i, j, l, k)
PEInts = PEInts + (FourIndInts(i, j, k, l)) - FourIndInts(i, j, l, k)
end do
end do
end do
end do
end if
if (tOnePartOrbEnMax .or. tOneElIntMax) then
do i = NoOcc + 1, NoOrbs
MaxTerm = 0.0_dp
MaxTerm = TMAT2DRot(i, i)
if (tOnePartOrbEnMax) then
do j = 1, NoOcc
MaxTerm = MaxTerm + (2 * FourIndInts(i, j, i, j)) - FourIndInts(i, j, j, i)
end do
MaxTerm = MaxTerm - EpsilonMin
MaxTerm = MaxTerm**OrbEnMaxAlpha
end if
PotEnergy = PotEnergy + MaxTerm
end do
end if
if (tHijSqrdMin) then
HijSqrdPotEnergy = 0.0_dp
do i = NoOcc + 1, NoOrbs
do j = NoOcc + 1, NoOrbs
if (j > i) then
PotEnergy = PotEnergy + (TMAT2DRot(i, j)**2)
HijSqrdPotEnergy = HijSqrdPotEnergy + (TMAT2DRot(i, j)**2)
end if
end do
end do
end if
if (tVirtCoulombMax) then
ERPotEnergy = 0.0_dp
ijOccVirtPotEnergy = 0.0_dp
do i = 1, NoOrbs
if (i <= NoOcc) then
do j = NoOcc + 1, NoOrbs
ijOccVirtPotEnergy = ijOccVirtPotEnergy + FourIndInts(i, j, i, j)
end do
end if
if (i > NoOcc) then
ERPotEnergy = ERPotEnergy + FourIndInts(i, i, i, i)
do j = NoOcc + 1, NoOrbs
if (j <= i) cycle
PotEnergy = PotEnergy + FourIndInts(i, j, i, j)
TwoEInts = TwoEInts + FourIndInts(i, j, i, j)
end do
end if
end do
end if
if (tHFSingDoubExcMax) then
do i = 1, NoOcc
do j = 1, NoOcc
do k = NoOcc + 1, NoOrbs
do l = NoOcc + 1, NoOrbs
PotEnergy = PotEnergy + (FourIndInts(i, j, k, l)**2)
end do
! Sing excitations <ij|ik> where i and j are occ, k virt.
PotEnergy = PotEnergy + (FourIndInts(i, j, i, k)**2)
end do
end do
end do
end if
end subroutine CalcPotentials