RefillUMATandTMAT2D Subroutine

    subroutine RefillUMATandTMAT2D()

        integer :: l, k, j, i, a, b, g, d, c, nBasis2, ierr
        integer(TagIntType) :: TMAT2DPartTag
        real(dp) :: NewTMAT
        real(dp), allocatable :: TMAT2DPart(:, :)
#ifdef CMPLX_
        call stop_all('RefillUMATandTMAT2D', 'Rotating orbitals not implemented for complex orbitals.')
#endif

        if (tStoreSpinOrbs) then
            allocate(TMAT2DPart((nBasis - NoFrozenVirt), nBasis), stat=ierr)
            call LogMemAlloc('TMAT2DPart', (nBasis - NoFrozenVirt) * nBasis, 8, 'RefillUMAT', TMAT2DPartTag, ierr)
            if (NoDumpTruncs > 1) then
                allocate(TMAT2DNew((nBasis - NoFrozenVirt), (nBasis - NoFrozenVirt)), stat=ierr)
                call LogMemAlloc('TMAT2DNew', (nBasis - NoFrozenVirt)**2, 8, 'RefillUMAT', TMAT2DNewTag, ierr)
                TMAT2DNew(:, :) = 0.0_dp
            end if
        else
            allocate(TMAT2DPart((nBasis - (NoFrozenVirt * 2)), nBasis), stat=ierr)
            call LogMemAlloc('TMAT2DPart', (nBasis - (NoFrozenVirt * 2)) * nBasis, 8, 'RefillUMAT', TMAT2DPartTag, ierr)
            if (NoDumpTruncs > 1) then
                allocate(TMAT2DNew((nBasis - NoFrozenVirt), (nBasis - NoFrozenVirt)), stat=ierr)
                call LogMemAlloc('TMAT2DNew', (nBasis - NoFrozenVirt)**2, 8, 'RefillUMAT', TMAT2DNewTag, ierr)
                TMAT2DNew(:, :) = 0.0_dp
            end if
        end if
        TMAT2DPart(:, :) = 0.0_dp

        RefillUMAT_Time%timer_name = 'RefillUMATandTMAT'
        call set_timer(RefillUMAT_Time, 30)

        do i = 1, nBasis
            write(stdout, *) SymLabelList2_rot(i), SymLabelList3_rot(i)
        end do

        ! Make the UMAT elements the four index integrals. These are calculated
        ! by transforming the HF orbitals using the coefficients that have been
        ! found.
        if (NoDumpTruncs <= 1) then
            do l = 1, (NoOrbs - (NoFrozenVirt))
                if (tTurnStoreSpinOff) then
                    d = CEILING(real(SymLabelList3_rot(l), dp) / 2.0_dp)
                else
                    d = SymLabelList3_rot(l)
                end if
                do k = 1, (NoOrbs - (NoFrozenVirt))

                    if (tTurnStoreSpinOff) then
                        g = CEILING(real(SymLabelList3_rot(k), dp) / 2.0_dp)
                    else
                        g = SymLabelList3_rot(k)
                    end if

                    do j = 1, (NoOrbs - (NoFrozenVirt))

                        if (tTurnStoreSpinOff) then
                            b = CEILING(real(SymLabelList3_rot(j), dp) / 2.0_dp)
                        else
                            b = SymLabelList3_rot(j)
                        end if
                        do i = 1, (NoOrbs - (NoFrozenVirt))

                            if (tTurnStoreSpinOff) then
                                a = CEILING(real(SymLabelList3_rot(i), dp) / 2.0_dp)
                            else
                                a = SymLabelList3_rot(i)
                            end if

                            if (tUseMP2VarDenMat .or. tFindCINatOrbs .or. tReadInCoeff) then
                                UMAT(UMatInd(a, b, g, d)) = (FourIndInts(i, k, j, l))
                            else
                                UMAT(UMatInd(a, b, g, d)) = (FourIndInts(i, j, k, l))
                            end if
                        end do
                    end do
                end do
            end do
        end if

        do a = 1, nBasis
            do k = 1, NoRotOrbs
                i = SymLabelList3_rot(k)
                NewTMAT = 0.0_dp
                do b = 1, NoOrbs
                    d = SymLabelList2_rot(b)
                    if (tStoreSpinOrbs) then
                        NewTMAT = NewTMAT + (CoeffT1(b, k) * real(TMAT2D(d, a), dp))
                    else
                        NewTMAT = NewTMAT + (CoeffT1(b, k) * real(TMAT2D(2 * d, a), dp))
                    end if
                end do
                if (tStoreSpinOrbs) then
                    TMAT2DPart(i, a) = NewTMAT
                else
                    TMAT2DPart(2 * i, a) = NewTMAT
                    TMAT2DPart(2 * i - 1, a) = NewTMAT
                end if
            end do
        end do

        if (tStoreSpinOrbs) then
            nBasis2 = nBasis - NoFrozenVirt
        else
            nBasis2 = nBasis - (NoFrozenVirt * 2)
        end if
        do k = 1, nBasis2
            do l = 1, NoRotOrbs
                j = SymLabelList3_rot(l)
                NewTMAT = 0.0_dp
                do a = 1, NoOrbs
                    c = SymLabelList2_rot(a)
                    if (tStoreSpinOrbs) then
                        NewTMAT = NewTMAT + (CoeffT1(a, l) * TMAT2DPart(k, c))
                    else
                        NewTMAT = NewTMAT + (CoeffT1(a, l) * TMAT2DPart(k, 2 * c))
                    end if
                end do
                if (tStoreSpinOrbs) then
                    if (NoDumpTruncs > 1) then
                        TMAT2DNew(k, j) = NewTMAT
                    else
                        TMAT2D(k, j) = (NewTMAT)
                    end if
                else
                    if (NoDumpTruncs > 1) then
                        TMAT2DNew(k, 2 * j) = NewTMAT
                        TMAT2DNew(k, 2 * j - 1) = NewTMAT
                    else
                        TMAT2D(k, 2 * j) = (NewTMAT)
                        TMAT2D(k, 2 * j - 1) = (NewTMAT)
                    end if
                end if
            end do
        end do

        deallocate(TMAT2DPart)
        call LogMemDeAlloc('RefillUMAT', TMAT2DPartTag)

        if (tROHistSingExc) call WriteSingHisttofile()

        call set_timer(RefillUMAT_Time, 30)

        if (tTurnStoreSpinOff) then
            tStoreSpinOrbs = .false.
            NoOrbs = nBasis / 2
        end if

        write(stdout, '(A, I5, A)') ' Printing the new ROFCIDUMP file for a truncation of ', NoFrozenVirt, ' orbitals.'
        if (tROFciDump .and. (NoDumpTruncs > 1)) then
            call PrintRepeatROFCIDUMP()
        else if (tROFciDUmp) then
            call PrintROFCIDUMP()
        end if

    end subroutine RefillUMATandTMAT2D