Configuration-interaction relativistic Hartree-Fock model

During the past decades, benefiting from the worldwide development of rare isotope beam facilities, a multitude of novel nuclear phenomena have been discovered, such as halo structures [1, 2], new magic numbers [3−5], and the so-called “island of inversion” [6, 7]. These discoveries have significantly extended the interests of nuclear physicists from stable nuclides to unstable ones. Meanwhile, these also present substantial challenges to our understanding of nuclear force and the methodology employed to address complex many-body nuclear system.

Configuration interaction shell model (CISM) has been one of the most powerful methods to understand the fundamental properties of nuclear structure, especially the sophisticated low-lying spectra [8]. Considering that the dimensionality of the configuration space grows exponentially with the increase of mass number, the CISM calculations are typically constrained to a truncated Hilbert space, namely the model space, where the core is frozen and only the degrees of freedom associated with valence nucleons are considered. As a result, an effective Hamiltonian should be introduced to take into account the configurations excluded from the model space. Practically, the success of shell model essentially depends on the choice of the configuration interactions [9−11].

Phenomenologically, the effective Hamiltonian, tailored to a specific model space, can be obtained by fitting the selected experimental data, e.g., the USD family of Hamiltonian for the $ sd $ shell [12, 13]. The fully microscopic derivation of the effective Hamiltonian based on the realistic forces remains a long-standing challenge for shell model calculations [11, 14−16]. Following time-dependent or time-independent perturbation theory [17, 18], the $ \widehat{Q} $-box resummation method, also known as the folded diagram approach [19, 20], has been widely used to derive the shell model Hamiltonian. It begins with the G-matrix or low-momentum nucleon-nucleon (NN) potential $ V_{\text{low-}k} $ to treat the strong short-range correlations induced by the realistic NN force [21, 22], then the $ \widehat{Q} $-box corrections and the folded terms are evaluated to account for the core polarization and other important correlations [23, 24].

Shell model calculations with such effective Hamiltonians present a satisfactory description of the binding energies and low-lying excitations of the nuclei with a few valence particles [25−27], while the agreement with experimental data deteriorates as more valence nucleons are added [15, 28, 29]. In fact, the main deficiency of such interactions lies in the monopole part [30], which is associated with the mean field of the specific nucleus. To improve the accuracy of theoretical descriptions, an empirical mass dependence of two-body interaction matrix elements (TBMEs) is typically adopted both for microscopic and phenomenological shell model calculations [9, 12, 31, 32]. Notably, the effects that account for the mass dependence of the TBMEs are nontrivial. On the one hand, the nuclear mean field changes consistently with the addition of valence nucleons, which indicates that the single-particle wave functions are mass dependent [9, 12]. By comparing the G-matrix obtained for different nuclei, Hosaka et al. pointed out that the G-matrix elements also exhibit a simple mass dependence [33]. Besides that, an important missing part of such monopole interactions is thought to be contributed by the three-body forces [30, 34, 35]. Considering the difficulty in handling three-body Hamiltonians for shell model diagonalization, the dominant contributions of the three-body forces are typically incorporated into a density-dependent two-body term [32, 36, 37].

In contrast to the CISM, with a given inter-particle interaction, the nuclear energy density functional (EDF) or the Hartree-Fock (HF) model can be applied to the study of nuclei across nearly the entire nuclide chart [38−40]. In the EDF or HF framework, the ground state properties of the specific nucleus are determined by a self-consistent mean field or HF potential, which is associated with the one-body density and can be derived through a variational procedure. Correspondingly, the essential physics encapsulated within the mean field or HF potential is manifested in the single-particle structure.

On the other hand, the empirical inter-particle interactions used in the EDF or HF model are obtained by fitting the properties of nuclei and nuclear matter at the saturation point [41−43]. These interactions can be considered as a kind of phenomenological G-matrix [38, 44], and their monopole and multipole parts are typically density-dependent to account for the nuclear in-medium effects [45, 46]. Notably, recent intensive works related to unstable nuclei have highlighted the importance of non-central parts of nuclear force, namely the spin-orbit and tensor forces, in the shell evolution [47−49], the isotopic shifts of charge radii [50−52], and the spin-isospin excitations [53−55], etc. Compared to the non-relativistic EDF, the spin-orbit force is naturally included in the relativistic scheme. Furthermore, the tensor force can be explicitly obtained via the exchange (Fock) terms in the relativistic Hartree-Fock (RHF) model [49, 56−62].

With several fitted parameters, the empirical EDFs or HF models and their extensions have achieved satisfactory descriptions of various nuclear phenomena from light to superheavy nuclei [63−68]. Naturally, one may expect that the density-dependent empirical interactions and the self-consistent single-particle basis given by the EDF or HF model could provide a reliable and general input for the shell model [69−77], where the effective Hamiltonians suffer from the issues of locality and mass dependence. In fact, Brown and the collaborators suggested that the mass-dependent monopole terms in the CISM calculations can be constrained by the EDF results, which contain the three-body, density-dependent and rearrangement contributions [75]. At the earlier time, Waroquier and the collaborators proposed a shell model approach based on the self-consistent Hartree-Fock single-particle basis with a universal Skyrme-type force [69]. The calculated two-particle (hole) excited spectra of doubly magic nuclei are in agreement with experimental results, with some adjusted or empirical single-particle energies. Their studies have shown that the effective TBMEs, derived from the Skyrme force SKE, and the calculated low-lying spectra strongly depend on the strength of the density-dependent zero-range force.

In this paper, we aim to integrate the strengths of the CISM and RHF theory to achieve a more robust method, namely the configuration-interaction relativistic Hartree-Fock (CI-RHF) model. On the one hand, we describe the nuclear system with an empirical but general Lagrangian, in which the spin-orbit and tensor forces are naturally included, and the nuclear in-medium effects are evaluated phenomenologically via the density dependencies of the meson-nucleon coupling strengths [56−58, 62, 78]. On the other hand, the RHF calculations provide an optimized single-particle basis to construct the model space, and in such a way, the energy of the frozen core can be calculated consistently. Following the spirit of the CISM, we derive the effective Hamiltonian self-consistently using the folded-diagram approach, and then diagonalize it in the model space to obtain the ground state and low-lying excitation properties. It is expected that the CI-RHF model can be applied to study the properties of a wide range of nuclei, without introducing additional parameters besides those well-defined in the phenomenological Lagrangian.

The paper is organized as follows. In Sec. 2, we present the general formalism of the CI-RHF model, and the calculation details of the effective Hamiltonian are provided in Appendix A. In Sec. 3, we first examine the convergence of intermediate-state excitations in calculating the effective TBMEs and illustrate the core-polarization effects. Subsequently, taking $ ^{16} {\rm{O}}$, $ ^{40} {\rm{Ca}}$ and $ ^{48} {\rm{Ca}}$ as examples, we investigate the beyond-mean-field corrections to the doubly magic nuclei. Furthermore, we obtain the properties of both ground states and low-lying excited states for neon isotopes with different model space, demonstrating the self-consistency of the CI-RHF model. Finally, a brief summary is given in Sec. 4.

Based on the meson propagated diagram of nuclear force [79], the Lagrangian of nuclear systems can be established considering the degrees of freedom associated with the nucleon field (ψ), the isoscalar scalar σ- and vector ω-meson fields, the isovector vector ρ- and pseudo-scalar π-meson fields, and the photon field (A) [56, 58, 80],

where $ \mathscr{L}_{I} $ contains the interactions between the nucleon and meson (photon) fields, and $ \mathscr L_\phi $ ($ \phi = \psi, \sigma, \omega, \rho, \pi $, and A) represents the free parts of nucleon and meson (photon) fields. Starting from the above Lagrangian, the Hamiltonian of nuclear systems can be derived from the Legendre transformation, and further the energy functional can be obtained following the standard procedure as detailed in Refs. [56, 62, 78, 80].

Restricted to the Hartree-Fock level, nucleons are treated as independent particles moving in an averaged field. Considering that collective correlations are significantly hindered by the shell structure, such an independent particle picture provides a satisfactory description for the ground states of closed-shell nuclei. Practically, the self-consistent RHF equation can be derived from a variation of the energy functional

where $ \varepsilon_{i} $ represents the single-particle energy, which excludes the rest mass of nucleon. The one-body kinetic energy $ T_{ii'} $ and the antisymmetric two-body interaction matrix elements $ \bar{V}_{iji'j'} $ can be obtained as follows,

where $ \Gamma_{\phi} $ and $ D_{\phi} $ denote the interaction vertex and meson (photon) propagators, respectively [58, 78, 80].

To further take into account the collective correlations explicitly, large-scale configuration mixing is introduced in the shell model approach. Starting from a complete set of single-particle basis, such as the self-consistent HF basis defined by the solutions of Eq. (2), the model space for the shell model, namely the P-space, can be constructed on top of a frozen core $ |\text{core}\rangle $,

where $ N_{\text{val}} $ is the number of valence nucleons, and $ |\text{core}\rangle $ corresponds to the wave function of a selected closed-shell core. It shall be stressed that the creation operation ($ c^\dagger $) is restricted within the valence space, and the outside space above it is unoccupied. Correspondingly, one can define the operator Q as the complement of P, namely $ Q=1-P $. In general, the nuclear system in a full Hilbert space can be described with the “bare” Hamiltonian $ H = H_{0} + V $, where $ H_{0} $ is the unperturbed Hamiltonian, and V contains the residual interactions. However, the effective Hamiltonian $ H^{\text{eff}} $ has to be considered when addressing the eigenvalue problem in a truncated model space,

There are several widely used methods to derive the effective interaction, such as the Krenciglowa-Kuo (KK) [81] and the Lee-Suzuki (LS) methods [82], which are actually the implementations of the folded diagram approach. It is worth mentioning that, in these methods, the model space is assumed to be degenerate, which is a relatively strong limitation, especially for nuclei far from the β-stability line. Recently, extended versions of the Krenciglowa-Kuo (EKK) and Lee-Suzuki (ELS) methods have been proposed [83, 84], which can be applied to the non-degenerate case. In the EKK scheme, the effective Hamiltonian can be obtained by the following iterations,

where the index n denotes the step of the iterations, and the Bloch-Horowitz Hamiltonian satisfies $ H^{\text{BH}}(E) = PH_{0}P + \widehat{Q}(E) $. The so-called $ \widehat{Q} $-box is defined as,

Notably, Eq. (6) can be interpreted as a Taylor series expansion of $ H^{\text{eff}} $ around the starting energy $ E_{0} $. Therefore, although the $ \widehat{Q} $-box and the folded terms, specifically the sum on the right side of Eq. (6), depend on the $ E_{0} $, the total effective Hamiltonian $ H^{\text{eff}} $ is independent of the starting energy.

For a single-shell valence space, where the unperturbed Hamiltonian can be considered degenerate, the effective interaction $ V^{\text{eff}} $ in the KK scheme can be obtained from Eq. (6) by assuming $ E_{0}P=PH_{0}P $,

In practical calculations, one cannot evaluate the $ \widehat{Q} $-box exactly. However, it can be further expressed as an expansion around $ H_{0} $:

For better understanding, Fig. 1 shows the valence-linked $ \widehat Q $-box diagrams, which contain the first-order $ \widehat{Q} $-box in V, and the second-order core-polarization diagrams induced by one and two particle-hole (p-h) excitations, namely the $ V^{\text{3p1h}} $ and $ V^{\text{4p2h}} $ corrections. At an early stage, the G-matrix augmented with the second-order core-polarization corrections served as a satisfactory effective interaction for shell model calculations [25, 26]. After the initial success, Barrett and Kirson found that large third-order terms tend to cancel the second-order contribution, indicating a lack of order-by-order convergence in this perturbation expansion [85]. However, when the summation of the folded diagrams extends to arbitrary orders as illustrated in Eq. (8), the differences between the effective interactions containing the $ \widehat{Q} $-box correction up to second order and those up to third order are small [14].

Figure 1.  Valence-linked $ \widehat{Q}$-box diagrams contributing to the effective two-body Hamiltonian.

In this paper, we consider the $ \widehat{Q} $-box corrections up to the second order, in which the contributions of ladder diagrams are not included to avoid double-counting [38, 86]. The folded terms with $ k\ge10 $ are terminated since the derivatives divided by the factorials are sufficiently small. Moreover, the core and single-particle energies are calculated consistently and all the corrections to the monopole part are neglected to keep the description of ground states of closed-shell nuclei almost unchanged. It is noteworthy that the rearrangement terms, induced by the density dependencies of the meson-nucleon coupling strengths, should be taken into account in deriving the effective Hamiltonian. The detailed derivations are presented in Appendix A.

It is worth emphasizing that all the CI-RHF calculations are carried out in a self-consistent manner. For a typical calculation, we outline the successive steps as follows. Firstly, starting from an empirical but general Lagrangian (1), spherical relativistic Hartree-Fock (RHF) calculations can be initially performed for specific nucleus. Then, aiming at interested observables, one can construct a suitable model space (4) by referring to the single-particle structure. Based on the RHF single-particle basis and the given Lagrangian, the core and single-particle energies are calculated consistently, then the effective Hamiltonian (6) can be derived using the folded-diagram approach. Finally, by diagonalizing the effective Hamiltonian in the model space, namely solving Eq. (5), one can obtain the properties of both the ground state and the low-lying states of the nucleus.

As the first applications of the CI-RHF model, we focus on the ground state and low-lying excitation properties of light and medium-mass nuclei. In the subsequent calculations, we utilize different Lagrangians, specifically PKA1 [57] and PKO2 [59]. The PKO2 contains the degrees of freedom associated with the σ-scalar (σ-S), ω-vector (ω-V), ρ-vector (ρ-V) and photon A-vector (A-V) couplings, while PKA1 further considers π-pseudo-vector (π-PV) and ρ-tensor (ρ-T) couplings. The single-particle basis is provided by the spherical RHF calculations. Furthermore, the valence space for even-even nuclei $ ^{18} {\rm{O}}$ and $ ^{18\sim 28} {\rm{Ne}}$ is restricted to the $ sd $ shell, composed of the $ 1d_{5/2} $, $ 2s_{1/2} $ and $ 1d_{3/2} $ orbits. In the current work, the excitations of the Dirac sea are negligible due to the large energy denominator. Additionally, the BIGSTICK code is used for the diagonalization of shell model Hamiltonians [87].

In general, one needs about 10 oscillator quanta excitations to obtain fully converged effective TBMEs [88−90]. However, in this work, the effective TBMEs are calculated with the self-consistent Hartree-Fock single-particle basis. Before realistic applications, it is necessary to examine the convergence of intermediate-state excitations, namely, to determine the energy truncation of particle states (p).

Using the Lagrangian PKA1 and considering the simple nucleus $ ^{18} {\rm{O}}$ as an example, we calculate the “bare” TBMEs, namely the first-order $ \widehat{Q} $-box diagram, and the effective TBMEs with different energy truncations. As shown in Fig. 2, the TBMEs are indicated by the red dots. The vertical coordinate corresponds to the results obtained with energy truncation of $ \varepsilon \leqslant 80 $ MeV, while the horizontal coordinates represent the “bare” TBMEs and the effective TBMEs with energy truncations of $ \varepsilon \leqslant 20 $, $ 40 $, and $ 60 $ MeV, for Figs. 2(a-d), respectively. It is evident that the corrections induced by the core-polarization effects are significant, and the vast majority of these effects are captured when $ \varepsilon \leqslant 20 $ MeV, as illustrated in Figs. 2(a) and 2(b). Furthermore, it can be observed that the calculated effective TBMEs tend to converge when $ \varepsilon \leqslant 40 $ MeV, and they remain almost unchanged from $ \varepsilon \leqslant 60 $ MeV to $ \varepsilon \leqslant 80 $ MeV. Therefore, considering the practical computational limitations, we consider all possible particle-hole excitations with particle states energy truncation of $ \varepsilon \leqslant 80 $ MeV for light and medium-mass nuclei. Additionally, it should be noted that the effective TBMEs given by the Lagrangian PKO2 are also converged under such truncation.

Figure 2.  (color online) Comparison of two-body interaction matrix elements $V_{abcd;J}$ (MeV) of $^{18}{\rm{O}}$ obtained with different particle state energy truncations using the Lagrangian PKA1. The diagonal line indicates the equivalence of TBMEs within different energy truncations. See text for more explanations.

The core-polarization corrections were initially calculated by Bertsch for $ ^{18} {\rm{O}}$ and $ ^{42} {\rm{Sc}}$ [91], and subsequently their effects were quantitatively evaluated in shell model calculations by Kuo and Brown [25]. It has been shown that the dominant contributions of the core-polarization diagrams originate from the quadruple component, providing a microscopic interpretation for the empirical long-range quadruple force [23]. Additionally, Shurpin and the collaborators found that the core-polarization corrections are significantly influenced by the single-particle basis [88]. As the first attempt, starting from the density-dependent Lagrangians PKA1 and PKO2, Fig. 3 illustrates the yrast band structure of $ ^{18} {\rm{O}}$ calculated with the self-consistent Hartree-Fock basis. In contrast to the calculations using “bare” PKA1 and PKO2 interactions, which present rather compressed spectra, fairly good agreement with the experimental data is achieved after considering the core-polarization corrections.

Figure 3.  (color online) The calculated yrast bands of $^{18}{\rm{O}}$ based on Lagrangians PKA1 and PKO2, compared with experimental results from [92]. Core-polarization corrections are not included in the case $V^{\text{bare}}$. Note that the energy of the ground state ($0^{+}_{1}$) is set to zero.

In fact, our results are similar to those obtained using microscopic G-matrix with a harmonic-oscillator basis [14, 25], but not with the Brueckner Hartree-Fock basis [88]. The latter typically yields reduced core-polarization corrections, due to the small radii for hole states and large ones for particle states, resulting in compressed spectra. Recently, a self-consistent Dirac Brueckner-Hartree-Fock approach has been implemented [86], providing a better description of the binding energy and charge radii of $ ^{16} {\rm{O}}$, which could improve this issue.

As pointed out by Brown and Kuo [23], the core-polarization corrections primarily arise from the long-range quadrupole component of neutron-proton isoscalar interactions. Consequently, for the stable nucleus $ ^{18} {\rm{O}}$, it can be observed that the core-polarization corrections, which are associated with the quadrupole properties of nuclei, are nearly identical for different Lagrangians PKA1 and PKO2.

However, the spectra obtained from the “bare” PKA1 and PKO2 exhibit notable differences. Specifically, these differences are determined by the isovector two-body matrix elements $ V^{\text{bare}}_{aaaa;J} $, where index a denotes the neutron orbit $ \nu 1d_{5/2} $. As an illustration, we extracted the contributions from various meson couplings to these TBMEs, as shown in Table 1. It is observed that the σ-scalar and ω-vector couplings provide crucial contributions to these TBMEs as described by PKO2, but less for PKA1. This is due the fact that the ρ-tensor coupling in PKA1 presents notable contributions, which results in distinct differences in the spectra obtained from “bare” PKA1 and PKO2.

	J	V$^{\text{bare}}_{aaaa;J}$	σ-S+ω-V	ρ-V	ρ-T	π-PV
PKA1	0	-1.92	-0.68	0.92	-2.13	-0.03
	2	-1.13	-0.82	0.24	-0.48	-0.07
	4	-0.63	-0.56	0.11	-0.21	0.03
PKO2	0	-1.35	-2.48	1.13	0.00	0.00
	2	-1.02	-1.32	0.29	0.00	0.00
	4	-0.66	-0.79	0.13	0.00	0.00

Table 1.  The calculated “bare” two-body interaction matrix elements $V^{\text{bare}}_{aaaa;J}$ (MeV) of $^{18}{\rm{O}}$ obtained using the Lagrangians PKA1 and PKO2. Here the index a represents the neutron orbit $\nu1d_{5/2}$. The contributions from different meson couplings are also presented. See the text for more explanations.

Due to the presence of the significant shell gaps, the ground state properties of doubly magic nuclei are frequently employed to constrain phenomenological Lagrangians within the mean field framework, such as PKA1 and PKO2. After considering the core polarization corrections and configuration mixing, the beyond-mean-field effects can be contained explicitly in the CI-RHF calculations. To elucidate the distinctions between the RHF and CI-RHF calculations, we present the calculated ground state binding energies and charge radii of doubly magic nuclei $ ^{16} {\rm{O}}$, $ ^{40} {\rm{Ca}}$ and $ ^{48} {\rm{Ca}}$ in Table 2. It should be noted that in the CI-RHF calculations, the valence space for $ ^{16} {\rm{O}}$ consists of the p shell and $ sd $ shell, while for $ ^{40} {\rm{Ca}}$, it includes the $ 2s_{1/2} $, $ 1d_{3/2} $, $ 1f_{7/2} $, and $ 2p_{3/2} $ orbitals. Additionally, we adopt the $ pf $ shell as valence space for $ ^{48} {\rm{Ca}}$.

	Binding Energy (MeV)				Charge radii (fm)
	$^{16}{\rm{O}}$	$^{40}{\rm{Ca}}$	$^{48}{\rm{Ca}}$		$^{16}{\rm{O}}$	$^{40}{\rm{Ca}}$	$^{48}{\rm{Ca}}$
Exp.	127.62	342.05	416.00		2.699	3.478	3.477
PKA1(HF)	127.16	341.48	416.04		2.802	3.530	3.495
PKA1(CI)	128.41	343.08	416.95		2.860	3.535	3.495
PKO2(HF)	127.09	340.70	415.25		2.673	3.434	3.457
PKO2(CI)	127.60	341.83	415.85		2.715	3.439	3.457

Table 2.  The calculated binding energies and charge radii of $^{16}{\rm{O}}$, $^{40}{\rm{Ca}}$ and $^{48}{\rm{Ca}}$ using the RHF and CI-RHF model with Lagrangians PKA1 and PKO2, compared to the experimental values (taken from [93, 94]).

It is observed that both the binding energies and charge radii calculated using the CI-RHF model (denoted as CI in Table 2) are larger than those given by the RHF calculations (denoted as HF in Table 2). For PKO2, the RHF results underestimate both the binding energies and the charge radii, as compared to the experimental values. However, after explicitly incorporating beyond-mean-field effects, the CI-RHF calculations improve the descriptions of $ ^{16} {\rm{O}}$, $ ^{40} {\rm{Ca}}$, and $ ^{48} {\rm{Ca}}$. Conversely, for the Lagrangian PKA1, the CI-RHF results for these nuclei are slightly overestimated, whereas the RHF results align more closely with experimental data. Even though, the beyond-mean-field corrections to the ground state bulk properties of doubly magic nuclei are relatively small, which indicates that the phenomenological Lagrangians PKA1 and PKO2 can serve as reasonable inputs for the CI-RHF calculations.

In general, with the addition of more valence particles, it is necessary to consider the mass dependence of the harmonic oscillator basis [12] and the G-matrix [33], as well as the effects of the three-body force [35]. For instance, an empirical mass dependence $ (18/A)^{0.3} $ is assumed for the USD Hamiltonians [13]. Within the CI-RHF framework, it is expected that the initial Hartree-Fock calculations can capture the effect of changing mean field, while the complicated many-body effects can be evaluated through the density-dependence of the meson-nucleon coupling strength [56, 57, 59]. Consequently, for the nuclei with more than two valence nucleons, no additional parameters are introduced in the CI-RHF calculations.

The neon isotopes span the entire neutron $ sd $ shell, provide ideal examples for studying low-lying spectra and their isotopic evolution. In the following, we first restrict the valence space to the $ sd $ shell and present the results for the even-even isotopes $ ^{18\sim 28} {\rm{Ne}}$.

As shown in Table 3, the binding energies of neon isotopes are described reasonably by the CI-RHF calculations, despite a slight overestimation by the Lagrangian PKA1. These results are consistent with the behavior observed for doubly magic nuclei, where the Lagrangian PKA1 similarly overestimates the binding energies, while the Lagrangian PKO2 exhibits better agreement with experimental data. On the other hand, the CI-RHF calculations provide a reasonable description of the isotopic shift of charge radii for neon isotopes, despite the values being systematically underestimated by PKO2, as illustrated in Fig. 4. Considering that the π-pseudo-vector and ρ-tensor couplings can enhance the quadruple correlations and give rise to more notable deformation, as indicated by previous studies [62, 78], it is not difficult to understand that PKA1 yields larger charge radii than PKO2 for the ground states of neon isotopes. Moreover, it can be observed that $ ^{20} {\rm{Ne}}$ presents the largest charge radius. This is understandable because $ ^{20} {\rm{Ne}}$ is strongly deformed, leading to a significant occupation of protons (occupancy number 0.62 and 0.42 for PKA1 and PKO2, respectively) in the more extended orbits $ \pi2s_{1/2} $ and $ \pi1d_{3/2} $. Notably, the charge radius of $ ^{24} {\rm{Ne}}$ is significantly smaller compared to that of neighboring isotopes, which may serve as an indicator of the neutron subshell $ N = 14 $.

	$^{18}{\rm{Ne}}$	$^{20}{\rm{Ne}}$	$^{22}{\rm{Ne}}$	$^{24}{\rm{Ne}}$	$^{26}{\rm{Ne}}$	$^{28}{\rm{Ne}}$
Exp.	132.14	160.65	177.77	191.84	201.55	206.86
PKA1	133.07	163.45	181.04	194.25	204.80	209.83
PKO2	131.87	160.67	179.99	193.99	201.80	206.56

Table 3.  The calculated binding energies (in MeV) for the even-even $^{18\sim28}{\rm{Ne}}$ isotopes using the CI-RHF model with Lagrangians PKA1 and PKO2, compared with the experimental values (taken from [94]).

Figure 4.  (color online) The calculated charge radii $R_{\text{c}}$ of neon isotopes using the CI-RHF model with Lagrangians PKA1 and PKO2, compared with experimental data. The dashed line represents the trend of experimental values (taken from [93]).

Besides the ground state properties, the CI-RHF model also aims at the low-lying excitations using the same Lagrangian. Taking PKA1 as the candidate, Fig. 5 further shows the low-lying spectra of even-even isotopes $ ^{18\sim 28} {\rm{Ne}}$. The reduced electric quadruple transition probabilities from the initial state $ \Psi_{i} $ to the final state $ \Psi_{f} $ are calculated as follows,

Figure 5.  (color online) Comparison of low-lying spectra of neon isotopes between experimental data [92, 95] and theoretical results obtained using the CI-RHF model with Lagrangian PKA1. The available reduced transition probabilities $ B(E2) $ (in e²fm⁴) within and between bands are also shown. Note that the energies of the ground states ($ 0^{+}_{1} $) are renormalized to zero.

It should be noted that $ e_{\tau_z} $ denotes the bare charge for protons or neutrons in the full space, while effective charges must be introduced in a truncated shell model space to account for the 2$ \hbar\omega $ core excitations. In fact, a reasonable quantitative description of $ B(E2) $ would require the inclusion of higher-order corrections to effective charge [11, 96]. Considering the limitation of practical calculations, in the current work, we use empirical effective charges of 1.29 and 0.49 for protons and neutrons, respectively, which are averages over orbits and masses in the $ sd $ shell obtained through a least-squares fit [97].

As illustrated in Fig. 5, it is clear that the calculated low-lying spectra are in good agreement with the experimental data for both the yrast and non-yrast bands. Specifically, the excited energies for $ ^{20,22} {\rm{Ne}}$ are slightly overestimated by the present calculation, while it yields a relatively compressed spectrum for $ ^{18} {\rm{Ne}}$ and even-even isotopes $ ^{24\sim 28} {\rm{Ne}}$. The trend of the reduced transition probabilities $ B(E2; 2^{+}_{1} \rightarrow 0^{+}_{1}) $ is consistent with the experimental results, despite the value for $ ^{28} {\rm{Ne}}$ being overestimated. In fact, earlier experimental measurement yielded a $ B(E2) $ with larger central value (53.8$ \pm $27.2 e²fm⁴) for $ ^{28} {\rm{Ne}}$ [98], as given by Monte Carlo shell model calculation [99]. Such inconsistency remains unexplained. On the other hand, it should be noted that the effective charges are orbit-dependent and can be influenced by the single-particle wave functions of the specific nuclide. As one moves away from stable nuclei toward neutron-rich side, the effective charges typically decrease [100]. Consequently, aiming at reliable analysis and prediction for transition properties, it is necessary to calculate the effective charges with the self-consistent Hartree-Fock single-particle wave functions in the future.

Furthermore, it can be observed that both $ ^{24} {\rm{Ne}}$ and $ ^{26} {\rm{Ne}}$ exhibit higher excited state energies $ E(2^{+}_{1}) $ compared to their neighboring isotopes. In general, large shell gaps result in relatively high excited state energies $ E(2^{+}_{1}) $ and low quadrupole collectivity. Along the isotopic chain, the $ B(E2; 2^{+}_{1} \rightarrow 0^{+}_{1}) $ values increase from $ ^{18} {\rm{Ne}}$ to $ ^{20} {\rm{Ne}}$ and $ ^{22} {\rm{Ne}}$ but are significantly reduced at $ N=14 $ and $ N=16 $. $ ^{20} {\rm{Ne}}$ and $ ^{22} {\rm{Ne}}$ are strongly deformed [92, 101], with a similar energy ratio $ E(4^{+}_{1})/E(2^{+}_{1}) $ of approximately 2.7, both experimentally and theoretically. However, the ratios $ E(4^{+}_{1})/E(2^{+}_{1}) $ for $ ^{24} {\rm{Ne}}$ and $ ^{26} {\rm{Ne}}$ are close to 2.0, indicating near pure collectivity of vibration. Therefore, combined with the isotopic evolution of charge radii shown in the lower panel of Fig. 4, a notable subshell closure at $ N=14 $ is indicated in neon isotopes, as well as a less pronounced one at $ N=16 $.

As discussed above, the CI-RHF calculations provide a reasonable description of even-even $ ^{18\sim28} {\rm{Ne}}$ nuclei when the valence space is restricted to the $ sd $ shell. However, extensive studies on atomic masses [102, 103], charge radii [104], and excitation spectra [105−107] have revealed that the neutron magicity $ N=20 $ disappears for nuclei in the “island of inversion”, indicating the intrusion of the $ pf $ shell. Consequently, it becomes necessary to obtain a multi-shell Hamiltonian to account for the cross-shell excitations. For instance, the SDPF-M effective Hamiltonian was constructed to explain the anomalous properties of neutron-rich nuclei around $ N=20 $ [99], with the model space comprising the $ sd $ shell and additional $ 1f_{7/2} $ and $ 2p_{3/2} $ orbits. More recently, taking the full $ sd $ shell and neutron $ pf $ shell as valence space, the effective Hamiltonian SDPF-U-MIX has been proposed [108], to study a broad range of nuclear properties. Notably, experimental guidance on single-particle energies is typically lacking for nuclei with extreme $ N/Z $ ratio, and it is challenging to determine the effective TBMEs in such a large model space by fitting experimental data.

In the CI-RHF framework, starting from the empirical but general Lagrangian, it is convenient to derive the multi-shell Hamiltonian by incorporating the EKK method [83, 84], since no additional parameters are needed to be readjusted. Taking the neutron-rich nuclei $ ^{30} {\rm{Ne}}$ and $ ^{32} {\rm{Ne}}$ as examples, we performed the $ sd $-$ pf $ cross-shell calculations with the effective Hamiltonian derived from the Lagrangian PKA1, in comparison to the results given by the effective Hamiltonian SDPF-U-MIX. In order to reduce the dimension of configurations, the CI-RHF calculations with PKA1 include only the neutron $ 1f_{7/2} $ and $ 2p_{3/2} $ orbits in $ pf $ shell as active orbits, in addition to the $ sd $ shell for proton and neutron. For comparison, the left panel of Fig. 6 also show the theoretical results of $ ^{30} {\rm{Ne}}$ calculated with PKA1 and USDB Hamiltonians [13], for which the model space is restricted to the $ sd $ shell.

Figure 6.  (color online) Cross-shell calculations of the excitation spectra for $^{30}{\rm{Ne}}$ and $^{32}{\rm{Ne}}$, compared with experimental data [92, 106, 109]. The left panel shows theoretical results with the valence space restricted to the $sd$ shell. Note that the energies of the ground states ($0^{+}_{1}$) have been renormalized to zero. See the text for more explanations.

It is observed that the excited energies $ E(2^{+}_{1}) $ of $ ^{30} {\rm{Ne}}$ are significantly overestimated with the neutron closed-shell configuration for both the PKA1 and USDB Hamiltonians. After considering the neutron cross-shell excitations, the CI-RHF calculations using the same Lagrangian PKA1 successfully reproduced the low-lying spectra, yielding results comparable to those obtained with the SDPF-U-MIX Hamiltonian [108]. It is noteworthy that the effective single-particle energy of the neutron $ 2p_{3/2} $ orbit, obtained with Lagrangian PKA1, is approximately 1 MeV lower than that of the $ 1f_{7/2} $ orbit, in contrast to the conventional spherical shell ordering. In fact, recent cross sections measurements for Ne isotopes also suggest that the valence neutrons in $ ^{31}\text{Ne} $ occupy low-l orbits [110, 111], leading to a halo structure in this nucleus [112]. It would be intriguing to investigate the mechanism driving shell evolution and halo formation around the $ N=20 $ “island of inversion” in future studies.

In this paper, we have proposed a method to consistently construct the configuration interactions for shell model calculations based on relativistic Hartree-Fock (RHF) theory, which is referred as the CI-RHF model. In the CI-RHF framework, initial RHF calculations capture the effect of changing mean field, and the density-dependent meson-nucleon coupling allows for a consistent treatment of complex many-body effects. Consequently, the core and single-particle energies, as well as the effective two-body matrix elements (TBMEs), can be derived from the same Lagrangian for each nuclide, which provides a robust solution for determining the effective Hamiltonian for any given model space. Based on that, the CI-RHF model is applicable for studying the properties of a wide range of nuclei with a few well-defined parameters in the phenomenological Lagrangian.

Taking the simple nucleus $ ^{18} {\rm{O}}$ as an example, we have checked the convergence of intermediate-state excitations in evaluating the effective interactions and studied the core-polarization effects, which are crucial for reproducing the structure of the yrast band. Furthermore, we have shown that the CI-RHF model provide a satisfactory description of both the bulk properties and low-lying excited spectra for doubly magic nuclei and even-even nuclei $ ^{18\sim 28} {\rm{Ne}}$. The isotopic evolution characteristics of charge radii and low-lying spectra indicate a notable subshell at $ N=14 $ in the middle of the $ sd $ shell for neon isotopes. Finally, by taking into account neutron cross-shell excitations from $ sd $ shell orbits to the $ 1f_{7/2} $ and $ 2p_{3/2} $ orbits, the CI-RHF calculations successfully reproduced the low-lying spectra for more neutron-rich nuclei $ ^{30} {\rm{Ne}}$ and $ ^{32} {\rm{Ne}}$ without readjusting any parameters. These promising results demonstrate that the CI-RHF model can provide equally reasonable descriptions of both stable and unstable nuclei using the same Lagrangian, which is crucial for making reliable predictions for nuclei near the proton or neutron dripline.

It should be noted that the CI-RHF calculations remain limited due to computational restrictions on the diagonalization of the shell model Hamiltonians, similar as met in the conventional CISM calculations. It would be interesting to diagonalize the shell model Hamiltonians in the nonorthogonal deformed basis by incorporating the Generator Coordinate Method [113−115], which could help us reach regions beyond the limit of the CISM and achieve a unified description for nuclei throughout the entire nuclide chart.

To account for the complicated in-medium effects in nuclear many-body system, the density dependence is typically introduced into phenomenological Lagrangians or nuclear forces. For instance, the meson-nucleon coupling strengths depend on the one-body density in the density-dependent relativistic Hartree-Fock model [56]. In the following, we consider a general density-dependent Hamiltonian to describe the many-body system in a full Hilbert space:

where the kinetic energy matrices $ T_{ii'} $ and antisymmetric two-body interaction matrices $ \bar{V}_{iji'j'} $ are expressed as Eq. (3) in the relativistic framework. Considering a general many-body wave function $ |\Psi\rangle $, the one-body and two-body nucleon density matrices can be obtained as follows:

Then one may express the energy functional E with respect to the one-body and two-body density matrices as a Taylor series expansion around the starting one-body density $ \rho^{0}\equiv\{\rho^{0}_{11},\rho^{0}_{12},\cdots\} $:

where $ \bar{V}^{'(kl)} $ represents the first-order derivative with respect to the one-body nucleon density matrix $ \rho_{kl} $. It is clear that, besides the normal terms, namely the linear terms of the density matrices $ \rho_{ii'} $ or $ \rho_{iji'j'} $, there are higher-order terms induced by the rearrangement effects in the energy functional (13). It should be emphasized that, before the variation of energy functional, one cannot evaluate the rearrangement effects exactly, which include all the terms up to infinite order.

In practical calculations, to accelerate the convergence of the Taylor series expansion, the one-body density of the Hartree-Fock ground state is adopted as the starting point of the expansion, namely $ \rho^{0} = \rho^{\text{HF}} $. On the other hand, for all nonlinear terms, one can approximate the one-body or two-body nucleon density matrices with $ \rho^{0}_{ij} $ or $ \rho^{0}_{iji'j'} $ to ensure the linear density-dependence of the energy functional. Consequently, the energy functional is approximated as follows:

where the zero-order rearrangement term $ E^{\text{re}}_{0} $ and kinetic energies containing the rearrangement term corrections $ T^{\text{re}}_{ii'} $ are defined as:

It is clear that, in this case, only the monopole part of the Hamiltonian is modified by the rearrangement terms. In the future, one may evaluate the rearrangement effects more precisely after the variation of the energy functional.

The energy eigenvalue equation can be obtained from a variation of the energy functional (14) with respect to the coefficients $ C_{n} $, namely the probability amplitude of a specific basis $ \Psi_{n} $. Considering that the shell model calculations are typically restricted to a truncated model space (4), the effective Hamiltonian must be introduced to account for the configurations beyond the model space,

where the one-body and two-body nucleon density matrices from the initial state $ \Psi_{n} $ to final state $ \Psi_{m} $ are defined as:

It is clear that the input for shell model calculations is the effective Hamiltonian matrix elements $ H^{\text{eff}}_{mn} $, which include the core and single-particle parts as well as the residual two-body interaction matrix elements. For conventional shell model, the core energy $ E^{\text{eff}}_{\text{c}} $ is taken from experimental data, and the single-particle energies $ \varepsilon^{\text{eff}}_{ii'} $ are typically determined by fitting or experimental data.

In this work, starting from a general Lagrangian (1), both the core and single-particle energies, along with the rearrangement terms, are calculated consistently based on the Hartree-Fock single-particle basis:

where the $ N_{\text{c}} $ represents the number of nucleons in the core. It should be noted that the given Lagrangian (1) is obtained by fitting the binding energies of several doubly magic nuclei at the Hartree-Fock level, which means higher corrections to the monopole part of energy have already been absorbed into the phenomenological Lagrangian. Consequently, the higher-order corrections to the core energy, such as the contributions of ladder diagrams, are neglected to avoid double-counting. Moreover, the $ \widehat{S} $-box corrections to single-particle energies are also not considered to keep the description of ground state properties of closed-shell nuclei almost unchanged. For the two-body part, we need to calculate the “bare” TBMEs and the second-order core-polarization corrections, as well as the folded terms. Similarly, the corrections to the monopole part of TBMEs are removed.

In the following, we focus on the calculation details of effective interactions (8). Imposing spherical symmetry, the single-particle states ψ are specified by a set of quantum numbers $ \alpha \equiv (a, m_{a}) \equiv (\tau_{a}, n_{a}, l_{a}, j_{a}, m_{a}) $. In the coupled representation, the particle-particle matrix elements can be obtained as follows:

where $ \bar{V}_{\alpha\beta\gamma\delta} $ is the uncoupled TBMEs, and $ \hat{J}\equiv\sqrt{2J+1} $. Noted that the notation used in this paper is standard [116], including the Clebsch-Gordan coefficients and $ 6j $ symbols.

Furthermore, the $ V^{\text{3p1h}} $ correction in the $ \widehat{Q} $-box consists of the following four terms:

where $ E_{0} $ represents the starting energy in $ \widehat{Q} $-box, and the factors $ f_{3p1h} $, corresponding to different intermediate-states, are expressed as:

Moreover, the particle-hole matrix elements in Eq. (20) can be obtained as follows:

Additionally, the different energy denominators in Eq. (20) can both be simplified as $ \varepsilon_{h} - \varepsilon_{p} $, when the model space is degenerate, namely $ E_{0}=2\varepsilon_{a}=2\varepsilon_{b}=2\varepsilon_{c}=2\varepsilon_{d} $.

The $ V^{\text{4p2h}} $ correction in the $ \widehat{Q} $-box can be evaluated in the coupled representation:

Similarly, in a degenerate model space, the energy denominator can be further simplified as $ \varepsilon_{h}+\varepsilon_{h'} - \dfrac{1}{2}(\varepsilon_{a}+\varepsilon_{b}+\varepsilon_{c}+\varepsilon_{d}) $. Moreover, it should be pointed out that, in this work, we removed the contributions from the isoscalar two-particle cross-shell excitations, which could lead to divergent binding energies with the given Lagrangian (1). Additionally, the center-of-mass correction is evaluated with the projection-before-variation method [117].