SUBROUTINE SetupTMAT2(nBASISFRZ, iSS, iSize)
! In:
! nBasisFRZ: number of active basis functions (orbitals).
! iSS: ratio of nBasisFRZ to the number of active spatial orbitals.
! iSize: number of elements in TMat2/TMatSym2.
! Initial allocation of TMat2D2 or TMatSym2 (if using symmetry-compressed
! storage of the <i|h|j> integrals) for post-freezing.
! See also notes in SetupTMat.
IMPLICIT NONE
integer Nirrep, nBasisfrz, iSS, nBi, i, basirrep, t, ierr, iState, nStateIrrep
integer iSize, iunit
character(*), parameter :: thisroutine = 'SetupTMAT2'
! If this is a CPMD k-point calculation, then we're operating
! under Abelian symmetry: can use George's memory efficient
! TMAT.
if (tCPMD) tCPMDSymTMat = tKP
if (tVASP) tCPMDSymTMat = .true.
IF (tCPMDSymTMat) THEN
! Set up info for indexing scheme (see TMatInd for full description).
Nirrep = NSYMLABELS
nBi = nBasisFRZ / iSS
iSize = 0
if (associated(SymLabelIntsCum2)) then
deallocate(SymLabelIntsCum2)
call LogMemDealloc(thisroutine, tagSymLabelIntsCum2)
end if
if (allocated(StateSymMap2)) then
deallocate(StateSymMap2)
call LogMemDealloc(thisroutine, tagStateSymMap2)
end if
allocate(SymLabelIntsCum2(nIrrep))
call LogMemAlloc('SymLabelIntsCum2', nIrrep, 4, thisroutine, tagSymLabelIntsCum2)
allocate(SymLabelCountsCum2(nIrrep))
call LogMemAlloc('SymLabelCountsCum2', nIrrep, 4, thisroutine, tagSymLabelCountsCum2)
allocate(StateSymMap2(nBi))
call LogMemAlloc('StateSymMap2', nBi, 4, thisroutine, tagStateSymMap2)
SYMLABELINTSCUM2(1:Nirrep) = 0
SYMLABELCOUNTSCUM2(1:Nirrep) = 0
! Ensure no compiler warnings
basirrep = SYMLABELCOUNTS(2, 1)
do i = 1, Nirrep
!SYMLABELCOUNTS is now mbas only for the frozen orbitals
basirrep = SYMLABELCOUNTS(2, i)
iSize = iSize + (basirrep * (basirrep + 1)) / 2
! # of integrals in this symmetry class and all preceeding
! symmetry classes.
SYMLABELINTSCUM2(i) = iSize
! JSS: we no longer need SymLabelCountsCum, but it's a useful test.
! SymLabelCountsCum is the cumulative total # of basis functions
! in the preceeding symmetry classes.
IF (i == 1) THEN
SYMLABELCOUNTSCUM2(i) = 0
ELSE
DO t = 1, (i - 1)
SYMLABELCOUNTSCUM2(i) = SYMLABELCOUNTSCUM2(i) + SYMLABELCOUNTS(2, t)
end do
end if
! write(stdout,*) basirrep,SYMLABELINTSCUM(i),SYMLABELCOUNTSCUM(i)
! call neci_flush(stdout)
! JSS: Label states of symmetry i by the order in which they come.
nStateIrrep = 0
do iState = 1, nBi
if (SymClasses2(iState) == i) then
nStateIrrep = nStateIrrep + 1
StateSymMap2(iState) = nStateIrrep
end if
end do
end do
IF ((SYMLABELCOUNTSCUM2(Nirrep) + basirrep) /= nBI) THEN
iunit = get_free_unit()
open(iunit, file="SYMLABELCOUNTS", status="unknown")
DO i = 1, Nirrep
write(iunit, *) SYMLABELCOUNTS(2, i)
CALL neci_flush(iunit)
end do
call stop_all(thisroutine, 'Not all basis functions found while setting up TMAT2')
end if
!iSize=iSize+2
!This is to allow the index of '-1' in the array to give a zero value
!Refer to TMatSym(-1) for when <i|h|j> is zero by symmetry.
allocate(TMATSYM2(-1:iSize), STAT=ierr)
CALL LogMemAlloc('TMatSym2', iSize + 2, HElement_t_size * 8, thisroutine, tagTMATSYM2, ierr)
TMATSYM2 = (0.0_dp)
ELSE
if (tOneElecDiag) then
! In the UEG, the orbitals are eigenfunctions of KE operator, so TMAT is diagonal.
iSize = nBasisFRZ
allocate(TMAT2D2(nBasisFRZ, 1), STAT=ierr)
call LogMemAlloc('TMAT2D2', nBasisFRZ, HElement_t_size * 8, thisroutine, tagTMat2D2)
TMAT2D2 = (0.0_dp)
else
! Using a square array to hold <i|h|j> (incl. elements which are
! zero by symmetry).
iSize = nBasisFRZ * nBasisFRZ
allocate(TMAT2D2(nBasisFRZ, nBasisFRZ), STAT=ierr)
call LogMemAlloc('TMAT2D2', nBasisFRZ * nBasisFRZ, HElement_t_size * 8, thisroutine, tagTMat2D2)
TMAT2D2 = (0.0_dp)
end if
end if
END SUBROUTINE SetupTMat2