Can we relate approximate Hartree-Fock orbitals to true solutions in the basis-set limit? Assume that we've &quot;solved&quot; the RHF Hartree-Fock equations in some finite basis, that is, we have arrived at a self-consistent set of coefficients $c_{\mu{}i}$ such that each spatial molecular orbital is expanded as $\chi_i = \sum_{\mu=1}^K c_{\mu{}i}\varphi_i$. By construction, this set of orbitals is orthonormal, and forms a collection of eigenfunctions for the Fock operators in the given basis, so $f_i |\chi_i\rangle = \varepsilon_i |\chi_i\rangle$ where $\varepsilon_i$ is the $i$th orbital energy.
Question: what, if anything, can we say about these orbitals in terms of their relationship to the &quot;true&quot; orbitals in the basis set limit, denoted $\chi_i^\mathrm{BS}$? Here, I mean those orbitals that would be obtained if the calculation were to be extended to an (infinite) complete basis set. (I realise that the assumed existence of a single unique set of true orbitals may not be justified...)
For example, is it possible to estimate, bound, or otherwise describe the overlap $\langle \chi_i | \chi_i^\mathrm{BS} \rangle$ between an approximate solution and the true solution? How about matrix elements such as $\langle \chi_i | \mathcal{H} | \chi_i^\mathrm{BS} \rangle$, where $\mathcal{H}$ is the electronic Hamiltonian, or indeed any self-adjoint operator?
I suspect/hope that some kind of a statement towards the first part might be possible, given that the Roothan-Hall equations are effectively a collection of Galerkin equations. One of the standard hand-waving statements about solutions obtained using Galerkin methods is that they are in some sense projections of the true solution onto finite-dimensional subspaces, and have residual vectors that are &quot;orthogonal&quot; to that subspace(for a potentially complicated definition of &quot;orthogonal&quot;).
(I'm also interested in answers that would apply in any of the HF settings, i.e. UHF, ROHF, etc.)
In quantum chemistry numerical experimentation is the norm and rigorous mathematical proof the rare exception. As a result, not many rigorous results are know, not even for the good-old Hartree-Fock method. For a review of mathematical results I recommend you have a look at these references:

Claude Le Bris, Computational chemistry from the perspective of numerical analysis, Acta Numerica (2005), pp. 363–444
Pablo Echenique and J. L. Alonso, A mathematical and computational review of Hartree–Fock SCF methods in quantum chemistry, Molecular Physics, Vol. 105, Nos. 23–24, 10 December–20 December 2007, 3057–3098

In particular, section 3.3 of the first reference is relevant, but don't expect to find mathematically strong results. For example, here's a quote:

The choice of an AO basis for solving
a given problem mostly relies upon some practical know-how. The lack of
rigorous understanding is a pity, because the output of the calculations
(typically some molecular property) might be very sensitive to the choice
of the basis set. The only available measures of the quality of the basis
set are obtained, in the chemistry literature, by choosing test cases, i.e.,
reference systems, where the solution of the exact Schrödinger equation may
be computed, mostly through numerical computations and, when possible,
with the help of an analytic calculation.

From a practical point of view I don't see why occupied Hartree-Fock orbitals expanded in the usual Gaussian basis sets should not converge (in the integral norm sense, ie leaving aside the lack of cusp at the nucleus) to the exact solutions of the Hartree-Fock equations. I have little idea as to what happens to the virtual orbitals.