Dirac equation for positronium with Breit interaction
Daniel Baye and Ergash M. Tursunov
Phys. Rev. A 112, 022801 (2025) - Published 1 August, 2025
The two-body Dirac equation for positronium is studied with the Coulomb and Breit interactions. When adding an annihilation correction, its energies agree with experiment up to fourth order in powers of the fine structure constant if the Breit term is treated perturbatively. We analyze why this is not true if the Dirac equation is solved numerically. Symmetries of the 16-component Dirac spinors are first detailed. For each symmetry, coupled differential equations verified by their radial components are established in a self-adjoint form convenient for numerical calculations. Nondifferential relations between these components are derived. The coupled equations are first solved in the Coulomb case with the Lagrange-mesh method. With only four mesh points, accurate binding energies are obtained which match the best published results except for the and states. A detailed study of the radial spinor components of the ground states of para- and orthopositronium explains why the convergence of their energies stops at an accuracy of about atomic units. A curious property of the Lagrange-mesh method is emphasized. Analytical approximations of these radial components are derived and compared with numerical solutions. The obtained Coulomb-Dirac spinors provide an original variant of first-order perturbation treatment of the Breit operator. The Coulomb-Breit-Dirac coupled equations are then solved under the numerical conditions tested with the Coulomb case. The differences between these numerical energies and the correct perturbation results at fourth order are analyzed. The unphysical numerical values of energies are shown to be due to couplings with the negative continua.

