From Schrödinger to Ionisation Potential in Minutes: A Hands-On Mini Lesson with Hartree-Fock

Most people meet Hartree–Fock (HF) in theory class and never use it to compute a property they care about. In this short, practical lesson you’ll go from the electronic Schrödinger equation to an estimated ionisation potential (IP) using only the HOMO energy.

Under the Born–Oppenheimer approximation, we treat the nuclei as fixed at a given geometry, so the nuclear–nuclear repulsion energy is constant, and we solve only the electronic problem. Finally, we add the nuclear repulsion energy to the expectation value of the electronic Hamiltonian to obtain the total energy.

Using a Slater determinant as the wave function, the Pauli exclusion principle is automatically fulfilled. Minimising the energy with respect to variations of the spin orbitals in the wave function leads to the HF equation, which is a pseudo-eigenvalue equation because the Fock operator F depends on its own eigenfunctions.

According to Koopmans’ theorem, the orbital energy of an occupied orbital can be associated with the ionisation potential - i.e., the energy required to remove an electron from the corresponding orbital. The first IP is approximately the negative of the HOMO energy.

Below is a minimal example in Python using PySCF to calculate the first IP.

# Minimal demo (educational).
import numpy as np
from pyscf import gto, scf

# Water geometry (Å). Keep it simple for demo; full tutorial uses cc-pVDZ and ∆SCF.
mol = gto.Mole()
mol.atom = '''
O  0.000000  0.000000  0.000000
H  0.000000 -0.757000  0.587000
H  0.000000  0.757000  0.587000
'''
mol.basis = 'sto-3g'
mol.spin = 0
mol.charge = 0
mol.build()

mf = scf.RHF(mol).run()
ehf = mf.e_tot
mo_energies = mf.mo_energy   # in Hartree
homo_idx = mol.nelectron//2 - 1 # assuming a closed-shell system
eps_homo_ha = mo_energies[homo_idx]
eps_homo_ev = eps_homo_ha * 27.211386245988  # Ha -> eV
ip_koop_ev = -eps_homo_ev

print(f"RHF total energy (Ha): {ehf:.6f}")
print(f"HOMO energy (Ha / eV): {eps_homo_ha:.6f} / {eps_homo_ev:.3f}")
print(f"Koopmans IP estimate (eV): {ip_koop_ev:.3f}")

If you’re interested in getting:

• a more detailed explanation of HF theory,
• a clear discussion of Koopmans’ theorem and its limitations, and
• a ready-to-run Jupyter Notebook with an installation guide,

please check out my **Hartree–Fock Quickstart** tutorial on Gumroad.