subroutine FillCoeffT1
use RotateOrbsData, only: CoeffT1, SymLabelList3_rot, SymOrbs_rot, SymOrbs_rotTag, &
TruncEval, NoRotOrbs, EvaluesTrunc, EvaluesTruncTag
use LoggingData, only: tTruncRODump, tTruncDumpbyVal
integer :: l, k, i, j, NoRotAlphBet, io1, io2, ierr
character(len=*), parameter :: t_r = 'FillCoeffT1'
character(len=5) :: Label
character(len=20) :: LabelFull
real(dp) :: OccEnergies(1:NoRotOrbs)
FillCoeff_Time%timer_name = 'FillCoeff'
call set_timer(FillCoeff_Time, 30)
write(stdout, *) 'NatOrbMat'
do i = 1, nBasis
do j = 1, nBasis
write(stdout, '(F10.6)', advance='no') NatOrbMat(j, i)
end do
write(stdout, *) ''
end do
if (tTruncRODump) then
if (tTruncDumpbyVal) then
NoFrozenVirt = 0
if (tStoreSpinOrbs) then
do i = SpatOrbs, 1, -1
if (Evalues(i) > TruncEval) exit
if (Evalues(i + SpatOrbs) > TruncEval) exit
NoFrozenVirt = NoFrozenVirt + 2
end do
if (NoFrozenVirt >= (NoOrbs - NEl)) call stop_all(t_r, 'Freezing all virtual orbitals.')
else
do i = SpatOrbs, 1, -1
if (Evalues(i) > TruncEval) exit
NoFrozenVirt = NoFrozenVirt + 1
end do
if (NoFrozenVirt >= (SpatOrbs - (NEl / 2))) call stop_all(t_r, 'Freezing all virtual orbitals.')
end if
NoRotOrbs = NoOrbs - NoFrozenVirt
end if
allocate(SymOrbs_rot(NoOrbs), stat=ierr)
call LogMemAlloc('SymOrbs_rot', NoOrbs, 4, t_r, SymOrbs_rotTag, ierr)
SymOrbs_rot(:) = 0
allocate(EvaluesTrunc(NoOrbs - NoFrozenVirt), stat=ierr)
call LogMemAlloc('EvaluesTrunc', NoOrbs - NoFrozenVirt, 4, t_r, EvaluesTruncTag, ierr)
EvaluesTrunc(:) = 0.0_dp
if (tStoreSpinOrbs) then
NoRotAlphBet = SpatOrbs - (NoFrozenVirt / 2)
else
NoRotAlphBet = NoOrbs - NoFrozenVirt
end if
if (tStoreSpinOrbs) then
! As we reorder these so that they are truncated, we also need
! to pair up symmetries.
call CalcOccEnergies(OccEnergies)
! First nOccBeta, then nOccAlpha.
do i = 1, (2 * nOccBeta), 2
k = 1
do while (abs(OccEnergies(k)) < 1.0e-10_dp)
k = k + 2
end do
do j = 1, (2 * nOccBeta), 2
if ((OccEnergies(j) < OccEnergies(k)) .and. (abs(OccEnergies(j)) > 1.0e-10_dp)) k = j
end do
l = ceiling(real(k, dp) / 2.0_dp)
CoeffT1(:, i) = NatOrbMat(:, l)
EvaluesTrunc(i) = Evalues(l)
SymOrbs_rot(i) = SymOrbs_rotTemp(l)
OccEnergies(k) = 0.0_dp
end do
do i = 2, (2 * nOccAlpha), 2
k = 2
do while (abs(OccEnergies(k)) < 1.0e-10_dp)
k = k + 2
end do
do j = 2, (2 * nOccAlpha), 2
if ((OccEnergies(j) < OccEnergies(k)) .and. (abs(OccEnergies(j)) > 1.0e-10_dp)) k = j
end do
l = (k / 2) + SpatOrbs
CoeffT1(:, i) = NatOrbMat(:, l)
EvaluesTrunc(i) = Evalues(l)
SymOrbs_rot(i) = SymOrbs_rotTemp(l)
OccEnergies(k) = 0.0_dp
end do
! Need to fill coeffT1 so that it goes alpha beta alpha beta.
k = (2 * nOccBeta) + 1
do i = nOccBeta + 1, NoRotAlphBet
CoeffT1(:, k) = NatOrbMat(:, i)
EvaluesTrunc(k) = Evalues(i)
SymOrbs_rot(k) = SymOrbs_rotTemp(i)
k = k + 2
end do
k = (2 * nOccAlpha) + 2
do i = SpatOrbs + nOccAlpha + 1, SpatOrbs + NoRotAlphBet
CoeffT1(:, k) = NatOrbMat(:, i)
SymOrbs_rot(k) = SymOrbs_rotTemp(i)
EvaluesTrunc(k) = Evalues(i)
k = k + 2
end do
else
! Order occupied in terms of energy again - this makes sure
! freezing etc doesn't get screwed up.
call CalcOccEnergies(OccEnergies)
! OccEnergies has the orbital energies as they are ordered
! currently - need to put NatOrbMat into CoeffT1 so that this
! goes from lowest energy to highest.
do i = 1, NEl / 2
k = 1
do while (abs(OccEnergies(k)) < 1.0e-10_dp)
k = k + 1
end do
do j = 1, NEl / 2
if ((OccEnergies(j) < OccEnergies(k)) .and. (abs(OccEnergies(j)) > 1.0e-10_dp)) k = j
end do
CoeffT1(:, i) = NatOrbMat(:, k)
EvaluesTrunc(i) = Evalues(k)
SymOrbs_rot(i) = SymOrbs_rotTemp(k)
OccEnergies(k) = 0.0_dp
end do
do i = (NEl / 2) + 1, NoRotAlphBet
CoeffT1(:, i) = NatOrbMat(:, i)
EvaluesTrunc(i) = Evalues(i)
SymOrbs_rot(i) = SymOrbs_rotTemp(i)
end do
end if
else
do i = 1, NoOrbs
CoeffT1(:, i) = NatOrbMat(:, i)
end do
end if
if (tTruncRODump) then
Label = ''
LabelFull = ''
write(Label, '(I5)') NoFrozenVirt
LabelFull = 'EVALUES-TRUNC-'//adjustl(Label)
io1 = get_free_unit()
open(io1, file=LabelFull, status='unknown')
if (tStoreSpinOrbs) then
write(io1, *) NoOrbs - NoFrozenVirt
do i = 1, NoOrbs - NoFrozenVirt, 2
write(io1, '(I5,ES20.10,I5,A5,I5,ES20.10,I5)') i, EvaluesTrunc(i), SymOrbs_rot(i), &
' * ', i + 1, EvaluesTrunc(i + 1), SymOrbs_rot(i + 1)
end do
else
write(io1, *) NoOrbs - NoFrozenVirt
do i = 1, NoOrbs - NoFrozenVirt
write(io1, '(ES20.10,I5)') EvaluesTrunc(i), SymOrbs_rot(i)
end do
end if
close(io1)
else
io2 = get_free_unit()
open(io2, file='EVALUES', status='unknown')
write(io2, *) NoOrbs
if (tStoreSpinOrbs) then
k = 0
do i = 1, NoOrbs, 2
k = k + 1
if (tTruncRODump) then
write(io2, '(2I5,ES20.10,I5,A5,I5,ES20.10,I5)') (NoOrbs - i + 1), i, Evalues(k), SymOrbs_rot(i), &
' * ', i + 1, Evalues(k + SpatOrbs), SymOrbs_rot(i + 1)
else
write(io2, '(2I5,ES20.10,I5,A5,I5,ES20.10,I5)') (NoOrbs - i + 1), i, Evalues(k), &
int(G1(SymLabelList3_rot(k))%Sym%S, 4), ' * ', &
i + 1, Evalues(k + SpatOrbs), int(G1(SymLabelList3_rot(k + SpatOrbs))%Sym%S, 4)
end if
end do
else
do i = 1, SpatOrbs
write(io2, '(3I5,ES20.10)') i, NoOrbs - i + 1, (NoOrbs - i + 1) * 2, Evalues(i)
end do
end if
close(io2)
end if
call HistNatOrbEvalues()
call halt_timer(FillCoeff_Time)
end subroutine FillCoeffT1