Auxiliary-Basis Expansion and Projection of Density
An auxiliary basis represents a scalar electron density by a finite set of coefficients attached to the atoms. The auxiliary functions considered here are centered on atoms and are distinct from the atomic-orbital basis. In quantum chemistry, such functions are commonly introduced for density fitting and related resolution-of-identity methods, but the term auxiliary basis does not tie the expansion to one particular algorithm.
This chapter first defines the density coefficients and the real-space \(L^2\) error between two densities. The auxiliary-basis overlap matrix then provides the exact coefficient-space expression for this error and determines the least-squares projection. Finite systems are treated first. For periodic systems, cell-shift-indexed overlap blocks are folded into the unit-cell overlap matrix used by the periodic density coefficients. Alternative projection metrics are discussed separately.
The target is the nonnegative scalar electron density \(n(\boldsymbol r)\) defined in Kohn–Sham Density Functional Theory, not the signed electric charge density. Hartree atomic units are used unless stated otherwise.
Auxiliary Basis
For a finite system, let
label one auxiliary-basis function. Here \(i\) is its atomic center, \(n\) distinguishes radial functions with the same angular momentum, and \((\ell,m)\) are angular indices. A typical function centered at \(\boldsymbol R_i\) is
The spherical-harmonic basis follows the same convention-neutral notation as Localized Atomic-Orbital Basis. The radial functions need not be analytic, normalized, or mutually orthogonal. Functions on different atoms overlap as well, so the complete multicenter auxiliary basis is generally nonorthogonal and may contain near-linear dependencies.
A density represented in this basis is
where the scalars \(d_\lambda\) are the density coefficients, collected in the column vector \(\mathsf d=[d_\lambda]\). Their values depend on the normalization and mixing of the auxiliary-basis functions. They are not point samples of the density, atomic populations, or observables by themselves. The reconstructed field \(n_{\mathsf d}(\boldsymbol r)\) is the physical object represented by the coefficients.
The term basis is used in the usual quantum-chemical sense and does not assert exact completeness or strict linear independence. When the functions are redundant, the represented density may be unique even though its coefficient vector is not.
L2 Density Error and the Auxiliary-Basis Overlap Matrix
For any two square-integrable scalar densities \(n(\boldsymbol r)\) and \(n'(\boldsymbol r)\) on a domain \(\mathcal D\), define their squared \(L^2\) density error by
The unsquared \(L^2\) density error is the square root of this quantity. For a finite system, \(\mathcal D=\mathbb R^3\) will be left implicit below.
Define the elements of the auxiliary-basis overlap matrix \(\mathsf S^\chi=[S^\chi_{\lambda\lambda'}]\) by
This is the auxiliary-basis counterpart of the AO overlap matrix \(\mathsf S=[S_{\mu\nu}]\). The matrix \(\mathsf S^\chi\) is Hermitian and positive semidefinite, and it is positive definite when the retained auxiliary-basis functions are linearly independent.
Consider two densities represented in the same auxiliary basis,
With \(\Delta\mathsf d\equiv\mathsf d-\mathsf d'\), their squared \(L^2\) density error is
For a real auxiliary basis and real coefficients, the Hermitian conjugate reduces to a transpose. Thus the Euclidean coefficient error \(\lVert\mathsf d-\mathsf d'\rVert_2^2\) is generally not the \(L^2\) error between the represented densities. The auxiliary-basis overlap matrix supplies the required metric. If the basis is linearly dependent, distinct coefficient vectors can represent the same density; their difference then lies in the null space of \(\mathsf S^\chi\) and has zero \(L^2\) density error.
Cholesky Form of the Density Error
When the retained auxiliary basis is linearly independent, \(\mathsf S^\chi\) is positive definite and admits the upper-triangular Cholesky factorization
The density error can therefore be evaluated as
This norm form is manifestly nonnegative and avoids a small negative result that can arise from roundoff when the quadratic form is evaluated directly. Cholesky factorization does not, however, turn a singular or indefinite overlap matrix into a positive-definite one. If the auxiliary basis has linear dependencies, the null modes must first be removed, or a rank-revealing factorization must be used on the retained subspace. Adding a diagonal shift \(\delta\mathsf I\) merely to make Cholesky succeed changes the density error by \(\delta\lVert\mathsf d-\mathsf d'\rVert_2^2\) and therefore changes the metric.
Finite-System L2 Projection
Given a target density \(n(\boldsymbol r)\), its \(L^2\) projection onto the auxiliary basis chooses \(\mathsf d\) to minimize
Define the density basis integrals
Let \(\mathsf b=[b_\lambda]\) denote the basis integrals of the density.
The objective can then be written as
Stationarity gives the normal equations
This distinction is essential in a nonorthogonal basis. The entries \(b_\lambda\) are covariant overlaps of the density with individual basis functions; they equal the expansion coefficients only when the auxiliary basis is orthonormal. In general, the normal equations define the coefficient vector and must be solved as a linear system. If the basis is redundant, the coefficients require a specified treatment of the null space even when the projected density itself is unique.
At the solution, the fitted residual is orthogonal to every retained auxiliary-basis function:
Solving the Normal Equations
The overlap matrix is a Gram matrix. For any coefficient vector \(\mathsf x\),
Its rank is therefore the dimension of the function space spanned by the retained auxiliary functions. Linear independence makes \(\mathsf S^\chi\) full rank and positive definite; an exact linear dependency makes it singular. A matrix may nevertheless be full rank but ill-conditioned when one positive eigenvalue is much smaller than the largest one. Full rank guarantees a unique coefficient vector, whereas good conditioning determines whether that vector is numerically stable.
When \(\mathsf S^\chi\) is positive definite, write its upper-triangular Cholesky factorization as
The normal equations are then solved without forming an inverse:
The intermediate vector satisfies \(\mathsf c=\mathsf U^\chi\mathsf d=(\mathsf U^{\chi\dagger})^{-1}\mathsf b\) and is the coefficient vector in the Euclidean coordinates of the Cholesky density-error norm defined above. Factoring an \(N\times N\) dense overlap costs \(O(N^3)\); after the factor is available, solving for each additional right-hand side costs \(O(N^2)\).
For a rank-deficient positive-semidefinite overlap, let
In exact arithmetic the Moore–Penrose pseudoinverse retains every positive eigenvalue and discards the exact zero modes. In finite precision, a stated relative threshold \(\tau\) defines the retained set \(I_\tau=\{a:\lambda_a>\tau\lambda_{\max}\}\), the effective numerical rank \(r_\tau=|I_\tau|\), and
This is the minimum-Euclidean-norm least-squares coefficient vector after the discarded directions have been treated as a numerical null space. Exact null directions represent no density and may be removed without changing the projected field. Discarding a small but genuinely positive direction instead changes the retained projection space, so the threshold is part of the numerical definition rather than a harmless implementation detail.
For a full-rank positive-definite overlap, the spectral condition number is
It is infinite for an exactly singular matrix; for a pseudoinverse calculation, the corresponding retained-space condition number uses the smallest retained eigenvalue. A large condition number means that small errors in the overlap or basis integrals can be strongly amplified in \(\mathsf d\), even when Cholesky factorization still succeeds. Near-null directions can consequently produce large mutually cancelling coefficients while changing the represented density much less. A numerical audit should therefore state the eigenvalue threshold, effective rank or nullity, \(\lambda_{\min}/\lambda_{\max}\) or \(\kappa_2\), and any negative eigenvalues separately. An exact overlap Gram matrix has no negative eigenvalues, so a significant negative value indicates numerical error rather than ordinary rank deficiency.
The same solution rules apply to the periodic equations below after replacing \((\mathsf S^\chi,\mathsf b)\) by the folded unit-cell quantities \((\widetilde{\mathsf S}^\chi,\widetilde{\mathsf b})\). The individual cell-shift blocks are not positive-semidefinite matrices and must be folded before either Cholesky factorization or rank analysis.
Projection from an LCAO Density Matrix
For the LCAO convention of Localized Atomic-Orbital Basis, the scalar electron density of a finite system is
Only the diagonal spin blocks enter because the scalar density is the spin trace of the local spin-density matrix. Define the three-center overlap
Substitution into the density-projection definition gives
The density coefficients are then obtained from
The AO overlap matrix \(\mathbf S\) does not appear in this contraction. The density matrix \(\mathbf D\) defined in Localized Atomic-Orbital Basis directly reconstructs the real-space density. By contrast, \(\mathbf S\mathbf D\mathbf S\) contains the covariant matrix elements of the density operator and must not be substituted for \(\mathbf D\) in the pair-density expansion.
Electron Number
Let
Let \(\mathsf q=[q_\lambda]\) denote the corresponding column vector.
The electron number represented by the fitted density is
An unconstrained \(L^2\) projection does not in general guarantee \(N_{\mathrm e}[n_{\mathsf d}]=N_{\mathrm e}\). Exact electron-number conservation can instead be imposed with a Lagrange multiplier. For a real auxiliary basis, the constrained normal equations are
This is a different projection: enforcing the integral redistributes error among the auxiliary-basis channels, and its minimized squared \(L^2\) error can only equal or exceed the unconstrained minimum. Electron-number conservation should therefore be specified as an additional constraint rather than assumed from the meaning of density projection.
Periodic Density Expansion
Now let the system be periodic with lattice \(\mathcal L\), unit cell \(\Omega\), and atoms \(i\) in that cell. The lattice and cell conventions are those of Periodic Systems. The auxiliary-basis index \(\lambda=(i,n,\ell,m)\) now labels a function associated with an atom in the unit cell. Its translate into cell \(\boldsymbol R\) is
Translation symmetry assigns the same coefficient \(d_\lambda\) to every translated copy. The periodic represented density is therefore
Define the auxiliary-basis function folded into the unit cell by
Throughout this chapter, a tilde marks a quantity obtained by periodic folding. The folded function is lattice periodic, and the density expansion becomes
Each auxiliary-basis function is assumed to have a finite radial cutoff and to be exactly zero beyond it. Consequently, only finitely many translated copies contribute at any point in the unit cell.
Cell-Shift-Indexed and Folded Overlap Matrices
In direct correspondence with the periodic AO overlap matrix, define the cell-shift-indexed auxiliary-basis overlap blocks by displacing the right function relative to the left:
Let \(\mathsf S^\chi(\boldsymbol R)\) collect the elements \(S^\chi_{\lambda\lambda'}(\boldsymbol R)\). It has the same displacement convention as the AO overlap block \(\mathsf S(\boldsymbol R)\). These localized, cell-shift-indexed blocks are the overlap data directly evaluated and typically stored by a periodic localized-basis calculation. Their collection is the unfolded auxiliary-basis overlap. It obeys
The overlap matrix acting on the coefficient vector of one unit cell is instead the folded auxiliary-basis overlap matrix \(\widetilde{\mathsf S}^\chi= [\widetilde S^\chi_{\lambda\lambda'}]\), whose elements are
Expanding both folded functions and using \(\boldsymbol R=\boldsymbol R_2-\boldsymbol R_1\) gives
In the second line, the translated integration cells generated by the sum over \(\boldsymbol R_1\) tile all space for each fixed relative translation \(\boldsymbol R\). The folded matrix is therefore Hermitian even though an individual \(\mathsf S^\chi(\boldsymbol R)\) block need not be.
Equivalently, if a Bloch overlap matrix is introduced with the same phase convention as the AO basis,
then the folded overlap is its \(\Gamma\)-point value:
L2 Density Error per Cell
A nonzero periodic density is not square integrable over all space. The physically relevant \(L^2\) density error is therefore defined per reference cell:
For two periodic densities represented by coefficient vectors \(\mathsf d\) and \(\mathsf d'\), this error is
Using the unfolded overlap blocks, the same result is
Thus the real-space error between periodic represented densities can be computed either with the folded matrix or directly from the lattice sum of the cell-shift-indexed overlap blocks.
The individual blocks \(\mathsf S^\chi(\boldsymbol R)\) are generally neither Hermitian nor positive definite, so they cannot be Cholesky-factorized independently. The lattice sum must be formed first. If the folded density basis is linearly independent, define
The per-unit-cell density error then becomes
Equivalently, \(\widetilde{\mathsf U}^\chi= \operatorname{chol}_{\mathrm{upper}}[\mathsf S^\chi(\boldsymbol k=\boldsymbol0)]\): the cell-periodic density expansion occupies only the \(\Gamma\) sector. A more general coefficient field with nonzero wave-vector components would instead require a separate factorization of each Hermitian \(\mathsf S^\chi(\boldsymbol k)\). In numerical evaluation, the \(\mathsf S^\chi(-\boldsymbol R)= \mathsf S^\chi(\boldsymbol R)^\dagger\) pairs should be combined consistently, or the folded matrix should be Hermitian-symmetrized before its Cholesky factorization. The same null-space caveat as in the finite system applies.
Periodic L2 Projection
The periodic \(L^2\) projection of a target density \(n(\boldsymbol r)\) minimizes
Define the periodic density basis integrals by
Let \(\widetilde{\mathsf b}=[\widetilde b_\lambda]\) denote the folded basis integrals of the density.
The objective is
so the coefficients per unit cell satisfy
The folded and unfolded density basis integrals are equivalent:
Periodicity of \(n(\boldsymbol r)\) turns the sum over translated copies of \(\chi_\lambda\) in one cell into the all-space integral of one reference copy. The all-space expression is finite because \(\chi_{\lambda\boldsymbol 0}\) has finite support. It does not treat the periodic density itself as an all-space square-integrable function. Only the auxiliary-basis functions are folded; the already-periodic scalar density is not folded a second time.
Electron Number per Cell
The represented electron number per cell is
where
Let \(\widetilde{\mathsf q}=[\widetilde q_\lambda]\) denote the folded electron-number vector.
The same optional electron-number constraint used for a finite system applies with \(\widetilde{\mathsf q}\) and the electron count per cell.
Projection from a Periodic LCAO Density Matrix
The periodic density-matrix blocks in Localized Atomic-Orbital Basis use \(\boldsymbol R\) for the cell of the unstarred row orbital relative to the starred column orbital. Taking the spin trace of the density reconstructed there gives
Using the unfolded density-projection expression, define the localized three-center integrals
Then
Here \(\boldsymbol R_0\) locates the starred orbital relative to the reference auxiliary-basis function, while \(\boldsymbol R\) is the displacement already carried by the density-matrix block. Localization or numerical screening reduces the formally infinite sums to the triples with overlapping spatial support. The final coefficients are obtained by solving
The same projection can be evaluated directly from the Bloch density matrices \(\mathbf D(\boldsymbol k)\). Define
At the continuum level,
The real-space and \(\boldsymbol k\)-space expressions agree when \(\mathbf D(\boldsymbol R)\) and \(\mathbf D(\boldsymbol k)\) follow the continuous Brillouin-zone transform and use the same orbital convention.
Metrics Other Than Real-Space Overlap
The preceding formulas use the ordinary real-space overlap because it gives the orthogonal projection that minimizes the squared \(L^2\) density error. More generally, a positive Hermitian kernel \(K\) defines
The corresponding overlap-matrix elements and projection-vector components are
Let \(\mathsf M^{\chi,K}=[M^{\chi,K}_{\lambda\lambda'}]\) and \(\mathsf b^K=[b_\lambda^K]\). The projection satisfies
The real-space overlap metric is obtained from \(K(\boldsymbol r,\boldsymbol r')=\delta(\boldsymbol r-\boldsymbol r')\), for which \(\mathsf M^{\chi,K}=\mathsf S^\chi\). The Coulomb metric uses
and emphasizes errors in electrostatic quantities rather than the local squared density error. In a periodic system, a Coulomb metric additionally requires a precise treatment of the long-range kernel, neutrality, and the zero-wave-vector component. A truncated Coulomb kernel is one possible definition. The metric must therefore accompany the coefficients as part of their kernel metadata; periodic Coulomb fitting is not defined by the overlap-matrix folding formulas alone.
References
- X. Ren et al., “Resolution-of-identity approach to Hartree–Fock, hybrid density functionals, RPA, MP2 and GW with numeric atom-centered orbital basis functions,” New J. Phys. 14, 053020 (2012), doi:10.1088/1367-2630/14/5/053020.
- K. R. Briling, A. Fabrizio, and C. Corminboeuf, “Impact of quantum-chemical metrics on the machine learning prediction of electron density,” J. Chem. Phys. 155, 024107 (2021), doi:10.1063/5.0055393.
- A. M. Lewis et al., “Learning electron densities in the condensed phase,” J. Chem. Theory Comput. 17, 7203 (2021), doi:10.1021/acs.jctc.1c00576.