Classically in Hartree-Fock implementations we define the wave-function in term of a basis set:

$$\phi_i = \sum_{\alpha} C_{\alpha i}\chi_\alpha$$

If I understand correctly the problem is reduced to find a set of these coefficients which minimise the total energy. Again, if my understanding is correct these coefficients are initially given in the basis set. If you look at STO-3G for example, there are always two numbers: the coefficient and the exponent.

What I am confused about is: In the algorithm we have to get a "first guess" of the density matrix. The density matrix is defined as:

$$D_{uv} = \sum_{i} C_{ui}C_{vi}$$

Why not just take the one given in the basis set as an initial guess?

  • 1
    $\begingroup$ No, the initial guesses for the basis set coefficients are not part of the basis set. Gaussian basis sets usually include two parameters per Gaussian primitive. One is the Gaussian exponent and the other is the expansion coefficient for the current basis function. This expansion coefficient is fixed, and is used to combine multiple Gaussian primitives into a single basis function. $\endgroup$
    – Hayden S
    Jun 19, 2022 at 21:57
  • $\begingroup$ Ok that's what I understood after I saw that the matrix C is a square matrix (which doesn't make sense with my original understanding). The expansion coefficients (written "C" in many books...) for each gaussian primitive function are not the same as the one in the first equation of the question, which are called basis set coefficients... $\endgroup$
    – Okano
    Jun 19, 2022 at 22:30

2 Answers 2


The question is confusing two things: the linear combination of atomic orbitals (LCAO) scheme where the molecular orbitals are expanded as a sum of the atomic orbitals \begin{equation} \psi_{i\sigma} ({\bf r}) = \sum_{\alpha} C_{\alpha i} \chi_\alpha ({\bf r}) \label{eq:lcao} \end{equation} and the actual definition for the atomic basis functions, which in the case of contracted Gaussian basis functions is \begin{equation} \chi_\alpha ({\bf r}) = r^{l_\alpha} \left[\sum_{k} d_{k \alpha} \exp(-\zeta_kr^2)\right] Y_{l_\alpha}^{m_\alpha}(\hat{\bf r}) \end{equation} where I have assumed a pure spherical basis for ease of notation. The LCAO scheme is not limited to Gaussian basis sets; the linear expansion ansatz is used in every basis set method, regardless of the nature of the basis functions (e.g. Gaussian or Slater-type orbital, numerical atomic orbital, plane waves, finite elements, etc). I have recently written two open access review articles on these topics, see Int. J. Quantum Chem. 119, e25968 (2019) for discussion on various kinds of basis functions and Molecules 25, 1218 (2020) for the self-consistent field algorithm.

It is true that the contracted functions in many basis sets such as the cc-pVXZ basis sets represent occupied atomic orbitals in Hartree-Fock calculations. However, this is not true for all basis sets. As programs tend to aim for generality, they typically aim to implement initial guesses that are not basis set dependent, such as the superposition of atomic densities - which is equivalent to occupying the contracted functions in certain basis sets but also works for basis sets where the individual basis functions do not look like occupied orbitals. Refer to J. Chem. Theory Comput. 15, 1593 (2019) for further discussion on initial guesses.

I would like to note that unlike stated in the comments, the matrix ${\bf C}$ is not square in general; this is unfortunately still a very common misconception in computational chemistry. There are often linear dependencies in the basis set (e.g. a basis function on one atom can be also well-described by the sum of basis functions on the neighboring atoms), and removing them removes some columns from ${\bf C}$. Stereotypical cases of this happening are when you use diffuse functions that are needed to describe long-range interactions and weakly bound electrons, or when two nuclei are close together. However, even such pathological linear dependencies can be cured and Coulomb barriers to nuclear fusion modeled when suitable numerical algorithms are chosen to pick a suitable sub-basis (see the links). The choice of the linearly independent basis is also discussed in the overview on SCF theory linked above.


Actually, it does sometimes. If you remember that HF scheme is to solve the pseudo-eigenvalue-eigenvector problem, the equation is FC=SCE where F is the Fock matrix and the C is the molecular orbital coefficient as shown in your equation. For a crude assumption, the initial Fock matrix can be the summation of the kinetic T and external potential V which is often referred core-hamiltonian H (all operators here are represented as the atomic orbital basis). For H, the only thing you need is the atomic orbital basis set information, and the initial C can be guessed with this H matrix. In your comment "matrix C is a square matrix", I think this is because the program deleted the unoccupied coefficient information for the sake of calculating the density matrix D. When making D, the occupation number is also multiplied there and, since the unoccupied occupation number is 0, coefficient matrix often treated like that.


You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .