Study of the Hyperfine Structure of the Low-Lying Ca II and Ra I-II Levels: Applying the MCDHF Models Developed for Ba I-II

1. Introduction

Alkaline-earth metals from Mg to Ra, particularly in their singly ionized forms (e.g., Ca II, Ra II), are ideal for high-precision calculations and measurements in astrophysics and metrology due to their simple electronic structure [1]. Neutral radium (Ra I) is also crucial for beyond-Standard-Model (BSM) physics.

In astrophysics, Ca II infrared triplet lines (8498, 8542, 8665 Å) serve as key calcium abundance indicators in stellar atmospheres, detectable even at extremely low metallicities $\left(\right[\mathrm{Fe}/\mathrm{H}]\le -3.0)$ [2], enabling studies of metal-poor stars [3] and galactic chemical evolution [4]. Ra II UV-Vis lines (3814, 4682, 6705 Å) probe r-process nucleosynthesis, as Ra is produced exclusively by neutron capture in neutron-rich environments where the s-process ends at bismuth (Z=83) [5,6]. Its likely production in neutron star mergers makes Ra detection in kilonovae spectra a direct indicator of such events [7].

In fundamental physics, Ca II and Ra II are crucial for frequency standards and precision measurements. Quadrupole moments in metastable atomic ions introduce systematic energy shifts in such systems [8]. High theoretical precision $(<1\%)$ is required for meaningful comparisons with parity non-conservation (PNC) experiments [9]. While heavier atoms like Ra offer more favorable PNC conditions, lighter elements such as Ca II serve as essential test cases for validating theoretical methods before tackling heavier systems.

Neutral radium (Ra I) offers unique opportunities for fundamental physics tests. The search for a permanent electron electric dipole moment (EDM) is among the most stringent BSM tests, with polar molecules like RaF serving as exceptionally sensitive probes due to their large internal electric fields and Ra’s heavy, octupole-deformed nucleus [10]. Precise ab initio calculations for Ra are essential, providing wavefunctions needed to interpret molecular EDM measurements and constrain $CP$-violating interactions. Multiconfiguration Dirac-Hartree-Fock (MCDHF) calculations for Ra I address this by capturing critical enhancement mechanisms, notably the [Rn]7s7p ³$P_1^o$ and [Rn]$7s6d{}^{3}D_2$ levels, separated by only $5.41$ ${\mathrm{cm}}^{-1}$. Combined with Ra high atomic number $(Z=88)$ and nuclear octupole deformation, this yields P- and T-reversal-odd effect enhancements up to five orders of magnitude greater than in other systems [11,12]. Electron correlation effects, which dominate the atomic EDM signal, are systematically accounted for in combination of MCDHF and relativistic configuration interaction (MCDHF+RCI) calculations. Ra II, with its 8230 Å clock transition and narrow linewidth, is also a candidate for optical atomic clocks. The low-cost diode lasers required for Ra II, compared to expensive UV lasers, make it attractive for both laboratory and transportable applications [13].

Accounting for the hyperfine structure (HFS) of the levels is essential. In astrophysics and atomic clock applications, HFS broadening and the resulting asymmetries in line profiles introduce inaccuracies in line positions and strengths. For BSM research, HFS constants enable extraction of electric field gradients to calibrate molecular theory against high-resolution RaF spectroscopy at ISOLDE (CERN), where 14 excited states have been characterized with $99.6\%$ theory-experiment agreement and laser-coolable transitions identified [10]. Accurate ab initio atomic wavefunctions and correlation energies could enhance RaF sensitivity to e-EDM searches, where experimental data are scarce.

For Ca II and Ra I-II, the HFS of certain levels remain experimentally unstudied, including the [Rn]$6d{}^{2}D_{3/2,5/2}$ levels in ^223,225Ra II and the [Rn]$7s6d{}^{3}D_{1,2,3}$ and [Rn]$7s7p{}^{3}P_2^o$ levels in Ra I. While extensive MCDHF calculations are available for Mg II [14], Ca II [8], Sr II, and Ba II [15], comprehensive calculations for Ra II remain limited.

In this work, we present theoretical calculations of magnetic dipole (A) and electric quadrupole (B) hyperfine structure constants for low-lying levels of Ca II and Ra I-II. Using the MCDHF method with Relativistic Configuration Interaction (RCI) within the GRASP framework [16], we extend our previously validated computational strategies from Sr II and Ba I-II [15] to the lighter Ca II and the heavier Ra I-II systems, employing for the first time a natural orbital basis transformation [17] approach for HFS calculations in heavy Ra I-II. Particular attention is given to the Bohr-Weisskopf (BW) effect, which is systematically accounted for in our calculations. By providing theoretical values for levels where experimental data are unavailable, and demonstrating agreement with existing measurements within our uncertainty estimates, we enable refined tests of atomic structure and nuclear effects. These results offer essential theoretical input for astrophysical abundance studies via Ca II infrared triplet lines, precision metrology using Ra II clock transitions, and searches for physics beyond the Standard Model through Ra I molecular and atomic electric dipole moment investigations.

This paper is organized as follows. Section 2 describes the MCDHF and RCI methods for calculating HFS constants and provides details about the transformation to natural orbitals. Section 3 presents Ca II calculations to validate the theoretical method from the earlier Sr II/Ba II study [15], followed by Ra II and Ra I HFS studies. Section 4 concludes with an outlook.

2. Method and Computational Details

2.1. Multiconfiguration Dirac–Hartree–Fock Method

The Multiconfiguration Dirac-Hartree-Fock (MCDHF) method, described in [18], is a fully relativistic approach to atomic structure calculations. It is particularly well suited for the study of medium- and high-Z elements as well as highly charged ions. Within this framework, an atomic state characterized by the quantum numbers $\Gamma$, J, $M_J$, and parity $\pi$ is represented by an Atomic State Function (ASF), denoted by $\Psi (\Gamma \pi J{M}_J)$. The ASF is expressed as a linear combination of Configuration State Functions (CSFs) $\Phi \left({\gamma }_i\pi J{M}_J\right)$,

$$\Psi (\Gamma \pi J{M}_J)={\displaystyle \sum _{i=1}^{N_{\mathrm{CSFs}}}}{c}_i^{\Gamma \pi J}\Phi \left({\gamma }_i\pi J{M}_J\right),$$

(1)

where the label ${\gamma }_i$ specifies the electronic configuration, the angular momentum coupling scheme, and any additional quantum numbers required to uniquely identify the CSF. The expansion coefficients $\left\{c_i^{\Gamma \pi J}\right\}$, together with the corresponding atomic-state energy E, are obtained by solving the secular equation

$$\mathit{H}\mathit{c}=E\mathit{c},$$

(2)

where $\mathit{H}$ denotes the Hamiltonian matrix.

Applying the variational principle to the total energy functional leads to the MCDHF equations, a set of coupled integro-differential equations that are solved self-consistently for the radial components of the one-electron spin-orbitals.

Subsequent Relativistic Configuration Interaction (RCI) calculations are performed by solving Eq. (2) using an enlarged CSF expansion, if desired, while keeping the orbitals obtained from the preceding MCDHF optimization fixed. At this stage, additional relativistic contributions may be incorporated into the Hamiltonian, including the Breit interaction and quantum electrodynamic (QED) corrections, thereby improving the accuracy of the calculated atomic properties.

2.2. Natural Orbitals

Due to numerical convergence challenges, and to reduce computation time, the MCDHF method employs a layer-by-layer (LBL) [19] approach in which only the newly introduced orbitals for the layer considered are optimized while the remaining are kept frozen. Therefore, the frozen orbitals, especially the spectroscopic reference orbitals, cannot react to the addition of the newly introduced layer. Natural orbitals (NO) provide an alternative orbital basis, first introduced in quantum chemistry, that addresses this limitation. They modify the orbitals of the reference states in a way that partly captures the response to the new layers. The NOs are obtained by diagonalizing the density matrix, which provides a more natural description of the electron correlation. One of the first applications to atomic physics was by Lindgren et al. [20] in order to determine hyperfine interaction in alkali atoms in many-body perturbation calculations. Schiffmann et al. [17] developed the NO formalism in the MCDHF framework and applied it in the hyperfine structure calculation of Na I [21], exhibiting results in better agreement with experiment than using the traditional LBL approach. Transformation to NOs also showed to considerably improve the hyperfine interaction constants and transition rates in other systems, e.g., In I [22].

For each value of $\kappa$, the relativistic angular quantum number, the density matrix is built such as

$${\rho }^{\kappa }={\rho }_{n{n}^{\prime }}^{\kappa },$$

(3)

with

$${\rho }_{n{n}^{\prime }}^{\kappa }=\sum _{ij}{c}_i^{★}{v}_{n{n}^{\prime }}^{ij,\kappa }{c}_j,$$

(4)

where $v_{n{n}^{\prime }\kappa }^{ij}$ are the spin-angular coefficients. The NOs $\tilde{\varphi }={\tilde{\varphi }}_{n\kappa }$ are obtained from the eigenvectors ${\mathit{U}}^{\kappa }$ of the density matrix diagonalizing it

$$\begin{array}{cc}\hfill {\mathit{U}}^{\kappa †}{\rho }^{\kappa }{\mathit{U}}^{\kappa }& ={\tilde{\rho }}^{\kappa },\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill {\tilde{\varphi }}^{\kappa }& ={\varphi }^{\kappa }{\mathit{U}}^{\kappa },\hfill \end{array}$$

(6)

where the eigenvalues are the occupation numbers ${\eta }_{n^{\prime }}^{\kappa }$ of the NOs

$$\sum _{\kappa }\sum _{n^{\prime }}{\eta }_{n^{\prime }}^{\kappa }=N_e,0\le {\eta }_{n^{\prime }}^{\kappa }\le 1,$$

(7)

with $N_e$ the number of electrons. The large and the small components of the radial function of the NOs are therefore

$$\begin{array}{cc}\hfill {\tilde{P}}_{n^{\prime }\kappa }\left(r\right)& =\sum _n{u}_{n{n}^{\prime }}^{\kappa }{P}_{n\kappa }\left(r\right),\hfill \\ \hfill {\tilde{Q}}_{n^{\prime }\kappa }\left(r\right)& =\sum _n{u}_{n{n}^{\prime }}^{\kappa }{Q}_{n\kappa }\left(r\right),\hfill \end{array}$$

(8)

with the coefficients $u_{n{n}^{\prime }}^{\kappa }$ the matrix elements of the unitary matrix ${\mathit{U}}^{\kappa }$ in eq. (5) and $P_{n\kappa }\left(r\right)$, $Q_{n\kappa }\left(r\right)$ the large and the small components of the radial function in the MCDHF equations [18].

2.3. Hyperfine Structure

The hyperfine structure (HFS) of an atomic fine-structure level arises from the interaction between the atomic electrons and the electromagnetic multipole moments of the nucleus. This interaction can be described through the multipole expansion

$$H_{\mathrm{HFS}}=\sum _{k\ge 1}{\mathit{T}}^{\left(k\right)}·{\mathit{M}}^{\left(k\right)},$$

(9)

where ${\mathit{T}}^{\left(k\right)}$ and ${\mathit{M}}^{\left(k\right)}$ are spherical tensor operators of rank k acting in the electronic and nuclear spaces, respectively. The terms corresponding to $k=1$ and $k=2$ represent the magnetic dipole (M1) and electric quadrupole (E2) interactions.

For an atom containing $N_e$ electrons, the electronic tensor operators are expressed as

$$\begin{array}{cc}\hfill {\mathit{T}}^{\left(1\right)}& ={\displaystyle \sum _{j=1}^{N_e}}{\mathit{t}}^{\left(1\right)}\left(j\right)={\displaystyle \sum _{j=1}^{N_e}}-i\sqrt{2}\alpha r_j^{-2}\left(\alpha j·{\mathit{C}}^{\left(1\right)}({\theta }_j,{\phi }_j)\right),\hfill \end{array}$$

(10)

$$\begin{array}{cc}\hfill {\mathit{T}}^{\left(2\right)}& =\sum {j=1}^{N_e}{\mathit{t}}^{\left(2\right)}\left(j\right)={\displaystyle \sum _{j=1}^{N_e}}-r_j^{-3}{\mathit{C}}^{\left(2\right)}({\theta }_j,{\phi }_j),\hfill \end{array}$$

(11)

where $i^2=-1$, $\alpha$ is the fine-structure constant, and $\alpha$ denotes the vector of Dirac matrices related to the Pauli matrices.

The nuclear magnetic dipole moment ${\mu }_I$ and electric quadrupole moment Q are defined through the matrix elements of the nuclear tensor operators,

$$\begin{array}{cc}\hfill {\mu }_I& =\langle II\left|M^{\left(1\right)}\right|II\rangle ,\hfill \end{array}$$

(12)

$$\begin{array}{cc}\hfill Q& =2\langle II\left|M^{\left(2\right)}\right|II\rangle ,\hfill \end{array}$$

(13)

where $|II\rangle$ denotes the nuclear state with maximal spin projection, $M_I=I$ [23]. The operators ${\mathit{M}}^{\left(1\right)}$ and ${\mathit{M}}^{\left(2\right)}$ correspond to the magnetic dipole and electric quadrupole moments of the nucleus, respectively.

Since the hyperfine interaction is generally much weaker than the fine-structure splitting, it can be treated within first-order perturbation theory. The corresponding hyperfine energy shift is given by

$$E_{\mathrm{HFS}}=\langle \Gamma JIF{M}_F|H_{\mathrm{HFS}}|\Gamma JIF{M}_F\rangle =E_{\mathrm{HFS}}^{\mathrm{M}1}+E_{\mathrm{HFS}}^{\mathrm{E}2},$$

(14)

where $|\Gamma JIF{M}_F\rangle$ is the zeroth-order wave function of the coupled electron-nucleus system [19], with F and $M_F$ denoting the total hyperfine angular momentum and its projection on the quantization axis, respectively.

It is convenient to express the hyperfine energy shifts in terms of the magnetic dipole (A) and electric quadrupole (B) hyperfine constants:

$$\begin{array}{cc}\hfill E_{\mathrm{HFS}}^{\mathrm{M}1}& ={\displaystyle \frac{1}{2}}{A}_{\Gamma J}K,\hfill \end{array}$$

(15)

$$\begin{array}{cc}\hfill E_{\mathrm{HFS}}^{\mathrm{E}2}& =B_{\Gamma J}{\displaystyle \frac{\frac{3}{4}K(K+1)-I(I+1)J(J+1)}{2I(2I-1)J(2J-1)}},\hfill \end{array}$$

(16)

where

$$\begin{array}{cc}\hfill A_{\Gamma J}& ={\displaystyle \frac{{\mu }_I}{I}}{\displaystyle \frac{1}{\sqrt{J(J+1)(2J+1)}}}\langle \Gamma J\left|\right|{\mathit{T}}^{\left(1\right)}\left|\right|\Gamma J\rangle ,\hfill \end{array}$$

(17)

$$\begin{array}{cc}\hfill B_{\Gamma J}& =2Q\sqrt{{\displaystyle \frac{J(2J-1)}{(J+1)(2J+1)(2J+3)}}}\langle \Gamma J\left|\right|{\mathit{T}}^{\left(2\right)}\left|\right|\Gamma J\rangle ,\hfill \end{array}$$

(18)

and

$$K=F(F+1)-J(J+1)-I(I+1).$$

(19)

The reduced matrix elements $\langle \Gamma J\left|\right|{\mathit{O}}^{\left(k\right)}\left|\right|\Gamma J\rangle$, with ${\mathit{O}}^{\left(k\right)}={\mathit{T}}^{\left(1\right)}$ or ${\mathit{T}}^{\left(2\right)}$, are evaluated using the atomic state functions (ASFs) defined in Eq. (1).

2.4. Computational Details

The MCDHF calculations and the evaluation of the hyperfine structure (HFS) constants were performed using GRASP2018 [24,25]. To improve computational efficiency, the recently developed GRASPG package by Si et al. [26], which exploits Configuration State Function Generators (CSFGs), was employed for the RCI calculations.

Following the methodology introduced by Papoulia et al. [27], the final HFS constants are reported as recommended values derived from the largest RCI calculations obtained with different orbital optimization strategies. The results from these calculations are treated as independent estimates, allowing the determination of both a mean value and a standard deviation for each constant. This procedure provides a robust prediction of the HFS constants together with a quantitative uncertainty assessment, which is particularly valuable in the absence of experimental data. The recommended values are expressed as $\mu \pm 2\sigma$, corresponding to a $95\%$ confidence interval.

3. Results and Discussion

3.1. Ca II

Calcium (Ca) has 26 known isotopes, of which 5 are stable [28]. ⁴³Ca has an abundance of only 0.135% but is the only stable isotope that has a hyperfine structure due to its non-zero nuclear spin $(I=7/2)$. The levels of interest in ⁴³Ca II are the even [Ar]$4s{}^{2}S_{1/2}$ and [Ar]$3d{}^{2}D_{3/2,5/2}$ as well as the odd [Ar]$4p{}^{2}P_{1/2,3/2}^o$.

Due to its relatively simple electronic structure, ⁴³Ca II has already been studied in depth theoretically by Silverans et al. [29] and Yu et al. [30] using the Relativistic Many-Body Perturbation Theory (RMBPT) approach, Itano [8] with the MCDHF method, Chakraborty et al. [31] with the Relativistic Coupled Cluster (RCC) method, and Safronova et al. [32] using the relavistic all-order method where all Single, Double and partial Triple (SDpT) excitations of the wavefunctions are included to all order of perturbation theory. Goble et al. [9] and Nortershauser et al. [1] performed optical spectroscopy experiments using hollow cathode discharge and fast ion beam collinear laser spectroscopy, whereas Garcia Ruiz et al. [33] adopted bunched-beam collinear laser spectroscopy.

The MCDHF method used in this paper relies on the AS [19] and LBL [25] approaches, which enable the optimization of only a subset of orbitals at each calculation step. For example, for the ground ²$S_{1/2}$ level, an initial Dirac-Hartree-Fock (DHF) layer is considered. This layer is built using the $\{4s,3p\}$ spectroscopic orbitals in nonrelativistic notation, where the highest n quantum number is quoted for each l-symmetry, and no electron correlation is included at this stage. The ASs are then expanded with correlation orbitals up to $l_{\max }=4$ (g-orbitals) to directly apply the Ba/Sr model [15] and achieve high precision in the calculations [34]. For each AS given in Table 1, only the newly introduced orbitals are optimized, while the remaining orbitals are kept frozen. The mixing coefficients of the ASFs are then determined, and Quantum Electrodynamics (QED) corrections are included at the RCI step. Finally, the HFS constants are computed.

Following Bieroń et al. [35], the impact of different classes of orbital substitutions has been investigated. Using an DHF orbital basis, more correlations are accumulated by allowing single and double (SD) substitutions from deeper orbitals employing Single References (SR) for each level using the fourth correlation layer (AS4) given in Table 1.

In Table 2, the individual contribution of single (S) orbital substitution to the magnetic dipole HFS constant of the ground level, $A{(}^2{S}_{1/2})$, is presented. Figure 1 gives a more visual representation of these contributions.

While the contribution of the $3p$ orbital is negligible, +0.08 MHz, the s orbital substitutions dominate, with the $3s$ orbital contributing $-81.54$ MHz. This highlights the critical role of s single substitutions in the HFS. To align with the Ba/Sr methodology [15], the $2p$ and $2s$ orbitals are included in the optimization step, even though the $2s$ substitution contributes less than a percent. Consequently, $1s$ correlations are excluded from the optimization.

As shown in Table 3 and Figure 2, the contributions for the excited states are similar, with the $3p$ contribution being negligible compared to the $3s,2p$ and $2s$ ones. The contributions for the B hyperfine constants indicate that the most important contribution comes from $2p$ correlations, which is mainly due to the orbital polarization or Sternheimer effect [23]. The trends above are also observed in the ²$D_{5/2}$ and the odd P levels.

Thanks to this analysis, the optimization of the radial orbitals can be aimed for relevant correlations for the hyperfine structure. The GRASP package [24] is considering spectroscopic orbitals in two ways. Let us take the example of the ground $1s^{2}2s^{2}2p^{6}3s^{2}3p^{6}4s{}^{2}S_{1/2}$ (or [Ar]4s ²$S_{1/2}$) level. The core orbitals are $1s^{2}2s^{2}2p^{6}3s^{2}3p^6$, while the valence orbital is $4s$. We can therefore build our MCDHF optimizations on core-valence (CV) correlations. To directly apply the MCDHF+RCI model developed for ⁸⁷Sr II and ¹³⁷Ba II [15], we used the equivalent CV single and restricted double (SrD) orbital optimization strategy with $nl\ge 3s$ substitutions. Here, "restricted" means at most one hole per active core orbital, i.e., allowing substitutions from $3s^2$ and $3p^6$ to $4s$ with at most one substitution from each core orbital to the valence one [36]. A second optimization strategy (MCDHF CV SrD $nl\ge 2s$) taking care of the $2s$ and $2p$ core orbitals that are important to correctly describe the HFS for the excited levels was also considered. With the ASs described in Table 1, seven layers of correlation orbitals were considered in order to achieve a convergence in our results. The latter is obtained when the HFS constants do not vary within some a priori accepted variation interval. Results for the HFS constants at each optimization step are given in Table 5.

As expected, CV correlations are crucial to describe the hyperfine interaction. The MCDHF CV SrD $nl\ge 3s$ shows a stable value achieved for the fourth active set (AS4). The ground level already agrees with experiment to within $\approx 7$ MHz (0.5%). This agreement is observed for each HFS constant except the $A{(}^2{D}_{5/2})$ value with a persistent $\approx 9$ MHz difference (223%). This behavior is well known from the Sr/Ba study [15]. To try to overcome this discrepancy, the second orbital set, MCDHF CV SrD $nl\ge 2s$ is considered. The agreement worsens for all HFS constants but is improved for the $A{(}^2{D}_{5/2})$ constant that still shows a $\approx 3$ MHz $(\approx 89\%)$ difference.

Following MCDHF calculations, deeper core-core (CC) correlations are added in the RCI step. The first RCI CC SD $nl\ge 3s$ model allows SD substitutions from $3p$ and $3s$ orbitals to the valence orbital, with no substitutions between the core orbitals in order to save computational time by decreasing the number of CSFs by a factor of two. As we want to account for as much correlation as possible, a second RCI strategy based on CC SD $nl\in \{3p,3s,2p,2s\}$ substitutions is considered. Adding correlations from the $1s$ core orbital was shown to be negligible, with less than a percent change in the values for the ground level. Therefore, they were also omitted. Using the MCDHF CV SrD nl$\ge 3s$ optimization, the extended RCI CC SD $nl\in \{3p,3s,2p,2s\}$ model modifies the electric quadrupole A hyperfine constants by $\approx 5\%$ for the ground level. The largest change appears for $A{(}^2{D}_{5/2})$ which varies by $\approx 60\%$. For the magnetic dipole B constant, we observe shifts of $\approx 1$MHz for the D levels, while $B{(}^2{P}_{3/2}^o)$ remains unchanged. The second MCDHF CV SrD $nl\ge 2s$ optimization yielded identical results with the same RCI approach.

The uncertainty in the theoretical HFS values can be evaluated through multiple approaches. One method involves treating the largest RCI calculations from distinct MCDHF optimizations as independent values. A mean value and standard deviation are then calculated from this set. The final result is expressed as $\mu \pm 2\sigma$, which corresponds to a 95% confidence interval.

However, this uncertainty estimation could only be applied to the layer-by-layer (LBL) calculation. The presented MCDHF+RCI methods were originally developed for barium (Ba) [15], which has a larger electronic core including d orbitals. When directly applied to lighter elements (such as calcium, with lower atomic number Z), the largest MCDHF+RCI method encountered limitations. Specifically, for calculations using the Natural Orbital (NO) basis, single (S) substitutions from all orbitals were included. In this case, both RCI calculations produced identical results. As a consequence, no uncertainty estimate could be obtained. This highlights a limitation of the model when applied to lighter elements.

Recent studies by Sahoo et al. [14] and Atalay [22] suggested a different approach for the use of a NO basis for HFS calculations. Compared to Nezosi et al. [15], where the biggest optimized MCDHF orbitals resulting from traditional single and restricted double (SrD) substitutions are transformed to NO using the rdensity program [17], here the NO basis is obtained via S substitutions (with at most one substitution from all core orbitals) from a multireference (MR) for the S, D and P levels studied in this paper. The largest MCDHF optimization calculation following this approach is then transformed to natural orbitals and then used in RCI calculations. This new approach is appealing as it can relax all orbitals in a single, compact calculation, as only a total of $4,196$ CSFs was needed compared to the $4,858$ CSFs for the MCDHF CV $nl\ge 2s$ optimization only for the ²$S_{1/2}$ HFS calculations at AS7. This approach is especially suitable for high-Z, as the number of CSFs grows rapidly with the number of active orbitals.

At first glance, equations (17,18) show that $A\propto \langle r^{-2}\rangle$ and $B\propto \langle r^{-3}\rangle$. Johnson [37] showed, in the Pauli approximation where the small $Q_{n\kappa }\left(r\right)$ and the large $P_{n\kappa }\left(r\right)$ can be linked by

$$Q_{n\kappa }\left(r\right)\approx -{\displaystyle \frac{1}{2c}}\left({\displaystyle \frac{d}{dr}}+{\displaystyle \frac{\kappa }{r}}\right)P_{n\kappa }\left(r\right),$$

(20)

that $\langle r^{-2}\rangle$ is linked to $\langle r^{-3}\rangle$, indicating that the A and B HFS constants are both $\propto \langle r^{-3}\rangle$ in the Pauli approximation.

Any modification in $\langle r\rangle$ caused by the transformation to NO could result in a modification in $\langle r^{-3}\rangle$. Table 4 compares the change in these expectation values for each d-subshell for a layer-by-layer (LBL) calculation and after transformation to NO. The spectroscopic $3d$ orbital is highlighted. Even the smallest change in $\langle r\rangle$ induces variation in $\langle r^{-2}\rangle$ as it can be seen for the $3d$ valence orbital

$$\begin{array}{cc}\hfill {\langle r\rangle }_{3d}^{\mathrm{LBL}}=2.5227a_0& \to {\langle r\rangle }_{3d}^{\mathrm{NO}}=2.3033a_0,\hfill \\ \hfill {\langle r^{-3}\rangle }_{3d}^{\mathrm{LBL}}=0.5937a_0^{-3}& \to {\langle r^{-3}\rangle }_{3d}^{\mathrm{NO}}=0.6947a_0^{-3},\hfill \end{array}$$

(21)

inducing a direct change in $A{(}^2{D}_{3/2})$ and $A{(}^2{D}_{5/2})$

$$\begin{array}{cc}\hfill A{{(}^2{D}_{3/2})}^{\mathrm{LBL}}=-51\mathrm{MHz}& \to A{{(}^2{D}_{3/2})}^{\mathrm{NO}}=-49\mathrm{MHz},\hfill \\ \hfill A{{(}^2{D}_{5/2})}^{\mathrm{LBL}}=-12\mathrm{MHz}& \to A{{(}^2{D}_{5/2})}^{\mathrm{NO}}=-11\mathrm{MHz}.\hfill \end{array}$$

(22)

The orbital contraction is more visible on Figure 3 between 3 and $6a_0$. The top left graph shows the large component $P_{3d}\left(r\right)$ of the valence orbital as a function of r, for the LBL calculation (in blue) and the transformed NO (in dashed red) obtained from the S substitutions MR calculation. The bottom-left graph shows the absolute difference between these two orbital bases. There are multiple zones where the NO valence orbital is more contracted than the LBL one. On the right-hand side, the $[3;6]a_0$ region is zoomed in where the contraction appears clearly with the red arrow. Even this small contraction (about $3\%$ in this region) can cause a difference in the A HFS constants. The same behavior applies to the valence orbitals of ground and the P excited levels, as well as in the B HFS constants.

Using the NO basis, we repeated all RCI calculations (labeled NO in Table 5). Results for both optimizations differ by only a few MHz. We recommend the final NO values as it accounts for orbital contraction and shows better agreement with experimental values, except for $A{(}^2{P}_{1/2}^o)$.

For an accurate theoretical determination of the HFS constants, one should account for the Bohr-Weisskopf (BW) effect [38], a correction of the finite nuclear magnetization distribution of the nucleus, giving a contribution to the HFS. Chakraborty et al. [31] studied the BW correction effect on the calculated HFS constants in ⁴³Ca II, with limited correction (at most 0.30% for the ²$D_{5/2}$ level). The BW corrected A HFS constants for the ²$S_{1/2}$ and ²$D_{3/2,5/2}$ levels are given in Table 6 but could not correct the fact that the values calculated with the LBL orbital basis are closer to the experiment than the ones calculated with the NO orbitals.

Table 5. Magnetic dipole A and electric quadrupole B hyperfine structure constants (in MHz) of ⁴³$\mathrm{Ca}\mathrm{II}$$(I=7/2,\mu =-1.31833{\mu }_N,Q=-0.0473\mathrm{b})$ [39], compared to avaialble theoretical and experimental data. The recommended values are given in the form $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^★: Number of CSFs obtained after reduction using the rcsfginteract_csfg code. ^a: Th. Silverans et al. (RMBPT) [29]; ^b: Th. Yu et al. (RMBPT) [30]; ^c: Th. Itano (MCDHF) [8]; ^d: Th. Chakraborty et al. (RCC) [31]; ^e: Th. Safronova et al. (SDpT) [32]; ^f: Expt. Goble et al. [9]; ^g: Expt. Nortershauser et al. [1]; ^h: Expt. Garcia Ruiz et al. [33]; ⁱ: Expt. Silverans et al. [29]

Step	2${\mathit{S}}_{1/2}$		2${\mathit{D}}_{3/2}$			2${\mathit{D}}_{5/2}$			2${\mathit{P}}_{1/2}^{\mathit{o}}$		2${\mathit{P}}_{3/2}^{\mathit{o}}$
	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 3s$
DHF	1	−606.08	1	−34.31	−2.67	1	−14.61	−3.77	1	−108.90	1	−21.11	−4.93
MCDHF (+AS1)	37	−740.61	120	−45.94	−1.92	130	−11.57	−2.74	53	−133.06	74	−25.64	−6.75
MCDHF (+AS2)	152	−774.30	433	−48.78	−2.67	496	−12.32	−3.79	181	−142.25	278	−32.92	−7.19
MCDHF (+AS3)	380	−802.39	944	−50.54	−2.65	1,100	−12.12	−3.77	423	−145.46	688	−34.14	−7.29
MCDHF (+AS4)	722	−796.91	1,653	−50.63	−2.78	1,942	−12.52	−3.94	781	−146.26	1,300	−33.58	−7.33
MCDHF (+AS5)	1,178	−800.83	2,560	−50.81	−2.72	3,022	−12.20	−3.87	1,255	−147.96	2,114	−33.77	−7.36
MCDHF (+AS6)	1,748	−800.26	3,665	−50.93	−2.76	4,340	−12.33	−3.92	1,845	−147.82	3,130	−33.76	−7.35
MCDHF (+AS7)	2,432	−801.87	4,968	−50.91	−2.75	5,896	−12.26	−3.90	2,551	−146.87	4,348	−33.79	−7.35
RCI CC SD $nl\ge 3s$ (+AS7)	9,407★	−737.72	27,462★	−44.88	−2.75	38,541★	−9.50	−3.91	9,825★	−125.47	23,013★	−28.59	−6.33
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7)	18,806★	−763.60	54,916★	−44.63	−4.64	77,074★	−4.94	−4.00	19,642★	−131.02	46,018★	−30.50	−6.68
RCI CC SD $nl\ge 3s$ (+AS7, NO)	9,407★	−764.91	27,462★	−46.83	−2.85	38,541★	−10.47	−4.06	9,825★	−126.58	23,013★	−28.81	−6.35
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7, NO)	18,806★	−810.51	54,916★	−48.34	−3.48	77,074★	−5.23	−4.95	19,642★	−133.17	46,018★	−31.30	−6.76
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 2s$
DHF	1	−606.08	1	−34.31	−2.67	1	−14.61	−3.77	1	−108.90	1	−21.11	−4.93
MCDHF (+AS1)	72	−745.89	238	−46.16	−2.05	258	−9.52	−2.90	104	−135.35	146	−25.82	−6.89
MCDHF (+AS2)	301	−782.91	863	−49.68	−3.11	989	−8.55	−4.42	359	−145.30	553	−33.02	−7.37
MCDHF (+AS3)	756	−841.46	1,884	−52.18	−3.19	2,196	−7.10	−4.54	842	−155.70	1,372	−35.74	−7.89
MCDHF (+AS4)	1,439	−845.50	3,301	−53.03	−3.47	3,879	−7.34	−4.94	1,557	−157.58	2,595	−36.59	−8.00
MCDHF (+AS5)	2,350	−854.23	5,114	−53.21	−3.36	6,038	−7.08	−4.79	2,504	−158.81	4,222	−37.50	−8.02
MCDHF (+AS6)	3,489	−854.63	7,323	−53.33	−3.41	8,673	−7.26	−4.86	3,682	−159.24	6,253	−37.52	−8.00
MCDHF (+AS7)	4,856	−854.43	9,928	−53.33	−3.39	11,784	−7.18	−4.81	5,094	−159.55	8,688	−37.51	−8.02
RCI CC SD $nl\ge 3s$ (+AS7)	11,831★	−782.79	32,422★	−46.82	−3.29	44,429★	−5.30	−4.70	12,368★	−135.67	27,353★	−31.58	−6.87
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7)	18,806★	−760.96	54,916★	−44.42	−3.25	77,074★	−-4.61	−3.98	19,642★	−130.84	46,018★	−30.35	−6.61
RCI CC SD $nl\ge 3s$ (+AS7, NO)	11,831★	−815.92	32,422★	−49.24	−3.42	44,429★	−5.29	−4.88	12,368★	−136.23	27,353★	−32.12	−6.64
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7, NO)	18,806★	−810.50	54,916★	−48.34	−3.48	77,074★	−5.23	−4.95	19,642★	−133.17	46,018★	−31.30	−6.76
This work (LBL)		−763.16 ± 2.85		−44.51 ± 0.32	−2.81 ± 0.02		−4.91 ± 0.12	−3.99 ± 0.03		−130.93 ± 0.25		−30.43 ± 0.21	−5.73 ± 0.10
This work (NO), recommended		−811.38		−48.32	−3.00		−5.24	−4.27		−133.17		−31.30	−5.83
Th.a		−819								−148		−30.9	-6.7
Th.b		−805.35		−47.82	−3.0		−3.55	−4.4		−143.07		−30.49	−6.7
Th.c				−47.27	−2.94		−4.84	−4.18
Th.d		−806 ± 3		−47.8 ± 0.3	−3.06 ± 0.07		−3.30 ± 0.04	−4.24 ± 0.09
Th.e		−801.28		−47.63	−2.87		−4.24	−4.08		−144.96		−30.34	−6.64
Expt.f		−797.5 ± 2.4								−158.0 ± 3.3		−29.7 ± 1.6
Expt.g				−47.3 ± 0.2	−3.7 ± 1.9		−3.8 ± 0.6	−3.9 ± 6.0		−145.4 ± 0.1		−31.0 ± 0.2	−6.9 ± 1.7
Expt.h		−806.87 ± 0.42										−31.10 ± 0.30	−4.2 ± 1.3
Expt.i		−805 ± 2								−145.5 ± 1.0		−31.9 ± 0.2	−6.7 ± 1.4

Table 5. Magnetic dipole A and electric quadrupole B hyperfine structure constants (in MHz) of ⁴³$\mathrm{Ca}\mathrm{II}$$(I=7/2,\mu =-1.31833{\mu }_N,Q=-0.0473\mathrm{b})$ [39], compared to avaialble theoretical and experimental data. The recommended values are given in the form $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^★: Number of CSFs obtained after reduction using the rcsfginteract_csfg code. ^a: Th. Silverans et al. (RMBPT) [29]; ^b: Th. Yu et al. (RMBPT) [30]; ^c: Th. Itano (MCDHF) [8]; ^d: Th. Chakraborty et al. (RCC) [31]; ^e: Th. Safronova et al. (SDpT) [32]; ^f: Expt. Goble et al. [9]; ^g: Expt. Nortershauser et al. [1]; ^h: Expt. Garcia Ruiz et al. [33]; ⁱ: Expt. Silverans et al. [29]

Step	2${\mathit{S}}_{1/2}$		2${\mathit{D}}_{3/2}$			2${\mathit{D}}_{5/2}$			2${\mathit{P}}_{1/2}^{\mathit{o}}$		2${\mathit{P}}_{3/2}^{\mathit{o}}$
	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 3s$
DHF	1	−606.08	1	−34.31	−2.67	1	−14.61	−3.77	1	−108.90	1	−21.11	−4.93
MCDHF (+AS1)	37	−740.61	120	−45.94	−1.92	130	−11.57	−2.74	53	−133.06	74	−25.64	−6.75
MCDHF (+AS2)	152	−774.30	433	−48.78	−2.67	496	−12.32	−3.79	181	−142.25	278	−32.92	−7.19
MCDHF (+AS3)	380	−802.39	944	−50.54	−2.65	1,100	−12.12	−3.77	423	−145.46	688	−34.14	−7.29
MCDHF (+AS4)	722	−796.91	1,653	−50.63	−2.78	1,942	−12.52	−3.94	781	−146.26	1,300	−33.58	−7.33
MCDHF (+AS5)	1,178	−800.83	2,560	−50.81	−2.72	3,022	−12.20	−3.87	1,255	−147.96	2,114	−33.77	−7.36
MCDHF (+AS6)	1,748	−800.26	3,665	−50.93	−2.76	4,340	−12.33	−3.92	1,845	−147.82	3,130	−33.76	−7.35
MCDHF (+AS7)	2,432	−801.87	4,968	−50.91	−2.75	5,896	−12.26	−3.90	2,551	−146.87	4,348	−33.79	−7.35
RCI CC SD $nl\ge 3s$ (+AS7)	9,407★	−737.72	27,462★	−44.88	−2.75	38,541★	−9.50	−3.91	9,825★	−125.47	23,013★	−28.59	−6.33
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7)	18,806★	−763.60	54,916★	−44.63	−4.64	77,074★	−4.94	−4.00	19,642★	−131.02	46,018★	−30.50	−6.68
RCI CC SD $nl\ge 3s$ (+AS7, NO)	9,407★	−764.91	27,462★	−46.83	−2.85	38,541★	−10.47	−4.06	9,825★	−126.58	23,013★	−28.81	−6.35
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7, NO)	18,806★	−810.51	54,916★	−48.34	−3.48	77,074★	−5.23	−4.95	19,642★	−133.17	46,018★	−31.30	−6.76
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 2s$
DHF	1	−606.08	1	−34.31	−2.67	1	−14.61	−3.77	1	−108.90	1	−21.11	−4.93
MCDHF (+AS1)	72	−745.89	238	−46.16	−2.05	258	−9.52	−2.90	104	−135.35	146	−25.82	−6.89
MCDHF (+AS2)	301	−782.91	863	−49.68	−3.11	989	−8.55	−4.42	359	−145.30	553	−33.02	−7.37
MCDHF (+AS3)	756	−841.46	1,884	−52.18	−3.19	2,196	−7.10	−4.54	842	−155.70	1,372	−35.74	−7.89
MCDHF (+AS4)	1,439	−845.50	3,301	−53.03	−3.47	3,879	−7.34	−4.94	1,557	−157.58	2,595	−36.59	−8.00
MCDHF (+AS5)	2,350	−854.23	5,114	−53.21	−3.36	6,038	−7.08	−4.79	2,504	−158.81	4,222	−37.50	−8.02
MCDHF (+AS6)	3,489	−854.63	7,323	−53.33	−3.41	8,673	−7.26	−4.86	3,682	−159.24	6,253	−37.52	−8.00
MCDHF (+AS7)	4,856	−854.43	9,928	−53.33	−3.39	11,784	−7.18	−4.81	5,094	−159.55	8,688	−37.51	−8.02
RCI CC SD $nl\ge 3s$ (+AS7)	11,831★	−782.79	32,422★	−46.82	−3.29	44,429★	−5.30	−4.70	12,368★	−135.67	27,353★	−31.58	−6.87
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7)	18,806★	−760.96	54,916★	−44.42	−3.25	77,074★	−-4.61	−3.98	19,642★	−130.84	46,018★	−30.35	−6.61
RCI CC SD $nl\ge 3s$ (+AS7, NO)	11,831★	−815.92	32,422★	−49.24	−3.42	44,429★	−5.29	−4.88	12,368★	−136.23	27,353★	−32.12	−6.64
RCI CC SD $nl\in \{3p,3s,2p,2s\}$ (+AS7, NO)	18,806★	−810.50	54,916★	−48.34	−3.48	77,074★	−5.23	−4.95	19,642★	−133.17	46,018★	−31.30	−6.76
This work (LBL)		−763.16 ± 2.85		−44.51 ± 0.32	−2.81 ± 0.02		−4.91 ± 0.12	−3.99 ± 0.03		−130.93 ± 0.25		−30.43 ± 0.21	−5.73 ± 0.10
This work (NO), recommended		−811.38		−48.32	−3.00		−5.24	−4.27		−133.17		−31.30	−5.83
Th.a		−819								−148		−30.9	-6.7
Th.b		−805.35		−47.82	−3.0		−3.55	−4.4		−143.07		−30.49	−6.7
Th.c				−47.27	−2.94		−4.84	−4.18
Th.d		−806 ± 3		−47.8 ± 0.3	−3.06 ± 0.07		−3.30 ± 0.04	−4.24 ± 0.09
Th.e		−801.28		−47.63	−2.87		−4.24	−4.08		−144.96		−30.34	−6.64
Expt.f		−797.5 ± 2.4								−158.0 ± 3.3		−29.7 ± 1.6
Expt.g				−47.3 ± 0.2	−3.7 ± 1.9		−3.8 ± 0.6	−3.9 ± 6.0		−145.4 ± 0.1		−31.0 ± 0.2	−6.9 ± 1.7
Expt.h		−806.87 ± 0.42										−31.10 ± 0.30	−4.2 ± 1.3
Expt.i		−805 ± 2								−145.5 ± 1.0		−31.9 ± 0.2	−6.7 ± 1.4

Table 6. Bohr-Weisskopf (BW) corrected magnetic dipole A hyperfine structure constants (in MHz) of ⁴³$\mathrm{Ca}\mathrm{II}$$(I=7/2,\mu =-1.31833{\mu }_N$ [40], $Q=-0.0473\mathrm{b}$ [39]), using the extracted relative BW corrections of Chakraborty et al. [31], compared to available theoretical and experimental data. The recommended values are given in the form of $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^a: Th. Silverans et al. (RMBPT) [29]; ^b: Th. Yu et al. (RMBPT) [30]; ^c: Th. Itano (MCDHF) [8]; ^d: Th. Chakraborty et al. (RCC) [31]; ^e: Th. Safronova et al. (SDpT) [32]; ^f: Expt. Goble et al. [9]; ^g: Expt. Nortershauser et al. [1]; ^h: Expt. Garcia Ruiz et al. [33]; ⁱ: Expt. Silverans et al. [29]

Reference	2${\mathit{S}}_{1/2}$	2${\mathit{D}}_{3/2}$	2${\mathit{D}}_{5/2}$
	43Ca II
This work (LBL)	−762.29 ± 3.73	−44.51 ± 0.34	−4.89 ± 0.11
This work (NO)	−810.49	−48.30	−5.22
Th.a	−819
Th.b	−805.35	−47.82	−3.55
Th.c		−47.27	−4.84
Th.d	−806 $\pm$3	−47.8 ± 0.3	−3.30 ± 0.04
Th.e	−801.28	−47.63	−4.24
Expt.f	−795.5 ± 2.4
Expt.g		−47.3 ± 0.2	−3.8 ± 0.6
Expt.h	−806.87 ± 0.42
Expt.i	−805 ± 2

Table 6. Bohr-Weisskopf (BW) corrected magnetic dipole A hyperfine structure constants (in MHz) of ⁴³$\mathrm{Ca}\mathrm{II}$$(I=7/2,\mu =-1.31833{\mu }_N$ [40], $Q=-0.0473\mathrm{b}$ [39]), using the extracted relative BW corrections of Chakraborty et al. [31], compared to available theoretical and experimental data. The recommended values are given in the form of $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^a: Th. Silverans et al. (RMBPT) [29]; ^b: Th. Yu et al. (RMBPT) [30]; ^c: Th. Itano (MCDHF) [8]; ^d: Th. Chakraborty et al. (RCC) [31]; ^e: Th. Safronova et al. (SDpT) [32]; ^f: Expt. Goble et al. [9]; ^g: Expt. Nortershauser et al. [1]; ^h: Expt. Garcia Ruiz et al. [33]; ⁱ: Expt. Silverans et al. [29]

Reference	2${\mathit{S}}_{1/2}$	2${\mathit{D}}_{3/2}$	2${\mathit{D}}_{5/2}$
	43Ca II
This work (LBL)	−762.29 ± 3.73	−44.51 ± 0.34	−4.89 ± 0.11
This work (NO)	−810.49	−48.30	−5.22
Th.a	−819
Th.b	−805.35	−47.82	−3.55
Th.c		−47.27	−4.84
Th.d	−806 $\pm$3	−47.8 ± 0.3	−3.30 ± 0.04
Th.e	−801.28	−47.63	−4.24
Expt.f	−795.5 ± 2.4
Expt.g		−47.3 ± 0.2	−3.8 ± 0.6
Expt.h	−806.87 ± 0.42
Expt.i	−805 ± 2

3.2. Ra II

Radium (Ra) is the heaviest alkaline earth metal, with 34 known radioactive isotopes. The most stable, ²²⁶Ra, has a half-life of 1,600 years [28]. The isotopes ²²³Ra and ²²⁵Ra, with nuclear spin $I=3/2$ and half-lives of $11.44$ and $14.82$ days respectively, only exist in trace amounts. Their short half-lives make experiments challenging, so theoretical hyperfine structure calculations are needed to provide experimental bounds. The levels of interest are the even [Rn]7s ²$S_{1/2}$ and [Rn]6d ²$D_{3/2,5/2}$ levels, as well as the odd [Rn]7p ²$P_{1/2,3/2}^o$ levels, for which only a few experimental values are available. Using collinear fast-beam laser spectroscopy, Wendt et al. [41] determined the HFS of both $J=1/2$ levels, while Neu et al. [42] studied both the ground ²$S_{1/2}$ and ²$P_{3/2}^o$ levels. Theoretical values were obtained using RCC [43][44], Atomic Many-body Perturbation Theory in the Screened Coulomb Interaction (AMPSCI) [45] or the relativistic all-order SDpT method [46]. This is the first paper to compute HFS constants using the MCDHF and RCI method and to benchmark the new GRASPG [26] package for high-Z ions.

With the AS given in Table 7, the code limited the investigation of orbital substitutions to those no higher than $3d$. Table 8 and Figure 4 show that, for the ground level, the highest contribution comes from $6s$ substitutions. $5d$ and $5p$ substitutions can be added to the optimization model in order to follow the Ba/Sr [15] strategies.

The obtained values for ²²³Ra II for both optimization strategies are given in Table 10. The first MCDHF CV SrD $nl\ge 6s$ gives electric quadrupole A HFS constants that are near the literature values with at most $14\%$ difference for $A{(}^2{D}_{5/2})$. However, the magnetic dipole B constants need more correlations to be accurately determined, as shown by the discrepancies for the D levels $(\approx 50\%)$. However, when an experimental value exists, the calculated value agrees with only $\approx 9\%$ difference for the $B{(}^2{P}_{3/2}^o)$ constant. Adding deeper correlations with the second MCDHF CV SrD $nl\ge 5p$ deteriorates the agreements with the experimental values [41,42] and the other theoretical values [43,44,45] for all the A HFS constants except $A{(}^2{D}_{5/2})$. The B HFS constants are $40\%$ closer to the calculated theoretical values of Sahoo et al. [43] (RCC). The agreement is preserved for the $B{(}^2{P}_{3/2}^o)$ constant $(\approx 9\%)$ but is now too large compared to the experiment value. We are still convinced that adding correlation orbitals with higher l encapsulates the orbital polarization that was neglected in the first MCDHF CV SrD $nl\ge 6s$ strategy and is important for a reliable determination of the B hyperfine constant [23].

After the MCDHF calculations, higher-order core-core (CC) correlations are added using RCI. Like in the Ca II case, two models are considered for each optimization strategy by allowing single and double (SD) substitutions. For all A HFS constants and for each orbital basis, both RCI models (RCI CC SD $nl\ge 6s$ and RCI CC SD $nl\in \{6p,6s,5p,5s\}$) lower the positive values and increase the negative values. However, the added CC correlations are not enough to have a significantly better agreement for the B HFS constants for the D levels while using the MCDHF CV SrD $nl\ge 6s$ orbital basis. The agreement even worsened for the ²$P_{3/2}^o$ level with nearly a 100 MHz difference with the experimental value. Overall, the best agreement is obtained with the MCDHF CV SrD $nl\ge 5p$ + RCI CC SD $nl\in \{6p,6s,5p,5s\}$ strategy. The agreement between the calculated value and RCC value [43] is still satisfactory for the $A{(}^2{D}_{5/2})$ constant.

We recall that, for Sr (Z = 38) and Ba (Z = 56), the $A{(}^2{D}_{5/2})$ HFS constant was problematic, with an opposite sign between calculated and literature values [15]. For Ra (Z = 88), the HFS constant does not show the same behavior. To understand this more clearly, we again move to the non-relativistic formalism, where the magnetic dipole HFS constant A is given as a sum of orbital, spin-dipole, and Fermi contact terms [47]. Table 9 shows the contribution of the different nonrelativistic hyperfine factors on the total A constant. For the Hartree-Fock (HF) calculation, the contact term is zero. Adding core-valence (CV) correlations by opening $6p$ and $6s$ for substitutions to one layer of correlation orbitals $\{7s,7p,7d\}$, we see a positive interference for the ²$D_{3/2}$ level and a negative one for the ²$D_{5/2}$ level. Compared with Sr and Ba the negative interference between the different terms is less pronounced for Ra making it less sensitive to neglected higher-order correlation effects from triple (T) and quadruple (Q) substitutions that are the likely causes for the lack of agreement with the experiment for the former elements.

In Ca II, the recommended values for this work were obtained using an NO basis. Following the same approach for Ra II, RCI calculations were repeated, with results presented in Table 10. While the LBL orbital basis remains suitable for the ground level, the NO basis provides better agreement for all excited levels, with an outstanding $<1\%$ difference with the experimental value of $A{(}^2{P}_{3/2}^o)$. Nevertheless, the uncertainty calculation for $A{(}^2{S}_{1/2})$ encompasses both the NO result and the experimental value.

Table 10. Magnetic dipole A and electric quadrupole B hyperfine structure constants (in MHz) of ²²³$\mathrm{Ra}\mathrm{II}$$(I=3/2,\mu =+0.2692{\mu }_N,Q=1.22\mathrm{b}$ [48]), compared to available theoretical and experimental data. The recommended values are given in the form of $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^★: Number of CSFs obtained after reduction using the rcsfginteract_csfg code; ^a: Th. Sahoo et al. (RCC) [43]; ^b: Th. Cserveny et al. (AMPSCI) [45]; ^c: Th. Li et al. (RCC) [44]; ^d: Th. Pal et al. (SDpT) [46];^e: Expt. Wendt et al. [41]; ^f: Expt. Neu et al. [42]

Step	2${\mathit{S}}_{1/2}$		2${\mathit{D}}_{3/2}$			2${\mathit{D}}_{5/2}$			2${\mathit{P}}_{1/2}^{\mathit{o}}$		2${\mathit{P}}_{3/2}^{\mathit{o}}$
	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 6s$
DHF	1	2,811	1	55	272	1	20	288	1	490	1	39	562
MCDHF (+AS1)	70	3,335	156	76	129	184	−35	157	90	620	150	51	732
MCDHF (+AS2)	252	3,408	505	72	194	601	−29	251	289	618	493	61	773
MCDHF (+AS3)	548	3,514	1,052	73	185	1,256	−30	235	604	642	1,038	53	789
MCDHF (+AS4)	958	3,487	1,797	72	185	2,149	−28	238	1,035	640	1,785	50	787
MCDHF (+AS5)	1,482	3,492	2,740	72	184	3,280	−28	235	1,582	642	2,734	52	791
MCDHF (+AS6)	2,120	3,480	3,881	72	184	4,649	−28	235	2,245	638	3,885	52	787
MCDHF (+AS7)	2,872	3,481	5,220	72	185	6,256	−28	236	3,024	638	5,238	52	787
RCI CC SD $nl\ge 6s$ (+AS7)	11,391★	3,222	29,379★	65	234	41,692★	−21	295	11,861★	570	28,513★	47	706
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7)	22,774★	3,349	58,750★	66	247	83,376★	−23	379	23,714★	589	57,018★	51	754
RCI CC SD $nl\ge 6s$ (+AS7, NO)	11,391★	3,336	29,379★	67	231	41,692★	−26	297	11,861★	606	28,513★	50	747
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7, NO)	22,774★	3,477	58,750★	69	317	83,376★	−26	405	23,714★	627	57,018★	55	806
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 5p$
DHF	1	2,811	1	55	272	1	20	288	1	490	1	39	562
MCDHF (+AS1)	184	3,424	440	80	190	526	−34	234	236	645	404	50	766
MCDHF (+AS2)	680	3,685	1,415	81	328	1,705	−27	418	775	655	1,347	62	836
MCDHF (+AS3)	1,490	3,877	2,944	85	327	3,558	−26	406	1,630	730	2,846	58	922
MCDHF (+AS4)	2,614	3,850	5,027	84	338	6,085	−24	523	2,801	737	4,801	57	932
MCDHF (+AS5)	4,052	3,877	7,664	84	336	9,286	−24	419	4,288	747	7,512	59	942
MCDHF (+AS6)	5,804	3,867	10,855	84	338	13,161	−24	423	6,091	743	10,679	58	937
MCDHF (+AS7)	7,870	3,871	14,600	84	338	17,710	−24	423	8,210	745	14,402	58	940
RCI CC SD $nl\ge 6s$ (+AS7)	16,389★	3,522	38,759★	74	366	53,146★	−18	455	17,047★	650	37,677★	51	823
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7)	25,714★	3,583	64,343★	76	365	90,280★	−21	454	26,668★	632	62,303★	54	808
RCI CC SD $nl\ge 6s$ (+AS7, NO)	16,389★	3,583	38,759★	81	396	53,146★	−22	494	17,047★	678	37,677★	54	869
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7, NO)	25,714★	3,709	64,343★	81	393	90,280★	−23	492	26,668★	669	62,303★	57	864
This work (LBL)		3,466 ± 331		71 ± 14	306 ± 167		−22 ± 3	417 ± 106		611 ± 61		53 ± 4	781 ± 76
This work (NO), recommended		3,593 ± 328		75 ± 17	355 ± 107		−25 ± 4	449 ± 123		648 ± 59		56 ± 3	835 ± 82
Th.a		3,567.26		77.08	383.88		−23.90	477.09
Th.b		3,424 ± 71		95 ± 20			$-29$± 10			675 ± 13		57.7 ± 60
Th.c		3,426 ± 82		79.7 ± 1.7			$-21.8$± 2.3			660 ± 11.4		56.5 ± 1.2
Th.d		3,450		79.56			$-24.08$			671.5		54.40
Expt.e		3,398.3 ± 2.9								667.1 ± 2.1
Expt.f		3,404 ± 1.9										56.5 ± 8	864.8 ± 1.9

Table 10. Magnetic dipole A and electric quadrupole B hyperfine structure constants (in MHz) of ²²³$\mathrm{Ra}\mathrm{II}$$(I=3/2,\mu =+0.2692{\mu }_N,Q=1.22\mathrm{b}$ [48]), compared to available theoretical and experimental data. The recommended values are given in the form of $\mu \pm 2\sigma$, exhibiting a $95\%$ confidence level. ^★: Number of CSFs obtained after reduction using the rcsfginteract_csfg code; ^a: Th. Sahoo et al. (RCC) [43]; ^b: Th. Cserveny et al. (AMPSCI) [45]; ^c: Th. Li et al. (RCC) [44]; ^d: Th. Pal et al. (SDpT) [46];^e: Expt. Wendt et al. [41]; ^f: Expt. Neu et al. [42]

Step	2${\mathit{S}}_{1/2}$		2${\mathit{D}}_{3/2}$			2${\mathit{D}}_{5/2}$			2${\mathit{P}}_{1/2}^{\mathit{o}}$		2${\mathit{P}}_{3/2}^{\mathit{o}}$
	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	${\mathbf{N}}_{\mathbf{CSFs}}$	$\mathit{A}$	$\mathit{B}$
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 6s$
DHF	1	2,811	1	55	272	1	20	288	1	490	1	39	562
MCDHF (+AS1)	70	3,335	156	76	129	184	−35	157	90	620	150	51	732
MCDHF (+AS2)	252	3,408	505	72	194	601	−29	251	289	618	493	61	773
MCDHF (+AS3)	548	3,514	1,052	73	185	1,256	−30	235	604	642	1,038	53	789
MCDHF (+AS4)	958	3,487	1,797	72	185	2,149	−28	238	1,035	640	1,785	50	787
MCDHF (+AS5)	1,482	3,492	2,740	72	184	3,280	−28	235	1,582	642	2,734	52	791
MCDHF (+AS6)	2,120	3,480	3,881	72	184	4,649	−28	235	2,245	638	3,885	52	787
MCDHF (+AS7)	2,872	3,481	5,220	72	185	6,256	−28	236	3,024	638	5,238	52	787
RCI CC SD $nl\ge 6s$ (+AS7)	11,391★	3,222	29,379★	65	234	41,692★	−21	295	11,861★	570	28,513★	47	706
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7)	22,774★	3,349	58,750★	66	247	83,376★	−23	379	23,714★	589	57,018★	51	754
RCI CC SD $nl\ge 6s$ (+AS7, NO)	11,391★	3,336	29,379★	67	231	41,692★	−26	297	11,861★	606	28,513★	50	747
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7, NO)	22,774★	3,477	58,750★	69	317	83,376★	−26	405	23,714★	627	57,018★	55	806
$\mathrm{Optimization}\mathrm{CV}\mathrm{SrD}nl\ge 5p$
DHF	1	2,811	1	55	272	1	20	288	1	490	1	39	562
MCDHF (+AS1)	184	3,424	440	80	190	526	−34	234	236	645	404	50	766
MCDHF (+AS2)	680	3,685	1,415	81	328	1,705	−27	418	775	655	1,347	62	836
MCDHF (+AS3)	1,490	3,877	2,944	85	327	3,558	−26	406	1,630	730	2,846	58	922
MCDHF (+AS4)	2,614	3,850	5,027	84	338	6,085	−24	523	2,801	737	4,801	57	932
MCDHF (+AS5)	4,052	3,877	7,664	84	336	9,286	−24	419	4,288	747	7,512	59	942
MCDHF (+AS6)	5,804	3,867	10,855	84	338	13,161	−24	423	6,091	743	10,679	58	937
MCDHF (+AS7)	7,870	3,871	14,600	84	338	17,710	−24	423	8,210	745	14,402	58	940
RCI CC SD $nl\ge 6s$ (+AS7)	16,389★	3,522	38,759★	74	366	53,146★	−18	455	17,047★	650	37,677★	51	823
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7)	25,714★	3,583	64,343★	76	365	90,280★	−21	454	26,668★	632	62,303★	54	808
RCI CC SD $nl\ge 6s$ (+AS7, NO)	16,389★	3,583	38,759★	81	396	53,146★	−22	494	17,047★	678	37,677★	54	869
RCI CC SD $nl\in \{6p,6s,5p,5s\}$ (+AS7, NO)	25,714★	3,709	64,343★	81	393	90,280★	−23	492	26,668★	669	62,303★	57	864
This work (LBL)		3,466 ± 331		71 ± 14	306 ± 167		−22 ± 3	417 ± 106		611 ± 61		53 ± 4	781 ± 76
This work (NO), recommended		3,593 ± 328		75 ± 17	355 ± 107		−25 ± 4	449 ± 123		648 ± 59		56 ± 3	835 ± 82
Th.a		3,567.26		77.08	383.88		−23.90	477.09
Th.b		3,424 ± 71		95 ± 20			$-29$± 10			675 ± 13		57.7 ± 60
Th.c		3,426 ± 82		79.7 ± 1.7			$-21.8$± 2.3			660 ± 11.4		56.5 ± 1.2
Th.d		3,450		79.56			$-24.08$			671.5		54.40
Expt.e		3,398.3 ± 2.9								667.1 ± 2.1
Expt.f		3,404 ± 1.9										56.5 ± 8	864.8 ± 1.9

The Bohr-Weisskopf (BW) effect could affect significantly our HFS A constants in such heavy element as radium [44,45]. Table 11 reports and compares our BW corrected values, as determined by using the extracted relative BW corrections of Cserveny et al. [45], with the available theoretical and experimental data. These corrections are more important for the S and D states than for the P states [45]. It can be seen that now our BW corrected NO values are closer to experiment than our LBL ones as expected.

Nuclear physicists are also interested in ²²⁵Ra. As $I{(}^{225}\mathrm{Ra})=1/2$ strictly implies that $Q{(}^{225}\mathrm{Ra})=0$ using the Wigner-Eckart theorem [23], only the A hyperfine structure constant can be evaluated for this isotope. Using ²²³Ra II hyperfine constants in Table 10, the ²²⁵Ra II constants are determined thanks to the HFS electronic constants [27]

$$\begin{array}{cc}\hfill A_{\mathrm{el}.}& =A\left({\displaystyle \frac{I}{{\mu }_I}}\right)\hfill \end{array}.$$

(23)

Table 11 reports the BW corrected values for the A HFS constants of ²²⁵Ra II, obtained through the use of the extracted relative BW corrections of Cserveny et al. [45], compared with the data of the literature. It shows the overall better agreement of the HFS constants calculated with the transformed NO basis with the notable exception of the ²$S_{1/2}$ ground level where our LBL value is closer to experiment. This exception remains unexplained as it was not the case for ²²³$\mathrm{Ra}\mathrm{II}$ and as these values were calculated using Eq. 23.

3.3. Ra I

Neutral radium (Ra I) has two valence electrons. Only a few number of levels have had their HFS experimentally studied. The odd [Rn]7s7p ³$P_1^o$ have been studied by Lynch et al. [50] at the ISOLDE facility and [Rn]7s7p ³$P_2^o,{}^{1}P_1^o$ have been studied by Wendt et al. [41]. The HFS of the even [Rn]$7s6d{}^{1}D_2,{}^{3}D_{1,2,3}$ levels are only known theoretically. Hassouneh et al. [51] and Bieroń et al. [52] studied them using the MCDHF method, while Dzuba et al. [53] used a Relativistic Hartree-Fock + Configuration Interaction (RHF+CI) method to determine the HFS of ²¹³Ra I $(I=1/2,\mu =0.6133{\mu }_N)$ Their values were used to compute the HFS constants for ²²⁵Ra II using eq. (23). The latter two were computationally limited by the number of CSFs or by using a pseudo-relativistic Schrödinger Hamiltonian that could not be enough to study the high relativistic effects of electrons in heavy atoms.

Nezosi et al. [15] Ba I HFS study highlighted that a multireference (MR) is needed for heavy neutral elements such as Ra I while using the MCDHF method. As the main goal of this paper is mainly to transpose Sr II/Ba I-II MCDHF strategies, the configuration that interacts the most with the reference is be added for both studied parities. The even [Rn]$7s6d{}^{1}D_2,{}^{3}D_{1,2,3}$ levels were augmented by a [Rn]$6d^2$ configuration, while the odd [Rn]7s7p ³$P_{1,2}^o,{}^{1}P_1^o$ levels were augmented by a [Rn]$7p6d$ configuration. Seven correlation layers were considered in Table 12 for the optimization of the MR formed by the two parity configurations.

The first calculations were performed in the layer-by-layer (LBL) approach and the results are displayed in Table 13. Due to the two valence-electron structure of Ra I, an additional optimization strategy could be considered based on valence-valence (VV) correlations. Computations shows that the MCDHF VV SrD optimization strategy is not sufficient to take into account the core polarization needed for an accurate determination of the hyperfine structure, with differences ranging from $\approx 40$MHz for $B{(}^3{D}_1)$ to $\approx 260$ MHz for $B{(}^3{P}_2^o)$. Adding core-valence (CV) correlations in the MCDHF VV SrD + CV SrD $nl\ge 6s$ strategy adds the important polarization of the core caused by the valence electrons with already satisfactory agreement with literature values for the A constants. The cost to include these important correlations is to increase the number of CSFs by nearly 70 times for the even levels, add 100 times for the odd levels. One final layer of optimization in MCDHF VV SrD + CV SrD $nl\ge 5d$ was added to match the Ba I [15] model. It again shows that the more correlations, the better, especially for the B hyperfine constants, when asymmetrical d ($5d$) core orbital substitutions are allowed. This behavior was already observed in Ba and reinforces the idea that substitutions from orbitals with a higher l take into account the orbital polarization that was previously neglected and are therefore important for a reliable determination of B [23]. This third layer of optimization increases the HFS constants for all the levels in absolute value.

Adding higher-order correlations with SD substitutions in a RCI step highly increases the number of CSFs to a total of 2 933 035 for the biggest MCDHF VV SrD + CV SrD $nl\ge 5d$ + RCI CC SD $nl\in \{6s,6p,5s,5p\}$ step. Therefore, even if T substitutions are known in Na I to correct CC correlations to the HFS [21], they could not be added in Ra I. Nevertheless, the largest model shows the best agreement for this approach, with only 1 MHz of discrepancy for $B{(}^3{D}_2)$, but can still under-perform for $B{(}^3{D}_2)$ or $B{(}^3{P}_{1,2}^o)$ with $\approx 60$ and $\approx 20\%$ difference with literature values, respectively.

Taking the LBL picture as a whole when determining the final values in the $\mu \pm 2\sigma$ format shows satisfactory results for $A{(}^3{D}_1)$ with only $2\%$ difference with the latest available MCDHF values [51]. Our calculations are in agreement with the available experimental value for $A{(}^3{P}_2^o)$ [41] with $\approx 2\%$ discrepency and outperforming Hassouneh et al. [51] MCDHF calculations for this level while giving an uncertainty estimation.

Moving to the NO basis allows us to get an overall better agreement than the LBL calculation, as it was already the case for ⁴³Ca II and ²²³Ra II. The final values with uncertainty estimates cover other theoretical calculations. When the NO basis degrades the values compared to the LBL basis, it is always only by a few percent.

Using the correction factor of Bieroń [52], our results for ²²³Ra I in Table 13 now include the BW effect. A rescaling of the A hyperfine constants using eq. (23) for ²²⁵Ra I is also done. The corrected results are shown in Table 14. The situation is less clear than in the case of Ra II. There is no evidence of a quasi systematic improvement of the agreement with the measurements when one transforms LBL orbitals to NO orbitals. Now, one can see from Table 14 that the BW corrected MCDHF-RCI calculations of both Hassouneh and Salah [51] and Bieroń [52] show better general agreement with experiment than our BW corrected values. There are still two exceptions: the ³$P_2^o$ level in ²²³Ra I and the ³$D_1$ level in ²²⁵Ra I. For these levels, our BW corrected NO and LBL HFS A constants are respectively in better agreement with the measurements of Wendt et al. [41] and Guest et al. [54]. It shows that the differences seen here in Ra I is more due to the correlation model adopted than the BW effect. However, this BW correction factor $(\approx 0.96)$ of Bieroń [52] applied to each levels is a limited approach as [31,44,45] showed the dependence of the BW effect to the level studied.

4. Conclusions

This work presented a comprehensive computational study of the hyperfine structure (HFS) of low-lying levels in Ca II and Ra I-II, employing the MCDHF and RCI method [8] within the GRASP framework [55] enhanced with Configuration State Function Generators (CSFGs) [26].

Building upon our established models for strontium (Sr) and barium (Ba) [15], we successfully extended this methodology to both a lighter alkaline-earth ion, calcium (Ca), and the heavy radioactive element radium (Ra). The implementation of a natural orbital (NO) basis for HFS calculations in heavy elements represented a significant methodological advancement, accounting for orbital relaxation and systematically improving accuracy. Furthermore, we introduced a robust statistical approach for uncertainty estimation, providing experimenters with reliable theoretical bounds for their measurements.

The calcium ion (Ca II) served as a crucial validation case, allowing us to refine our computational approach through detailed analysis of orbital substitution contributions. This investigation revealed the dominant electron correlations influencing HFS and demonstrated that an NO basis induces valence orbital contraction, consistent with observations in lighter elements [21]. Our final values for Ca II exhibited excellent agreement with experimental data while offering precise theoretical uncertainty estimates, thereby confirming the reliability of our method for subsequent applications to more complex systems. However, the application of the theoretical uncertainty determination approach within the new NO transformation framework yielded MCDHF+RCI calculations identical to those performed with the LBL orbital set, making impossible the estimation of uncertainty bars. This limitation was caused by the direct adoption of MCDHF models originally developed for Ba II [15], which incorporate a larger electronic core. Furthermore, the inclusion of the Bohr-Weisskopf (BW) effect using the extracted values of Chakraborty et al. [31] did not improve our results.

For singly ionized radium (Ra II), we presented the first MCDHF calculations of HFS constants for the [Rn]6d ²$D_{3/2,5/2}$ levels, supplementing existing data for the ground and excited P levels. Despite computational limitations that restricted our correlation analysis to the $3d$ core orbital, we identified $6s$ substitutions as providing the most significant contributions. The NO basis yielded the best agreement with experimental values for the [Rn]7p ²$P_{3/2}^o$ level, with discrepancies below one percent. Our recommended values, complete with uncertainty estimates, encompassed existing theoretical predictions and established essential reference data for future spectroscopic studies. However, unlike other one-valence atomic structures like Sr II and Ba II, the sign change between the final and the experimental value for $A{(}^2{D}_{5/2})$ did not appear and, instead, a good agreement was seen with other accurate calculations [43,44,45]. A succinct MCHF calculation for the [Rn]6d ²$D_{3/2,5/2}$ levels still showed that this constant might be sensible to interference in the CSFs expansion with cancellation between the non-relativistic orbital and Fermi contact terms.

The neutral radium (Ra I) atom was examined to provide HFS constants critical for determining electric field gradients and interpreting molecular spectra, particularly for RaF. Our calculations for ²²³Ra I produced electric quadrupole B constants that aligned well with literature values while incorporating uncertainty quantification. The [Rn]7s6d ³$D_2$ and [Rn]7s7p ³$P_1^o$ levels, which are separated by merely $5.41$ ${\mathrm{cm}}^{-1}$ and are of particular importance in BSM physics due to their enhanced sensitivity to parity and time-reversal (PT) violating effects, were also studied. The theoretical bounds we provided will guide future experimental and theoretical investigations in this domain.

This HFS investigation marked the first successful application of the newly developed GRASPG [26] code to high-Z elements, demonstrating its capacity to efficiently handle complex atomic systems. Future enhancements, including the integration of partitioned correlation function interaction (PCFI) [56], promise to further improve the precision of hyperfine structure calculations by directly targeting the most significant orbital correlations. Experimental validation of our theoretical models would enable more refined uncertainty estimates. While higher-order substitutions remain computationally intensive, with configuration state function counts already approaching three million for Ra I, natural orbitals emerge as a particularly promising approach. A dedicated study exploring various natural orbital methodologies for second-row elements is currently in progress, highlighting the potential of this approach for atomic structure calculations. In addition, it would be interesting to test the accuracy of our methods for higher excited states, as they may be needed in astrophysical applications.