```
subroutine FinalizeNewOrbs()
! At the end of the orbital rotation, have a set of coefficients CoeffT1 which transform
! the HF orbitals into a set of linear
! combinations ui which minimise |<ij|kl>|^2. This is the final subroutine after
! all iterations (but before the memory deallocation)
! that calculates the final 4 index integrals to be used in the NECI calculation.
use sym_mod, only: GenSymStatePairs
integer :: i, a, j
real(dp) :: TotGSConstraints, GSConstraint(TotNoConstraints)
HElement_t(dp) :: CoeffTemp(SpatOrbs, SpatOrbs)
! First need to do a final explicit orthonormalisation. The orbitals
! are very close to being orthonormal, but not exactly. Need to make
! sure they are exact orthonormal using Gram Schmit.
if (tStoreSpinOrbs .and. (.not. tMaxHLGap)) then
CoeffTemp(:, :) = 0.0_dp
do i = 1, SpatOrbs
do j = 1, SpatOrbs
CoeffTemp(i, j) = CoeffT1(2 * i, 2 * j)
end do
end do
call GRAMSCHMIDT_NECI(CoeffTemp, SpatOrbs)
CoeffT1(:, :) = 0.0_dp
do i = 1, SpatOrbs
do j = 1, SpatOrbs
CoeffT1(2 * i, 2 * j) = CoeffTemp(i, j)
CoeffT1((2 * i) - 1, (2 * j) - 1) = CoeffTemp(i, j)
end do
end do
else if (.not. tMaxHLGap) then
call GRAMSCHMIDT_NECI(CoeffT1, NoOrbs)
end if
! Put routine in here that takes this rotation matrix, CoeffT1, and forms raises it to the power of a small number, alpha.
! Changeing this number allows us to see the change in plateau level with various rotations.
! Write out some final results of interest, like values of the constraints, values of new coefficients.
write(stdout, *) 'The final transformation coefficients after gram schmidt orthonormalisation'
do i = 1, NoOrbs
do a = 1, NoOrbs
write(stdout, '(F10.4)', advance='no') CoeffT1(a, i)
end do
write(stdout, *) ''
end do
call WriteTransformMat()
call CalcConstraints(CoeffT1, GSConstraint, TotGSConstraints)
write(stdout, *) 'Final Potential Energy before orthogonalisation', PotEnergy
call Transform2ElInts()
! Use these final coefficients to find the FourIndInts(i,j,k,l).
! These are now the <ij|kl> integrals we now want to use instead of the HF UMat.
! New potential energy is calculated in this routine using the orthogonalised coefficients.
! Compare to that before this, to make sure the orthogonalisation hasn't shifted them back to a non-minimal place.
write(stdout, *) 'Final Potential Energy after orthogonalisation', PotEnergy
! Calculate the fock matrix, and print it out to see how much the off diagonal terms contribute.
! Also print out the sum of the diagonal elements to compare to the original value.
call CalcFOCKMatrix()
call RefillUMATandTMAT2D()
! UMat is the 4 index integral matrix (2 electron), whereas TMAT2D is the 2 index integral (1 el) matrix
! This is the keyword that tells the NECI calculation that the orbitals are not HF. It means that contributions to
! the energy from walkers on singly occupied determinants are included in the values printed.
! Making it true here allows us to go directly from a Rotation into a spawn if required.
tRotatedOrbs = .true.
call GENSymStatePairs(SpatOrbs, .false.)
end subroutine FinalizeNewOrbs
```

NECI was developed by Alavi Group

© 2023

Documentation generated by FORD