Orbital angular momentum eigenfunctions for fast and numerically stable evaluations of closed-form pseudopotential matrix elements.

The computation of s-type Gaussian pseudopotential matrix elements involving low powers of the distance from the pseudopotential center using Gaussian orbitals can be reduced to familiar integrals. They may be directly expressed as either simple three-center overlap integrals for even powers of the radial distance from the pseudopotential center or related to the three-center nuclear integrals of a Gaussian charge distribution for odd powers. Orbital angular momentum about each atom is added to these integrals by solid-harmonic differentiation with respect to its center. The solid-harmonic addition theorem allows all the integrals to be factored into products of invariant one-dimensional integrals involving the Gaussian exponents and angular factors that contain the azimuthal quantum numbers but are independent of all Gaussian exponents. Precomputing the angular factors allow looping over all Gaussian exponents about the three centers. The fact that solid harmonics are eigenstates of angular momentum removes the singularities seen in previous treatments of pseudopotential matrix elements.