TF-RSE-IonCF v3.2: A Unified Closed-Form Model for the Total Electron Binding Energy of All Neutral Atoms and Ions (Z = 1–92, q = 0 to Z–1)

1. Introduction

The total electron binding energy (TBE) — the energy required to remove all electrons from an atom or ion in its ground state — is a fundamental quantity in atomic physics, plasma modeling, astrophysics and nuclear engineering. In hot plasmas and astrophysical environments matter is highly ionized, so the ionic TBE (the sum of all successive ionization potentials of the remaining electrons) is required for every charge state, not merely for the neutral atom. Ab initio methods (Hartree–Fock, density functional theory, relativistic many-body perturbation theory) deliver high accuracy at a cost of seconds to minutes per ion, which is prohibitive for workflows that evaluate millions of charge states, such as the energy budgets of collisional-radiative plasma models, nucleosynthesis networks and radiation-damage cascades.

The model presented here unifies the neutral and ionic problems in a single closed form. Its neutral-anchor layer is TF-RSE v2, a phenomenological extension of the Thomas–Fermi (TF) statistical model that predicts the TBE of all 92 neutral atoms (Z = 1–92) with a mean relative error of 0.40%. TF-RSE v2 retains the TF power-law skeleton and adds three corrections: a relativistic regularization factor R(Z) of Dirac form, a trigonometric shell-oscillation factor S(Z), and a piecewise light-element correction L(Z). By itself it contains no charge-state dependence: it returns a single value per element, the neutral-atom TBE.

Extending such a model to ions is nontrivial. A naive ionic extrapolation — multiplying the neutral value by a smooth electron-fraction factor F(N/Z) — fails for two reasons. First, the relativistic correction of a highly charged ion differs from that of the neutral atom: as outer (weakly relativistic) electrons are removed, the average relativistic weight of the remaining electrons rises, an effect that a Z-only factor R(Z) cannot represent. Second, the neutral shell factor S(Z) is blind to the successive shell closures encountered as electrons are stripped: the successive ionization potentials display a pronounced sawtooth, jumping by up to an order of magnitude at each closed shell, and a smooth electron-fraction factor averages these jumps away. The naive extrapolation was measured at a mean error exceeding 10% over all charge states, with errors above 90% for deep inner-shell ionization of light elements.

TF-RSE-IonCF (Thomas–Fermi Relativistic-Shell-Extended Ionic Closed-Form model) resolves both problems by construction. Rather than perturbing the neutral value, the ionic TBE is decomposed into shell-resolved contributions anchored at four reference limits that are known to high accuracy: the point-nucleus Dirac energy at N = 1 (finite-nuclear-size and QED corrections neglected; +0.35% from the NIST value at Z = 92); a screened helium-like energy at N = 2; a screened neon-like energy at N = 10; and the TF-RSE v2 neutral value at N = Z. Power-law blends interpolate between the anchors. Two physical deficiencies are then closed by two targeted upgrades: (Upgrade-3) an effective nuclear charge Z_eff(Z, N) is introduced into the relativistic factor; and (Upgrade-2) the shell factor S(Z) is replaced by an electron-number-dependent sawtooth factor S(Z, N). A final upgrade (Upgrade-1/4) adds explicit shell-closure bump corrections and a global refit over five isoelectronic element sequences, turning the model from a proof of concept into a complete data product covering the whole periodic table. The result is a single formula system whose neutral limit (N = Z) is by construction the v2 model and whose hydrogenic limit (N = 1) is the point-nucleus Dirac value.

A separate issue, addressed in v3.1 and strengthened in v3.2, is physical consistency of the evaluated sequence itself. The closed-form blend reproduces the four anchors by construction and is accurate to 0.26% in the grid interior, but for 74 of the 92 elements the multiplicative shell corrections produce a slight overshoot above the neutral value — or a small oscillation — within a few charge states of the neutral boundary, so that a finite difference of the raw sequence would yield a small number of negative or non-monotonic ionization energies. Version 3.1 repaired these defects with a three-tier safeguard (hard cap; increasing-trend boundary-layer reconstruction; minimal-touch projection). Independent evaluation, including the referee’s analysis, subsequently exposed a residual defect of v3.1: the boundary-layer scan triggered only on non-monotonicity, so a near-neutral link could survive with a monotone but quantitatively unphysical implied ionization energy — most severely for the heaviest actinides, where the raw first ionization energy of uranium evaluates to +2588 eV (experiment 6.19 eV), and altogether for 19 elements whose implied first ionization energies fall outside any physically plausible window. Version 3.2 closes this defect by extending the Tier-2 trigger from a pure monotonicity test to a physical-plausibility window: a boundary link k (counted inward from the neutral anchor, k = 1 being the first ionization energy) is reconstructed whenever its implied ionization energy falls outside 0.5 eV ≤ IP_k ≤ 4·IP_base(period)·1.8^(k−1). The window upper bound exceeds every experimental first ionization energy in the periodic table (maximum 24.6 eV, He) by a factor of at least three in every period, so no physically real valence link is ever removed; the lower bound additionally removes degenerate near-zero links (the v3.1 Na N = 10 link of ≈1×10⁻¹¹ eV). After the extension, all 92 elements have strictly positive, monotone and trend-correct boundary ionization energies, at a modest and code-reproducible footprint (646 corrected points, 15.1%; Section 3.1.1).

The remainder of the paper is organized as follows. Section 2 presents the theoretical framework. Section 3 gives the complete formula system, the fitted parameters and the standard operating procedure. Section 4 reports validation on all 4278 (Z, N) reference points, the two-stage fitting strategy, the ablation study, two anti-leakage tests, and an uncertainty model. Section 5 discusses applications and limitations, and Section 6 concludes.

2. Theoretical Framework

2.1. Neutral-Anchor Layer: TF-RSE v2

The neutral-atom layer of the unified model is TF-RSE v2, whose accuracy is preserved exactly: at N = Z the ionic model returns the v2 neutral value by construction. TF-RSE v2 reads

E_TBE(Z) = a Z^b E_h R(Z) S(Z) L(Z). (1)

Here E_h = 13.605693 eV is the Rydberg energy (not the Hartree energy, which is 27.211386 eV), a = 2.0224 and b = 2.1414 are the TF amplitude and effective exponent, R(Z) the relativistic factor, S(Z) the shell factor and L(Z) the piecewise light-element correction. The relativistic factor of Dirac form is

R(Z) = [1 − (Z/137.036)²]^(−1/2), (2)

which grows from 1.000 at Z = 1 to 1.349 at Z = 92. The shell factor is

S(Z) = 1 + c sin(d Z^(1/3)) + e cos(d Z^(1/3)), (3)

with c = −0.2404, d = 0.9737, e = −0.2254 restoring periodic structure. The light-element correction is piecewise:

L(Z) = Α + Α e^(−γZ) [Z ≤ 8]; 1 [Z = 9]; 1.0149 [Z = 10]; 1 + A e^(−B(Z−10)) [10 < Z < 18]; 1 [Z ≥ 18], (4)

with Α = 0.9973, c = −0.9925, γ = 1.3240, A = 0.01955 and B = 0.06539. Table 1 collects the eleven parameters of the neutral layer.

Fitted on all 92 neutral atoms (Z = 1–92) this layer achieves a mean relative error of 0.40% (median 0.39%, maximum 1.53% at uranium), with 89 of 92 elements below 1% and 100% below 5%. Leave-one-out cross-validation yields 0.70% (0.49% excluding hydrogen), and a genuine out-of-sample test on the superheavy elements Z = 93–103, never used in fitting, yields 3.6%. All eleven parameters are kept frozen throughout the ionic extension: the neutral layer enters the model only as the exact boundary condition E_v2(Z).

2.2. Shell-Resolved Ionic Decomposition

The central structural idea is to anchor the ionic TBE at the charge states where it is known to high accuracy, and to blend between these anchors. Four anchors are used:

(i)

Hydrogenic anchor, N = 1. A single electron in the field of the bare point nucleus has the Dirac binding energy [7]

E_K1(Z) = Z² E_h D(Z), D(Z) = 2[1 − sqrt(1 − (Zα)²)]/(Zα)², (5)

where α = 1/137.035999084 is the fine-structure constant and D(Z) the Dirac factor (D → 1 for small Z). This is the exact point-nucleus Dirac value; finite-nuclear-size and QED corrections are not included. It deviates from the NIST value by +0.35% at Z = 92 (132.28 keV vs. 131.82 keV) — a deviation larger than the mean error of the full model — so the N = 1 anchor is “exact” only in the point-nucleus Dirac sense, not relative to experiment.

(ii)

Helium-like anchor, N = 2. E_K2(Z) = 2 (Z − s_K)² E_h D(Z), (6)

with s_K = 0.2399 fitted. Physically s_K is the Slater screening of one 1s electron by the other.

(iii)

Neon-like anchor, N = 10. E_K10(Z) = E_K2(Z) + 8 (Z − s_L)²/4 · E_h D(Z), s_L = s_L0 + s_L1 Z, (7)

with s_L0 = 3.6555 and s_L1 = 0.0204. The Z-dependent screening accounts for the growing penetration of the L shell.

• (iv) Neutral anchor, N = Z: E = E_v2(Z), the TF-RSE v2 value. For light elements the hydrogenic anchors (i)–(iii) can exceed the neutral value; each anchor is therefore clamped from above by E_v2(Z), which guarantees a smooth, exact reduction to the neutral model.

Between the anchors the TBE is blended by power laws. For 2 < N ≤ 10 the blend runs from E_K2 to E_K10 with exponent q_L; for N > 10 it runs from E_K10 to E_v2 with exponent q_M. Because the blends vanish at the anchor points, the model reproduces the anchors at N = 1, 2, 10 (to the anchor accuracy stated above) and at N = Z (to the v2 accuracy) for every element.

2.3. Upgrade-3: Effective Charge in the Relativistic Factor

Z_eff(Z, N) = Z [1 + r_1(1 − x) + r_2(1 − x)²], x = N/Z, (8)

with r_1 = −0.8000 and r_2 = 0.8280 fitted. The relativistic correction is applied as the ratio R(Z_eff)/R(Z), which equals unity at the neutral anchor (x = 1) and therefore leaves the published neutral accuracy untouched. The ablation study (Section 4.5) shows this single upgrade reduces the heavy-ion mean error from 2.67% to 1.41%.

We state plainly, in response to the referee’s question, that Z_eff(Z, N) is a purely phenomenological calibration parameter and carries no physical interpretation as a screened charge seen by an individual electron. Its role is limited to re-weighting the Dirac-form relativistic factor: as outer, weakly relativistic electrons are removed, the average relativistic weight of the remaining electron cloud rises, and Eq. (8) parameterizes this re-weighting empirically. The formal extrapolation Z_eff → 1.028 Z as N/Z → 0 therefore has no physical meaning and does not represent a charge exceeding the bare nuclear charge; moreover it lies strictly outside the model domain, because the relativistic ratio is evaluated only in the N > 10 branch (N = Z returns the neutral anchor earlier), where Z_eff/Z < 1 holds mathematically for every element. The hard constraint Z_eff ≤ Z in the released code is defensive — a no-op within the model domain — and is retained only as protection against out-of-domain use. A first-principles derivation of the re-weighting law (8) is an open problem. We have nevertheless tested the natural alternative of constraining Z_eff ≤ Z by construction, i.e. restricting Eq. (8) to the sign-definite family Z_eff = Z(1 − q_1 u − q_2 u^m) with q_i ≥ 0 and u = 1 − N/Z, for both m = 2 and m = 3. Refitted on the same validation grid (Z ≥ 11, 10 < N < Z), the best constrained form reaches only RMS 2.30%, versus 0.245% for the retained unconstrained quadratic with the defensive cap — a factor of ≈9 degradation. The reason is that the data require a non-monotone Z_eff(u): it dips to ≈0.81 Z near u ≈ 0.5 and rises back toward Z for deeply stripped ions, a shape no sign-definite monotone family can represent (the constrained optimum drives q_2 → 0 in both cases). The unconstrained quadratic with the defensive cap is therefore retained as the best phenomenological choice; the constraint test is fully reproducible with the script Supplementary_Script_ZeffTest.py in the supplementary archive.

2.4. Upgrade-2: Electron-Number Shell Sawtooth S(Z, N)

S_sh(N) = 1 + c_sh sin(d_sh N^(1/3)) + e_sh cos(d_sh N^(1/3)), (9)

M(y) = 1 + [S_sh(N)/S_sh(Z) − 1] env(y), env(y) = [sin(π y)]^ep, (10)

with c_sh = 0.0098, d_sh = 2.0531, e_sh = 0.0024, ep = 0.6376. The envelope env(y) vanishes at both ends of each blend interval (so the anchors stay exact) and peaks in between. The ablation shows this upgrade cuts the error at closed-shell electron counts (N = 10, 18, 36, 54, 86) from 1.10% to 0.81%.

2.5. Upgrade-1/4: Shell-Closure Bumps and Global Refit

B(Z, N) = Σ_k A_k (Z/m_k)^(g_b) exp[−((N − m_k)/w_b)²], m_k = 10, 18, 36, 54, 86, (11)

applied as a factor (1 + B(1 − y)) that vanishes at the N = Z anchor. The amplitudes A_k, the width w_b and the Z-scaling exponent g_b are fitted. Finally, all 18 parameters are refitted globally against the full 4278-point data product with the neutral layer frozen.

3. Formula System and Standard Operating Procedure

3.1. Core Equations

For integer Z in [1, 92] and integer N in [1, Z], the total electron binding energy is assembled as follows. By construction E(Z, Z) = E_v2(Z) exactly for every element (the neutral anchor, including Be and B where the N ≤ 10 blend would otherwise undershoot the neutral value). For N < Z, with x = N/Z, D(Z) the Dirac factor, E_v2(Z) the TF-RSE v2 neutral value, and the anchors E_K1, E_K2, E_K10 of Section 2.2 (each clamped above by E_v2):

E(Z, N) = E_K1, N = 1. (12a)

E(Z, N) = E_K2, N = 2. (12b)

E(Z, N) = [E_K2 + (E_K10 − E_K2) y^(q_L)] (1 + B(1 − y)) M(y), y = (N − 2)/(min(10, Z) − 2), 2 < N ≤ 10. (12c)

E(Z, N) = [E_K10 + (E_v2 − E_K10) y^(q_M)] [R(Z_eff)/R(Z)] (1 + B(1 − y)) M(y), y = (N − 10)/(Z − 10), N > 10. (12d)

All quantities are defined in Section 2; the complete ionic parameter set is listed in Table 2 and the neutral-layer parameters in Table 1. The result is obtained in eV and divided by 1000 for keV.

3.1.1. Physical-Consistency Safeguard (v3.2)

The closed-form evaluation (12a–12d) reproduces the four anchors by construction and is accurate to 0.26% mean error over the grid interior. On the raw sequence, however, 74 of the 92 elements exhibit a small non-monotonic defect: for 69 of them the multiplicative shell corrections (bump, sawtooth, relativistic ratio) push one or more of the last charge states strictly above the neutral anchor E_v2 — producing a negative binding energy for the last electron upon finite differencing — and two further elements (Li and Na) sit on E_v2 to within floating-point precision at degenerate anchors (310 raw defect links in total; four anchor equalities at Mg, Al, Si, P (N = 10) are created by the cap pinning of Tier 1; the deepest raw defect lies 10 charge states from neutrality, at Hg N = 70; the deepest point reached by the Tier-3 minimal-touch projection lies 15 charge states from neutrality, at Tm N = 54 in the lanthanide region). In addition — the defect pointed out by the referee and repaired in v3.2 — the raw blend can also produce near-neutral links that are formally monotone but whose implied ionization energies are physically impossible for valence electrons: 19 elements have an implied first ionization energy outside the plausibility window defined below, most severely uranium (+2588 eV vs. the experimental 6.19 eV). Version 3.2 repairs both classes of defect with a three-tier safeguard that acts at the source and is minimal by construction:

Tier 1 — hard physical cap. For every N < Z the raw value is capped at the neutral anchor:

E_raw(Z, N) ← min(E_raw(Z, N), E_v2(Z)). (13)

This enforces the first-principles requirement that the last electron binds positively (E(Z, Z) − E(Z, Z−1) > 0), since TBE(Z, N) < TBE(Z, Z) for all N < Z. The cap touches 505 raw points.

Tier 2 — boundary-layer reconstruction (N > Z − 5) with physical-plausibility window (extended in v3.2). A backward scan from N = Z − 1 to N = Z − 4 examines every boundary link k (counted inward from the neutral anchor, so k = 1 is the first ionization energy). The link is reconstructed whenever its implied ionization energy IP_k = E(Z, N+1) − E(Z, N) falls outside the window

0.5 eV ≤ IP_k ≤ 4 · IP_base(period) · 1.8^(k−1), (14)

i.e. whenever the link is non-monotone (IP_k ≤ 0), degenerate (0 < IP_k < 0.5 eV), or monotone but physically implausible (IP_k above the upper bound). The upper bound exceeds every experimental first ionization energy in the periodic table (maximum 24.6 eV, He) by a factor of at least three in every period, so no physically real valence link is ever removed; the real shell-crossing jumps of the deeper boundary links are preserved because the bound scales with the same 1.8^(k−1) trend (e.g. the Na second ionization energy of 47.3 eV falls below the k = 2 bound of 57.6 eV). Reconstruction replaces the value by

E(Z, N) ← E(Z, N+1) − max(IP_est, 0.5 eV), IP_est = IP_base(period) · 1.8^k, (15)

where IP_base = 10, 10, 8, 7, 8, 6, 5 eV for periods 1–7 (an order-of-magnitude estimate of the last-electron ionization energy) and k counts the corrected steps from the neutral anchor. The factor 1.8 > 1 encodes the experimental stripping trend — successive ionization energies INCREASE across the valence shell as screening decreases (e.g. Cu: 7.7, 20.3, 36.8, 57.4 eV; Dy: 5.9, 11.6, 22.9, 41.2 eV). The v3.2 window extension adds 64 reconstructed points relative to v3.1 (71 ceiling triggers, 1 floor trigger — the v3.1 degenerate Na N = 10 link of ≈1×10⁻¹¹ eV, now repaired — minus 8 points released elsewhere), confined to the boundary layer of the 20 elements Z = 8–12 and Z = 78–92. The reconstructed boundary points inherit the accuracy of the neutral anchor (0.40% mean); they are intended as positive, monotone, trend-correct estimates, not as quantitative valence ionization potentials (Section 5).

Tier 3 — minimal-touch projection (N ≤ Z − 5). Residual interior oscillation links are repaired by the smallest possible monotone projection:

E(Z, N) ← E(Z, N+1) − 0.5 eV (only where E(Z, N) ≥ E(Z, N+1)). (16)

No heuristic cascade is applied beyond what monotonicity strictly requires: each corrected oscillation-link point is moved to just below its upper neighbour, so the cumulative TBE is preserved to within the local oscillation amplitude (relative adjustment ≤ 0.25% at the largest core link; the clamp-degenerate anchor replacement at P (Z = 15, N = 10), which falls in the core region by the N ≤ Z − 5 convention, is a separate category disclosed in Table D2 with a 0.69% replacement). The artificial 0.5 eV step appears only in the differential at the affected link and is flagged in Section 5.

Anchor protection. The scan (14)–(16) never modifies the hydrogenic anchor N = 1 (outside the scan range) or the neutral anchor N = Z (returned directly). The helium-like (N = 2) and neon-like (N = 10) anchors are modified only at the five clamp-degenerate points where the anchor value is pinned to E_v2 by construction and would otherwise equal the neutral anchor: Li (Z = 3, N = 2) and Mg, Al, Si, P (Z = 12–15, N = 10). At these five points the safeguard value is MORE accurate than the clamped anchor (Section 4.4 and Appendix D); their values are unchanged from v3.1. At all other elements the N = 2 and N = 10 anchors remain exactly the closed-form values. One degenerate boundary case is disclosed for completeness: at Z = 1 the hydrogenic and neutral anchors coincide (N = 1 = Z); the code returns the neutral anchor E_v2(1) = 13.608109 eV, which deviates from the Dirac value 13.605874 eV by +0.0164%. The claims of exactness at N = 1 and N = Z therefore hold exactly (in the point-nucleus sense of Section 2.2) for Z ≥ 2.

Net effect on the full 4278-point grid (Appendix D), recomputed directly from the released v3.2 code: 646 points are corrected relative to the raw closed form (15.1%); Tier 2 performs 319 boundary-layer reconstructions (247 monotonicity triggers, 1 floor trigger, 71 ceiling triggers) and Tier 3 performs 327 minimal-touch core adjustments; 86 elements require at least one correction; the largest single adjustment relative to the raw closed-form value is 7.84 keV (U, N = 88; boundary layer of the heaviest actinides, where the raw blend underestimates the near-neutral drop by several keV), while the largest core-grid adjustment remains 1.14 keV (Lu, N = 63, deep-core oscillation link); the deepest correction lies 15 charge states from neutrality (Tm, N = 54; unchanged from v3.1); zero negative ionization energies, zero non-monotonic elements, and zero out-of-window first ionization energies remain (v3.1: 19 elements out of window); the first ionization energies of all 92 elements are now strictly positive, monotone and trend-correct (Si 8.0 eV, Fe 7.0 eV, U 5.0 eV vs. experiment 8.15, 7.90, 6.19 eV); against the 4278-point reference dataset the mean relative error is 0.259% for the raw closed form and 0.260% for the safeguarded output (core grid N ≤ Z−5: 0.230% → 0.238%; on the external 264-point spot-check of Section 4.9 the core mean changes only from 0.304% to 0.321%), and all 29 fitted parameters are unchanged.

3.2. Fitted Parameters

Table 2. The 18 fitted parameters of TF-RSE-IonCF v3.2 (fitted on 4278 NIST reference points; the 11 neutral v2 parameters are kept fixed). 

Parameter	Value	Role/region
s_K	+0.2399	K-shell screening constant (N = 2 anchor)
s_L0	+3.6555	L-shell screening, constant term (N = 10 anchor)
s_L1	+0.0204	L-shell screening, Z slope
q_L	+0.7401	blend exponent, 2 < N ≤ 10
q_M	+0.6431	blend exponent, N > 10
r_1	−0.8000	Z_eff linear coefficient (Upgrade-3)
r_2	+0.8280	Z_eff quadratic coefficient (Upgrade-3)
c_sh	+0.0098	sawtooth sine amplitude (Upgrade-2)
d_sh	+2.0531	sawtooth frequency (Upgrade-2)
e_sh	+0.0024	sawtooth cosine amplitude (Upgrade-2)
e_p	+0.6376	sawtooth envelope exponent
A_10	−0.1342	closure-bump amplitude, N = 10
A_18	+0.2318	closure-bump amplitude, N = 18
A_36	−0.0455	closure-bump amplitude, N = 36
A_54	+0.2260	closure-bump amplitude, N = 54
A_86	+0.4764	closure-bump amplitude, N = 86
w_b	+26.1947	closure-bump width
g_b	+0.7139	closure-bump Z-scaling exponent

Table 2. The 18 fitted parameters of TF-RSE-IonCF v3.2 (fitted on 4278 NIST reference points; the 11 neutral v2 parameters are kept fixed). 

Parameter	Value	Role/region
s_K	+0.2399	K-shell screening constant (N = 2 anchor)
s_L0	+3.6555	L-shell screening, constant term (N = 10 anchor)
s_L1	+0.0204	L-shell screening, Z slope
q_L	+0.7401	blend exponent, 2 < N ≤ 10
q_M	+0.6431	blend exponent, N > 10
r_1	−0.8000	Z_eff linear coefficient (Upgrade-3)
r_2	+0.8280	Z_eff quadratic coefficient (Upgrade-3)
c_sh	+0.0098	sawtooth sine amplitude (Upgrade-2)
d_sh	+2.0531	sawtooth frequency (Upgrade-2)
e_sh	+0.0024	sawtooth cosine amplitude (Upgrade-2)
e_p	+0.6376	sawtooth envelope exponent
A_10	−0.1342	closure-bump amplitude, N = 10
A_18	+0.2318	closure-bump amplitude, N = 18
A_36	−0.0455	closure-bump amplitude, N = 36
A_54	+0.2260	closure-bump amplitude, N = 54
A_86	+0.4764	closure-bump amplitude, N = 86
w_b	+26.1947	closure-bump width
g_b	+0.7139	closure-bump Z-scaling exponent

3.3. Standard Operating Procedure (SOP)

Step 1: input Z (1–92) and N (1–Z); compute x = N/Z and the Dirac factor D(Z). The released reference implementation accepts Z and N directly as command-line arguments (python Supplementary_Code.py Z [N]; with N omitted the neutral atom N = Z is returned; with no arguments a self-test demonstration is run).

Step 2: compute the neutral value E_v2(Z) with the TF-RSE v2 formula (Section 2.1).

Step 3: compute the anchors E_K1, E_K2, E_K10 (Section 2.2), clamping each at E_v2.

Step 4: select the branch by N; compute the blend variable y and the power-law base energy.

Step 5: for N > 10 compute Z_eff (with Z_eff ≤ Z) and the relativistic ratio R(Z_eff)/R(Z) (Section 2.3).

Step 6: compute the sawtooth modulation M(y) (Section 2.4) and the closure bumps B(Z, N) (Section 2.5).

Step 7: assemble E(Z, N) per Section 3.1; apply the Tier-1 cap E ≤ E_v2 (N < Z).

Step 8 (safeguard): evaluate the full charge-state sequence TBE(N), N = 1..Z; apply the backward scan (Tier 2 boundary layer with the plausibility window (14), Tier 3 core) per Section 3.1.1; return TBE(N).

A worked example (Fe at N = 13) is given in Appendix B; the full reference implementation is listed in Appendix A.

3.4. Computational Complexity

The core model uses only elementary arithmetic (powers, trigonometric and exponential functions) and evaluates in about one microsecond per ion. For the 86 corrected elements the safeguard requires the full charge-state sequence, O(Z) per element (~Z × 8 µs); all 4278 charge states are produced in approximately 70 milliseconds when pre-computed and cached. The three-tier design (cap → boundary reconstruction → minimal-touch projection) is non-iterative and deterministic.

4. Validation Results

4.1. Reference Data

The reference dataset is built from the NIST Atomic Spectra Database (version 5.12): for each element Z = 1–92 the complete set of successive ionization potentials IE(q) (charge q → q+1) was extracted and summed to obtain the total binding energy of every ion, TBE(Z, N) = Σ_{q=Z−N}^{Z−1} IE(q). This yields 4278 (Z, N) reference points, covering every element and every charge state q = 0 to Z−1. The neutral-column sums agree to within 0.05% with the reference values used in the TF-RSE v2 paper (Fe: 34.619 vs. 34.613 keV; U: 761.52 vs. 761.87 keV). These are experimental reference data, not v2 model output: the v2 model predictions themselves (34.743 keV for Fe, 773.46 keV for U) differ from experiment by the v2 model residuals (0.40% mean over Z = 1–92), and the present model reproduces those v2 predictions exactly at N = Z by construction.

4.2. Two-Stage Fitting Strategy: Why the Neutral Layer Is Frozen

Freeing the five global neutral parameters (a, b, c, d, e) and refitting all 23 free parameters jointly against the 4278 points reduces the overall mean error only marginally (0.259% → 0.240%) but degrades the neutral-atom accuracy from 0.40% to 0.57% — a 43% loss — because the 92 neutral points carry only 2.2% of the total weight. The neutral layer is the exact boundary condition of the ionic decomposition; its accuracy is more valuable than a marginal global gain. The 11 neutral parameters are therefore fitted once on the 92 neutral atoms and frozen; the 18 ionic parameters are then fitted on the 4278 points with the neutral layer fixed.

4.3. Statistical Summary

Table 3. Statistical summary of TF-RSE-IonCF v3.2 over all 4278 (Z, N) points. 

Quantity	Value
Reference points (Z, N)	4278 (92 elements, all charge states q = 0..Z−1)
Data source	NIST ASD v5.12 successive ionization potentials (summed)
Neutral layer (92 atoms, q = 0)	0.401% mean (frozen v2 parameters)
Mean relative error, core closed form (all 4278 points)	0.259%
Median relative error, core closed form	0.151%
Maximum relative error, core closed form	4.63% (Ne, N = 6)
Points with error < 1%/< 5%, core closed form	97.7%/100%
Mean relative error, v3.2 safeguarded output	0.260% (median 0.161%, maximum 4.33% at F, N = 5)
Points with error < 1%/< 5%, safeguarded output	97.8%/100%
Core grid only (N ≤ Z−5), raw → safeguarded	0.230% → 0.238%
Leave-one-element-out CV (mean)	0.263% (max 4.74%)
Five-sequence blind test (87 elements)	0.266% mean
Neutral-anchor consistency (N = Z)	exact (TF-RSE v2 recovered)
Hydrogenic anchor (N = 1)	point-nucleus Dirac value (+0.35% from NIST at Z = 92)
First ionization energies (v3.2 safeguard)	strictly positive, monotone, trend-correct for all 92 elements; estimates only (Section 5)
Parameters/reference points	18/4278 (+ 11 fixed neutral)

Table 3. Statistical summary of TF-RSE-IonCF v3.2 over all 4278 (Z, N) points. 

Quantity	Value
Reference points (Z, N)	4278 (92 elements, all charge states q = 0..Z−1)
Data source	NIST ASD v5.12 successive ionization potentials (summed)
Neutral layer (92 atoms, q = 0)	0.401% mean (frozen v2 parameters)
Mean relative error, core closed form (all 4278 points)	0.259%
Median relative error, core closed form	0.151%
Maximum relative error, core closed form	4.63% (Ne, N = 6)
Points with error < 1%/< 5%, core closed form	97.7%/100%
Mean relative error, v3.2 safeguarded output	0.260% (median 0.161%, maximum 4.33% at F, N = 5)
Points with error < 1%/< 5%, safeguarded output	97.8%/100%
Core grid only (N ≤ Z−5), raw → safeguarded	0.230% → 0.238%
Leave-one-element-out CV (mean)	0.263% (max 4.74%)
Five-sequence blind test (87 elements)	0.266% mean
Neutral-anchor consistency (N = Z)	exact (TF-RSE v2 recovered)
Hydrogenic anchor (N = 1)	point-nucleus Dirac value (+0.35% from NIST at Z = 92)
First ionization energies (v3.2 safeguard)	strictly positive, monotone, trend-correct for all 92 elements; estimates only (Section 5)
Parameters/reference points	18/4278 (+ 11 fixed neutral)

All statistics above were regenerated against the reference dataset (Dataset_S1.csv) with the released v3.2 code. Relative to v3.1 the safeguard modifies 64 additional points, all within or adjacent to the boundary layer whose values are disclosed as estimates (Section 5); the core grid (N ≤ Z−5) closed form and all Table 4 model values are unchanged, and the safeguarded full-grid mean error improves slightly relative to the v3.1 safeguarded output (0.264% → 0.260%), with the maximum error reduced from 4.63% to 4.33%.

4.4. Representative Validation Across the Table

Table 4 lists representative validation for seven elements. All values are the v3.2 model outputs (after the physical-consistency safeguard); boundary-layer points (Si N = 11, and the neutral-anchor rows) reflect the v3.2 reconstruction of Section 3.1.1 and are unchanged from v3.1. Model values are given to the same precision as the reference; the error column is computed from unrounded values.

Table 4. Representative validation for seven elements across charge states (keV), v3.2 model output. The N = Z row is the neutral anchor; the N = 1 row is the hydrogenic Dirac anchor. Si N = 11 is a boundary-layer reconstruction (raw closed-form value 7.904 keV, +0.94%; unchanged from v3.1).

Z	El.	N	q	Model (keV)	NIST (keV)	Err (%)
6	C	1	+5	0.490	0.490	+0.01
6	C	2	+4	0.903	0.882	+2.40
6	C	3	+3	0.945	0.947	−0.13
6	C	5	+1	1.005	1.019	−1.38
6	C	6	0	1.030	1.030	−0.01
14	Si	1	+13	2.674	2.673	+0.02
14	Si	2	+12	5.166	5.111	+1.07
14	Si	4	+10	6.077	6.111	−0.54
14	Si	7	+7	7.062	7.167	−1.46
14	Si	11	+3	7.834	7.831	+0.05
14	Si	14	0	7.883	7.889	−0.07
26	Fe	1	+25	9.282	9.278	+0.04
26	Fe	2	+24	18.22	18.11	+0.64
26	Fe	7	+19	27.07	27.16	−0.35
26	Fe	13	+13	32.65	32.58	+0.22
26	Fe	20	+6	34.41	34.34	+0.21
26	Fe	26	0	34.74	34.62	+0.36
36	Kr	1	+35	17.95	17.94	+0.07
36	Kr	2	+34	35.42	35.23	+0.53
36	Kr	9	+27	59.92	60.33	−0.67
36	Kr	18	+18	70.99	71.10	−0.16
36	Kr	27	+9	75.41	75.13	+0.37
36	Kr	36	0	76.22	75.87	+0.46
54	Xe	1	+53	41.35	41.30	+0.11
54	Xe	2	+52	81.96	81.57	+0.48
54	Xe	14	+40	164.0	163.9	+0.07
54	Xe	27	+27	191.1	192.0	−0.47
54	Xe	41	+13	201.3	200.8	+0.26
54	Xe	54	0	202.1	202.4	−0.16
74	W	1	+73	80.91	80.76	+0.19
74	W	2	+72	160.8	159.9	+0.52
74	W	19	+55	358.0	357.1	+0.26
74	W	37	+37	417.3	416.9	+0.09
74	W	56	+18	435.6	435.7	−0.03
74	W	74	0	436.5	439.0	−0.56
92	U	1	+91	132.3	131.8	+0.35
92	U	2	+90	263.2	261.4	+0.69
92	U	23	+69	626.4	624.0	+0.40
92	U	46	+46	723.4	724.4	−0.14
92	U	69	+23	758.5	756.2	+0.30
92	U	92	0	773.5	761.5	+1.57

Table 4. Representative validation for seven elements across charge states (keV), v3.2 model output. The N = Z row is the neutral anchor; the N = 1 row is the hydrogenic Dirac anchor. Si N = 11 is a boundary-layer reconstruction (raw closed-form value 7.904 keV, +0.94%; unchanged from v3.1).

Z	El.	N	q	Model (keV)	NIST (keV)	Err (%)
6	C	1	+5	0.490	0.490	+0.01
6	C	2	+4	0.903	0.882	+2.40
6	C	3	+3	0.945	0.947	−0.13
6	C	5	+1	1.005	1.019	−1.38
6	C	6	0	1.030	1.030	−0.01
14	Si	1	+13	2.674	2.673	+0.02
14	Si	2	+12	5.166	5.111	+1.07
14	Si	4	+10	6.077	6.111	−0.54
14	Si	7	+7	7.062	7.167	−1.46
14	Si	11	+3	7.834	7.831	+0.05
14	Si	14	0	7.883	7.889	−0.07
26	Fe	1	+25	9.282	9.278	+0.04
26	Fe	2	+24	18.22	18.11	+0.64
26	Fe	7	+19	27.07	27.16	−0.35
26	Fe	13	+13	32.65	32.58	+0.22
26	Fe	20	+6	34.41	34.34	+0.21
26	Fe	26	0	34.74	34.62	+0.36
36	Kr	1	+35	17.95	17.94	+0.07
36	Kr	2	+34	35.42	35.23	+0.53
36	Kr	9	+27	59.92	60.33	−0.67
36	Kr	18	+18	70.99	71.10	−0.16
36	Kr	27	+9	75.41	75.13	+0.37
36	Kr	36	0	76.22	75.87	+0.46
54	Xe	1	+53	41.35	41.30	+0.11
54	Xe	2	+52	81.96	81.57	+0.48
54	Xe	14	+40	164.0	163.9	+0.07
54	Xe	27	+27	191.1	192.0	−0.47
54	Xe	41	+13	201.3	200.8	+0.26
54	Xe	54	0	202.1	202.4	−0.16
74	W	1	+73	80.91	80.76	+0.19
74	W	2	+72	160.8	159.9	+0.52
74	W	19	+55	358.0	357.1	+0.26
74	W	37	+37	417.3	416.9	+0.09
74	W	56	+18	435.6	435.7	−0.03
74	W	74	0	436.5	439.0	−0.56
92	U	1	+91	132.3	131.8	+0.35
92	U	2	+90	263.2	261.4	+0.69
92	U	23	+69	626.4	624.0	+0.40
92	U	46	+46	723.4	724.4	−0.14
92	U	69	+23	758.5	756.2	+0.30
92	U	92	0	773.5	761.5	+1.57

4.5. Ablation Study: Contribution of Each Upgrade

Table 5. Ablation study. Mean relative error over all points, heavy ions (Z ≥ 40, q ≥ 4), and closed-shell electron counts. 

Stage	Configuration	All (%)	Heavy-ion (%)	Closed-shell (%)
A	Backbone only (shell-resolved blend)	2.27	2.67	2.26
B	A + Upgrade-3: Z_eff in R	1.27	1.41	1.10
C	B + Upgrade-2: shell sawtooth S(Z,N)	0.72	0.65	0.81
D	C + Upgrade-1/4: closure bumps (full)	0.26	0.19	0.22

Table 5. Ablation study. Mean relative error over all points, heavy ions (Z ≥ 40, q ≥ 4), and closed-shell electron counts. 

Stage	Configuration	All (%)	Heavy-ion (%)	Closed-shell (%)
A	Backbone only (shell-resolved blend)	2.27	2.67	2.26
B	A + Upgrade-3: Z_eff in R	1.27	1.41	1.10
C	B + Upgrade-2: shell sawtooth S(Z,N)	0.72	0.65	0.81
D	C + Upgrade-1/4: closure bumps (full)	0.26	0.19	0.22

4.6. Anti-Leakage Test: Five-Sequence Training, 87-Element Blind Test

All 18 parameters were fitted using only the five element sequences Fe (Z = 26), Kr (36), Xe (54), W (74) and U (92) — 282 points — and then applied blindly to the remaining 87 elements (3996 points) without any refitting. The training sequences achieve 0.21% mean error; the 87 unseen elements achieve 0.27% mean error (maximum 8.2% at a single light-element point), essentially identical to the full fit.

4.7. Leave-One-Element-Out Cross-Validation

Removing one entire element at a time, refitting the 18 parameters on the other 91 elements, and predicting the held-out element gives a mean error of 0.263% (median 0.156%, maximum 4.74%), with 100% of the 4278 held-out predictions below 5%. The largest residuals occur for the lightest elements (Ne, F, Na, O: 2–3% mean), where the statistical treatment is least adequate.

4.8. Uncertainty Quantification

Table 6. Recommended uncertainty U by Z region (95th percentile of residuals, rounded up). 

Z region	Points	Mean err (%)	95th pct (%)	Recommended U (%)
Z = 1–10	55	1.25	3.46	3.5
Z = 11–18	116	0.74	2.55	3.0
Z = 19–36	495	0.35	0.92	1.0
Z = 37–54	819	0.24	0.65	1.0
Z = 55–71	1071	0.21	0.56	1.0
Z = 72–92	1722	0.21	0.63	1.0

(Computed on the safeguarded v3.2 output against the reference dataset.).

Table 6. Recommended uncertainty U by Z region (95th percentile of residuals, rounded up). 

Z region	Points	Mean err (%)	95th pct (%)	Recommended U (%)
Z = 1–10	55	1.25	3.46	3.5
Z = 11–18	116	0.74	2.55	3.0
Z = 19–36	495	0.35	0.92	1.0
Z = 37–54	819	0.24	0.65	1.0
Z = 55–71	1071	0.21	0.56	1.0
Z = 72–92	1722	0.21	0.63	1.0

(Computed on the safeguarded v3.2 output against the reference dataset.).

4.9. Independent Spot-Check Against NIST Successive Ionization Potentials

As an independent guard against transcription errors in the fitted-parameter table, the released code was additionally spot-checked against freshly retrieved NIST ASD v5.12 successive-ionization data for seven elements spanning the full table (Na, Mg, Cu, Zn, Nb, Dy, Re; 264 charge states). Cumulative sums of the experimental ionization potentials were compared with the model output: mean relative error 0.36% (raw core 0.36%), median 0.20%, 95% of points below 1%, 100% below 5% — consistent with the global statistics of Table 3 given the deliberate overweighting of light elements (Na, Mg) in this spot-check. Near the neutral boundary the v3.2 reconstruction gives strictly positive, monotone, trend-correct ionization energies for every element (e.g. Cu last four reconstructed IPs 7.0, 12.6, 22.7, 35.1 eV vs. experiment 7.7, 20.3, 36.8, 57.4 eV; Dy 6.0, 10.8, 19.4, 35.0 eV vs. 5.9, 11.6, 22.9, 41.2 eV; U 5.0, 9.0, 16.2, 29.2 eV vs. an experimental first IP of 6.19 eV), in contrast to the v3.0 heuristic which produced decreasing sequences there, and in contrast to the v3.1 safeguard which left monotone-but-unphysical links unrepaired (U raw first IP +2588 eV).

5. Discussion

TF-RSE-IonCF occupies a niche that, to the authors knowledge, no other closed-form model fills: it predicts the total binding energy of any ion of any element Z = 1–92 in microseconds with a mean error of about 0.26%, while remaining exactly consistent with a validated 0.40% neutral-atom model at q = 0 and with the point-nucleus Dirac value at q = Z−1 (+0.35% from experiment at Z = 92). Ab initio ionic calculations reach ~0.1% at much higher cost; semi-empirical electron-binding formulas (FRDM2012) reach only 2–3% and lack both charge-state dependence and shell structure.

Applications include: (1) cumulative binding-energy budgets in collisional-radiative plasma modeling — the TBE enters the plasma energy content and the energy released on recombination; we emphasize, as the referee correctly notes, that the level-resolved ionization potentials required for collisional-radiative rate coefficients cannot be obtained from this model and must be taken from experimental databases or ab initio calculations; (2) astrophysical opacity and nucleosynthesis calculations; (3) nuclear mass evaluation (ionic electron-binding corrections in the extraction of bare nuclear masses from mass-spectrometric measurements of highly charged ions); (4) radiation-damage and ion-beam simulations.

Physical-consistency safeguard (v3.2). The three-tier safeguard guarantees positive ionization energies and strict monotonicity of TBE(N) over the full 4278-point grid, and the v3.2 plausibility window (14) additionally guarantees that every near-neutral link is physically plausible, while the statistical accuracy of the core grid (N ≤ Z−5) is essentially untouched (mean error 0.230% raw vs. 0.238% safeguarded). We state plainly, as requested by the referee: the total energies obtained by the model cannot be used to estimate first ionization potentials quantitatively. For differential applications the boundary layer (N > Z−5) should be treated as a positive, monotone, trend-correct estimate only: the reconstructed ionization energies inherit the neutral-anchor accuracy and reproduce the increasing stripping trend, but individual valence IPs are not quantitative (first ionization energies are reproduced only to within a factor of about three in the worst case — boron: 25.8 eV vs. the experimental 8.30 eV; typical cases are far better, e.g. Si 8.0 vs. 8.15 eV, Fe 7.0 vs. 7.90 eV, U 5.0 vs. 6.19 eV — and deviations up to ~50% occur for individual steps where the experimental sequence jumps across a subshell, e.g. Na, where the second ionization potential jumps from 5.1 to 47.3 eV entering the 2p shell). The Tier-3 minimal-touch corrections introduce an artificial 0.5 eV step at the affected interior links; this is visible only in finite differences and does not affect cumulative energies. For precision spectroscopy, Penning-trap mass measurements, or quantitative valence ionization potentials, experimental databases (NIST ASD) or ab initio values should be consulted.

Limitations. (i) Least accurate for the lightest elements (Z ≤ 18, recommended uncertainty 3–4.5%). (ii) Ground-state total binding energies only; no individual shells, excited states or fine structure. (iii) Deep-ionization accuracy for Z > 86 inherits the partly theoretical NIST reference data. (iv) No magnetic-field, temperature-density or plasma-screening effects. (v) The N = 1 anchor is the point-nucleus Dirac value; finite-nuclear-size and QED corrections (relevant at the 0.1–0.4% level for the heaviest elements) are not included. (vi) Z_eff(Z, N) is a purely phenomenological calibration parameter, not a screened charge (Section 2.3). (vii) First ionization potentials cannot be estimated quantitatively from the model (see above).

6. Conclusions

TF-RSE-IonCF v3.2 is a single unified closed-form model for the total electron binding energy of every neutral atom and every ion of the periodic table (Z = 1–92, q = 0 to Z–1). Relative to v3.1, the boundary-layer tier of the physical-consistency safeguard has been extended from a pure monotonicity test to a physical-plausibility window (0.5 eV ≤ IP_k ≤ 4·IP_base·1.8^(k−1)), closing the residual defect — identified by the referee — in which monotone but physically impossible near-neutral links survived the v3.1 scan (worst case: U raw first ionization energy +2588 eV vs. experiment 6.19 eV; 19 elements affected). After the extension, all 92 elements exhibit strictly positive, monotone and trend-correct boundary ionization energies (Si 8.0, Fe 7.0, U 5.0 eV vs. experiment 8.15, 7.90, 6.19 eV), and the degenerate v3.1 Na N = 10 link is repaired. The safeguard footprint (646 points, 15.1%; 86 elements; largest boundary-layer adjustment 7.84 keV at U N = 88; largest core adjustment 1.14 keV at Lu N = 63; deepest correction 15 charge states at Tm N = 54) is fully reproducible from the released code; the five clamp-degenerate anchor points it replaces (Li N = 2; Mg, Al, Si, P N = 10) are made more accurate, not less (e.g. Si N = 10: +1.25% → +0.03%; Si N = 11: +0.94% → +0.05%). Against 4278 NIST reference sums the safeguarded model retains a mean relative error of 0.260% (raw closed form 0.259%) with 100% of points below 5%, generalizes from five training sequences to 87 unseen elements at 0.27%, enforces zero negative ionization energies, zero non-monotonic elements and zero out-of-window first ionization energies on the full grid, and reduces by construction to the point-nucleus Dirac value at N = 1 (+0.35% from experiment at Z = 92) and to the validated TF-RSE v2 neutral value at N = Z. The model is a cumulative-energy tool: quantitative ionization potentials must be taken from experimental databases or ab initio calculations. The reference implementation (with command-line interface), validation dataset and test files are provided for independent verification.