Relativistic Hartree-Fock model for axial-symmetric nuclei with quadruple and octupole deformations

Over the decades, octupole deformation in nuclei has been an active research field in nuclear physics [1, 2]. It is widely acknowledged that the shape of ground state for most deformed nuclei is symmetric under space inversion and therefore the dominant intrinsic deformation is of quadruple character. Even though, there are a handful of nuclear systems, in which the reflection symmetry is not preserved, giving rise to an octupole deformation or pear-like shape. Microscopically, the occurrence of the octupole deformation is attributed to the coupling between two single-particle orbits near the Fermi surface, which differ by $ \Delta l $ = 3 and $ \Delta j $ = 3. Such nuclei are located close to proton and neutron numbers of 34, 56, 88 and neutron number of 134 [1].

In recent years, the occurrences of octupole deformation have been confirmed in various regions of nuclear chart, such as $ ^{224} {\rm{Ra}}$ [3], $ ^{144} {\rm{Ba}}$ [4] and $ ^{146} {\rm{Ba}}$ [5], $ ^{228} {\rm{Th}}$ [6], and $ ^{96} {\rm{Zr}}$ [7]. More interestingly, a multitude of nuclear novel phenomena related to the octupole deformation were also discovered, such as the alternating-parity rotational bands, low-lying $ 1^- $ and $ 3^- $ states, and the $ E3 $ transitions [1]. In parallel, considerable theoretical efforts have been devoted to understanding these colorful phenomena, including the macroscopic-microscopic models [8−10], the self-consistent mean-field models based on nuclear density functional theory [11−14], the interacting boson model [15−19], the geometrical collective models [20, 21], and the cluster models [22, 23].

As one of the representatives, the relativistic mean field (RMF) theory based on the meson exchange diagram of nuclear force [24], also referred as the covariant density functional theory (CDFT), has achieved great success in exploring the structural properties of nuclei which spread over almost the whole nuclear chart [25−31]. In particular, due to the covariant representation, the RMF theory owns the advantages of providing the natural interpretation on the strong spin-orbit couplings [32, 33] and the origin of the pseudo-spin symmetry [34, 35]. However, in the RMF theory, the Fock terms, the inseparable part of the meson exchange diagram of nuclear force, are dropped for the sake of simplicity. As a result, the important degrees of freedom associated with the π- and ρ-tensor couplings, that work only through the Fock terms, are missing within the RMF scheme.

In the past decades, the relativistic Hartree-Fock (RHF) theory and its extension, which incorporate the density-dependent meson-nucleon coupling strengths [36−38], have achieved comparable accuracy as popular mean-field models in describing various nuclear structure properties, using the RHF Lagrangians PKOi(i = 1,2,3) [36, 39] and PKA1 [37]. Because of the Fock terms, notable improvements have been achieved in the self-consistent description of nuclear structure properties, such as shell evolutions [39−41], symmetry energy [42, 43], new magicity [44, 45], pseudo-spin symmetry [46, 47], and spin and isospin excitations [48, 49]. In particular, the tensor force, an important ingredient of nuclear force, has been considered naturally via the Fock terms [50, 51], e.g., by the π- and ρ-tensor couplings. Not only limited within the RHF scheme, it was revealed that the tensor force can play an important role in nuclear shell evolution [39, 41], symmetry energy [50, 52] and nuclear excitations [53, 54].

Recently, using the Dirac Woods-Saxon (DWS) basis [55], both the RHF and relativistic Hartree-Fock-Bogoliubov (RHFB) theories [36−38] were extended for quadruple-deformed nuclei with axial symmetry, giving the D-RHF and D-RHFB models [56, 57], respectively. It shall be stressed that the Dirac equations become integro partial-differential ones, which are solved by expanding the Dirac spinors on the spherical DWS basis [56, 57]. In fact, the DWS basis [55] owns the advantages of providing reasonable asymptotic behaviors for the wave functions of unstable nuclei and a microscopic insight into deformed nuclei. Taking $ ^{20} {\rm{Ne}}$ as an example, it was illustrated that the tensor force carried by the π-coupling plays an essential role in determining the deformed single-particle structures [56]. Furthermore, the D-RHFB model with PKA1 reproduces well both the even-parity ground state and halo structure of $ ^{11} {\rm{Be}}$, from which a microscopic picture of deformed halo was indicated in terms of the DWS basis [58].

For a long time, it remains an interested topic to understand the occurrence of stable octupole deformation in nuclear chart, in particular for the role of nuclear tensor force. As encouraged by the successes of both D-RHF and D-RHFB models, we are motivated to extend the RHF model by incorporating both quadruple and octupole deformations, giving the O-RHF model in this work. Similarly, the Dirac equation will be solved by expanding the Dirac spinor on the spherical DWS basis, regarding its advantages. Due to the complicated Fock terms, the axial symmetry is still imposed for octupole-deformed nuclei. Meanwhile, the pairing correlations are treated within the BCS scheme by using the central part of the finite-range Gogny force D1S [59] as the pairing force.

The contents are organized as follows. In Sec. 2, the general formalism of the O-RHF model is presented by using the spherical DWS basis. Sec. 3 presents the results and discussions, including the space truncation, the convergence check, and the description of the nucleus $ ^{144} {\rm{Ba}}$, where the effects of the Fock terms and the tensor force were emphasized. Finally, a summary is given in Sec. 4.

This section will provide a brief introduction of the general formalism of the relativistic Hartree-Fock (RHF) theory. In order to provide readers with a comprehensive understanding, some details of the RHF energy functional and Dirac equations for the nuclei with the reflection asymmetry and axial symmetry will be also presented, by using the spherical DWS basis.

Under the meson-exchange diagram of nuclear force, the Lagrangian of nuclear systems, the starting point of the RHF theory, can be obtained by considering the degrees of freedom associated with nucleon (ψ), the isoscalar σ- and ω-mesons, the isovector ρ- and π-mesons, and the photon field ($ A^\mu $). Starting from the general Lagrangian as detailed in Ref. [56], the Hamiltonian can be derived via the Legendre transformation as,

where the kinetic energy term (T) and the potential energy one ($ V_\phi $) read as,

Here x and $ {\boldsymbol{x}}' $ are used to denote the space coordinate vectors, and ψ standing for the Dirac spinor of nucleon. In the potential energy terms $ V_\phi $, ϕ represents the two-body interaction channels, including the scalar (σ-S), the vector (ω-V, ρ-V and A-V), the tensor (ρ-T), the vector-tensor (ρ-VT), and the pseudo-vector (π-PV) couplings. Accordingly, the interaction vertex $ \Gamma_{\phi}(x,x') $ are of the following form [56, 60],

In the above expressions, the symbol $ \vec{\tau} $ represents the isospin vector with $ \tau_3 $ for the projection, $ x = (t, {\boldsymbol{x}}) $ with the bold type for space vectors, and M and $ m_\phi $ for the masses of nucleon and mesons, respectively. After neglecting the retardation effects, namely ignoring the time component of the four-momentum carried by the mesons and photon, the propagators $ D_{\phi}({\boldsymbol{x}}, {\boldsymbol{x}}') $ in the potential terms $ V_\phi $ read as,

where ϕ represents the σ-, ω-, ρ- and π-meson fields.

Restricted on the level of the mean field approach, the no-sea approximation is introduced as usual, which is amount to neglect the contributions from the negative energy states [60]. Thus, the nucleon field operator ψ in the Hamiltonian can be quantized as,

where $ x = (t,{\boldsymbol{x}}) $, the index α denotes the positive-energy solutions of Dirac equation, and $ \psi_\alpha({\boldsymbol{x}}) $, $ \varepsilon_\alpha $ and $ c_\alpha $ ($ c_\alpha^\dagger $) represents the single-particle wave functions, energies and related annihilation (creation) operators, respectively. Subsequently, nuclear energy functional E can be obtained from the expectation of the Hamiltonian with respect to the Hartree-Fock ground state $ \left|{\Phi_0}\right> $ [56, 60],

where $ \left|{-}\right> $ represents the vacuum and A is nuclear mass number. For the two-body interaction $ V_\phi $, the above expectation leads to two types of contributions, namely the direct Hartree and the exchange Fock terms. If only the Hartree terms are taken into account, this leads to the so-called RMF theory. If both Hartree and Fock terms are considered, it gives the RHF theory.

Taking the variational of the energy functional E, an integro-differential Dirac equation can be obtained as,

where the single-particle Hamiltonian h consists of three parts, i.e., the kinetic term $ h^{\text{kin}} $, the local potential term $ h^D $, and the non-local terms $ h^E $ contributed by the Fock terms. The details can be found in Ref. [56]. It shall be commented that the integral Dirac equation is hard to be solved directly in coordinate space, especially when considering nuclear deformations.

In the popular RHF Lagrangians PKOi ($ i = 1,2,3) $ and PKA1, the density dependencies are introduced for the meson-nucleon coupling strengths $ g_\phi $ ($ \phi = \sigma, \omega^\mu, \vec{\rho}^\mu $) and $ f_{\phi'} $ ($ \phi' = \vec{\rho}^\mu, \vec{\pi} $) in equation (4), which are taken as functions of nucleon density $ \rho_b = \bar\psi\gamma^0\psi $ [36, 37]. It is noteworthy that the density dependencies of the coupling strengths lead to an additional contribution to the self-energy $ \Sigma_0 $, i.e., the rearrangement terms $ \Sigma_R $ [36−38], which should be taken into account to preserve the energy-momentum conservation [61]. Together with the Fock terms, the density dependencies of the coupling strengths remarkably increase the numerical complexity, when the deformation degree of freedom is involved [56, 57].

As limited by complicated Fock terms, we restrict ourselves under the axial symmetry for the deformed nuclei with quadruple and octupole deformations. Due to the reflection asymmetry, the parity does not remain as a good quantum number, and the projection m of the total angular momentum remains a good one. In order to label the single-particle states, the symbol $ \alpha = (\nu m) $ is introduced, by using the index ν to denote the orbits with different energies in a given m-block. In the following, we will derive the RHF energy functional (7) and solve the integro-differential equation (8) by expanding the single-particle wave functions on the spherical DWS basis [55].

In expanding the wave functions, both positive and negative energy states in the spherical DWS basis shall be considered. This does not conflict with the no-sea approximation, which is considered in quantizing the nucleon field operator ψ [see Eq. (6)] by neglecting the contributions from the Dirac sea, e.g., to nucleon densities. Conversely, considering the negative energy states of the spherical DWS basis is demanded by the mathematical completeness. It was proved that the DWS states with negative energies can be essential for correctly describing nuclear matter [62] and nuclear structure [57]. In terms of the DWS basis, the expansion of the wave function of the orbit $ \alpha = (\nu m) $ reads as,

where $ \psi_{am} $ represents the wave functions of the spherical DWS basis states denoted by $ a = (n\kappa) $, where n represents the principle number and $ \kappa = \pm (j+1/2) $ with $ j = l\mp1/2 $. For the sake of convenience, the expansion coefficients $ C_{a,\alpha} $ are set as real number. In order to provide compact expressions, the sum over the principal number n in above equation is absorbed by $ \psi_{ \nu \kappa m} $ as,

with $ {{\cal{G}}}_{\alpha\kappa} = \sum\nolimits_n C_{n\kappa,\alpha} G_{n\kappa} $ and $ {{\cal{F}}}_{\alpha\kappa} = \sum\nolimits_n C_{n\kappa,\alpha} F_{n\kappa} $, and $ \Omega_{\kappa m} $ (also referred as $ \Omega_{j m}^l $) is the spherical spinor [63]. It shall be noted that the sum in the expansion (9) contains the DWS state $ \left|{n\kappa m}\right> $ with both even- and odd-parity, due to the reflection asymmetry.

For the propagators (5), the following expansion in the spherical coordinate space ($ r,\vartheta,\varphi $) [56, 60, 63] is used to derive the energy functional,

where $ {\bf{\Omega}} = (\vartheta,\varphi) $, and the index $ \lambda_y $ denotes the expansion terms of the Yukawa propagator, with $ R_{\lambda_y\lambda_y} $ representing the radial part. For details, please refer to Ref. [56].

Due to the density dependence, the coupling strengths in the interaction vertex (4) are expanded in a series of spherical harmonic functions $ Y_{\lambda_p0} $,

where $ g_\phi $ represents the coupling strengths $ g_\sigma $, $ g_\omega $, $ g_\rho $ and $ f_\pi $ in Eq. (4), and the index $ \lambda_p $ denotes the expansion terms of the density-dependent parameters in the effective Lagrangians. In contrast to the D-RHF and D-RHFB model [56, 57] which preserve the reflection symmetry, here $ \lambda_p $ does not remain as only even integers due to the octupole deformation.

In this work, we concentrate on the RHF Lagrangians PKOi ($ i = 1,2,3 $) [36, 39]. For PKO2, that share the same degrees of freedom as the popular RMF models, it contains the σ-S, ω-V, ρ-V and A-V couplings. In addition to that, PKO1 and PKO3 take the π-PV coupling into account. Thus, the energy function E can be expressed as,

where the kinetic energy $ E^{\text{kin}} = \left<{\Phi_0}\right| T\left|{\Phi_0}\right> $, and $ E_\phi^D $ and $ E_\phi^E $ are respectively the Hartree and Fock terms of the potential energies, namely $ E_\phi = \left<{\Phi_0}\right|V_\phi\left|{\Phi_0}\right> $ with $ \phi = \sigma $-S, ω-V, ρ-V, π-PV and A-V.

With the expansions (9), (11) and (12), one can derives the energy functional terms $ E^{\text{kin}} $, $ E_{\phi}^{D} $ and $ E_{\phi}^E $ of Eq. (13) in great details. The detailed expressions for these terms are not provided in this paper, as their forms are identical to those given in equations (37), (40) and (41) of Ref. [56], with the exception of the following points. As a consequence of the reflection asymmetry, the sums with respect to $ \lambda_d $ in Eq. (40) of Ref. [56] will be released with both even and odd numbers. For the Fock terms, the restrictions on the parity index, π and $ \pi' $, in Eq. (41) and Appendix 3 of Ref. [56], are not valid for octupole-deformed nuclei. In order to avoid misunderstanding, the symbols $ \bar{{\mathscr D}}_{\kappa_1m_1;\kappa_2m_2}^{\lambda \bar\mu} $ and $ \bar{{\cal{Q}}}_{\kappa_1m_1;\kappa_2m_2}^{\lambda\bar\mu\sigma} $ appearing in the non-local self-energies, namely Eqs. (A5) and (A7) in Ref. [56], are redefined for octupole-deformed nuclei as,

where $ \pi_1 = 0 $ and $ 1 $ for $ \kappa_1 $ quantity representing the blocks with the even- and odd-parity, respectively.

For the open-shell nuclei, the pairing correlations are considered within the BCS scheme, and the central part of the finite-range Gogny force D1S [59] is adopted as the pairing force, similar as we did in the D-RHF model [56]. While, different from the D-RHF model where only the contribution of the main component $ J = 0 $ is included [56], the full contributions from all the J-components are considered in this work.

Since the Dirac spinors are expanded on the spherical DWS basis, the variation of the RHF energy functional (13) with respect to the expansion coefficient $ C_{a,\alpha} $ leads to a series of eigenvalue equations as,

where the symbol $ \widehat C_\alpha $ is a column matrix composed of the coefficients $ C_{a,\alpha} $. By diagonalizing the Hamiltonian $ H_{m} $ for given m-block, the eigenvalue, i.e., the single-particle energy $ \varepsilon_\alpha $ can be obtained, as well as the eigenvector $ \widehat C_\alpha $ for deformed single-particle orbit $ \alpha = (\nu m) $.

Similar as the single-particle Hamiltonian in Eq. (8), the matrix $ H_{m} $ in Eq. (16) consists of three parts, i.e., the kinetic $ H^{\text{kin}} $, local $ H^D $ and non-local $ H^E $ terms,

where the subscript $ (m) $ is omitted. These terms have the same form as the ones in Eqs. (54), (55) and (56) of Ref. [56]. Consistently, $ \lambda_d $ in Eq. (55) of Ref. [56] is also released with both the even and odd numbers, and the restrictions on the parity index π in Eq. (56) of Ref. [56] do not apply to octupole-deformed nuclei.

In order to provide some standard references for future applications of the O-RHF model, we firstly test the space truncations of the spherical DWS basis and the cutoff on $ \lambda_p $ in Eq. (12), using $ ^{144} {\rm{Ba}}$ as an example, which is known to exhibit octupole deformation [4]. In determining the spherical DWS basis, it is in general precise enough to set the spherical box size as $ 20 $ fm with a radial mesh step of $ 0.1 $ fm, and the continuum states with both positive and negative energies are discretized within such a spherical box. Furthermore, with determined space truncations and the cutoff of $ \lambda_p $, we analyze the structural properties of $ ^{144} {\rm{Ba}}$, using the RHF Lagrangians PKOi ($ i = 1,2,3 $) [36, 39] and RMF one DD-ME2 [64]. Particular efforts will be devoted to the quadruple and octupole deformation effects, and the role of the tensor force component carried by the π-PV coupling.

In the O-RHF calculations, two independent truncations shall be subjected to rigorous examination. These are the cutoff on the quantities $ (n\kappa) $ of the spherical DWS basis [Eq. (9)] and the expansion term $ \lambda_p $ of the density-dependent coupling strengths (12). It is noteworthy that the expansion terms $ \lambda_y $ in the propagators (11) are truncated automatically by the selected DWS basis states and expansion terms $ \lambda_p $ of coupling strengths, as detailed in Eq. (A2) of Ref. [56]. Due to the reflection asymmetry, it is necessary to incorporate both odd and even $ \lambda_p $ terms in expanding the coupling strengths (12). This differs from the D-RHF and D-RHFB models [56, 57], in which only even $ \lambda_p $ terms are needed for quadruple deformed nuclei with the reflection symmetry. In general, it is sufficient to set $ \lambda_p = 0, 1, 2, \cdots, 8 $ for the majority of nuclei.

Practically, the maximum value of m, designated as $ m_{\max} $, depends on the nucleus under consideration. Let's set $ K_m $ to denote the number of κ-blocks included in expanding the Dirac spinors with $ m_{\max} $. Consistently, the $ m_{\max} $ and $ K_m $ together determine the maximum value of $ |\kappa| $, i.e., $ \kappa_{\max} = m_{\max} + K_{m} - 1/2 $. Thus, for an arbitrary Dirac spinor $ \psi_{\nu m} $, the κ-quantities in the expansion (9) read as $ \kappa = \pm (m + 1/2), \pm(m + 3/2), \cdots, \pm\kappa_{\max} $, including the DWS states with both even- and odd-parity due to the reflection asymmetry. In principle, the $ m_{\text{max}} $ value can be determined referring to the conventional shell-model picture. For $ ^{144}_{\ 88} {\rm{Ba}}$$ _{56} $, the $ m_{\max} $ value is determined as $ 15/2 $ for both neutron and proton orbits, in accordance with the magic number $ 184 $ as predicted in the conventional shell-model picture. For such large $ m_{\max} $ value, it is accurate enough to set $ K_{m} = 4 $. As a comment, for the cases with fairly large deformation, some high-j orbits may penetrate the well-known major shells, and a rigorous test shall be conducted, e.g., by increasing the $ m_{\max} $ value.

In a manner analogous to the D-RHFB model [57], the maximum values of the principal number n for each κ-block in the expansion (9) are determined by the energy cutoff $ E_{\pm}^C $. This cutoff is defined by the single-particle energy ε of the spherical DWS state, with the sign $ + $ or $ - $ indicating the positive or negative energy states, respectively. Specifically, the states with positive/negative energies, that $ E_{+}^C+ M>\varepsilon> E_{-}^C-M $, are considered in the expansion (9). It is important to exercise caution when testing the energy cutoff $ E_\pm^C $.

As aforementioned, the Fock terms increase largely the numerical complexity. As a result, the O-RHF calculations become extremely time consuming. Similar as in developing the D-RHF and D-RHFB models, we applied the GPU parallel speedup technology to calculate the nonlocal Fock mean fields and the pairing matrix elements. This notably decreases the computing time, which makes the extension to superheavy region becoming expectable. Even though, for a single calculation of the nuclide $ ^{144} {\rm{Ba}}$ with $ m_{\max} = 15/2 $, it still costs about 19.1 and 51.9 hours respectively for PKO2 and PKO3, by fully using eight GPU units (Tesla A100-40G) and eighty CPU cores. When the $ m_{\max} $ value increases, the time cost grows up by times. Fortunately, it will not increase much when enlarging the $ E_\pm^C $ values.

Figures 1a-1b and 1c-1d present the tests of the energy cutoff $ E_\pm^C $ for the total energy E (MeV) and the deformation $ (\beta_2,\beta_3) $ of $ ^{144} {\rm{Ba}}$, respectively. The results were obtained by employing the RHF Lagrangians PKO2 and PKO3, and the RMF one DD-ME2. Figures 1a and 1c illustrate the convergence with respect to $ E_{+}^C $, with $ E_{-}^C $ set as $ 0 $ MeV, and Figs. 1b and 1d demonstrate the convergence with respect to $ E_{-}^C $, in which $ E_{+}^C $ is fixed as $ 350 $ MeV. It is clear that both the total energy E [Fig. 1a] and the deformation $ (\beta_2, \beta_3) $ [Fig. 1c] are converged when $ E_{+}^C > 200 $ MeV for $ ^{144} {\rm{Ba}}$. When changing the $ E_-^C $ value from zero to $ -50 $ MeV, namely considering the continuum negative energy states of the DWS basis in the expansion (9), the total energy E remains almost unchanged for DD-ME2 and PKO3, while PKO2 shows a slight but visual change, as illustrated in Fig. 1b.

Figure 1.  (Color Online) Total energy E (MeV) [plots a and b] and deformations $(\beta_2, \beta_3)$ [plots c and d] for $^{144}$Ba with respect to the positive ($+$) and negative ($-$) energy cutoff $E_{\pm}^C$ (MeV) in expanding the Dirac spinors $\psi_{\nu m}$. The results are calculated by DD-ME2, PKO2 and PKO3.

In contrast to the total energy E, the deformations $ (\beta_2, \beta_3) $ appear to be more sensitive to the negative energy cutoff $ E_-^C $. As illustrated from Figs. 1c to 1d, an extended truncation of the negative energy states of the DWS basis is required to give converged deformations, especially for the octupole deformation $ \beta_3 $ given by the RHF Lagrangians PKO2 and PKO3. Despite these evident alterations, both deformations $ \beta_2 $ and $ \beta_3 $ converge rapidly with respect to the negative $ E_{-}^C $ values, namely $ E_-^C<-100 $ MeV.

In order to understand the convergence of the octupole deformation $ \beta_3 $ as illustrated in Figs. 1c and 1d, Fig. 2 shows the proton single-particle spectra of $ ^{144} {\rm{Ba}}$ calculated by DD-ME2 and PKO2 with $ E_-^C = 0 $ and $ -50 $ MeV, where the Fermi levels are denoted by $ E_{\text{F}} $, and $ m_\nu $ for the deformed single-particle orbits. For the RHF Lagrangian PKO3, the results are not shown due to similar description as PKO2.

Figure 2.  (color online) Proton single-particle spectra of $^{144}$Ba given by DD-ME2 (plot a) and PKO2 (plot b) with negative energy cutoff $E_{-}^C = 0$ MeV and $-50$ MeV, in which the positive energy cutoff is fixed as $E_{+}^C = 350$ MeV. The ultra-thick bars represent the occupation probabilities of the orbits $m_\nu$, with the index ν representing the νth state in m-block, and $E_F$ denotes the Fermi levels.

Consistent with relatively small alterations from Figs. 1c to 1d, the proton single-particle spectra given by DD-ME2 remain almost unchanged when the $ E_-^C $ value varying from zero to $ -50 $ MeV, see Fig. 2a. Conversely, as illustrated in Fig. 2b, the PKO2 results exhibit remarkable alterations. In fact, similar systematics were observed for the neutron single-particle spectra, which are not shown. Namely, the results given by DD-ME2 remain almost unchanged, while the PKO2 results manifest rather distinct alterations, when the $ E_-^C $ value changes from zero to $ -50 $ MeV.

In this work, the deformed single-particle orbits are expanded on the spherical DWS basis. This provides us an insight into the microscopic properties of the octupole nucleus $ ^{144} {\rm{Ba}}$, as well as the obvious alterations in Fig. 2b. Table 1 lists the proportions (%) of the main DWS waves for the proton valence orbits $ 3/2_8 $ and $ 1/2_{13} $ of $ ^{144} {\rm{Ba}}$, which are all the DWS basis states with positive energies. In contrast to the DD-ME2 results, the proportions of the main spherical DWS waves given by PKO2 undergo a notable redistribution from $ E_-^C = 0 $ to $ E_-^C = -50 $ MeV, especially for the $ 1h_{11/2} $ and $ 2d_{5/2} $ waves, which are essential for the occurrence of the octupole deformation. This indicates that the completeness regarding the DWS basis states with negative energies can be essential for correctly describing the octupole deformation.

orbits	$E_{-}^C$	$2d_{5/2}$	$1g_{7/2}$	$1g_{9/2}$	$1f_{5/2}$	$1h_{11/2}$
$3/2_{8}$	0	5.9	65.5	5.4	2.4	11.5
	−50	6.6	59.5	2.6	6.9	12.6
$1/2_{13}$	0	21.8	24.4	7.7	0.7	27.1
	−50	25.8	20.1	2.6	2.1	30.1
$3/2_{8}$	0	2.1	43.8	25.3	4.0	2.4
	−50	6.5	53.8	6.2	3.5	20.2
$1/2_{13}$	0	2.9	23.5	37.7	3.3	5.0
	−50	24.7	2.8	8.8	1.1	30.1

Table 1.  Proportions (%) of the main expansion components of proton orbits $3/2_8$ and $1/2_{13}$, These results are calculated with DD-ME2 (upper panel) and PKO2 (lower panel) by selecting the negative energy cutoff $E_{-}^C$, respectively as 0 MeV and -50 MeV, with the positive energy cutoff $E_{+}^C = 350$ MeV.

In order to promise the completeness of the expansion (9), it is necessary to consider the spherical DWS states with negative energies, although the relevant contributions are rather small. As emphasized above, the mixing of the high-j waves like $ h_{11/2} $ is essential for giving rise to the octupole deformation. In contrast to the low-j waves, fewer high-j waves with negative energies are bound in the DWS basis, which corresponds to $ E_-^C = 0 $ MeV. Thus, when the high-j waves become essential, it is not enough by considering only bound negative-energy states ($ E_-^C = 0 $) to guarantee the completeness of the expansion (9). Therefore, the notable alterations of the deformation ($ \beta_2, \beta_3 $) from Figs. 1c to 1d are nothing but an indication of the improved completeness. In particular, the Fock terms get more two-body correlations involved, in contrast to the Hartree terms [57]. As a result, the RHF calculations with PKO2 and PKO3 appear to be more sensitive to the negative energy cutoff $ E_-^C $ than the RMF one, as illustrated in Figs. 1 and 2.

With the given space truncations, namely $ E_+^C = 350 $ MeV, $ E_-^C = -150 $ MeV, $ m_{\max} = 15/2 $ and $ K_{m} = 4 $, we perform the O-RHF calculations for the octupole nucleus $ ^{144} {\rm{Ba}}$, using the RHF Lagrangians PKOi ($ i = 1,2,3 $) [36, 39] and the RMF one DD-ME2 [64]. Table 2 lists the total energy E (MeV), the deformations $ (\beta_2, \beta_3) $ and charge radius $ r_c $ (fm) for the ground state (g.s.), the quadruple local minimum (q.m.) and the case with the spherical shape (sph.), respectively.

	Cases	E	$(\beta_2, \beta_3)$	$r_c$	$\Delta E$
Exp.		$-1190.22$	(0.18, 0.17)	4.9236
DD-ME2	g.s.	−1187.95	(0.23, 0.09)	5.0340	−1.51
	q.m.	−1186.44	(0.20, 0.00)	5.0013	−5.77
	sph.	−1180.87	(0.00, 0.00)	4.9013
PKO2	g.s.	−1189.47	(0.22, 0.14)	5.0234	−3.69
	q.m.	−1185.78	(0.19, 0.00)	4.9971	−3.28
	sph.	−1182.50	(0.00, 0.00)	4.9828
PKO1	g.s.	−1192.54	(0.24, 0.14)	5.0525	−3.15
	q.m.	−1189.39	(0.21, 0.00)	5.0224	−4.27
	sph.	−1185.12	(0.00, 0.00)	5.0047
PKO3	g.s.	−1191.77	(0.23, 0.13)	5.0365	−2.69
	q.m.	−1189.08	(0.21, 0.00)	5.0066	−4.61
	sph.	−1184.47	(0.00, 0.00)	4.9878

Table 2.  Total energy E (MeV), quadruple and octupole deformations $(\beta_2, \beta_3)$, and charge radius $r_c$ (fm) of $^{144}$Ba, calculated by PKOi ($i = 1,2,3$) and DD-ME2 for the cases of the octupole g.s., the quadruple local minima (q.m.) and the spherical shape (sph.). The experimental data from Refs. [4, 65, 66] are shown as the references. The last column $\Delta E$ (MeV) is the energy differences between the cases of the g.s. and the q.m., and the ones between the q.m. and sph. cases, respectively.

As illustrated in Table 2, $ ^{144} {\rm{Ba}}$ is predicted to be deeper bound from the sph. to the q.m. cases, and the total energy E becomes more closely aligned with the experimental value [65]. Further considering the octupole deformation, all the selected Lagrangians yield the octupole g.s. for $ ^{144} {\rm{Ba}}$, and PKO2 shows the best agreement with the experimental data, including the total energy E, the deformations $ (\beta_2, \beta_3) $ and charge radius $ r_c $. It is noteworthy that a value of $ \beta_3 = 0.17(_{-6}^{+4}) $ was derived experimentally with $ \beta_2 = 0.18 $ [4]. As illustrated in Table 2, the $ \beta_3 $ value given by the RMF model DD-ME2 is even less than the lower limit $ 0.11 $, similar as previous mean-field calculations [8, 67−69]. In contrast, the values of $ \beta_3 $ given by the RHF models PKOi $ (i = 1,2,3) $ are all within the uncertainty of the experimental value, which deserves to be explored with further detail.

As shown in Table 2, it is commonly seen that the quadruple deformations $ \beta_2 $ for the g.s. slightly increase from the q.m. cases. Thus, the energy differences between the cases of the g.s. and q.m., and the cases of the q.m. and sph., namely the $ \Delta E $ values (MeV) in the last column of Table 2 can be approximately regarded as the effects of the octupole and quadruple deformations, respectively. It is observed that the RHF Lagrangians PKOi ($ i = 1,2,3 $) produce more gains in binding caused by the octupole deformation than the RMF one DD-ME2. Roughly speaking, this can be attributed to the effects of the Fock terms, which get more two-body correlations into account than the Hartree terms. Furthermore, as referred to the value of $ \sim 0.5 $ MeV given by the point coupling RMF model PC-PK1 [13], as well as the one of $ \sim1.8 $ MeV given by the projected HFB calculations [12], the selected RHF models also predict more pronounced effects of octupole deformation.

On the other hand, it is somewhat beyond our expectation that PKO1 and PKO3, which contain the degree of freedom associated with the π-PV coupling, yield weaker octupole enhancements than PKO2. It is worthwhile to mention that PKO2 does not contain the π-PV coupling and PKO3 carries stronger π-PV coupling than PKO1. It is observed that the $ \Delta E $ values between the g.s. and the q.m. cases decrease in sequence from PKO2 to PKO1 and further to PKO3, and vice versa for the values of $ \Delta E $ between the q.m. and the sph. cases. This suggests that the π-coupling does not favor the octupole deformation, but rather enhances the quadruple deformation. Indeed, PKO1 and PKO3 have been shown to produce an enhanced quadruple deformation effect in comparison to PKO2. For instance, in the well-deformed nucleus $ ^{20} {\rm{Ne}}$, the tensor force component carried by the π-PV coupling has been demonstrated to enhance the deformation effects [56]. Moreover, the ρ-tensor coupling, due to the nature of Lorentz tensor coupling [57], can also enhance the deformation effects, as evidenced in $ ^{24} {\rm{Mg}}$ [57] and $ ^{11} {\rm{Be}}$ [58, 70]. This may be applicable for nuclei with octupole deformation, which calls for the future implementation of the ρ-tensor coupling in the O-RHF model. A qualitative understanding of the systematics of octupole and quadruple enhancements can be obtained by combining these with the single-particle spectra of $ ^{144} {\rm{Ba}}$.

Figures 3a and 3b show the neutron and proton spectra of $ ^{144} {\rm{Ba}}$ given by PKO2, respectively. It shall be commented that the single-particle spectra given by PKO2 are not significantly different from the ones given by PKO1 and PKO3, which are not shown here. As illustrated in Fig. 3, from q.m. to the octupole g.s., both neutron and proton valence orbits (marked in blue color) become much deeper bound, and the induced large shell gaps (denoted by double arrows) stabilize the octupole deformation for the g.s. of $ ^{144} {\rm{Ba}}$. In order to further understand the effects of octupole deformation, Table 3 shows the proportions (%) of the main DWS waves for the neutron and proton valence orbits. It is seen that the intrusions of both high-lying neutron $ 1i_{13/2} $ and proton $ 1h_{11/2} $ waves play a key role in giving rise to the octupole deformation, which is commonly supported by the selected models in this work.

Figure 3.  (color online) Neutron (left) and proton (right) spectra of $^{144}$Ba for the sph., q.m. and octupole g.s. given by PKO2. The ultrathick bars represent the occupation probabilities of the orbits $m_\nu$ and $E_F$ denotes the Fermi levels.

Neutron		$1i_{13/2}$	$3p_{3/2}$	$2f_{5/2}$	$2f_{7/2}$	$1h_{9/2}$	$1h_{11/2}$
$3/2_{12}$	q.m.	0.0	2.3	16.5	11.7	62.4	4.4
	g.s.	24.5	2.9	12.1	22.5	11.5	14.1
$1/2_{18}$	q.m.	0.0	23.9	0.8	42.7	18.8	8.3
	g.s.	23.9	2.4	15.7	12.4	22.0	3.1
$1/2_{17}$	q.m.	0.0	6.5	27.4	4.0	48.2	3.1
	g.s.	3.5	15.6	15.7	5.6	11.3	14.6
Proton		$1h_{11/2}$	$3s_{1/2}$	$2d_{3/2}$	$2d_{5/2}$	$1g_{7/2}$	$1g_{9/2}$
$3/2_{8}$	q.m.	0.0	−	4.1	3.0	86.9	3.8
	g.s.	20.2	−	1.9	6.6	53.5	5.9
$1/2_{13}$	q.m.	0.0	14.1	4.4	55.7	14.2	9.1
	g.s.	29.8	7.1	5.7	25.3	2.3	8.6
$1/2_{12}$	q.m.	0.0	2.9	18.8	2.2	70.2	2.6
	g.s.	3.4	3.3	21.4	0.0	49.6	3.3

Table 3.  Proportions (%) of the main DWS waves in the neutron orbits $3/2_{12}$, $1/2_{18}$ and $1/2_{17}$ (upper panel) and the proton orbits $3/2_8$, $1/2_{13}$ and $1/2_{12}$(lower panel) of $^{144}$Ba. These results are calculated with the RHF Lagrangian PKO2.

Taking proton as an example, it is clear that the couplings between the $ 2d_{5/2} $ and the intruding $ 1h_{11/2} $ waves, which fulfill the conditions $ \Delta l = 3 $ and $ \Delta j = 3 $, play an essential role in giving a stable octupole deformation. At the q.m., the orbits $ (3/2_8, 1/2_{12}) $ and $ 1/2_{13} $ are dominated by the $ 1g_{7/2} $ and $ 2d_{5/2} $ waves (in bold types), respectively, as illustrated in the lower panel of Table 3. Consistently, when the octupole deformation is involved, more evident couplings between the $ 1h_{11/2} $ and $ 2d_{5/2} $ waves are observed for the orbit $ 1/2_{13} $ than the ones $ (3/2_8, 1/2_{12}) $. This explains partly the octupole enhancements, especially for the proton orbit $ 1/2_{13} $, which becomes much more deeper bound from the q.m. to the octupole g.s. as illustrated in Fig. 3b. For neutron, the couplings between the $ 2f_{7/2} $ and the intruding $ 1i_{13/2} $ waves are significant for giving rise to the octupole deformation, as illustrated in the upper panel of Table 3. In particular, the proportions of both neutron and proton pseudo-spin partners, respectively ($ 2f_{7/2}, 1h_{9/2} $) and ($ 2d_{5/2}, 1g_{7/2} $), undergo significant redistributions with the penetrations of the high-j waves $ 1i_{13/2} $ and $ 1h_{11/2} $. This can favor the emergence of octupole deformation in $ ^{144} {\rm{Ba}}$.

As demonstrated by the $ \Delta E $ values given by PKOi $ (i = 1,2,3) $ in Table 2, the π-PV coupling, which contributes only via the Fock terms, seems to present an effect against the octupole deformation, but enhance the effect of the quadruple deformation for $ ^{144} {\rm{Ba}}$. Combined Fig. 3 and Table 3, it is not difficult to understand the role of the π-PV coupling, in which the tensor force components lead to the repulsive couplings between the $ j_> = l+1/2 $ ($ j_< = l-1/2 $) and $ j_> $ ($ j_< $) waves, but attractive ones between the $ j_< $ and $ j_> $ waves [71]. At the q.m. as illustrated in Fig. 3 and Table 3, the valence neutron ($ 3/2_{12}, 1/2_{17} $) and proton $ (3/2_8, 1/2_{12}) $ orbits are dominated by the $ j_< $ waves, i.e., the $ 1h_{9/2} $ and $ 1g_{7/2} $ waves, respectively. In contrast, as the results of strong spin-orbit couplings, the core of $ ^{144} {\rm{Ba}}$ corresponds to $ N = 82 $ and $ Z = 50 $, both representing the nature of $ j_> $ on average. Thus, due to the attractive tensor force effects carried by the π-PV coupling, PKO1 and PKO3 present a stronger quadruple enhancement than PKO2, similar as revealed in Ref. [56]. Notice that the RMF Lagrangian DD-ME2, which do not contain either the Fock terms or the π-PV coupling, present a stronger quadruple enhancement than the RHF models PKOi ($ i = 1,2,3 $). This can be attributed to the fact that the modeling of nuclear binding, mainly the interplay between strong attractive σ-S and repulsive ω-V couplings, has been changed with the presence of the Fock terms, particularly by the much enhanced isovector ρ-V coupling, as compared to the RMF models like DD-ME2 [36, 47].

As mentioned above, the intrusion of the $ j_> $ waves, namely neutron $ 1i_{13/2} $ and proton $ 1h_{11/2} $ waves, is essential for the emergence of the octupole deformation in $ ^{144} {\rm{Ba}}$. However, due to the nature of tensor force, this results in repulsive tensor couplings between the core of $ ^{144} {\rm{Ba}}$ and the intruding waves, because both represent the nature of $ j_> $. Meanwhile, as illustrated in Table 3, the $ j_< $ waves in both neutron and proton valence orbits are averagely reduced from the q.m. to the octupole g.s., which further weakens the attractive tensor couplings between the core and the $ j_< $ waves of valence nucleons. Thus, it is not difficult to understand that the π-PV coupling, mainly due to its tensor force components, manifests an effect against the octuple deformation for $ ^{144} {\rm{Ba}}$, as revealed by the $ \Delta E $ values in Table 2. This may be applicable for the regions with $ N/Z = 34 $, $ 56 $, $ 88 $ and $ N = 136 $, in which the couplings $ (1g_{9/2}, 2p_{3/2}) $, $ (1h_{11/2}, 2d_{5/2}) $, $ (1i_{13/2}, 2f_{7/2}) $ and $ (1j_{15/2}, 2g_{9/2}) $ play an essential role in giving rise to the octupole deformation, respectively. It is due to the fact that both intruding waves $ (1g_{9/2}, 1h_{11/2}, 1i_{13/2}, 1j_{15/2}) $ and the cores ($ N = 34, 56, 88,136 $ and/or $ Z = 34, 56, 88 $) carry the nature of $ j_> $, between which the tensor force couplings are repulsive. This also deserves a special focus in future applications of the O-RHF model.

In this work, the relativistic Hartree-Fock (RHF) theory is extended to accommodate the octupole-deformed nuclei, giving the O-RHF models. In terms of the spherical Dirac Woods-Saxon (DWS) basis, the general formalism of the O-RHF model are introduced, and the space truncations of the DWS basis are tested as well. Taking $ ^{144} {\rm{Ba}}$ as an example, the mechanism related to the occurrence of the octupole deformation is analyzed with a special focus on the Fock terms and the tensor force effects carried by the π-PV coupling.

It is found that the selected RHF models can reproduce the octupole deformation of $ ^{144} {\rm{Ba}}$ within the uncertainty of the experimental measurements, in contrast to the other popular mean field models. Moreover, the intrusion of the high-lying neutron $ 1i_{13/2} $ and proton $ 1h_{11/2} $ waves, which gives rise to the octupole deformation in $ ^{144} {\rm{Ba}}$, is notably enhanced by the Fock terms. This leads to an enhanced effect of octupole deformation. In particular, due to the repulsive tensor couplings between the core of $ ^{144} {\rm{Ba}}$ and the intruding $ j_> $ waves ($ 1i_{13/2} $ and $ 1h_{11/2} $), the tensor force components carried by the π-PV coupling present an effect against the octupole deformation for $ ^{144} {\rm{Ba}}$, which deserves special attentions in future study of octupole nuclei.