The Hartree-Fock equations derived

Units are SI, with energies in eV, and no metric signature enters: everything here is non-relativistic. Write x = (r, σ), ∫dx for the spin sum and the d³r integral, h = −(ħ²/2m)∇² − Ze²/(4πε₀r) and v(x, y) = e²/(4πε₀|r − r′|). Topic 10.05 derived the energy of the determinant of occupied spin-orbitals φᵢ, quoted here: E = Σᵢ⟨i|h|i⟩ + ½Σᵢⱼ(⟨ij|v|ij⟩ − ⟨ij|v|ji⟩). That formula assumed ⟨i|j⟩ = δᵢⱼ, so the orbitals cannot be varied freely: N² conditions must hold throughout. Attach a multiplier to each, L = E − Σᵢⱼ λⱼᵢ(⟨i|j⟩ − δᵢⱼ). Since ⟨i|j⟩* = ⟨j|i⟩ the conditions come in conjugate pairs, and L is real for every variation only if λⱼᵢ* = λᵢⱼ: the multiplier matrix is Hermitian. A complex orbital has two independent real parts, or equivalently φₖ and φₖ*; varying φₖ* is enough, because the equation from φₖ is its complex conjugate.

Change φₖ* by δφₖ* and keep first order. The one-body sum gives ⟨δk|h|k⟩. In the pair sum φₖ* appears twice, as i = k and as j = k, and because v(x, y) = v(y, x) the two pieces are equal: they cancel the ½ and leave Σⱼ(⟨δk j|v|k j⟩ − ⟨δk j|v|j k⟩). The constraint term gives Σⱼ λⱼₖ⟨δk|j⟩. The coefficient of the arbitrary δφₖ*(x) must vanish, so hφₖ + Jφₖ − Kφₖ = Σⱼ λⱼₖφⱼ, with J(x) = Σⱼ∫dy |φⱼ(y)|² v(x, y) multiplying and (Kf)(x) = Σⱼ φⱼ(x)∫dy φⱼ*(y) v(x, y) f(y) integrating. Dimensions: |φⱼ|² is in m⁻³, v in joules and d³r′ in m³, so J is an energy, and K returns an energy times f. One operator, F = h + J − K, now serves every orbital, and it is Hermitian. That is the derivation: no mean field was assumed, it came out of the variation.

Fφₖ = Σⱼ λⱼₖφⱼ is not yet an eigenvalue equation. Replace the occupied set by φ′ₖ = Σⱼ φⱼUⱼₖ with U unitary. The determinant is linear in each column, so Φ′ = (det U)Φ with |det U| = 1: the same state up to a phase. J and K depend on the orbitals only through Σⱼ φⱼ(x)φⱼ*(y), which U leaves alone, so F is unchanged, while the multipliers become U†λU. A Hermitian matrix is diagonalised by a unitary one, and that choice gives the canonical equations Fφᵢ = εᵢφᵢ, derived. What the derivation fixes is the occupied subspace, not the orbitals in it. Two spin-up orbitals of a symmetric molecule, one even and one odd, rotated by 45° become one orbital on the left and one on the right, with the same density, energy and many-electron wavefunction. The canonical basis is singled out only later, when its εᵢ acquire a meaning in Koopmans' theorem, topic 10.07.

J is local: the Coulomb potential of the whole charge cloud of all N occupied orbitals, the acted-on one included, blind to spin. K is non-local, since (Kf)(r) needs f at every r′, and its spin sum gives zero unless φⱼ and f share a spin. On an occupied orbital the j = k terms are the same integral, ∫dy |φₖ(y)|² v(x, y) φₖ(x), so (Jₖ − Kₖ)φₖ = 0, derived: each electron feels the other N − 1 and none of itself. For one electron that leaves F = h, so the one-particle limit is the Schrödinger equation itself and Hartree-Fock is exact there, as the Easy example checks against hydrogen. Drop K and that lone electron meets its own cloud: on hydrogen's 1s orbital the quoted self-energy (5/8)e²/(4πε₀a₀) = 17.0 eV outweighs its 13.6 eV binding, and the electron is pushed off. Hartree's 1928 equations removed the self-term by giving each orbital the potential of the others only. The operator then differs from orbital to orbital, so orthogonality is not built in, and the exchange that antisymmetry brings is missing.

F depends on the orbitals it is meant to produce, so Fφ = εφ is nonlinear. It is solved by iteration: guess orbitals, build J and K, diagonalise F in a basis (Roothaan's matrix form), occupy the N lowest, rebuild, and stop when the density settles. Without interaction one pass suffices; with it, convergence is not guaranteed. Take two electrons on two sites with hopping t and on-site repulsion U, sharing one restricted orbital with site weights p₁, p₂ and polarisation x = p₁ − p₂. Each spin feels only its partner's density on its site (K cancels its own), so F has diagonal Up₁ and Up₂ and off-diagonal −t, and its lower eigenvector gives x′ = −Ux/√(U²x² + 4t²), derived. At the solution x = 0 the error is multiplied by −U/(2t) per build: above U = 2t it grows, flipping sign, until the iteration alternates between two states, neither a solution. Cancès and Le Bris proved in 2000 that, under mild conditions, plain iteration either converges or ends in such a cycle (quoted). Mixing a fraction α of each new density makes the factor 1 − α(1 + U/(2t)), which must lie between −1 and 1.

Self-consistency means δL = 0 and nothing more. The equations have many solutions: occupying an excited orbital instead of the lowest, keeping a symmetry that the lowest determinant breaks (restricted against unrestricted, topic 10.09), or sitting at a saddle where the energy still falls along some orbital rotation. Only the second variation, a stability analysis, shows a converged answer to be a local minimum, and nothing shows it to be the global one. Even the best determinant leaves out correlation. Topic 10.08 computes closed-shell helium at −77.87 eV, quoted here, against the measured −79.0 eV. The Hartree-Fock limit for H⁻, also quoted, is −13.28 eV, 0.33 eV above a hydrogen atom and a free electron, so the method cannot bind an ion that is measured bound by 0.754 eV. The equations are non-relativistic as well: fine structure, the Breit interaction and the Lamb shift lie outside them.