These calculations are all carried out in the framework of generalized Kohn-Sham theory ([5], and Chapter 7 of Ref. [3]). Details of the LDA, LSD, RLDA, and ScRLDA formalisms are described in independent sections below; here we focus on matters that are generic to all these approximations.

We utilize the central-field approximation, with conventional labelling of principal and angular momentum quantum numbers of electronic orbitals. We limit our calculations to the ground-state electronic configurations of the first 92 neutral atoms and singly-charged cations of the periodic table; the specific configurations used are described below. In cases of partially filled electronic subshells, fractional occupancies are assigned to orbitals with different azimuthal quantum number, *m*, to accomplish a spherical averaging of the charge distribution. In the case of RLDA, this extends to averaging over subshells with the same orbital angular momentum *ℓ* but different values of total angular momentum *j*. Thus, for example, if there were 2 electrons in a *p* shell, we assign an electron population of 4/3 to the *p*_{3/2} states and 2/3 to the *p*_{1/2} states.

The results presented here derive from four codes that were written independently. These were found to give results of good mutual consistency, provided that the numerical approximations within each code were varied until a very high degree of convergence was obtained within each code. The original authors of the four codes are, in alphabetical order

- Sverre Froyen
- Donald Hamann
- Eric Shirley
- Ilia Tupitsyn and Svetlana Kotochigova

These codes were used by us with the permission of their authors. However, all codes required modification to obtain numerical convergence to the target accuracy of 1 microHartree in total energy. These modifications were not subject to review by the original authors. Some of these codes circulate relatively freely within the electronic structure community, so we must caution readers that a given available version of any one of these codes need not yield results identical to those presented here. Our purpose in this study was to use robust tested tools to accomplish a specific task, not to provide a relative ranking of various codes. For this reason, we do not give details of individual code performance beyond what is necessary to describe the uncertainties in the results, and we refer to each code by a numerical label between 1 and 4, chosen arbitrarily.

The local-density approximation (LDA) requires that the exchange-correlation potential be given as a function of the electron density at a given point in space. For this part of our study, we use the form of the exchange-correlation potential given by Vosko, Wilk, and Nusair (VWN) [4]. The form is a fit to the Ceperley-Alder electron gas study [6]. The VWN functional reproduces the random-phase-approximation (RPA) results for a uniform electron gas in the high-density limit, it reproduces the spin-stiffness constant calculated in the RPA in the paramagnetic limit of a uniform electron gas, and it is uniformly differentiable as a function of the electron gas parameter, *r _{s}*. It is also in standard use, or available as an option, in many electronic structure codes, and thereby provides a convenient reference potential for checking the accuracy of numerical calculations.

We now summarize the form of the VWN functional. The exchange-correlation energy per electron is separated into two parts, an exchange term and a correlation term.

In the RLDA and ScRLDA calculations, we use the relativistic corrections to the energy-density functional proposed by MacDonald and Vosko [7].

Suitable choice of a radial grid is key to obtaining accurate numerical solutions of the integro-differential equations of density-functional theory. The codes in our test suite make different choices for the radial grid. Two codes make perhaps the simplest choice, an exponentially increasing grid

$$r_n = r_{\rm min} \left({{r_{\rm max}}\over{r_{\rm min}}}\right)^{n/N} $$

(Eq. 22)

with three parameters: the minimum radius, *r*_{min}, the maximum radius, *r*_{max}, and the number of intervals, *N*. The application of the exponential grid to the atomic Schrödinger equation has been discussed by Desclaux [8]. For one code we used *N* = 15788, *r*_{min} = 1/(160 *Z*), and *r*_{max} = 50. (All distances are in atomic units.) Another code used *N* ≤ 8000, *r*_{min} = 10^{-6}/*Z*, and *r*_{max} = 800 *Z*^{-1/2}; in this case, the energies were extrapolated to *n* → ∞ using an *N*^{-2} or *N*^{-4} dependence, depending on the quantity in question.

Another code chooses a grid which is nearly linear near the origin, and exponentially increasing at large *r*,

$$r_n = a ({\rm e}^{b (n-1)} - 1) ~ ,$$

(Eq. 23)

which again is determined by three parameters, *a*, *b*, and *N*. This grid includes the origin explicitly as *r*_{0}. In this case, we took *a* = 4.34 10^{-6}/*Z*, *b* = 0.002 304, and *r*_{max} = 50, leading to *N* = 7058 for H, increasing to *N* = 9021 for U, and to *r*_{1} = 10^{-7} for H, decreasing to 1.1 10^{-9} for U.

A fourth code uses a change of variable technique:

$$\rho = \ln r. $$

(Eq. 24)

A uniform grid is taken in the transformed variable from *ρ*(*r*_{min}) to *ρ*(*r*_{max}) where the parameters are taken to be *r*_{min} = 0.01 e^{-4}/*Z*, for atomic number *Z*, and *r*_{max} = 50. The number of points increased from *N* = 2113 for H to *N* = 2837 for U. The density of points chosen in the latter two codes, linear near the origin and exponentially increasing at large *r*, is similar to that suggested from theoretical considerations [9].

The calculations used the 1986 CODATA recommended value for α^{-1}, which was 137.0359895(61), where the digits in parentheses are the one-standard-deviation uncertainty in the last digits of the given value. The value of α^{-1} has changed in more recent compilations of fundamental constants; see the Fundamental Physical Constants homepage, where values are kept current. In test calculations using the 2006 updated values, the total energy of U (i.e., Uranium) changes by about 461 microhartrees, and the 1s Kohn-Sham eigenvalue by 108 microhartrees in the RLDA approximation, with the shift in the direction of smaller binding energies. This is consistent with the recommended value for α^{-1} having become slightly larger, implying that the recommended value for alpha has become slightly smaller, thereby reducing relativistic effects that contribute to stronger overall binding.

Created October 15, 2015, Updated June 2, 2021