The Hartree-Fock method

• Prof. Paul-Henri Heenen, Physique Nucléaire Théorique et Physique Mathématique, Université Libre de Bruxelles

• Prof. Michel R. Godefroid, Service de Chimie quantique et Photophysique - Université Libre de Bruxelles

The Hartree-Fock (HF) method is a variational method that provides the wave function of a many-body system assumed to be in the form of a Slater determinant for fermions and of a product wave function for bosons. It treats correctly the statistics of the many-body system, antisymmetry for fermions and symmetry for bosons under the exchange of particles. The variational parameters of the method are the single-particle wave functions composing the many-body wave function. We will focus the present article on the Hartree-Fock method for fermionic systems.

Let a many-body system be described by a non-relativistic Hamiltonian composed of a one-body term, denoted \(t\), representing the kinetic energy and possibly a central potential \(V\) like the Coulomb attractive potential between the electrons and the nucleus in an atom, and a two-body interaction \(v\): \[\tag{1} \hat{H}= \sum_{\alpha \beta} t_{\alpha \beta} a^{\dagger}_{\alpha} a_{\beta} +\frac{1}{2}\sum_{\alpha \beta\gamma \delta} v_{\alpha \beta \gamma \delta} a^{\dagger}_\alpha a^{\dagger}_{\beta} a_{\delta} a_{\gamma} \; , \] where the indices \( \alpha \beta \gamma \delta \) label the single particle states in a complete orthonormal basis. The one and two-body matrix elements are denoted by \( t_{\alpha \beta}= (T_{\alpha \beta}+V_{\alpha \beta})=\langle \alpha |T +V|\beta \rangle \) and \( v_{\alpha \beta \gamma \delta}=\langle \alpha \beta|v|\gamma \delta \rangle \). Three-body terms can be also included requiring only straightforward changes in the equations. This system is characterized by its one-body density whose matrix elements \(\rho_{\beta \alpha }\) are given by\[\tag{2} \rho_{\beta \alpha } =\langle \Phi | a^{\dagger}_{\alpha } a_{\beta}|\Phi \rangle \; \] For a system of \(A\) fermions (note: This notation is natural in nuclear physics, since \(A=N+Z\), where \(N\) and \(Z\) are the numbers of neutrons and protons, respectively. In atomic physics however, \(A\) has nothing to do with the so-called mass number \(A\) and must be understood in the present context as the number of electrons \(N_e\)), a Slater determinant composed of \(A\) orbitals \({\{\vert \phi_\alpha \rangle = a^{\dagger}_{\alpha} |0 \rangle ; \alpha = 1,\ldots,A \}}\) chosen in the complete basis can be written in second quantization: \[ \tag{3} \left |\Phi \right \rangle = \prod_{\alpha=1}^A a^{\dagger}_{\alpha} \left | 0 \right \rangle \; . \] The density operator associated with a Slater determinant has a very simple form\[ \tag{4} \hat {\rho} =\sum _{\alpha=1}^A | \alpha \rangle \langle \alpha |= \sum _{\alpha=1}^\infty n_{\alpha} | \alpha \rangle \langle \alpha | \; , \] which shows that the operator \(\hat {\rho}\) has \(A\) eigenvalues equal to 1, all the others being zero. Note also the idempotence property of the density operator\[ \tag{5} \hat {\rho}^2 =\hat {\rho} \; , \] demonstrating that the operator \(\hat{\rho}\) is a projector. One can show that if a density matrix has the property that it is equal to its square, the associated wave function is a Slater determinant (Ring & Schuck, 2000).

The exact ground state corresponding to the Hamiltonian (1) can in principle be determined by using a linear combination of all possible Slater determinants constructed for A fermions with the single particle states of the basis. However, this combination should already be infinite for a two-fermion system and a truncation scheme has to be introduced. The simplest approximation is a single Slater determinant. The ground state energy might be poorly approximated if an arbitrary single-particle basis were chosen. The Hartree-Fock (HF) approximation enables one to determine the best—in the meaning of giving the lowest energy—set of single particle states that is optimized for each Hamiltonian and for a given number of particles. Since a single Slater determinant corresponds to non-interacting particles, the method is also called the independent-particle model. Note that the reduction of the many-body wave function to a single determinant inevitably breaks some symmetries of the original Hamiltonian. One cannot for instance construct a Slater determinant that is invariant by translation except for the trivial case of single-particle states that are plane waves. The HF single-particle states are obtained by a linear transformation of the orthonormal basis with creation operators \(a^{\dagger}_{\alpha}\) to a new basis with operators \(b^{\dagger}_{i}\). The matrix element of the Hamiltonian for a Slater determinant composed of orbitals labeled by an index \(i\) can be calculated explicitly in terms of \( \rho \) (Ring & Schuck, 2000) and (Blaizot & Ripka, 1986). \[ \tag{6} E^{HF} [\rho] = \left \langle \Phi | \hat{H} | \Phi \right \rangle = \sum_{ij} \left ( T_{ij} + V_{ij} \right ) \rho_{ji} + \frac{1}{2} \sum_{ijkl} \rho_{ki} \bar{v}_{ijkl} \rho_{lj} \] where \(\bar{v}_{ijkl}\) is the antisymmetrized 2-body matrix element \((v_{ijkl}-v_{ijlk})\). The HF wave function is the Slater determinant for which the expectation value of the Hamiltonian is stationary with respect to variations in the single particle wave functions \(\delta\phi_i\) \[ \tag{7} \delta E = \delta \left \langle \Phi \vert \hat{H} \vert \Phi \right \rangle = 0 \; . \] subject to orthonormality constraints \(\langle \phi_i \vert \phi_j \rangle = \delta_{ij} \) that are imposed to keep a simple expression of the expectation value (Froese Fischer, 1977). Varying the energy functional (6), with respect to the single particle wave functions and introducing a matrix of Lagrange multipliers \(\epsilon_{i,j}\) associated with the orthonormality constraints, one obtains the \(A\) Hartree-Fock integro-differential equations in the configuration space (Ring & Schuck, 2000), \[ \tag{8} \begin{array}{lcl} \epsilon_k \phi_k \left( \vec{r} \right) & = & \left( -\frac{\hbar^2}{2m}\Delta + V \left( \vec{r} \right) \right) \phi_k \left( \vec{r} \right) + \left( \int d^3 \vec{r}^{\prime} v( \vec{r} - \vec{r}^{\prime} ) \sum_{j=1}^A | \phi_j \left( \vec{r}^{\prime} \right) |^2 \right) \phi_k \left( \vec{r} \right) \\ & - & \sum_{j=1}^A \left( \int d^3 \vec{r}^{\prime} v (\vec{r} - \vec{r}^{\prime} ) \phi_j \left( \vec{r}^{\prime} \right)^* \phi_k \left( \vec{r^{\prime}} \right) \right) \phi_j (\vec{r}) - \sum_{j=1, j\neq k}^{A} \epsilon_{k,j} \phi_j (\vec{r}) \; . \end{array} \] The diagonal element of the matrix of Lagrange multipliers, \(\epsilon_{k,k}\) has been rewritten \(\epsilon_{k}\) for simplicity. Each of these equations has a form similar to a Schrödinger equation for the single-particle states. The second term, on the right-hand side, called the Hartree potential, is the average potential: \[\tag{9} U(\vec{r}) = \int d^3 \vec{r}^{\prime} v \left( \vec{r} - \vec{r}^{\prime} \right) \sum_{j=1}^A \left | \phi_j \left( \vec{r}^{\prime} \right ) \right |^2 \] which has the simple interpretation of the potential generated by the density distribution of the particles. The third term of equation (8) called the Fock term is the exchange potential. Both terms define the mean field. The last term arises from the orthogonality constraints. The HF equations constitute a set of coupled integro-differential equations and are therefore not trivial to solve. These equations form a self-consistent problem in the sense that the wave functions determine the mean-field which in turn determines the wave functions. In practice, these equations are solved using an iterative procedure. The HF Hamiltonian is a one-body operator: \[ \tag{10} \hat{h}^{HF} = \sum_{lk} h_{lk}^{HF} a^{\dagger}_{l} a_{k} \] that is defined by equation (8). Since there is an equivalence between a Slater determinant and a density matrix that is a projector, the HF hamiltonian can equivalently be determined by a variation of the energy with respect to the matrix elements of \(\rho\): \[ \tag{11} h_{lk}^{HF} = \frac{\partial E^{HF}[\rho]}{\partial \rho_{kl}} = t_{lk}+ \sum_{ij=1}^A (\bar{v}_{likj}\rho_{ji}) \] From the symmetry of the interaction, the HF operator is hermitian. One can easily show that the Hartree-Fock Hamiltonian and the density operator commute and can be diagonalized simultaneously. However, this requires that one can define an orthonormal set of eigenvectors of \(h ^{HF}\), which is possible only if the non-diagonal elements of the matrix of Lagrange multipliers \(\epsilon\) vanish. One can then find a set of orbitals that diagonalizes the Hartree-Fock Hamiltonian : \[ \tag{12} h^{HF} | i \rangle = \epsilon_i | i \rangle \; , \] defining a single-particle basis with corresponding single-particle energies \(\{\epsilon_i \}\): \[\tag{13} h_{ji}^{HF}=\epsilon_i \delta_{ji} \; . \] Equations (12) are sometimes referred to as the canonical form of the Hartree-Fock equations. In the HF approximation, the many-body Slater determinant is built on the wave functions corresponding to the \(A\) lowest eigenvalues of \(\hat{h}^ {HF} \) to which the eigenvalues \(1\) of the density operator are attributed. These states are the occupied \((o)\) or hole \((h)\) states while all the others corresponding to the eigenvalues \(0\) of the density matrix are virtual \((v)\) (also called unoccupied) or particle \((p)\) states. The diagonal form of the density matrix implies that the matrix elements of the HF Hamiltonian vanish between virtual and occupied states: \[\tag{14} \langle v \vert h^{HF} \vert o \rangle = 0 \; . \] Combining (14), with the completeness relation, \(\sum_i \vert i \rangle \langle i \vert = 1\), it is easy to show that the HF operator only produces occupied orbitals when acting on an occupied orbital: \[\tag{15} h^{HF} \vert o \rangle = \sum_{o'} \vert o' \rangle \langle o' \vert h^{HF} \vert o \rangle \; . \] The energy given by eq.(6) can be rewritten in the HF basis. Taking into account that the density matrix is diagonal in this basis, with eigenvalues equal to \(1\) and \(0\), one has: \[ \tag{16} E^{HF} = \sum_{i=1}^A t_{ii} + \frac{1}{2} \sum_{i,j=1}^A \bar{v}_{ij,ij} \] One can also use the single particle energies: \[\tag{17} \epsilon_i = t_{ii} + \sum_{j=1}^A \left ( \bar{v}_{ijij} \right ) \] to rewrite the total energy in the form: \[ \tag{18} E^{HF} = \sum_{i=1}^A \epsilon_{i} - \frac{1}{2} \sum_{i,j=1}^A \bar{v}_{ij,ij} \] One sees that the total energy is not equal to the sum of the single-particle energies. Indeed, these energies include a term generated by the two-body interaction of a given particle with all the others. When the single-particle energies are added, these interactions are counted twice, leading to the second term in eqn. (18). The two forms of the energy calculated according to eqns. (16) and (18), are equal only for the solution of the HF equations and this equality constitutes a very stringent test of convergence. The two-body interaction depends for some systems on the local density of particles. This is often the case in nuclear physics, where this density dependence can be justified by the elimination of the very repulsive core of most bare nucleon-nucleon interaction (for a recent account on this problem and the clarification of some misconceptions, see Ref. (Bogner & Furnstahl & Schwenk). The Skyrme forces, today more often labeled Skyrme energy density functionals, are a popular example of such density-dependent interactions (Bender & Heenen & Reinhard, 2003). In this case, an extra term appears in the HF hamiltonian defined in equation (11) that is given by \(\frac{1}{2} \sum_{ijpq}\langle pq \vert \frac{\partial v}{\partial \rho_{kl}} \vert ij \rangle \rho_{ip}\rho_{jq}\). The HF equations (8) have been formulated without any assumptions on the symmetry properties of the single-particle wave functions. They are the basic equations of the so-called unrestricted Hartree-Fock (UHF) formalism, a terminology that is not used in nuclear physics but is common in atomic and molecular physics (Shavitt & Bartlett, 2009). The only variational quantity of the HF method is the total energy. There is therefore a priori no reason to impose to the Hartree-Fock wave function the symmetries of the exact Hamiltonian if breaking symmetries lowers the HF total energy. A state of broken symmetries does not carry the quantum numbers of the eigenstates of the Hamiltonian (Blaizot & Ripka, 1985, Ring & Schuck, 2000). The interest of symmetry breaking is that it allows to incorporate many-body correlations without loosing the simple independent-particles picture. It has led in nuclear physics to the very powerful concept of deformed nuclei that enables to describe in an economic way many experimental data (Hamamoto & Mottelson, 2011). Symmetries of the exact state are usually imposed during the energy optimization in atomic (Froese Fischer, 1977) and molecular (Helgaker, Jørgensen and Olsen, 2000) HF calculations. In this formalism known as the restricted Hartree-Fock (RHF) theory, the total non-relativistic wave function used in the self-consistent field variational process is an eigenfunction of the total \(\bf{S}^2\) and projected \(S_z\) spins, of the total \(\bf{L}^2\) and projected \(L_z\) angular momenta for atoms, and is required to transform as an irreducible representation (IR) of the appropriate point group for molecules. This symmetry-adaptation is accomplished by:

1. requiring the atomic or molecular \(m_s = \pm 1/2\) spin-orbitals to have the same spatial parts,
2. the spatial orbitals to transform according to the IR of SO(3) for atoms and of the molecular point group for molecules,
3. to write if necessary the wave function as a Configuration State Function (ie. a symmetry-adapted linear combination of Slater determinants), rather than a single Slater determinant

By assuming a spherical symmetry, the atomic RHF equations are rewritten in spherical coordinates and reduced to a system of coupled radial equations, one for each \(nl\)-subshell, independently of the \(2(2l+1)\) projection states \((m_l m_s)\) of the orbital angular momentum \(l\) and spin \(s\) of the particles. This RHF has its natural extension in the relativistic Dirac-Fock scheme (Grant, 2007) for which the orbital energies are degenerate for each value of \(j\) as a function of \(m\). The Hartree-Fock wave function has some interesting properties, as first illustrated by Koopmans' theorem. Making the assumption that the mean-field is unchanged by the addition or the removal of a single fermion to a system with an even number of particles, the wave-function of the odd system is given by: \[ \tilde{\Phi}_{o} = b_o | \Phi ^{HF} \rangle \] for the removal of a particle. The only change \(\delta \rho_o\) in the density matrix is then the removal of the contribution coming from the occupied state \(o\). The energy of the odd system is then: \[\tag{19} E^{A-1}_o = E[\rho -\delta \rho_o] \] \[\tag{20} E^{A-1}_o = E^A - \sum_{ij}h_{ij} (\delta \rho_o)_{ij} +\frac{1}{2} \sum_{ijkl} \frac{\delta ^2 E}{\delta \rho_{ij}\delta \rho_{kl}} \delta (\rho_o)_{ij}\delta (\rho_o)_{kl} \] Expressed in the HF basis, the second term of this equation is the HF single particle energy of the orbital $o$ and the third term vanishes because of antisymmetry. Note that extra term appears for density-dependent interactions. One concludes from this expression that the energy required to remove a particle from the state \(o\) is equal to \(\epsilon_o\), ie. the corresponding HF single-particle energy. The derivation of Koopmans' theorem supposes that one can neglect the rearrangement of the mean field due to the removal or the addition of a particle. Its validity is limited, especially in nuclear physics but it has been widely used for estimating ionization potentials or electron affinities in atomic and molecular systems. Brillouin's theorem is another property of the Hartree-Fock solution that can explain its relatively high quality in atomic and molecular electronic structure calculations. It implies that there is no first-order mixing of the Hartree-Fock solution with states obtained from single substitutions of the type \[ \tilde{\Phi}_{o \rightarrow v} = b^{\dagger}_v b_o | \Phi ^{HF} \rangle \; , \] demonstrating that the HF method partially takes into account the particle-hole part of the interaction. One can show that the off-diagonal Lagrange multipliers \(\epsilon_{ij}\) can be eliminated when the fermion system corresponds to a closed shell system in which all j-shells are fully occupied. In practice, this case is the only one in nuclear physics for which the HF approximation is valid (see below). In atomic and molecular physics, the symmetry restrictions imposed in the RHF formalism to describe open-shell configurations for electronic systems with point-group symmetry have some unavoidable consequences (Nesbet, 2005). With this respect, it is worthwhile to stress that the off-diagonal Lagrange multipliers cannot always be eliminated, as it has been shown in the very first applications of the HF approximation for the Lithium atomic ground state, yet described by a single Slater determinant (Slater, 1960). However, replacing the differential equations by systems of non-linear equations and generalized eigenvalue problems, projection operators can be applied to solve this issue (Froese Fischer, 2011). The use of symmetry-adapted N-electron functions in the RHF approximation also requires some adaptation of Brillouin's and Koopmans' theorems (Froese Fischer, 1977). The Hartree-Fock method has some intrinsic limitations, mostly due to the basic assumption that particles move independently in some average potential produced by all the particles. In nuclear physics, many-body correlation effects are captured through the symmetry breaking of the UHF solution in the sense that a Slater determinant built on deformed single-particle wave functions can be expanded as a function of Slater determinants in a spherical basis and includes very excited spherical particle-hole excitations. This approach has the advantage of preserving the simple picture of independent particles. The rotational symmetry can then be restored afterwards by projecting the deformed mean-field wave function on good total momentum, leading to a projection after variation method. For atoms and molecules however, the variation after projection method is more often adopted and electron correlation is strictly defined as the difference between the RHF energy and the exact eigenvalue of the non-relativistic Schrödinger equation for the many-body system. Note however that some electron correlation is implicitly included in the HF approximation through the use of antisymmetric wave function that prevents two electrons having the same spin projection to occupy the same space region. This effect is known as the Fermi correlation.

The Hartree-Fock wave function can be obtained by solving numerically the radial integro-differential HF equations (Froese Fischer, 1977). One can also use algebraic approaches. They consist in expanding the one-fermion orbitals in some suitable set of analytical basis functions, reaching the numerical Hartree-Fock limit provided the basis is large enough and complete with respect to square-integrable functions. For molecular polyatomic systems, the molecular orbitals are expanded in a set of atomic orbitals whose expansion coefficients are the variational parameters. In this scheme, the Hartree-Fock equations are reformulated as the Roothaan-Hall self-consistent field equations (Helgaker, Jørgensen and Olsen, 2000).

In addition to Fermi correlation recognized in the HF model, the instantaneous correlation in electron motions due to their mutual repulsion that is neglected in the average field picture of the HF model is often crucial to get an accurate description of the electronic phenomena in atoms, molecules, and solids. For atomic and molecular systems, the HF solution often accounts for more then 99% of the total electronic energy, and has an overlap of 95% with wave functions obtained using more sophisticated methods. The 5% left in the latter case are difficult to capture and can dramatically affect other observables than energies such as isotope shifts, hyperfine structures and transition probabilities (Froese Fischer, Brage and Jönsson, 1997).

The Hartree-Fock method is often applied to get an approximate description of excited states that are not the lowest of their symmetry (Froese Fischer, 1977). In this case, one determines a stationary energy through the selection of the orbital solution having the desired number of radial nodes (Froese Fischer, Brage and Jönsson, 1997). Although the accuracy required today in the description of nuclei, atoms and molecules cannot be obtained by the HF method, the latter is often used as the starting point of several more elaborated methods (see for instance (Shavitt and Bartlett, 2009; Ring and Schuck, 2000)). Amongst the latter, let us quote in particular:

1. the configuration interaction (CI) method that diagonalizes the total Hamiltonian in the configuration space enlarged by including multiple excitations \(( b^{\dagger}_{v'} b^{\dagger}_v \ldots b_{o'} b_o \ldots ) | \Phi ^{HF} \rangle \) from the single Slater determinant or reference configuration state function, without further orbital optimization,
2. the multiconfiguration approach (MCHF/MCSCF) (Froese Fischer, 1977, Helgaker, Jørgensen and Olsen, 2000, Grant, 2007), optimizing both the orbital set and the CI mixing coefficients, particularly adapted to describe the electron correlation due to near degeneracies that make the single-configuration HF model inappropriate,
3. the generator coordinate method, mainly used in nuclear physics, has some similarity with the two previous methods. The total Hamiltonian is diagonalized in a basis defined by mean-field wave functions generated by a constraint on a collective variable, like the quadrupole moment of a nucleus (Ring and Schuck, 2000 and Bender & Heenen & Reinhard, 2003). This method allows also to restore symmetries that could be broken at the mean-field level of approximation.
4. the many-body perturbation (MBPT) and Coupled-Cluster theories (Shavitt and Bartlett, 2009), also existing in the relativistic version (Lindgren, 2011),
5. the HF-BCS and Hartee-Fock-Bogoliubov methods that enable to introduce the correlations due to supraconductivity,
6. the random phase approximation (RPA) that has been developed to describe collective excitations that are a coherent superposition of single particle excitations,
7. the density functional theory that leads to equations that have a form similar to mean-field equations but that incorporates correlations not included in an HF scheme. It has been shown to be very successful in condensed matter physics but also in computational chemistry,
8. the non relativistic framework presented in this article has been extended to incorporate the relativistic Dirac-Fock method for nuclei, atoms and molecules and the Relativistic Mean Field method for nuclei that extend the non relativistic framework of the present article.

References

• Ring, P and Schuck, P (2000). The Nuclear Many-Body Problem Springer-Verlag, Berlin, Heidlberg. ISBN 3-540-09820-8.

• Blaizot, J-P and Ripka, G (1986). Quantum Theory of Finite Systems MIT Press, Cambridge. ISBN 0262022141.

• Froese Fischer, C (1977). The Hartree-Fock Method for Atoms: A Numerical Approach John Wiley and Sons, New York. ISBN 047125990X.

• Bogner, SK; Furnstahl, RJ and Schwenk, A. (2010). Progress in Particle and Nuclear Physics 65 94: .

• Bender, M; Heenen, P-H and Reinhard, P-G (2003). Review of Modern Physics 75 121: .

• Shavitt, I and Bartlett, R-J (2009). Many-Body Methods in Chemistry and Physics: MBPT and Coupled-Cluster Theory Cambridge Molecular Science Cambridge University Press, Cambridge. ISBN 9780521818322

• Hamamoto(2011). Shape deformations in atomic nuclei Scholarpedia : 13. http://www.scholarpedia.org/article/Shape_deformations_in_atomic_nuclei

• Helgaker, T; Jørgensen, P and Olsen, J (2000). Molecular Electronic-Structure Theory Wiley, Chichester. ISBN 0471967556.

• Grant, I-P (2007). Relativistic Quantum Theory of Atoms and Molecules Springer Series on Atomic, Optical, and Plasma Physics, Volume 40 Springer, New York. ISBN 978-0-387-34671-7

• Nesbet, R-K (2005). Variational Principles and Methods in Theoretical Physics and Chemistry Cambridge University Press, Cambridge. ISBN 0521675758

• Slater, J-C (1960). Quantum Theory of Atomic Structure, Vol. 2 Cambridge University Press, New York.

• Froese Fischer , C (2011). Computational Physics Communications 182 1315: .

• Froese Fischer, C; Brage, T and Jönsson, P (1997). Computational Atomic Structure: An MCHF Approach, Institute of Physics, Bristol. ISBN 0750304669

• Lindgren, I (2011). Relativistic Many-Body Theory: A New Field-Theoretical Approach Springer Series on Atomic, Optical and Plasma Physics, Volume 63 Springer, Berlin. ISBN 1441983082

Sponsored by: Dr. Nicolas Alamanos, CEA, IRFU; F-91191 Gif-sur-Yvette, France, France

Reviewed by: Dr. Thomas DUGUET, CEA/DSM/IRFU/SPhN Orme des Merisiers F-91191 Gif sur Yvette Cedex

Reviewed by: Dr. Jean-Paul Blaizot, IPhT- Saclay

Accepted on: 2012-09-16 15:26:12 GMT