GetTMatEl Function

    elemental function GetTMatEl(i, j) result(ret)

        ! Return the one electron integral <i|h|j>
        !
        ! In: i,j - Spin orbitals

        integer, intent(in) :: i, j
        HElement_t(dp) :: ret
#ifdef CMPLX_
        HElement_t(dp) :: t
#endif

        if (tCPMDSymTMat) then
            ! TMat is Hermitian, rather than symmetric.
            ! Only the upper diagonal of each symmetry block is stored
            if (j >= i) then
                ret = TMatSym(TMatInd(i, j))
            else
                ! Work around a bug in gfortran's parser: it doesn't like
                ! doing conjg(TMatSym).
#ifdef CMPLX_
                t = TMatSym(TmatInd(i, j))
                ret = conjg(t)
#else
                ret = TMatSym(TmatInd(i, j))
#endif
            end if
        else
            if (tOneElecDiag) then
                if (i /= j) then
                    ret = 0.0_dp
                else
                    ret = TMat2D(i, 1)
                end if
            else
                ret = TMat2D(i, j)
            end if
        end if
    end function GetTMatEl