The tensor part of the Skyrme energy density functional. I. Spherical nuclei

TL;DR

This study systematically investigates how the J^2 tensor term influences spherical nuclei properties within the Skyrme energy functional, revealing limitations in current functional flexibility to incorporate large tensor effects.

Contribution

It introduces a comprehensive set of 36 Skyrme parameterizations exploring tensor term variations and assesses their impact on nuclear observables.

Findings

Tensor terms significantly affect spin-orbit splittings and spectra.

Current Skyrme functionals lack flexibility for large tensor contributions.

Tensor effects do not improve agreement with experimental data.

Abstract

We perform a systematic study of the impact of the J^2 tensor term in the Skyrme energy functional on properties of spherical nuclei. In the Skyrme energy functional, the tensor terms originate both from zero-range central and tensor forces. We build a set of 36 parameterizations, which covers a wide range of the parameter space of the isoscalar and isovector tensor term coupling constants, with a fit protocol very similar to that of the successful SLy parameterizations. We analyze the impact of the tensor terms on a large variety of observables in spherical mean-field calculations, such as the spin-orbit splittings and single-particle spectra of doubly-magic nuclei, the evolution of spin-orbit splittings along chains of semi-magic nuclei, mass residuals of spherical nuclei, and known anomalies of charge radii. Our main conclusion is that the currently used central and spin-orbit parts of the Skyrme energy density functional are not flexible enough to allow for the presence of large tensor terms.

Figures

Equations

$$\mathcal{E}={\mathcal{E}}_{\text{kin}}+{\mathcal{E}}_{\text{Skyrme}}+{\mathcal{E}}_{\text{pairing}}+{\mathcal{E}}_{\text{Coulomb}}+{\mathcal{E}}_{\text{corr}}.$$

$$v^{\text{Skyrme}}=v^{\text{c}}+v^{\text{t}}+v^{\text{LS}}.$$

$${\mathcal{E}}_{\text{Skyrme}}=\int d^{3}r{\mathcal{H}}^{\text{Skyrme}}(\mathbf{r}),$$

$${\rho }_q(\mathbf{r}\sigma ,{\mathbf{r}}^{\prime }{\sigma }^{\prime })$$

$${\rho }_q(\mathbf{r},{\mathbf{r}}^{\prime })$$

$${\mathbf{s}}_q(\mathbf{r},{\mathbf{r}}^{\prime })$$

$${\rho }_q(\mathbf{r})$$

$${\mathbf{s}}_q(\mathbf{r})$$

$${\tau }_q(\mathbf{r})$$

$$T_{q,\mu }(\mathbf{r})$$

$${\mathbf{j}}_q(\mathbf{r})$$

$$J_{q,\mu \nu }(\mathbf{r})$$

$$F_{q,\mu }(\mathbf{r})$$

$$J_{\mu \nu }(\mathbf{r})=\frac{1}{3}{\delta }_{\mu \nu }{J}^{(0)}(\mathbf{r})+\frac{1}{2}\sum _{\kappa =x}^z{ϵ}_{\mu \nu \kappa }{J}_{\kappa }^{(1)}(\mathbf{r})+J_{\mu \nu }^{(2)}(\mathbf{r}),$$

$$J^{(0)}(\mathbf{r})$$

$$J_{\kappa }^{(1)}(\mathbf{r})$$

$$J_{\mu \nu }^{(2)}(\mathbf{r})$$

$${\rho }_0(\mathbf{r})$$

$${\rho }_1(\mathbf{r})$$

$$v^{\text{c}}(\mathbf{R},\mathbf{r})$$

$$\mathbf{r}$$

$$\mathbf{R}$$

$$\hat{\mathbf{k}}=-\frac{i}{2}({\nabla }_1-{\nabla }_2)$$

$${\hat{P}}_{\sigma }=\frac{1}{2}(1+{\hat{\sigma }}_1\cdot {\hat{\sigma }}_2).$$

$${\mathcal{H}}^{\text{c}}(\mathbf{r})$$

$${\mathcal{H}}_t^{\text{c,even}}$$

$${\mathcal{H}}_t^{\text{c,odd}}$$

$$v^{\text{LS}}(\mathbf{r})=i{W}_0({\hat{\sigma }}_1+{\hat{\sigma }}_2)\cdot {\hat{\mathbf{k}}}^{\prime }\times \delta (\mathbf{r})\hat{\mathbf{k}}$$

$${\mathcal{H}}^{\text{LS}}(\mathbf{r})$$

$${\mathcal{H}}_t^{\text{LS,even}}$$

$${\mathcal{H}}_t^{\text{LS,odd}}$$

$${\hat{S}}_{12}=3({\hat{\sigma }}_1\cdot {\mathbf{e}}_r)({\hat{\sigma }}_2\cdot {\mathbf{e}}_r)-{\hat{\sigma }}_1\cdot {\hat{\sigma }}_2,$$

$$v^{\text{t}}(\mathbf{r})$$

$${\mathcal{H}}^{\text{t}}(\mathbf{r})$$

$${\mathcal{H}}_t^{\text{t,even}}$$

$${\mathcal{H}}_t^{\text{t,odd}}$$

$${\mathcal{E}}_{\text{Skyrme}}$$

Full text

I Introduction

The strong nuclear spin-orbit interaction in nuclei is responsible for the observed magic numbers in heavy nuclei Goeppert Mayer (1948); Haxel et al. (1949); Feenberg and Hammack (1949); Goeppert Mayer (1949). While a simple spin-orbit interaction allows for the qualitative description of the global features of shell structure, the available data suggest that single-particle energies evolve with neutron and proton number in a manner that cannot be related to the geometrical growth of the single-particle potential with $N$ and $Z$. Many anomalies of shell structure have been identified that do not fit into simple experimental systematics, and that challenge any global model of nuclear structure.

The evolution of shell structure with $N$ and $Z$ as a feature of self-consistent mean-field models has been known for long. To quote the pioneering study of shell structure in a self-consistent model performed by Beiner et al. Beiner et al. (1975), the “most striking effect is the appearance of $N=16$, 34 and 56 as neutron magic numbers for unstable nuclei, together with a weakening of the shell closure at $N=20$ and 28”. Various mechanisms that modify the appearance of gaps in the single-particle spectra have been discussed in detail in the literature. The two most prominent ones that were worked out by Dobaczewski et al. in Ref. Dobaczewski et al. (1994), however, play mainly a role for weakly-bound exotic nuclei far from stability, as they are directly or indirectly related to the physics of loosely bound single-particle states, namely that the enhancement of the diffuseness of neutron density distribution reduces the spin-orbit coupling in neutron-rich nuclei on the one hand, and the interaction between bound orbitals and the continuum results in a quenching of shell effects in light and medium systems on the other hand. The former effect was also extensively discussed in the framework of relativistic models by Lalazissis et al. Lalazissis et al. (1998a, b), while the latter triggered a number of studies that discussed the potential relevance of this so-called “Boguliubov enhanced shell quenching” to explain the abundance pattern from the astrophysical $r$-process of nucleosynthesis Chen et al. (1995); Dobaczewski et al. (1995); Pearson et al. (1996); Pfeiffer et al. (1997).

These two effects take place in neutron-rich nuclei. In proton-rich nuclei, the Coulomb barrier suppresses both the diffuseness of the proton density and the coupling of bound proton states to the continuum. But the Coulomb interaction itself can also modify the shell structure: for super-heavy nuclei, it begins to destabilize the nucleus as a whole. Mean-field models predict that it amplifies the shell oscillations of the densities for incomplete filled oscillator shells, which leads to strong variations of the density profile that feed back onto the single-particle spectra Dechargé et al. (1999); Bender et al. (1999).

Interestingly, most theoretical papers about the evolution of shell structure from the last decade have speculated about new effects that mainly affect neutron shells in nuclei far from stability in the anticipation of the rare-isotope physics that might become accessible with the next generation of experimental facilities. The known anomalies, some of which have been known for a long time, and many more have been identified recently, concern also proton shells and already appear sufficiently close to stability that “exotic phenomena can be ruled out for their explanation” in most cases, to paraphrase the authors of Ref. Langanke et al. (2003). By contrast, this suggests that there exists a mechanism that induces a strong evolution of single-particle spectra already in stable nuclei that has been overlooked for long.

There is a prominent ingredient of the nucleon-nucleon interaction that has been ignored for decades in virtually all global nuclear structure models for medium and heavy nuclei, be it macroscopic-microscopic approaches or self-consistent mean-field methods. It is only very recently, that the systematic discrepancies between model predictions and experiment have triggered a renaissance of the tensor force in the description of finite medium- and heavy-mass nuclei.

The tensor force is a crucial and necessary ingredient of the bare nucleon-nucleon interaction Wiringa et al. (1995); Machleidt (2001), and consequently is contained in all ab-initio approaches that are available for light, mainly $p$-shell nuclei Pieper and Wiringa (2001); Navrátil and Ormand (2003). One of the first experimental signatures of the tensor force was the small, but finite quadrupole moment of the deuteron. In a boson-exchange picture of the bare nucleon-nucleon interaction, the tensor force originates from the exchange of pseudoscalar pions, which have both central and tensor couplings, see for example section 2.3 in Ref. Eisenberg and Greiner (1972) or appendix 13A of Ref. Nilsson and Ragnarsson (1995). In a nuclear many-body system, the bare tensor force induces a strong correlation between the spatial and spin orientations in the two-body density matrix. For two nucleons with parallel spins, the tensor force energetically favors the configuration where the distance vector is aligned with the spins, while for anti-parallel spins the tensor force prefers when the distance vector is perpendicular to the spins, see the discussion of Fig. 13 in Ref. Neff and Feldmeier (2003) and of Fig. 3 in Ref. Roth et al. (2004). The authors of these papers also demonstrate very nicely the well-known fact Bethe (1968); Negele (1970) that in an approach that starts from the bare nucleon-nucleon interaction, nuclei are not bound without taking into account the two-body correlations induced by the tensor force.

The role of the tensor force, however, manifests itself differently in self-consistent mean-field models, otherwise called energy density functional (EDF) methods, the tool of choice for medium and heavy nuclei. The latter methods use an independent-particle state as a reference state to express the energy of the correlated nuclear ground state. Thus, correlations are not explicitly present in the higher-order density matrices of the reference state, but rather included under the form of a more elaborate functional of the (local and nearly local parts of the) one-body density matrix of that reference state. In such a scheme, most of the effect of the bare tensor force on the binding energy is integrated out through the renormalization of the coupling constants associated with a central effective vertex, in a similar fashion as the tensor part of the bare interaction is renormalized into the central one when going from the bare nucleon-nucleon force to a Brueckner $G$ matrix. The tensor terms of the EDF relate to a residual tensor vertex, that gives nothing but a correction to the spin-orbit splittings, which for light $p$-shell nuclei might be of the same order as the contribution from the genuine spin-orbit force. The interplay of spin-orbit and tensor forces in the mean field of medium and heavy nuclei was explored in Refs. Scheerbaum (1976); Goodman and Borysowicz (1978); Zheng and Zamick (1991), where the particular role of spin-unsaturated shells was pointed out.

There are two widely used effective interactions for non-relativistic self-consistent mean-field models Bender et al. (2003a), the zero-range non-local Skyrme interaction Skyrme (1956, 1958a); Bell and Skyrme (1956); Skyrme (1958b) on the one hand and the finite-range Gogny force Gogny (1975); Dechargé and Gogny (1980) on the other hand.

In fact, the effective zero-range non-local interaction proposed by Skyrme in 1956 Skyrme (1956, 1958a); Bell and Skyrme (1956); Skyrme (1958b) already contained a zero-range tensor force. The first applications of Skyrme’s interaction in self-consistent mean-field models that became available around 1970, however, neglected the tensor force, and the simplified effective Skyrme interaction used in the seminal paper by Vautherin and Brink Vautherin and Brink (1972) soon became the standard Skyrme interaction that was used in most applications ever since. Until very recently, there was only very little exploratory work on Skyrme’s tensor force. In their early study, Stancu, Brink and Flocard Stancu et al. (1977), who added the tensor force perturbatively to the SIII parameterization, pointed out that some spin-orbit splittings in magic nuclei can be improved with a tensor force. A complete fit including the terms from the tensor force that contribute in spherical nuclei was attempted by Tondeur Tondeur (1983), with the relevant coupling constants of the spin-orbit and tensor terms adjusted to selected spin-orbit splittings in 16O, 48Ca and 208Pb. Another complete fit of a generalized Skyrme interaction including a tensor force was performed by Liu et al. Liu et al. (1991), but the authors did not investigate the effect of the tensor force in detail, nor was the resulting parameterization ever used in the literature thereafter.

Similarly, the seminal paper by Gogny Gogny (1975) on the evaluation of matrix elements of a finite-range force of Gaussian shape in an harmonic oscillator basis contains the expressions for a finite-range tensor force, which, however, was omitted in the parameterizations of Gogny’s force adjusted by the Bruyères-le-Châtel group Dechargé and Gogny (1980). It were Onishi and Negele Onishi and Negele (1978) who first published an effective interaction that combined a Gaussian two-body central force, a finite-range tensor force with a zero-range spin-orbit force and a zero-range non-local three-body force, which, however, also fell into oblivion.

The role of the tensor force is slightly different in Skyrme and Gogny interactions. In the Gogny force, the contributions from the central and tensor parts remain explicitly distinct, although, of course, this does not prevent a certain entanglement of their physical effects. In the context of Skyrme’s functional, however, the contribution of a zero-range tensor force to the spherical mean-field state of an even-even nucleus has exactly the same form as a particular exchange term from the non-local part of the central Skyrme force. When looking at spherical nuclei only, adding Skyrme’s tensor force simply allows one to decouple a term that is already provided by the central force. This indeed makes the effective-interaction-restricted functional more flexible, as the additional degrees of freedom from the tensor force remove an interdependence between the effective mass, the surface terms and the “tensor terms”. However, one must always keep in mind that both the central and tensor part of the effective vertex contribute to the so-called ${\mathbf{J}}^{2}_{t}$ “tensor” terms of the functional.111As we will outline below, and as was already pointed out in Ref. Beiner et al. (1975), this argument does not hold for deformed even-even nuclei or any situation where intrinsic time-reversal is broken, for example odd nuclei or dynamics. There, the tensor and non-local central parts of the effective Skyrme interaction give contributions to the mean-fields and the binding energy with different analytical expressions. This will be discussed in a companion article Bender et al. (2007).

In the context of relativistic mean-field models, the equivalent of the non-relativistic tensor force appears as the exchange term of effective fields with the quantum numbers of the pion, which by construction do not appear in the standard relativistic Hartree models. Only relativistic Hartree-Fock models contain this tensor force, with the first predictive parameterizations becoming available just recently Long et al. (2006a).

We also mention that there is a large body of work on the tensor force in the interacting shell model, see Ref. Fayache et al. (1997) for a review, that concentrates on a completely different aspect of the tensor force, namely its unique contribution to excitations with unnatural parity.

The recent interest in the effect of the tensor force in the context of self-consistent mean field models was triggered by the observed evolution of single-particle levels of one nucleon species in dependence of the number of the other nucleon species. Otsuka et al. Otsuka et al. (2005) proposed that at least part of the effect is caused by the proton-neutron tensor force from pion exchange. Many groups attempt now to explain known, but so far unresolved, anomalies of shell structure in terms of a tensor force. A particularly popular playground is the relative shift of the proton $1g_{7/2}$ and $1h_{11/2}$ levels in tin isotopes, which is interpreted as the reduction of the spin-orbit splittings of both levels with their respective partners with increasing neutron number Schiffer et al. (2004).

Otsuka et al. Otsuka et al. (2006) added a Gaussian tensor force, adjusted on the long-range part of a one-pion+$\rho$ exchange potential, to a standard Gogny force. After a consistent readjustment of the parameters of its central and spin-orbit parts, they were able to explain coherently the anomalous relative evolution of some single-particle levels without, however, being able to describe their absolute distance in energy. Dobaczewski Dobaczewski (2006) has pointed out that a perturbatively added tensor interaction with suitably chosen coupling constants in the Skyrme energy density functional does not only modify the evolution of shell structure, but does also improve the description of nuclear masses around magic nuclei. Brown et al. Brown et al. (2006) have fitted a Skyrme interaction with added zero-range tensor force with emphasis on the reproduction of single-particle spectra. While the authors appreciate the qualitatively correctly described evolution of relative level distances, they point out that the combination of zero-range spin-orbit and tensor forces does not and can not correctly describe the $\ell$-dependence of spin-orbit splittings. Colò et al. Colò et al. (2007), and Brink et al. Brink and Stancu (2007) have added Skyrme’s tensor force perturbatively to the existing standard parameterization SLy5 Chabanat et al. (1997, 1998), and to the SIII Beiner et al. (1975) one, respectively. They have investigated some single-particle energy differences: the $1h_{11/2}$ and $1g_{7/2}$ proton states in tin isotopes as well as $1i_{13/2}$ and $1h_{9/2}$ neutron states in $N=82$ isotones and propose similar parameters as in Ref. Brown et al. (2006). The effect of the tensor force on the centroid of the GT giant resonance is also estimated by Colò et al. using a sum-rule approach and found to be substantial. Long et al. Long et al. (2006b), demonstrate that the tensor force that emerges naturally in relativistic Hartree-Fock also improves the relative shifts of the proton $1g_{7/2}$ and $1h_{11/2}$ levels in tin isotopes.

The work on the tensor force published so far aims at an optimal single parameterization, that establishes a best fit to either the underlying bare tensor force Otsuka et al. (2006); Brown et al. (2006) or empirical data Tondeur (1983); Dobaczewski (2006); Colò et al. (2007). The published results, as well as our first exploratory studies, however, suggest that adding a tensor force to the existing mean-field models gives only a local improvement of the relative change of certain single-particle energies, but not necessarily a global improvement of single-particle spectra or other observables. In the framework of the Skyrme interaction, that we will employ throughout this work, there is also the already mentioned ambiguity that the contribution from the tensor force to spherical nuclei has the same structure as a term from the central force. In view of this situation, we will pursue a different strategy and investigate the effect of the tensor terms on a multitude of observables in nuclei though a set of Skyrme interactions with systematically varied coupling constants of the tensor terms.

The present study was motivated by the finding that the performance of the existing Skyrme-type effective interactions for masses and spectroscopic properties is limited by systematic deficiencies of the single-particle spectra Bender et al. (2006a, 2003b, b); Chatillon et al. (2006) that seem to be impossible to remove within the standard Skyrme interaction. The details of single-particle spectra were so far somewhat outside the focus of self-consistent mean-field methods, on the one hand as they do not correspond directly to empirical single-particle energies (we will come back to that below), and on the other hand because many of the observables that are usually calculated with self-consistent mean-field methods are not very sensitive to the exact placement of single-particle levels. By contrast, there is an enormous body of work that examines the infinite and semi-infinite nuclear matter properties of the effective interactions that are the analog of liquid-drop and droplet parameters in great detail. The reason is, of course, that the global trends over the whole chart of nuclei have to be understood before one can look into details. The last few years have seen an increasing demand on predictive power. Moreover, beyond-mean-field approaches of the projected generator coordinate method (GCM), or Bohr-Hamiltonian type, have become widely used tools to analyze and predict spectroscopic properties in medium and heavy nuclei, employing either Gogny or Skyrme interactions. The underlying single-particle spectra thus now deserve more attention, as many of the spectroscopic properties of interest turn out to be extremely sensitive to even subtle details of the single-particle spectra. As the tensor force is the most obvious missing piece in all standard mean-field interactions, it is the natural starting point for the systematic investigation of possible generalizations with the ultimate goal to improve the predictive power of the interactions for spectroscopy.

In the present paper, we will outline the formalism of a Skyrme interaction with added tensor force, describe the fit of the parameterizations, analyze the role of the tensor terms for single-particle spectra, masses and radii of spherical even-even nuclei. A second paper Bender et al. (2007) studies the surface and deformation properties of these Skyrme interactions for even-even nuclei, and future work will examine the stability of nuclear matter and the role of the time-odd terms from the tensor force in odd and rotating nuclei. Only deformed nuclei and, in particular, observables sensitive to the time-odd contributions, will possibly allow to distinguish clearly between the non-local central and tensor parts of the Skyrme force.

II The Skyrme interaction with tensor terms

II.1 The energy density functional

The usual ansatz for the Skyrme effective interaction Chabanat et al. (1997, 1998) leads to an energy density functional which can be written as the sum of a kinetic term, the Skyrme potential energy functional that models the effective strong interaction in the particle-hole channel, a pairing energy functional corresponding to a density-dependent contact pairing interaction, the Coulomb energy functional (calculated using the Slater approximation Slater (1951)) and correction terms to approximately remove the excitation energy from spurious motion caused by broken symmetries

[TABLE]

II.2 The Skyrme energy density functional

Throughout this work, we will use an effective Skyrme energy functional that corresponds to an antisymmetrized density-dependent two-body vertex in the particle-hole channel of the strong interaction, that can be decomposed into a central, spin-orbit and tensor contribution

[TABLE]

Other choices for the writing of the Skyrme energy functional are possible and have been made in the literature, which might affect the form of the effective interaction, its interpretation and the results obtained from it. We will come back to that in section II.4 below.

The Skyrme energy density functional is a functional of local densities and currents

[TABLE]

which has many technical advantages compared to finite-range forces such as the Gogny force. All exchange terms have the same structure as the direct terms, which greatly reduces the number of necessary integrations during a calculation.

II.2.1 Local densities and currents

Throughout this paper we will assume that we have pure proton and neutron states. The formal framework of the general case including proton-neutron mixing is discussed in Ref. Perlińska et al. (2004). Without making reference to any single-particle basis, we start from the density matrices of protons and neutrons in coordinate space Dobaczewski et al. (2000)

[TABLE]

where

[TABLE]

The Skyrme energy functional up to second order in derivatives that we will introduce below can be expressed in terms of seven local densities and currents Perlińska et al. (2004) that are defined as

[TABLE]

which are the density $\rho_{q}({\mathbf{r}})$, the kinetic density $\tau_{q}({\mathbf{r}})$, the current (vector) density ${\mathbf{j}}_{q}({\mathbf{r}})$, the spin (pseudovector) density ${\mathbf{s}}_{q}({\mathbf{r}})$, the spin kinetic (pseudovector) density ${\mathbf{T}}_{q}({\mathbf{r}})$, the spin-current (pseudotensor) density $J_{q,\mu\nu}({\mathbf{r}})$, and the tensor-kinetic (pseudovector) density ${\mathbf{F}}_{q}({\mathbf{r}})$. $\rho_{q}({\mathbf{r}})$, $\tau_{q}({\mathbf{r}})$ and $J_{q,\mu\nu}({\mathbf{r}})$ are time-even, while ${\mathbf{s}}_{q}({\mathbf{r}})$, ${\mathbf{T}}_{q}({\mathbf{r}})$, ${\mathbf{j}}_{q}({\mathbf{r}})$ and ${\mathbf{F}}_{q}({\mathbf{r}})$ are time-odd. For a detailed discussion of their symmetries see Ref. Dobaczewski et al. (2000). There are other local densities up to second order in derivatives that can be constructed, but when constructing an energy functional they either cannot be combined with others to terms with proper symmetries or they lead to terms that are not independent from the others Dobaczewski and Dudek (1996).

The cartesian spin-current pseudotensor density $J_{\mu\nu}$ can be decomposed into pseudoscalar, (anti-symmetric) vector and (symmetric) traceless pseudotensor parts, all of which have well-defined transformation properties under rotations

[TABLE]

where $\delta_{\mu\nu}$ is the Kronecker symbol and $\epsilon_{\mu\nu\kappa}$ the Levi-Civita tensor. The pseudoscalar, vector and pseudotensor parts expressed in terms of the cartesian tensor are given by

[TABLE]

The vector spin current density ${\mathbf{J}}^{(1)}({\mathbf{r}})\equiv{\mathbf{J}}({\mathbf{r}})$ is often called spin-orbit current, as it enters the spin-orbit energy density. 222Some authors call ${\mathbf{J}}({\mathbf{r}})$ spin density, which is ambiguous and confusing when discussing the complete energy density functional including terms that contain the time-odd ${\mathbf{s}}({\mathbf{r}})$.

For the formal discussion of the physical content of the Skyrme energy functional it is of advantage to recouple the proton and neutron densities to isoscalar and isovector densities, for example

[TABLE]

and similar for all others. As we assume pure proton and neutron states, only the $T_{z}=0$ component of the isovector density is non-zero, which we exploit to drop the index $T_{z}$ from the isovector densities $\rho_{1T_{z}}({\mathbf{r}})$ etc.

II.2.2 Skyrmes’s central force

We will use the standard density-dependent central Skyrme force

[TABLE]

where we use the shorthand notation

[TABLE]

while $\hat{{\mathbf{k}}}$ is the usual operator for relative momenta

[TABLE]

and $\hat{{\mathbf{k}}}^{\prime}$ its complex conjugated acting on the left. Finally, $\hat{P}_{\sigma}$ is the spin exchange operator that controls the relative strength of the $S=0$ and $S=1$ channels for a given term in the two-body interaction

[TABLE]

As said above, we restrict ourselves to a parameterization of the Skyrme energy functional as obtained from the average value of an effective two-body vertex in the reference Slater determinant. We decompose the isoscalar and isovector parts of the resulting energy density functional $\mathcal{H}^{\text{c}}$ into a part $\mathcal{H}_{t}^{\text{c,even}}$ that is composed entirely of time-even densities and currents, and a part $\mathcal{H}_{t}^{c,\text{odd}}$ that contains terms which are bilinear in time-odd densities and currents and vanishes in intrinsically time-reversal invariant systems

[TABLE]

Both $\mathcal{H}_{t}^{\text{c,even}}$ and $\mathcal{H}_{t}^{\text{c,odd}}$ are of course constructed such that they are time-even; they are given by Engel et al. (1975); Perlińska et al. (2004)

[TABLE]

where $A^{\rho}_{t}[\rho_{0}]$ and $A^{s}_{t}[\rho_{0}]$ are density dependent coupling constants that depend on the total (isoscalar) density. The detailed relations between the coupling constants of the functional and the central Skyrme force are given in appendix A. The notation reflects that two pairs of terms in $\mathcal{H}_{t}^{\text{c,even}}$ and $\mathcal{H}_{t}^{\text{c,odd}}$ are connected by the requirement of local gauge invariance of the Skyrme energy functional Dobaczewski and Dudek (1995).

II.2.3 A zero-range spin-orbit force

The spin-orbit force used with most standard Skyrme interactions

[TABLE]

is a special case of the one proposed by Bell and Skyrme Bell and Skyrme (1956); Skyrme (1958b). Again, the corresponding energy functional Engel et al. (1975); Perlińska et al. (2004) can be separated into a time-even and a time-odd term

[TABLE]

where

[TABLE]

which share the same coupling constant as again both terms are linked by the local gauge invariance of the energy functional. The relation between the $A^{\nabla\cdot J}_{t}$ and the one coupling constant of the two-body spin-orbit force $W_{0}$ is given in appendix A.

II.2.4 Skyrme’s tensor force

By convention, the tensor operator in the tensor force is constructed using the unit vectors in the direction of the relative coordinate ${\mathbf{e}}_{r}={\mathbf{r}}/|{\mathbf{r}}|$ and subtracting $\hat{\boldsymbol{\mathbf{\sigma}}}_{1}\cdot\hat{\boldsymbol{\mathbf{\sigma}}}_{2}$

[TABLE]

such that its mean value vanishes for a relative $S$ state, which decouples the central and tensor channels of the interaction. The operator $\hat{S}_{12}$ commutes with the total spin $[\hat{S}_{12},\hat{{\mathbf{S}}}^{2}]=0$, therefore it does not mix partial waves with different spin, i.e. spin singlet and spin triplet states. In particular, it does not act in spin singlet states at all, as $\hat{S}_{12}\hat{P}_{S=0}=0$ (see section 13.6 of Ref. Nilsson and Ragnarsson (1995)). As a consequence, there is no point in multiplying a tensor force with an exchange operator $(1+x_{t}\hat{P}_{\sigma})$ as done for the central force, as this will only lead to an overall rescaling of its strength.

The derivation of the general energy functional from a zero-range two-body tensor force is discussed in detail in Refs. Flocard (1975); Perlińska et al. (2004). We repeat here the details relevant for our discussion, starting from the two zero-range tensor forces proposed by Skyrme Skyrme (1956, 1958a)

[TABLE]

where ${\mathbf{r}}$, $\hat{{\mathbf{k}}}$ and $\hat{{\mathbf{k}}}^{\prime}$ are defined as above, Eqs. (II.2.2) and (12). The corresponding energy density functional can again be decomposed in a time-even and a time-odd part

[TABLE]

with Perlińska et al. (2004)

[TABLE]

where we already used the local gauge invariance of the energy functional Perlińska et al. (2004) for the expressions of the coupling constants. The actual expressions for the coupling constants expressed in terms of the two coupling constants $t_{e}$ and $t_{o}$ of the tensor forces are given in appendix A.

The “even” term proportional to $t_{e}$ in the two-body tensor force (20) mixes relative $S$ and $D$ waves, while the “odd” term proportional to $t_{o}$ mixes relative $P$ and $F$ waves. Thus, due to the fact that both act in spin-triplet states only, antisymmetrization implies that the former acts in isospin-singlet states (and hence contributes to the neutron-proton interaction only) and the latter in isospin-triplet states (contributing both to the like-particle and neutron-proton interactions). The central and spin-orbit interactions as we use them, however, do not contain $D$ or $F$ wave interactions. From this point of view, one might suspect a mismatch when combining the various interaction terms. From the point of view of the energy functional (II.2.4), however, all contributions from the zero-range tensor force are of the same second order in derivatives as the contributions from the non-local part of the central Skyrme force (II.2.2) and from the spin-orbit force (II.2.3).

In the time-even part of the energy functional $\mathcal{H}_{t}^{\text{t,even}}$, there appear three different combinations of the cartesian components of the spin current tensor. The term proportional to $B^{T}_{t}$ contains the symmetric combination $J_{\mu\nu}J_{\mu\nu}$ as it already appeared in the energy functional from the central Skyrme interaction (II.2.2), while the term proportional to $B^{F}_{t}$ contains two different terms, namely the antisymmetric combination $J_{\mu\nu}J_{\nu\mu}$ and the square of the trace of $J_{\nu\mu}$.

II.2.5 Combining central and tensor interactions

The Skyrme energy functional representing central, tensor, and spin-orbit interactions is given by

[TABLE]

This functional contains all possible bilinear terms up to second order in the derivatives that can be constructed from local densities and that are invariant under spatial and time inversion, rotations, and local gauge transformations Perlińska et al. (2004).

Some of the coupling constants are completely defined by the standard central Skyrme force, i.e. $C^{\rho}_{t}=A^{\rho}_{t}$, $C^{s}_{t}=A^{s}_{t}$, $C^{\tau}_{t}=A^{\tau}_{t}$, and $C^{\Delta\rho}_{t}=A^{\Delta\rho}_{t}$, two by the spin-orbit force, $C^{\nabla J}_{t}=A^{\nabla J}_{t}$, others by the tensor force, $C^{F}_{t}=B^{F}_{t}$ and $C^{\nabla s}_{t}=B^{\nabla s}_{t}$, while some are the sum of coupling constants from both central and tensor forces, $C^{T}_{t}=A^{T}_{t}+B^{T}_{t}$, and $C^{\Delta s}_{t}=A^{\Delta s}_{t}+B^{\Delta s}_{t}$.

The three terms bilinear in $J_{\mu\nu}$ can be recoupled into terms bilinear in its pseudoscalar, vector, and pseudotensor components $J^{(0)}$, $J^{(1)}$, and $J^{(2)}$, Eq. (8), which is prefered by some authors Perlińska et al. (2004)

[TABLE]

After combining (23) with the kinetic, Coulomb, pairing and other contributions from (1), the mean-field equations are obtained by standard functional derivative techniques from the total energy functional Bender et al. (2003a); Perlińska et al. (2004).

The complete Skyrme energy functional (23) has quite complicated a structure, and in the most general case leads to seven distinct mean fields in the single-particle Hamiltonian Perlińska et al. (2004). As already mentioned, we want to divide the examination of those terms that contain two derivatives and two Pauli matrices in the complete functional, i.e. those terms from the central Skyrme force that are often neglected and all the terms from the tensor Skyrme force, into three distinct steps: First, in the present paper, we enforce spherical symmetry which removes all time-odd densities and all but one out of the nine components of the spin current tensor $J_{\mu\nu}$ as will be outlined in the following section. A subsequent paper Bender et al. (2007) will discuss deformed even-even nuclei where the complete spin current tensor $J_{\mu\nu}$ is present, and future work will address the time-odd part of the energy functional (23).

II.3 The Skyrme energy functional in spherical symmetry

For the rest of this paper, we will concentrate on spherical nuclei, enforcing spherical symmetry of the $N$-body wave functions. As a consequence, the canonical single-particle wave functions $\Psi_{i}$ Dobaczewski et al. (1984) can be labeled by $j_{i}$, $\ell_{i}$ and $m_{i}$. The index $n_{i}$ labels the different states with same $j_{i}$ and $\ell_{i}$. The functions $\Psi_{i}$ separates into a radial part $\psi$ and an angular and spin part, represented by a tensor spherical harmonic $\Omega_{j\ell m}$

[TABLE]

Spherical symmetry also enforces that all magnetic substates of $\Psi_{nj\ell m}$ have the same occupation probability $v^{2}_{nj\ell m}\equiv v^{2}_{nj\ell}$ for all $-j\leq m\leq j$. For a static spherical state, all time-odd densities are zero ${\mathbf{s}}_{q}({\mathbf{r}})={\mathbf{T}}_{q}({\mathbf{r}})={\mathbf{j}}_{q}({\mathbf{r}})={\mathbf{F}}_{q}({\mathbf{r}})={\mathbf{0}}$, as are the corresponding mean fields in the single-particle Hamiltonian.

Enforcing spherical symmetry also greatly simplifies the spin-current tensor, both the pseudoscalar and pseudotensor parts of $J_{\mu\nu}$ vanish. From the vector spin-orbit current, only the radial component is non-zero, which is given by Vautherin and Brink (1972)

[TABLE]

so that there is only one out of the nine components of the spin-current tensor density that contributes in spherical nuclei. Unlike the total density $\rho$ and the kinetic density $\tau$, that are bulk properties of the nucleus and grow with the size of the nucleus, the spin-orbit current is a shell effect that shows strong fluctuations. Assume the two shells with same $n$ and $\ell$ which are split by the spin-orbit interaction, one coupled with the spin to $j=\ell+\tfrac{1}{2}$, the other to $j=\ell-\tfrac{1}{2}$. It is easy to verify that their contributions to $J_{q}(r)$ are equal but of opposite signs such that they cancel when (i) both shells are completely filled and (ii) their radial wave functions are identical $\psi_{n,\ell+1/2,\ell}=\psi_{n,\ell-1/2,\ell}$. Although the latter condition is never exactly fulfilled, this demonstrates that the spin-orbit current is not a bulk property, but a shell effect that strongly fluctuates with $N$ and $Z$. It nearly vanishes in so-called spin-saturated nuclei, where all spin-orbit partners are either completely occupied or empty, and it might be quite large when only the $j=\ell+1/2$ level out of one or even several pairs of spin-orbit partners is filled.

Altogether, the Skyrme part of the energy density functional in spherical nuclei is reduced to

[TABLE]

where we have introduced an effective coupling constant $C_{t}^{J}$ of the ${\mathbf{J}}_{t}^{2}$ tensor terms at sphericity, such that the corresponding contribution to the energy functional is given by

[TABLE]

The effective coupling constants can be separated back into contributions from the non-local central and tensor forces

[TABLE]

which are given by

[TABLE]

where we also give the expressions using the notation $T=3t_{e}$ and $U=3t_{o}$ employed in Flocard (1975); Stancu et al. (1977); Colò et al. (2007).

For the following discussion it will be also illuminating to recouple this expression to a representation that uses proton and neutron densities, where we use the notation introduced in Ref. Stancu et al. (1977)

[TABLE]

with

[TABLE]

The proton-neutron coupling constants $\alpha=\alpha_{C}+\alpha_{T}$ and $\beta=\beta_{C}+\beta_{T}$ can again be separated into contributions from central and tensor forces

[TABLE]

As could be expected, the isospin-singlet tensor force contributes only to the proton-neutron term, while the isospin-triplet tensor force contributes to both.

The spin-orbit potential of the neutrons is given by

[TABLE]

The expression for the protons is obtained exchanging the indices for protons and neutrons. In spherical symmetry, the tensor force gives a contribution to the spin-orbit potential, but does not alter the structure of the spin-orbit terms in the single-particle Hamiltonian as such. This will be different in the case of deformed mean fields Perlińska et al. (2004); Bender et al. (2007).

The dependence of the spin-orbit potential $W_{q}(r)$ on the spin-orbit current $J_{q}(r)$ through the tensor terms is the source of a potential instability. When the spin-orbit splitting becomes larger than the splitting of the centroids of single-particle states with different orbital angular momentum $\ell$, the reordering of levels might increase the number of spin-unsaturated levels, which increases the spin-orbit current $J_{n}$ and feeds back on the spin-orbit potential by increasing it even further, which ultimately leads to an unphysical shell structure. An example will be given in appendix B.

II.4 A brief history of tensor terms in the central Skyrme

energy functional

For the interpretation of the parameterizations we will describe below it is important to point out that within our choice of the effective Skyrme interaction as an antisymmetrized vertex the two coupling constants of the contribution from the central force to $\mathcal{H}^{T}$, Eq. (29), either represented through $A_{0}^{J}$, $A_{1}^{J}$ or through $\alpha_{C}$, $\beta_{C}$, are not independent from the coupling constants $A^{\tau}_{0}$, $A^{\tau}_{1}$, $A^{\Delta\rho}_{0}$, and $A^{\Delta\rho}_{1}$, that appear in Eq. (28). Through the expressions given in appendix A, all six of them are determined by the four coupling constants $t_{1}$, $x_{1}$, $t_{2}$, and $x_{2}$ from the central Skyrme force, Eq. (10). As a consequence, a tensor force is absolutely necessary to decouple the values of the $C_{t}^{J}$ from those of the $C^{\tau}_{t}$ and $C^{\Delta\rho}_{t}$, which determine the isoscalar and isovector effective masses and give the dominant contribution to the surface and surface asymmetry coefficients, respectively.

This interpretation of the Skyrme interaction is, however, far from being common practice and a source of confusion and potential inconsistencies in the literature. Many authors have used parameterizations of the central and spin-orbit Skyrme energy functional with coupling constants that in one way or the other do not exactly correspond to the functional obtained from Eqns. (10) and (16), which, depending on the point of view, can be seen as an approximation to or a generalization of the original Skyrme interaction. As the most popular modification concerns the tensor terms, a few comments on the subject are in order. Again, the practice goes back to the seminal paper by Vautherin and Brink Vautherin and Brink (1972), who state that “the contribution of this term to [the spin-orbit potential] is quite small. Since it is difficult to include such a term in the case of deformed nuclei, it has been neglected”. This choice was further motivated by the interpretation of the effective Skyrme interaction as a density-matrix expansion (DME) Negele (1970); Negele and Vautherin (1972, 1975); Campi and Bouyssy (1978). All early parameterizations as SI and SII Vautherin and Brink (1972), SIII-SVI Beiner et al. (1975), SkM Krivine et al. (1980) and SkM∗ Bartel et al. (1982) followed this example and did not contain the ${\mathbf{J}}^{2}$ terms. Beiner et al. Beiner et al. (1975) weakened the case for ${\mathbf{J}}^{2}$ terms further by pointing out that they might lead to unphysical single-particle spectra. During the 1980s and later, however, it became more popular to include them, for example in SkP Dobaczewski et al. (1984), the parameterizations T1-T9 by Tondeur et al. Tondeur et al. (1984), Eσ and Zσ by Friedrich and Reinhard Friedrich and Reinhard (1986). Some of the recent parameterizations come in pairs, where variants without and with ${\mathbf{J}}^{2}$ terms are fitted within the same fit protocol, for example (SLy4, SLy5) and (SLy6, SLy7) in Ref. Chabanat et al. (1998), or (SkO, SkO’) in Ref. Reinhard et al. (1999).

Interestingly, all but one parameterization of the central Skyrme interaction found in the literature set the coupling constants of the ${\mathbf{J}}^{2}$ terms either to their Skyrme force value (A) or strictly to zero. The exception is Ref. Tondeur (1983) by Tondeur, where an independent fit of the coupling constants of the ${\mathbf{J}}^{2}$ terms was attempted, making explicit reference to a DME interpretation of the energy functional.

Setting the coupling constants of a term to zero when one does not know how to adjust its parameters is of course an acceptable practise when permitted by the chosen framework. For Skyrme interactions fitted without the ${\mathbf{J}}^{2}$ terms, the situation becomes confusing when one looks at deformed nuclei and any situation that breaks time-reversal invariance. First of all, Galilean invariance of the energy functional dictates that the coupling constant of the ${\mathbf{s}}\cdot{\mathbf{T}}$ terms is also set to zero, as already indicated by the presentation of the energy functional in Eq. (23). Second, using a DME interpretation of the Skyrme energy functional in one place, but the interrelations from the two-body Skyrme force in all others is not entirely satisfactory. Many authors who drop the ${\mathbf{J}}^{2}$ terms rarely show scruples to keep most of the time-odd terms in the Skyrme energy functional (23) with coupling constants $A_{t}^{s}$ and $A_{t}^{\Delta s}$ from (A), although they are not at all constrained in the common fit protocols employing properties of even-even nuclei and spin-saturated nuclear matter. For a list of exceptions see Sect. II.A.2.d of Ref. Bender et al. (2003a). An alternative is to set up a hierarchy of terms, as it was attempted by Bonche, Flocard and Heenen in their mean-field and beyond codes, which set $A_{t}^{\Delta s}=0$ in addition to the coupling constant of the ${\mathbf{J}}^{2}$ terms, as all three terms have in common that they couple two Pauli matrices with two derivatives in different manners, see the footnote on page 129 of Bonche et al. (1987).

There are also inconsistent applications of parameterizations without ${\mathbf{J}}^{2}-{\mathbf{s}}\cdot{\mathbf{T}}$ terms to be found in the literature. For example, almost all applications of Skyrme interactions to the Landau parameters $g_{\ell}$ and $g_{\ell}^{\prime}$ and the properties of polarized nuclear matter, include the contribution from the ${\mathbf{s}}\cdot{\mathbf{T}}$ terms, although it should be dropped for parameterizations fitted without ${\mathbf{J}}^{2}$ terms. Similarly, most RPA and QRPA codes include them for simplicity, see the discussion in Refs. Engel et al. (1999); Bender et al. (2002); Terasaki et al. (2005).

As it is relevant for the subject of the present paper, we also mention another generalization of the Skyrme interaction that invokes the interpretation of the Skyrme energy functional in a DME framework. The spin-orbit force (16) fixes the isospin mix of the corresponding terms in the Skyrme energy functional (23) such that $A_{0}^{\nabla J}=3A_{1}^{\nabla J}$ (A). There are a few parameterizations as MSkA Sharma et al. (1995), SkI3 and SkI4 Reinhard and Flocard (1995), SkO and SkO’ Reinhard et al. (1999) and SLy10 Chabanat et al. (1998) that liberate the isospin degree of freedom in the spin-orbit functional. A DME interpretation of the energy functional is mandatory for this generalization. It is motivated by the better performance of standard relativistic mean-field models for the kink of the charge radii in Pb isotopes. Note that the standard RMF models are effective Hartree theories without exchange terms, and that the standard Lagrangians have very limited isovector degrees of freedom Bender et al. (2003a), both of which supress a strong isospin dependence of the spin-orbit interaction. It is interesting to note that the existing fits of Skyrme energy functionals with generalized spin-orbit interaction do not improve spin-orbit splittings Bender et al. (1999).

III The fits

III.1 General remarks

In order to study the effect of the ${\mathbf{J}}^{2}$ terms, we have built a set of 36 effective interactions that systematically cover the region of coupling constants $C^{J}_{0}$ and $C^{J}_{1}$ that give a reasonable description of finite nuclei in connection with the standard central and spin-orbit Skyrme forces. At variance with the perturbative approach used in Refs. Stancu et al. (1977); Colò et al. (2007), each of these parameterizations has been fitted separately, following a procedure nearly identical to that used for the construction of the SLy parameterizations Chabanat et al. (1997, 1998), so that we can keep the connection between the new fits with parameterizations that have been applied to a large variety of observables and phenomena. The Saclay-Lyon fit protocol focuses on the simultaneous reproduction of nuclear bulk properties such as binding energies and radii of finite nuclei and the empirical characteristics of infinite nuclear matter (i.e. symmetric and pure neutron matter). The latter establishes an important, though highly idealized, limiting case as it permits to confront the energy functional with calculations from first principles using the bare nucleon-nucleon force Akmal et al. (1998).

The region of effective coupling constants $(C^{J}_{0},C^{J}_{1})$ of the ${\mathbf{J}}^{2}$ terms acting in spherical nuclei, as defined in Eq. (28), that we will explore, is shown in Fig. 1. The parameterizations are labeled T$IJ$, where indices $I$ and $J$ refer to the proton-neutron ($\beta$) and like-particle ($\alpha$) coupling constants in Eq. (32) such that

[TABLE]

The corresponding values of $C^{J}_{t}$ can be obtained through Eq. (II.3) or from Fig. 1. On the one hand, we cover the positions of the most popular existing Skyrme interactions that take the ${\mathbf{J}}^{2}$ terms from the central force into account, which are SLy5 Chabanat et al. (1998), SkP Dobaczewski et al. (1984), Zσ Friedrich and Reinhard (1986), T6 Tondeur et al. (1984), SkO’ Reinhard et al. (1999) and BSk9 Goriely et al. (2005). On the other hand, among recent parameterizations including a tensor term, i.e. Skxta Brown et al. (2006), Skxtb Brown et al. (2006); Brown as well as those published by Colò et al. Colò et al. (2007) and Brink and Stancu Brink and Stancu (2007), most fall in a region of negative $C^{J}_{1}$ and vanishing $C^{J}_{0}$, that is to the lower left of Fig. 1. Parameterizations of this region, which also includes a part of the triangle advocated in the perturbative study of Stancu et al. Stancu et al. (1977), gave unsatisfactory results for many observables. Moreover, when attempting to fit parameterizations with large negative coupling constants, we sometimes obtained unrealistic single-particle spectra or even ran into the instabilities already mentioned and outlined in appendix B. Parameterizations further to the lower and upper right also have unrealistic deformations properties. The contribution from the ${\mathbf{J}}^{2}$ terms vanishes for T22, which will serve as the reference point. For the parameterizations T2$J$, only the proton-proton and neutron-neutron terms in $\mathcal{H}^{\text{t}}$ are non-zero ($\beta=0$), while for the parameterizations T$I$2, only the proton-neutron term in $\mathcal{H}^{\text{t}}$ contributes ($\alpha=0$). Note that the earlier parameterizations T6 and Zσ have a pure like-particle ${\mathbf{J}}^{2}$ terms as a consequence of the constraint $x_{1}=x_{2}=0$ employed for both (and most other early parameterizations of Skyrme’s interaction).

III.2 The fit protocol and procedure

The list of observables used to construct the cost function $\chi^{2}$ minimized during the fit (see Eq. (4.1) in Ref. Chabanat et al. (1997)) reads as follows: binding energies and charge radii of 40Ca, 48Ca, 56Ni, 90Zr, 132Sn and 208Pb; the binding energy of 100Sn; the spin-orbit splitting of the neutron $3p$ state in 208Pb; the empirical energy per particle and density at the saturation point of symmetric nuclear matter; and finally, the equation of state of neutron matter as predicted by Wiringa et al. Wiringa et al. (1995).

Furthermore, some properties of infinite nuclear matter are constrained through analytic relations between coupling constants in the same manner as they were in Refs. Chabanat et al. (1997, 1998): the incompressibility modulus $K_{\infty}$ is kept at 230 MeV, while the volume symmetry energy coefficient $a_{\tau}$ is set to 32 MeV. The isovector effective mass, expressed through the Thomas-Reiche-Kuhn sum rule enhancement factor $\kappa_{v}$, is taken such that $\kappa_{v}=0.25$.

When using a single density-dependent term in the central Skyrme force (10), the isoscalar effective mass $m^{\ast}_{0}$ cannot be chosen independently from the incompressibility modulus for a given exponent $\alpha$ of $\rho_{0}$. We follow here the prescription used for the SLy parameterizations Chabanat et al. (1997, 1998) and use $\alpha=1/6$, which leads to an isoscalar effective mass close to 0.7 in units of the bare nucleon mass for all T$IJ$ parameterizations. This value allows for a correct description of dynamical properties, as for example the energy of the giant quadrupole resonance Liu and Brown (1976). Using such a protocol we cannot reproduce the isovector effective mass consistent with recent ab-initio predictions Lesinski et al. (2006). Regarding the present exploratory study of the tensor terms this is not a critical limitation, in particular as the influence of this quantity on static properties of finite nuclei turns out to be small.

There are three modifications of the fit protocol compared to Chabanat et al. (1997, 1998). The obvious one is that the values for $C^{J}_{0}$ and $C^{J}_{1}$ are fixed beforehand as the parameters that will later on label and classify the fits. The second is that we have added the binding energies of 90Zr and 100Sn to the set of data. Indeed, we observed that the latter nucleus is usually significantly overbound when not included in the fit. The third is that we have dropped the constraint $x_{2}=-1$ that was imposed on the SLy parameterizations Chabanat et al. (1997, 1998) to ensure the stability of infinite homogeneous neutron matter against a transition into a ferromagnetic state. On the one hand, this stability criterion is completely determined by the coupling constants of the time-odd terms in the energy functional Bender et al. (2002), that we do not want to constrain here, accepting that the parameterizations might be of limited use beyond the present study. On the other hand, the tensor force brings many new contributions to the energy per particle of polarized nuclear matter that lead to a much more complex stability criterion. We postpone the entire discussion concerning the stability in polarized systems in the presence of a tensor force to future work that will also address finite-size instabilities Lesinski et al. (2006). It also has to be stressed that the actual stability criterion, as all properties of the time-odd part of the Skyrme energy functional, depends on the choices made for the interpretation of its coupling constants, i.e. antisymmetrized vertex or density functional Bender et al. (2002).

The properties of the finite nuclei entering the fit are computed using a Slater determinant without taking pairing into account. The cost function $\chi^{2}$ was minimized using a simulated annealing algorithm. The annealing schedule was an exponential one, with a characteristic time of 200 iterations (also referred to as “simulated quenching”) Thus, assuming a reasonably smooth cost function, we strive to obtain satisfactory convergence to its absolute minimum in a single run, allowing a systematic and straightforward production of a large series of forces. The coupling constants for all 36 parameterizations can be found in the Physical Review archive EPA .

Figure 2 displays the value of $\chi^{2}$ after minimization as a function of the recoupled coupling constants $\alpha$ and $\beta$. The first striking feature is the existence of a “valley” at $\beta=0$, i.e. a pure like-particle tensor term $\sim({\mathbf{J}}^{2}_{\mathrm{n}}+{\mathbf{J}}^{2}_{\mathrm{p}})$. The abrupt rise of $\chi^{2}$ around this value can be attributed to the term depending on nuclear binding energies, as sharp variations of energy residuals can be seen between neighboring magic nuclei with functionals of the T6$J$ series ($\beta=240$). For example, 48Ca and 90Zr tend to be significantly overbound in this case. We will come back later to discussing the implications for the quality of the functionals.

III.3 General properties of the fits

The coupling constants of the energy functional for spherical nuclei (28) obtained for T22 are very similar to those of SLy4, except for a slight readjustment coming from the inclusion of the binding energies of 90Zr and 100Sn in the fit as well as the abandoned constraint on $x_{2}$. With its value of $-0.945$, the $x_{2}$ obtained for T22 still stays close to the value $-1$ enforced for SLy4, which confirms that this is not too severe a constraint for parameterizations without effective ${\mathbf{J}}^{2}$ terms at sphericity. Increasing the effective tensor term coupling constants $C^{J}_{t}$, however, the values for $x_{2}$ start to deviate strongly from the region around $-1$, which is to a large extent due to the feedback from the contribution of the ${\mathbf{J}}^{2}$ terms to the surface and surface symmetry energy coefficients in the presence of constraints on isoscalar and isovector effective masses, all of which also depend on $x_{2}$. A more detailed discussion of the contribution of the ${\mathbf{J}}^{2}$ terms to the surface energy coefficients will be given elsewhere Bender et al. (2007).

From the constrained coupling constants $C^{J}_{0}$ and $C_{1}^{J}$, the respective contributions $B_{0}^{J}$ and $B_{1}^{J}$ from the tensor force can be deduced afterwards using the expressions given in Sect. II.3. Their values, shown in Fig. 3, are less regularly distributed, which is a consequence of the the non-linear interdependence of all coupling constants. Still, a general trend can be observed, such that all parameterizations are shifted towards the “south-west” compared to Fig. 1. In turn, this indicates that the contribution from the central Skyrme force always stays in the small region outlined by SkP, SLy5, Zσ, etc in Fig. 1, with values that range between $28$ and $104$ MeV fm5 for $A^{J}_{0}$ and between $38$ to $62$ MeV fm5 for $A^{J}_{1}$, respectively. This justifies a posteriori to use the tensor force as a motivation to decouple the ${\mathbf{J}}^{2}_{t}$ terms from the central part of the effective Skyrme vertex. We note in passing that all our parameterizations T$I$4 correspond to an almost pure proton-neutron or isospin-singlet tensor force, i.e. the term $\propto t_{e}$ in Eq. (20), as they are all located close to the $\alpha_{T}=0$ line.

We also find a particularly strong and systematic variation of the coupling constant $W_{0}$ of the spin-orbit force, which varies from $W_{0}=103.7$ MeV fm5 for T11 to $W_{0}=195.3$ MeV fm5 for T66, see Fig. 4. This variation is of course correlated to the strength of the tensor force. As already shown, the tensor force has the tendency to reduce the spin-orbit splittings in spin-unsaturated nuclei. To maintain a given spin-orbit splitting in such a nucleus, the spin-orbit coupling constant $W_{0}$ has to be increased.

IV Results and discussion

The calculations presented below include open-shell nuclei treated in the Hartree-Fock-Bogoliubov (HFB) framework. In the particle-particle channel, we use a zero-range interaction with a mixed surface/volume form factor (called DFTM pairing in Ref. Duguet et al. (2006)). The HFB equations were regularized with a cutoff at 60 MeV in the quasiparticle equivalent spectrum Bennaceur and Dobaczewski (2005). The pairing strength was adjusted in 120Sn with the particle-hole mean field calculated using the parameter set T33. The resulting strength was kept at the same value for all parameterizations, which is justified by the fact that the effective mass parameters are the same. Moreover, we thus avoid including, in the adjustment of the pairing strength, local effects linked with changes in details of the single-particle spectrum.

IV.1 Spin-orbit currents and potentials

As a first step in the analysis of the role of the tensor terms and their interplay with the spin-orbit interaction in spherical nuclei, we analyze the spin-orbit current density and its relative contribution to the spin-orbit potential. We choose the chain of nickel isotopes, $Z=28$, as it covers the largest number of spherical neutron shells and subshells ($N=20$, 28, 40 and 50) of any isotopic chain, two of which are spin-saturated ($N=20$ and 40), while the other two are not. Figure 5 displays the radial component of the neutron spin-orbit current ${\mathbf{J}}_{\mathrm{n}}$ for isotopes from the proton to the neutron drip-lines. The calculations are performed with T44, but the spin-orbit current is fairly independent from the parameterization. Starting from $N=20$, which corresponds to a completely filled and spin-saturated $sd$-shell, the next magic number at $N=28$ is reached by filling the $1f_{7/2}$ shell, which leads to the steeply rising bump in the plot of ${\mathbf{J}}_{\mathrm{n}}$ in the foreground, peaked around $r\simeq 3.5$ fm. Then, from $N=28$ to $N=40$ the rest of the $fp$ shell is filled, which first produces the small bump at small radii that corresponds to the filling of the $2p_{3/2}$ shell, but ultimately leads to a vanishing spin-orbit current when the $1f$ and $2p$ levels are completely filled for the $N=40$ isostope, visible as the deep valley in Fig. 5. Adding more neutrons, the filling of the $1g_{9/2}$ shell leads again to a strong neutron spin-orbit current at $N=50$. For the remaining isotopes up to the neutron drip line, the evolution of ${\mathbf{J}}_{\mathrm{n}}$ is slower with the filling of the $2d$ and $3s$ orbitals.

A few further comments are in order. First, the spin-orbit current clearly reflects the spatial probability distribution of the single-particle wave function in pairs of unsaturated spin-orbit partners. Within a given shell, the high-$\ell$ states contribute at the surface, represented by the solid line on the base of Fig. 5, while low-$\ell$ states contribute at the interior. The peak from the high-$\ell$ orbitals, however, is always located on the inside of the nuclear surface, as defined by the radius of half saturation density. Second, within a given shell, the largest contributions to the spin-orbit current density obviously come from the levels with largest $\ell$, as they have the largest degeneracy factors in (27), and because they do not have nodes, which leads to a single, sharply peaked contribution. Third, the spin-orbit current is not exactly zero for nominally “spin-saturated” nuclei, exemplified by the $N=20$ and $N=40$ isotopes in Fig. 5, as the radial single-particle wave functions are not exactly identical for all pairs of spin-orbit partners, which is a necessary requirement to obtain ${\mathbf{J}}_{\mathrm{n}}=0$ at all radii (Cf. the example of the $\nu$ $2d$ states in 132Sn in Fig. 16 below). Fourth, pairing and other correlations will always smooth the fluctuations of the spin-orbit current with nucleon numbers, as levels in the vicinity of the Fermi energy will never be completely filled or empty.

Next, we compare the contributions from the tensor terms and from the spin-orbit force to the spin-orbit potentials of protons and neutrons, Eq. (35). The contributions from the tensor force to the spin-orbit potential are proportional to the spin-orbit currents of protons and neutrons. For the Ni isotopes, the proton spin-orbit current is very similar to that of the neutrons at $N=28$ displayed in Fig. 5. For the parameterization T44 we use here as an example, we have contributions from both proton and neutron spin-orbit currents, which come with equal weights. Their combined contribution to the spin-orbit potential of the neutron $W_{n}$ might be as large as 4 MeV, see Fig. 6. This is more than a third of the maximum contribution from the spin-orbit force to $W_{n}$, see Fig. 7. The latter is proportional to a combination of the gradients of the proton and neutron densities, $2\nabla\rho_{n}(r)+\nabla\rho_{p}(r)$, see Eq. (35). As a consequence, it has a smooth behavior as a function of particle number, with slowly and monotonically varying width, depth and position. Only limited local variations can be seen on the interior due to small variations of the density profile originating from the successive filling of different orbits. Furthermore, one can easily verify that the contribution from the spin-orbit force is peaked at the surface of the nucleus (the solid line on the base plot). The strongest variation of the depth of this potential occurs just before the neutron drip line at $N=62$, where is becomes wider and shallower due to the development of a diffuse neutron skin, which reduces the gradient of the neutron density Dobaczewski et al. (1994); Lalazissis et al. (1998a, b).

Adding the contributions from the proton and neutron tensor terms to that from the spin-orbit force, the total neutron spin-orbit potential for neutrons in Ni isotopes is shown in Fig. 8. For the parameterization T44 used here (and most others in the sample of parameterizations used in this study) the dominating contributions from the spin-orbit and tensor forces to the spin-orbit potential are of opposite sign. For Ni isotopes, ${\mathbf{J}}_{p}$ is always quite large, while ${\mathbf{J}}_{n}$ varies as shown in Fig. 5. Notably, both are peaked inside of the surface. When examining the combined contribution from the spin-orbit and tensor forces to the spin-orbit potential (35), one must keep in mind that they are peaked at different radii. Moreover, the variation of tensor-term coupling constants among a set of parameterizations implies a rearrangement of the spin-orbit term strength, as will be discussed later. As a consequence, taking into account the tensor force modifies the width and localization of the spin-orbit potential $W_{q}(r)$ much more than it modifies its depth through the variation of the spin-orbit currents.

Our observations also confirm the finding of Otsuka et al. Otsuka et al. (2006) that the spin-orbit splittings might be more strongly modified by the tensor force than they are by neutron skins in neutron-rich nuclei through the reduction of the gradient of the density.

Figure 9 shows the spin-orbit potential of the protons for the chain of Ni isotopes. Here, the contribution from the spin-orbit force has a larger contribution coming from the gradient of the proton density that just grows with the mass number, without being subject to varying shell fluctuations. The same holds for the proton contribution from the tensor terms. Only the neutron contribution from the tensor terms varies rapidly, proportional to ${\mathbf{J}}_{n}$ displayed in Fig. 5, which has a very limited effect on the total spin-orbit potential, though.

With that, we can examine how the tensor terms affect the evolution of single-particle spectra. To that end, Fig. 10 shows the single-particle energies of protons and neutrons along the chain of Ni isotopes for the parameterization T22 with vanishing combined tensor terms, which will serve as a reference, while Fig. 11 shows the same for the parameterization T44 with proton-neutron and like-particle tensor terms of equal strength. For the latter, the variation of the neutron spin-orbit current with $N$ influences both neutron and proton single-particle spectra. The effect of the tensor terms is subtle, but clearly visible: for T22, the major change of the single-particle energies is their compression with increasing mass number, while for T44 the level distances oscillate on top of this background correlated to the neutron shell and sub-shell closures at $N=20$, 28, 40 and 50. As shown above, the neutron spin-orbit current vanishes for $N=20$, where it consequently has no effect on the spin-orbit potentials and splittings. By contrast, the neutron spin-orbit current is large for $N=28$ and 50, where its contribution to the spin-orbit potential reduces the splittings from the spin-orbit force.

The strong variation of the spin-orbit current with nucleon numbers is typical for light nuclei up to about mass 100. For heavier nuclei, its variation becomes much smaller. This is exemplified in Fig. 12 for the neutron spin-orbit current in the chain of Pb isotopes. There remain the fast fluctuations at small radii which we already saw for the Ni isotopes and that reflect the subsequent filling of low-$\ell$ levels with many nodes, but which have a very limited impact on the spin-orbit splittings when fed into the spin-orbit potential. The dominating peak of the spin-orbit current, just beneath the surface shows only small fluctuations, as the overlapping spin-orbit splittings of levels with different $\ell$ never give rise to a spin-saturated configuration in heavy nuclei.

Note that both the spin-orbit current ${\mathbf{J}}$ and the spin-orbit potential are exactly zero at $r=0$ as they are vectors with negative parity.

IV.2 Single-particle energies

As a next step, we analyze the modifications that the presence of ${\mathbf{J}}^{2}$ terms brings to single-particle energies in detail. Before we do so, a few general comments on the definition and interpretation of single-particle energies are in order. From an experimental point of view, empirical single-particle energies in a doubly-magic nucleus are determined as the separation energies between the even-even doubly magic nucleus and low-lying states in the adjacent odd-$A$ nuclei, i.e. they are differences of binding energies. In nuclear models, however, it is customary to discuss shell structure and single-particle energies in terms of the spectrum of eigenvalues $\epsilon_{i}$ of the Hartree-Fock mean-field Hamiltonian (in even-even nuclei), as we have done already in Figs. 10 and 11:

[TABLE]

In the nuclear EDF approach without pairing, the reference state is directly constructed as a Slater determinant of eigenstates of $\hat{h}$; hence, the corresponding eigenvalues are directly connected to the fundamental building blocks of the theory and reflect the mean field in the nucleus. The density of single-particle levels around the Fermi surface drives the magnitude of pairing correlations, the relative distance of single-particle levels at sphericity and their quantum numbers determine to a large extent the detailed structure of the deformation energy landscape which in turn, determines the collective spectroscopy. The spectroscopic properties of even-even nuclei, in particular when they exhibit shape coexistence, provide valuable benchmarks for the underlying single-particle spectrum Bender et al. (2006b). The link between the spectrum of single-particle energies on the one hand and the collective excitation spectrum on the other hand, however, always remains indirect.

On the other hand, “single-particle” states near the Fermi level of a magic nucleus can be observed by adding or removing a particle in one of these states, and thus correspond to the ground and excited states of the neighboring odd-mass nuclei. Assuming an infinitely stiff magic core, which is neither subject to any rearrangement or polarization, nor to any collective excitations following the addition (or removal) of a nucleon, the separation energies with the states in the odd-mass neighbors are equal to the single-particle energies as defined through (37). This highly idealized situation is modified by static Rutz et al. (1998) and dynamic Bernard and Nguyen Van Giai (1980); Litvinova and Ring (2006) correlations, often called “core polarization” (see chapter 7 of Ref. Bertsch (1972)) and “particle-vibration coupling” (see section 9.3.3 of Ref. Ring and Schuck (1980)) in the literature, that alter the separation energies. The main effect of the correlations is that they compress the spectrum, pulling down the levels from above the Fermi energy and pushing up those from below. The gross features, i.e. the ordering and relative placement of single-particle states, however, are more weakly affected by correlations. The particle-vibration coupling, however, is also responsible for the fractionization of the single-particle strength. When the latter is too large, the naive comparison between the calculated $\epsilon_{i}$ given by Eq. (37) and the energy of the lowest experimental state with the same quantum numbers is not even qualitatively meaningful anymore Brown et al. (2006).

We mention that a part of the static correlations originate from the non-vanishing time-odd densities in the mean-field ground-state of an odd-$A$ nucleus, that also cannot be truly spherical, so that the complete energy functional from Eq. (23) should be considered in a fully self-consistent calculation of the separation energies.

The effective single-particle energies that are used to characterize the underlying shell structure in the interacting shell model Caurier et al. (2005) have a slightly different meaning. Their definition usually renormalizes polarization and particle-vibration coupling effects around a doubly-magic nucleus whereas their evolution is discussed in terms of monopole shifts Duflo and Zuker (1999). A collection of effective single-particle energies and their evolution was collected by Grawe Grawe et al. (2005, 2006). Note that the SkX parameterization of the Skyrme energy functional by Brown and its variants Brown (1998); Brown et al. (2006) were constructed aiming at a description of effective single-particle energies along these lines.

It should be kept in mind that the obvious, coarse discrepancies between the calculated spectra of $\epsilon_{\mu}$ and the empirical single-particle energies are often larger than the uncertainties coming from the missing correlations, as long as one observes some elementary precautions. We took care to ensure that the states used in the analysis below were one-quasiparticle states weakly coupled to core phonons. First, we checked that the even-even nucleus of interest could be described as spherical, indicated by a sufficiently high-lying 2+ state. Second, we avoided all levels which were obviously correlated with the energies of 2+ states in the adjacent semi-magic series, as this indicates strong coupling with core excitations. Finally, we carefully examined states, lying above the 2+ energy and/or twice the pairing gap of adjacent semi-magic nuclei, in order to eliminate those more accurately described as an elementary core excitation coupled to one or more quasiparticles, which generally appear as a multiplet of states. We did not attempt to use energy centroids calculated with use of spectroscopic factors, as these are not systematically available. Indeed, our requirement is that if some collectivity is present, it should be similar among all nuclei considered, in order to be easily subtracted out. Empirical single-particle levels shown below are determined from the lowest states having given quantum numbers in an odd-mass nucleus.

IV.2.1 Spin-orbit splittings

The primary effect one expects from a tensor term is that it affects spin-orbit splittings by altering the strength of the spin-orbit field in spin-unsaturated nuclei, according to Eq. (35). One should remember, though, that the spin-orbit coupling itself is readjusted for each pair of coupling constants $C^{J}_{0}$, and $C^{J}_{1}$. The effect of this readjustment is generally opposite to that of the variation of the isoscalar tensor term coupling constant. It should thus be stressed that the effects described result from the balance between the variation of tensor and spin-orbit terms, which for most of our parameterizations pull into opposite directions.

Common wisdom states that the energy spacing between levels that are both above or both below the magic gap are not much affected by correlations, even when their absolute energy changes; hence it is common practice to confront only the spin-orbit splittings between pairs of particle or hole states with calculated single-particle energies from the spherical mean field. Figure 13 shows the relative error of single-particle splitting of such levels for doubly-magic nuclei throughout the chart of nuclei. The calculated values are typically 20 to $60\,\%$ larger than the experimental ones, with the exception of 16O, where the splittings of the neutron and proton $1p$ states are acceptably reproduced at least for the parameterizations T22, T24 and T42, i.e. those with the weakest tensor terms in the sample.

It is noteworthy that the calculated splittings depend much more sensitively on the tensor terms for light nuclei with spin-saturated shells (protons and neutrons in 16O, protons in 90Zr) than for the heavy doubly-magic 132Sn and 208Pb, which are quite robust against a variation of the tensor terms. The reason will become clear below.

IV.2.2 Connection between tensor and spin-orbit terms

The finding that our parameterizations systematically overestimate the spin-orbit splittings deserves an explanation. It was earlier already noted that all standard Skyrme interactions, including the SLy parameterizations that share our fit protocol, have an unresolved trend that overestimates the spin-orbit splittings in heavy nuclei Bender et al. (1999, 2003a); López-Quelle et al. (2000). Adding the tensor terms, however, further deteriorates the overall description of spin-orbit splittings, instead of improving it. It is particularly disturbing that the spin-orbit splitting of the $3p$ level in 208Pb that was used to constrain $W_{0}$ in the fit is overestimated by 30 to 40$\%$, which is larger than the relative tolerance of 20$\%$ included in the fit protocol. In fact, it turns out that the coupling constant $W_{0}$ of the spin-orbit force is more tightly constrained by the binding energies of light nuclei than by this or any other spin-orbit splitting. In the HF approach used during the fit, the structure of 40Ca, 48Ca, and 56Ni differs by the occupation of the neutron and proton $1f_{7/2}$ levels. First, we have to note that the terms in the energy functional that contain the spin-orbit current play an important role for the energy difference between 40Ca and 56Ni. The combined contribution from the tensor and spin-orbit terms varies from a near-zero value in the spin-saturated 40Ca to about $-60$ MeV in 56Ni for all our parameterizations, which is a large fraction of the $-142$ MeV difference in total binding energy between both nuclei. The $Z=40$ subshell and $Z=50$ shell are another example of abrupt variation of the spin-orbit current with the filling of the $1g_{9/2}$ level, which strongly affects the relative binding energy of $N=50$ isotones 90Zr and 100Sn. Second, the fit to phenomenological data can take advantage of the large relative variation of these terms to mock up missing physics in the energy functional that should contribute to the energy difference, but that is absent in it. The consequence will be a spurious increase of the spin-orbit and tensor term coupling constants. The resulting energy functional will correctly describe the mass difference, but not the physics of the spin-orbit and tensor terms.

In order to test the above interpretation, we performed a refit of selected T${IJ}$ parameterizations without taking into account the masses of 40Ca, 48Ca, 56Ni and 90Zr in the fit procedure. In the resulting parameterizations, the spin-orbit coefficient $W_{0}$ is typically $20\,\%$ lower than in the original ones. As a consequence, the empirical value for the spin-orbit splitting of the neutron $3p$ level in 208Pb is met well within tolerance, at the price of binding energy residuals in light nuclei being unacceptably large, i.e. 56Ni being underbound by 5 MeV while 40Ca and 90Zr are overbound by up to 10 MeV. While the global trend of the spin-orbit splittings shown in Fig. 13 is enormously improved with these fits, in particular for heavy nuclei, the overall agreement of the single-particle spectra with experiment is not, so that we had to discard these parameterizations. This finding hints at a deeply rooted deficiency of the Skyrme energy functional. The spin-orbit and, when present, tensor terms indeed do simulate missing physics of the energy functional at the price of unrealistic spin-orbit splittings. This also hints why perturbative studies, as those performed in Stancu et al. (1977); Colò et al. (2007) give much more promising results than what we will find below with our complete refits. We will discuss mass residuals in more detail in Sect. IV.3.1 below.

During the fit, the masses of light nuclei do not only compromise the spin-orbit splittings, they also establish a correlation between $W_{0}$ and $C^{J}_{0}$ in all our parameterizations. The combined spin-orbit and spin-current energy of a given spherical nucleus $(N,Z)$ is given by (keeping only the isoscalar part since we shall focus on the $N=Z$ nuclei 40Ca and 56Ni)

[TABLE]

with

[TABLE]

The difference of $E^{\text{spin}}_{0}$ between 56Ni and 40Ca

[TABLE]

turns out to be fairly independent from the parameterization. Averaged over all 36 parameterizations T$IJ$ used here, $\Delta E^{\text{spin}}$ has a value of $-58.991$ MeV with a standard deviation as small as 3.202 MeV, or $5.4\%$.

The integrals in Eqs. (39,40) are fairly independent from the actual parameterization. For a rough estimate, we can replace them in Eq. (38) by their average values. Plugged into Eq. (41) this yields

[TABLE]

Figure 14 compares the values of $C_{0}^{\nabla J}$ as obtained through (42) with the values for the actual parameterizations. The estimate works very well, which demonstrates that $C_{0}^{\nabla J}=-\tfrac{3}{4}W_{0}$ and $C^{J}_{0}$ are indeed correlated and cannot be varied independently within a high quality fit of the energy functional (28). As the combined strength of the spin-orbit and tensor terms in the energy functional is mainly determined by the mass difference of the two $N=Z$ nuclei 40Ca and 56Ni, the spin-orbit coupling constant $W_{0}$ depends more or less linearly on the isoscalar tensor coupling constant $C^{J}_{0}$, while for all practical purposes it is independent from the isovector one, see also Fig. 4 above.

IV.2.3 Splitting of high-$\ell$ states and the role of the radial

form factor

As stated above, it is common practice to confront only the spin-orbit splittings between pairs of particle or hole states with calculated single-particle energies from the spherical mean field. The spin-orbit splitting of intruder states is rarely examined. Figure 15 displays the relative deviation of the spin-orbit splittings of the intruder states with $\ell\geq 3$ that span across major shell closures and are thus given by the energy difference of a particle and a hole state. These splittings are not “safe”, i.e. they can be expected to be strongly decreased by polarization and correlation effects Rutz et al. (1998); Bernard and Nguyen Van Giai (1980); Litvinova and Ring (2006). To leave room for this effect, a mean-field calculation should overestimate the empirical spin-orbit splittings. We observe, however, that mean-field calculations done here give values that are quite close to the experimental ones, or even smaller for parameterizations with large positive isoscalar tensor coupling (cf. the evolution from T22 to T66).

This means that the spin-orbit splittings are not too large in general, as might be concluded from Fig. 13, but that there is a wrong trend of the splittings with $\ell$ with the strength of the spin-orbit potential establishing a compromise between the in-shell splittings of small $\ell$ orbits that are too large and the across-shell splittings of the intruders that are tentatively too small. In fact, the levels in Fig. 15 obviously have in common that their radial wave functions do not have nodes, while the levels in Fig. 13 have one or two nodes, with the notable exception of the $1p$ levels in 16O, for which we also find smaller deviations of the spin-orbit splittings than for the other levels in Fig. 13.

Underestimating the spin-orbit splittings of intruder levels has immediate and obvious consequences for the performance of an effective interaction, as this closes the magic gaps in the single-particle spectra and compromises the predictions for doubly-magic nuclei, as we will demonstrate in detail below. By contrast, the spin-orbit splittings of the low-$\ell$ states within the major shells have no obvious direct impact on bulk properties. Their deviation from empirical data is less dramatic, as the typical bulk observables discussed with mean-field approaches are not very sensitive to them. It is only in applications to spectroscopy that their deficiencies become evident. It is noteworthy that the parameterization T22 without effective tensor terms at sphericity provides a reasonable compromise between the tentatively underestimated splittings of the intruder levels shown in Fig. 15 and the tentatively overestimated splittings of the levels within major shells shown in Fig. 13 above, while for parameterizations with tensor terms this balance is lost.

There clearly is a proton-neutron staggering in Figs. 13 and 15, such that calculated proton splittings are relatively smaller than the neutron ones. The effect appears both when comparing proton and neutron levels with different $\ell$ in the same nucleus, and when comparing proton and neutron levels with the same $\ell$ in the same or different nuclei (see the $1h$ levels in 132Sn and 208Pb). The staggering for the intruder levels is even amplified for parameterizations with large proton-neutron tensor term, as T62, T64 or T66. The effect is particularly prominent for the heavy 132Sn and 208Pb with a large proton-to-neutron ratio $N/Z$, which might hint at unresolved isospin dependence of the spin-orbit interaction, although alternative explanations that involve how single-particle states in different shells should interact through tensor and spin-orbit forces are possible as well, see also the next paragraph.

Note that also the spin-orbit splittings of the low-$\ell$ levels shown in Fig. 13 exhibit a staggering, which is of smaller amplitude, though. It has been pointed out by Skalski Skalski (2001), that an exact treatment of the Coulomb exchange term (compared to the Slater approximation used here and nearly all existing literature) does indeed slightly increase the spin-orbit splittings of protons across major shells. This effect might give a clue to the staggering observed for the $N=Z$ nucleus 56Ni, but the magnitude of the effect reported in Skalski (2001) is too small to explain the large staggering we find for the heavier $N\neq Z$ nuclei.

Next, we use the example of 132Sn to demonstrate why the spin-orbit splittings of nodeless high-$\ell$ states are more sensitive to the tensor terms than low-$\ell$ states with one or several nodes, see Fig. 16. The lower panel shows the neutron spin-orbit potential in 132Sn for four different parameterizations, while the upper panel shows selected radial single-particle wave functions. The $\nu$ $1h_{11/2}$ and $\pi$ $1g_{9/2}$ levels give the main contribution to the neutron and proton spin-orbit currents in this nucleus, and consequently to the tensor contribution to the spin-orbit potential. Indeed, the largest differences between the spin-orbit potentials from the chosen parameterizations are caused by the varying contribution from the tensor terms and appear for the region between 3 and 6 fm, where the wave functions of the $1g$ and $1h$ states are peaked. This region corresponds to the inner flank of the spin-orbit potential well, while the outer flank is much less affected. While the $1g$ and $1h$ wave functions are peaked at the inner flank, the $2d$ orbitals have their node in this region. Consequently, the splittings of the $1g$ and $1h$ levels are strongly modified by the tensor terms, while those of the $2d$ orbitals are quite insensitive.

As a rule of thumb, the tensor contribution to the spin-orbit potential in doubly-magic nuclei comes mainly from the nodeless intruder states, which, when present, in turn mainly affect their own spin-orbit splittings, leaving the splittings of the low-$\ell$ states with one or more nodes nearly unchanged for reasons of geometrical overlap.

We note in passing that the slightly different radial wave functions of the $2d$ orbitals demonstrate nicely that their contribution to the spin-orbit current, Eq. (27), cannot completely cancel.

In fact, when regarding more specifically the evolution of the spin-orbit potential between the parameterizations T22 and T66, it is striking that for T66 it is essentially narrowed and its minimum slightly pushed towards larger radii, while its depth remains unaltered. Recalling that T66 shows a pathological behavior of too weak spin-orbit splitting of the intruder states, it appears that a correct $\ell$-dependence of spin-orbit splittings might require to modify the radial dependence of the spin-orbit potential such that it becomes wider towards smaller radii. This uncalled-for modification of the shape of the spin-orbit field has previously been put forward by Brown et al. Brown et al. (2006) as an argument for a negative like-particle ${\mathbf{J}}^{2}$ coupling constant $\alpha$. However, as will be discussed in paragraph IV.2.6 below, the evolution of single-particle levels along isotopic chains calls for $\alpha>0$, see also Brown et al. (2006). Additionally, as we will show in appendix B, large negative values of $\alpha$ pose the risk of instabilities towards the transition to states with unphysical shell structure.

IV.2.4 Single-particle spectra of doubly-magic nuclei

After we have examined the predictions for spin-orbit splittings, we will now turn to the overall quality of the single-particle spectra of doubly-magic nuclei. Figure 17 shows the single-particle spectrum of 132Sn. It is evident that as a consequence of the underestimated spin-orbit splittings of the intruder levels that we discussed in the last section, the spectrum is deteriorated for large positive isoscalar tensor term coupling constants $C^{J}_{0}$ (see T66), as, for example, a decrease of the spin-orbit splitting of the neutron $1h$ shell pushes the $1h_{11/2}$ further up, closing the $N=82$ gap. As a consequence, the presence of the tensor terms cannot remove the problem shared by all standard mean-field methods that always wrongly put the neutron $1h_{11/2}$ level above the $2d_{3/2}$ and $3s_{1/2}$ levels Bender et al. (2003a), which compromises the description of the entire mass region. For the same reason, the proton spectrum of 132Sn also excludes interactions with large positive $C^{J}_{0}$, which reduces the $Z=50$ gap between the $1g$ levels to unacceptable small values.

Figure 17 also shows the energy centroids of the $\nu$ $1h$ and $\pi$ $1g$ levels, defined as

[TABLE]

The position of the centroid is fairly independent from the parameterization. Assuming that the calculated energy of the centroid of an intruder state is more robust against corrections from core polarization and particle-vibration coupling that its spin-orbit splitting, we see that the $\nu\,1h$ centroid is clearly too high in energy by about 1 MeV. In combination with its tentatively too small spin-orbit splitting, see Fig. 15, this offers an explanation for the notorious wrong positioning of the $\nu$ $1h_{11/2}$, $2d_{3/2}$ and $3s_{1/2}$ levels in 132Sn Bender et al. (2003a). The near-degeneracy of the $\nu$ $2d_{3/2}$ and $3s_{1/2}$ levels is always well reproduced, while the $1h_{11/2}$ comes out much too high. As the $1h_{11/2}$ is the last occupied neutron level, self-consistency puts it close to the Fermi energy, which, in turn, pushes the $2d_{3/2}$ and $3s_{1/2}$ levels down in the spectrum.

The overall situation is similar for 208Pb, see Fig. 18. Again, the high-$\ell$ intruder states move too close to the $Z=82$ and $N=126$ gaps for large positive $C^{J}_{0}$. The effect is less obvious than for 132Sn as the intruders and their spin-orbit partners are further away from the gaps. Still, the level ordering and the size of the $Z=82$ gap become unacceptable for parameterizations with large tensor coupling constants. For strong tensor term coupling constants (both like-particle and proton-neutron), a $Z=92$ gap opens in the single-particle spectrum of the protons that is also frequently predicted by relativistic mean-field models Rutz et al. (1998); Bender et al. (1999) but absent in experiment Hauschild et al. (2001).

The single-particle spectra for the light doubly magic nuclei 40Ca (Fig. 19), 48Ca (Fig. 20), 56Ni (Fig. 21), 68Ni (Fig. 22) and 90Zr (Fig. 23), all have in common that the relative impact of the ${\mathbf{J}}^{2}$ terms on the ordering and relative distance of single-particle levels is even stronger than for the heavy nuclei discussed above. But not all of the strong dependence on the coupling constants of the ${\mathbf{J}}^{2}$ terms that we see in the figures is due to the actual contribution of the tensor terms to the spin-orbit potential. This is most obvious for 40Ca, where protons and neutrons are spin-saturated so that the ${\mathbf{J}}^{2}$ terms do not contribute to the spin-orbit potentials. Still, increasing their coupling constants increases the spin-orbit splittings, which manifests the readjustment of the spin-orbit force to a given set of $C^{J}_{0}$ and $C^{J}_{1}$ (see Fig. 4). The evolution of the spin-orbit splittings in 40Ca visible in Fig. 19 is the background which we have to keep in mind when discussing the impact of the tensor terms on nuclei with non-vanishing spin-orbit currents. Note that the spin-orbit coupling constant $W_{0}$ is correlated with isoscalar tensor coupling constant $C^{J}_{0}$, such that the single-particle spectra obtained with T24 and T42 are very similar, as they are for T26, T44 and T62.

For 48Ca, Fig. 20, the protons are still spin-saturated with vanishing proton spin-orbit current ${\mathbf{J}}_{\mathrm{p}}$, while for neutrons we have a large ${\mathbf{J}}_{\mathrm{n}}$. Depending on the nature of the tensor terms in the energy functional – i.e. like-particle or proton-neutron or a mixture of both – the spin-orbit current will either contribute to the spin-orbit potential of the neutrons or that of the protons or both, see Eq. (35). For the parameterizations with dominating like-particle ${\mathbf{J}}^{2}$ term, for example T24 and T26, the situation for the protons is the same as for 40Ca: there is no contribution from the tensor terms to the proton spin-orbit splittings, but compared to T22 the proton $Z=20$ gap is reduced through the readjustment of the spin-orbit force, leading to values that are too small. For the same parameterizations, the large contribution from ${\mathbf{J}}_{\mathrm{n}}$ to $W_{\mathrm{n}}$ opens up the $N=20$ gap to values that are tentatively too large, as it reduces the neutron spin-orbit splittings and thereby compensates, even overcompensates, the effect from the readjustment of the spin-orbit force. At the same time the $N=28$ gap is reduced. The opposite effect is seen for parameterizations with large proton-neutron tensor term, for example T42 or T62. For those, the proton spin-orbit splitting is reduced, opening up the $Z=20$ gap compared to T22, while the neutron spin-orbit splittings are increased by the background effect from the readjusted spin-orbit force.

For 56Ni, Fig. 21, we have large ${\mathbf{J}}_{\mathrm{n}}$ and ${\mathbf{J}}_{\mathrm{p}}$. In this $N=Z$ nucleus, the like-particle or proton-neutron parts of the tensor terms cannot be distinguished. The spectra depend only on the overall coupling constant of the isoscalar tensor term $C^{J}_{0}$, on the one hand directly through the contribution of the tensor terms to the spin-orbit potentials, and on the other hand through the background readjustment of $W_{0}$ that is correlated to $C^{J}_{0}$ as well. As already mentioned, results for T24 and T42 are very similar, as they are for T26, T44 and T62. All parameterizations have in common that the proton and neutron gaps at 28 are too small. The variation of the single-particle spectra among the parameterizations is smaller than for 40Ca, mainly because the tensor terms compensate the background drift from the readjustment of $W_{0}$.

The slightly neutron-rich 68Ni combines a spin-saturated sub-shell closure $N=40$ that gives a vanishing neutron spin-orbit current with the magic $Z=28$ that gives a strong proton spin-orbit current. The variation of the single-particle spectra in dependence of the coupling constants of the tensor terms is similar to those of 48Ca, with the roles of protons and neutrons exchanged.

The nucleus 90Zr combines the spin-saturated proton sub-shell closure $Z=40$ with the major neutron shell closure $N=50$. The high degeneracy of the occupied $\nu$ $1g_{9/2}$ level leads to a very strong neutron spin-orbit current, while the proton spin-orbit current is zero. Even in the absence of a tensor term contributing to their spin-orbit potential for parameterizations with pure like-particle tensor terms, the proton single-particle spectra are dramatically changed by the feedback effect from the readjusted spin-orbit force; see the evolution from T22 to T26. The $\pi$ $1g_{9/2}$ comes down, and closes the $Z=40$ sub-shell gap. For parameterizations with pure proton-neutron tensor term, one has the opposite effect, this time because the contribution from the tensor terms overcompensates the background effect from the spin-orbit force. The effect of the tensor terms on the neutron spin-orbit splittings is less dramatic, but still might be sizable.

We have to point out that the calculations displayed in Fig. 23 were performed without taking pairing into account, as the HFB scheme breaks down in the weak pairing regime of doubly magic nuclei. For some extreme (and unrealistic) parameterizations, however, the gaps disappear which, in turn, would lead to strong pairing correlations if the calculations were performed within the HFB scheme. This happens, for example, for neutrons in 90Zr when using T26 and T46. Interestingly, the pairing correlations for neutrons break the spin saturation, which leads to a substantial neutron spin-orbit current ${\mathbf{J}}_{\mathrm{n}}$. As these parameterizations use values of the like-particle coupling constant significantly larger than the neutron-proton one, ${\mathbf{J}}_{\mathrm{n}}$ feeds back onto the neutron spin-orbit potential only, Eq. (35). As the corresponding coupling constant $\alpha$ is positive for T26 and T46, the contribution from the tensor terms reduces the spin-orbit splittings, in particular those of the $1g_{9/2}$ and $1f_{5/2}$. As a result, this counteracts the reduction of the $N=40$ gap predicted by T26 and T46 in calculations without pairing.

IV.2.5 Evolution along isotopic chains: $np$ coupling

In the preceding sections, we have analyzed characteristics of the single-particle spectra for isolated doubly-magic nuclei. We found that larger tensor terms do not lead to an overall improvement of the single-particle spectra. However, we also argued that it might be essentially due to deficiencies of the central (and possibly spin-orbit) interactions and that it should not be used to discard the tensor terms as such. In any case, the results gathered so far on single-particle spectra of doubly-magic nuclei do not permit to narrow down a region of meaningful coupling constants of the tensor terms. The analysis must be complemented by looking at other observables. A better suited observable is provided by the evolution of spin-orbit splittings along an isotopic or isotonic chain, which ideally reflects the nucleon-number-dependent contribution from the ${\mathbf{J}}^{2}$ terms to the spin-orbit potentials. Unfortunately, safe experimental data for the evolution of spin-orbit partners are scarce; hence, one has to content oneself to the evolution of the energy distance of levels with different $\ell$, assuming that the effect is primarily caused by the evolution of the spin-orbit splittings of each level with its respective partner. A popular playground for such studies is the chain of Sn isotopes, where two such pairs of levels have gained attention; the $\pi$ $2d_{5/2}$ and $\pi$ $1g_{7/2}$ on the one hand, and the $\pi$ $1g_{7/2}$ and $\pi$ $1h_{11/2}$ on the other hand. Figure 24 shows these two sets of results for a selection of our parameterizations.

Experimentally, the $2d_{5/2}$ and $1g_{7/2}$ levels cross between $N=70$ and 72, such that the $2d_{5/2}$ provides the ground state of light odd-$A$ Sb isotopes, and $1g_{7/2}$ that of the heavy ones, see for example Ref. Shergur et al. (2005). The crossing as such is predicted by many mean-field interactions and most of the parameterizations of the Skyrme interaction we use here. It has also been studied in detail with the standard Gogny force (without any tensor term) using elaborate blocking calculations of the odd-$A$ nuclei Porquet et al. (2005). The crossing, however, is never predicted at the right neutron number, see Fig. 24. As we have learned above, we should not assume that the absolute distance of the two levels will be correctly described by any of our parameterizations (as the centroids of the $\ell$ shells will not have the proper distance and the spin-orbit splittings have a wrong $\ell$ dependence within a given shell). Hence, the neutron number where the crossing takes place cannot and should not be used as a quality criterion. What does characterize the tensor terms is the bend of the curves in Fig. 24, as ideally it reflects how the spin-orbit splittings of both levels change in the presence of the tensor terms. Similar caution has to be exercised in the analysis of the unusual relative evolution of the proton $1g_{7/2}$ and $1h_{11/2}$ levels that was brought to attention by Schieffer et al. Schiffer et al. (2004). Their spacing has been investigated in terms of the tensor force before Otsuka et al. (2005, 2006); Brown et al. (2006); Colò et al. (2007). Again, we pay attention to the qualitative nature of the bend without focusing too much on the precise value by which the splitting changes when going from $N\approx 58$ to $N=82$. Indeed, the matching of the lowest proton fragment with quantum number $1h_{11/2}$ seen experimentally with the corresponding empirical single-particle energy is unsafe because of the fractionization of the strength as discussed in Ref. Brown et al. (2006).

For both pairs of levels, the evolution of their distance can be attributed to the tensor coupling between the proton levels and neutrons filling the $1h_{11/2}$ level below the $N=82$ gap. Unfortunately, this introduces an additional source of uncertainty: as can be seen in Fig. 17, the ordering of the neutron levels in 132Sn is not properly reproduced by any of our parameterizations, with the $1h_{11/2}$ level being predicted above the $2d_{3/2}$ level, while it is the other way round in experiment. This means that in the calculations, the contribution from the $1h_{11/2}$ level to the neutron spin-orbit current builds up at larger $N$ than what can be expected in experiment. As a consequence, the prediction for the relative evolution of the levels might be shifted by up to four mass units to the right compared to experiment for both pairs of levels we examine here.

In the end, the trend of both splittings is best reproduced when using a positive value of the neutron-proton ${\mathbf{J}}_{\mathrm{n}}\cdot{\mathbf{J}}_{\mathrm{p}}$ coupling constant $\beta$ such that the filling of the neutron $1h_{11/2}$ shell decreases the spin-orbit splittings of the proton shells. The parameterizations from the T4$J$ and T6$J$ series indeed do reproduce the bend of empirical data, with, however, a clear shift in the neutron number where it occurs, as expected from the previous discussion. A value of $\beta=120$ MeV fm5, which corresponds to the series of T4$J$ parameterizations, matches its magnitude best (see for example T44).

A similar analysis can be performed for the proton $1f_{5/2}$ and $2p_{3/2}$ levels in the chain of Ni isotopes, see Fig. 25. This case is interesting as no distinctive feature can be observed in the empirical spectra, yet the standard parameterizations without tensor terms like T22 do not reproduce them. In fact, to keep the $1f_{5/2}$ and $2p_{3/2}$ at a constant distance, two competing effects have to cancel. First, the increasing diffuseness of the neutron density with increasing neutron number diminishes the proton spin-orbit splittings through its reduced gradient in the expression for the proton spin-orbit potential when going from $N=32$ to $N=40$. Second, the filling of the neutron $1f_{5/2}$ state reduces the neutron spin-orbit current which in turn increases the proton spin-orbit splittings for interactions with sizable proton-neutron tensor contribution to the proton spin-orbit potential when going from $N=32$ to $N=40$. The former effect can be clearly seen for parameterizations T2$J$ with vanishing proton-neutron tensor term, $\beta=0$. Again, parameterizations of the T4$J$ series seem to be the most appropriate to describe the evolution of these levels.

The evolution of single-particle levels is the tool of choice to determine the sign and magnitude of the proton-neutron tensor coupling constant. The value which we favor, as a result of our semi-qualitative analysis is $\beta=120$ MeV fm5. This value is only slightly larger than the value of 94 to 96 MeV fm5 advocated by Brown et al. in Ref. Brown et al. (2006), which was adjusted to theoretical level shifts in the chain of tin isotopes obtained from a $G$-matrix interaction. We can consider this as a reasonable agreement.

Let us defer the discussion of this value to the end of this section and study in the next paragraph the like-particle tensor-term coupling constant $\alpha$.

IV.2.6 Evolution along isotopic chains: $nn$ coupling

In order to narrow down an empirical value for the neutron-neutron tensor coupling constant, the ideal observable would be the evolution of neutron single-particle levels along an isotopic chain. Unfortunately, these are only accessible at the respective shell closures. We shall therefore compare neutron single-particle spectra of pairs of doubly-magic nuclei belonging to the same isotopic chain. Again, the necessity to extract pure single-particle effects calls for precautions. We choose pairs of particle or hole levels which are close enough in energy that their absolute spacing is not much affected by particle-vibration coupling. Of course, one also has to be careful if both states appear at relatively high excitation energy in the neighboring odd isotope because the fractionization of their strength could again interfere with the analysis. In the following, we choose pairs of orbitals which are as safe as possible.

To remove the uncertainties from the deficiencies of the central and spin-orbit parts of the effective interaction that we have identified above, we will look at a double difference, where, first, we construct the energy difference between the neutron $1d_{3/2}$ and $2s_{1/2}$ levels separately for 40Ca and 48Ca, and then compare the value of this difference in both nuclei

[TABLE]

Assuming that the problems from the central and spin-orbit forces discussed in Sects. IV.2.1 and IV.2.4 have the same effect in both nuclei, they will cancel out in $\delta^{\text{Ca}}$.

The interesting feature of this pair of states is that they are separated by more than 2 MeV in 40Ca, while they are nearly degenerate in 48Ca, see Figs. 19 and 20. Such a shift can only be reproduced with a positive (140-180 MeV fm5) value of $\alpha$, which decreases the splitting of the neutron $1d$ shell when the neutron $1f_{7/2}$ level is filled.

A similar analysis can be performed for the $1f_{5/2}$ and $2p_{1/2}$ neutron states in the Ni isotopes 56Ni and 68Ni

[TABLE]

Going from 56Ni to 68Ni, the neutron $1f_{5/2}$ level comes further down in energy than the $2p_{1/2}$ level for parameterizations without tensor terms (T22), see Figs. 21 and 22. The reason for this trend is the geometrical growth of the nucleus, which on the one hand lowers the centroid of the $1f$ levels in the widening potential well, and on the other hand pushes the spin-orbit field to larger radii, which has opposite effects on the splittings of $2p$ and $1f$ states. The like-particle tensor terms can compensate this trend through a reduction of the spin-orbit splitting of the $1f$ levels. The observed downward shift by $0.3$ MeV can be recovered with a value of $\alpha$ around 120 MeV fm5, see Fig. 26.

It is also gratifying to see that the analysis of Ca and Ni isotopes suggests nearly the same value for the like-particle tensor term coupling constant $\alpha$.

IV.3 Binding energies

Our ultimate goal, although far beyond the scope of the present paper, is the construction of a universal nuclear energy density functional that simultaneously describes bulk properties like masses and radii, giant resonances, and low-energy spectroscopy, such as quasiparticle configurations and collective rotational and vibrational states. To crosscheck how our findings on single-particle spectra and spin-orbit splittings translate into bulk properties, we will now analyze the evolution of mass residuals and charge radii along isotopic and isotonic chains. It has been repeatedly noted in the literature that the mass residuals from mean-field calculations show characteristic arches Dobaczewski et al. (1984); Friedrich and Reinhard (1986); Chabanat et al. (1998); Patyk et al. (1999); Bender et al. (2003a); Lunney et al. (2003); Dobaczewski et al. (2004); Bender et al. (2006a), where heavy mid-shell nuclei are usually underbound compared to the doubly magic ones that are located at the bottom of deep ravines. For light nuclei, the patterns are often less obvious. Part of this effect can be explained and removed taking large-amplitude correlations from collective shape degrees of freedom into account through suitable beyond-mean-field methods. In turn, this means that the mass residuals should leave room for the extra binding of mid-shell nuclei from correlations. However, it turns out that for typical effective interactions the amplitude of the arches is larger than what is brought by correlations Bender et al. (2006a). Furthermore, this effect seems not to be of the same size for isotopic and isotonic chains, which altogether hints at deficiencies of the current effective interactions.

Recently, Dobaczewski pointed out Dobaczewski (2006) that the strongly fluctuating contribution brought by the ${\mathbf{J}}^{2}$ terms to the total binding energy could remove at least some of the ravines found in the mass residuals around magic numbers. The hypothesis was motivated by calculations that evaluate the tensor terms either perturbatively, or self-consistently, using in this case an existing standard parameterization without tensor terms for the rest of the energy functional. Our set of refitted parameterizations with varied coupling constants of the tensor terms gives us a tool to check how much of the argument persists to a full fit.

IV.3.1 Semi-magic series

Figure 27 displays binding energy residuals along various isotopic and isotonic chains of semi-magic nuclei for a selection of our parameterizations: T22 is the reference with vanishing ${\mathbf{J}}^{2}$ terms at sphericity; T24 has a substantial like-particle coupling constant $\alpha$ and vanishing proton-neutron coupling constant $\beta$, which is similar to most of the published parameterizations which take the ${\mathbf{J}}^{2}$ terms from the central Skyrme force into account; T42 and T62 are parameterizations with substantial proton-neutron coupling constant $\beta$ and vanishing like-particle coupling constant; T44 has a mixture of like-particle and proton-neutron tensor terms that is close to what we found preferable for the evolution of spin-orbit splittings above; and T46 is a parameterization that gives the best root-mean-square residual of binding energies for spherical nuclei, as we will see below. Finally, T66 is a parameterization with large and equal proton-neutron and like-particle tensor-term coupling constants.

The tensor terms have opposite effects in light and heavy nuclei: The curves obtained with T22, the parameterization without ${\mathbf{J}}^{2}$ term contribution at sphericity, are relatively flat for the light isotopic and isotonic chains, but show very pronounced arches with an amplitude of 5 or even more MeV for the heavy Sn and Pb isotopic chains. By contrast, the most striking effect of the ${\mathbf{J}}^{2}$ terms is that they induce large fluctuations of the mass residuals in light nuclei, while they flatten the curves in the heavy ones.

The strong variation between the parameter sets for light nuclei are of course the direct consequence of the strong variation of the spin-orbit current ${\mathbf{J}}$ that enters the spin-orbit and tensor terms when going back and forth between nuclei where the configuration of at least one nucleon species is spin-saturated. The variations seen are a result of the modifications of tensor-term coupling constants and the associated readjustment of the spin-orbit strength $W_{0}$. For example, 48Ca is overbound with respect to 40Ca and 56Ni for parameterizations with a proton-neutron coupling constant $\beta>0$, while the like-particle coupling constant $\alpha$ has a more limited effect. Since only the neutron core is spin-unsaturated in this nucleus, this must be attributed to the increase in the readjusted spin-orbit strength $W_{0}$ (correlated with $C^{J}_{0}=\tfrac{1}{2}(\alpha+\beta)$) which dominates when $\beta$ is increased and $\alpha$ kept at zero, and counterbalances the effect of $\alpha$ when the latter varies. See the parameter sets T62 and T66 in Figures 27 and 28. The large overbinding of nuclei around 90Zr ($Z=40$, $N=50$) for parameterizations with large proton-neutron tensor coupling constant has the same origin. For a given parameterization and a given nucleus, the energy gain from the spin-orbit term seems to be almost always larger than the energy loss from the ${\mathbf{J}}^{2}$ one, see Fig. 28 for Ca isotopes and Fig. 29 for Sn isotopes. Of course other terms in the energy functional compensate for a part of the gain from the spin-orbit term, but the overall trends of the mass residuals suggest that the spin-orbit energy has a much larger contribution to the differences between the parameterizations visible in Fig. 27 than the ${\mathbf{J}}^{2}$ terms.

We have to note that the spin-orbit current does not completely vanish for the nominally proton and neutron spin-saturated 40Ca for parameterizations with large coupling constants of the ${\mathbf{J}}^{2}$ terms. For those, the gap at 20 is strongly (and nonphysically) reduced, see Fig. 19. The small gap at 20 does not suppress pairing correlations anymore in our HFB approach. The resulting scattering of particles from the $sd$ shell to the $fp$ shell breaks the spin-saturation, such that there is a finite, in some cases quite sizable, contribution from the spin-orbit term to the total binding energy. Owing to the compensation between all contributions, the total energy gain compared to a HF calculation without pairing is usually small and rests on the order of 200 keV for the parameterizations shown in Fig. 27.

It is also important to note that some of the light chains in Fig. 27 are sufficiently close to or even cross the $N=Z$ line that they are subject to the Wigner energy, which still lacks a satisfying explanation, not to mention a description in the framework of mean-field methods Satuła et al. (1997). The Wigner energy is not taken into account in our fits, while it turned out to be a crucial ingredient of any HFB Tondeur et al. (2000); Samyn et al. (2002); Goriely et al. (2003) or other mass formula. In fact, as shown in Fig. 14 of Ref. Bender et al. (2006a), the missing Wigner energy clearly sticks out from the mass residuals for SLy4 (which is very similar to T22) when they are plotted for isobaric chains. This local trend around $N=Z$ is, however, overlaced with a global trend with mass number, such that the missing Wigner energy cannot be spotted anymore when looking at the mass residuals for the isotopic chain of Ca isotopes, similar to what is seen for T22 in Fig. 27. Within our fit protocol, the correlation between the masses of 40Ca, 48Ca and 56Ni, that is brought by the spin-orbit force (see Sect. IV.2.2) does not tolerate a correction for the Wigner energy for standard central and spin-orbit Skyrme forces, as this will lead to an unacceptable underbinding of 48Ca. This, however, might change when the ${\mathbf{J}}^{2}$ terms are added. Indeed, Fig. 27 suggests that adding a phenomenological Wigner term around 40Ca and 56Ni to a parameter set like T44, which is consistent with the evolution of single-particle levels, would flatten the curves for the mass residuals in the Ca, Ni and $N=28$ chains. The mass residuals for the chain of oxygen isotopes that are not shown here would be improved in a similar manner. However, extreme caution should be exercised before jumping to premature conclusions, as the spin-orbit splittings and level distances in light nuclei are far from realistic for all our parameterizations; as a consequence it is difficult to judge if the room we find for the Wigner energy is fortuitous or indeed a feature of well-tuned ${\mathbf{J}}^{2}$ terms. Note that the HFB mass formulas that do include a correction for the Wigner energy side-by-side with the ${\mathbf{J}}^{2}$ terms from the central Skyrme force give satisfying mass residuals for light nuclei Tondeur et al. (2000); Samyn et al. (2002); Goriely et al. (2003), but have nuclear matter properties that are quite different from ours; cf. BSk1 and BSk6 with SLy4 in Table I of Ref. Reinhard et al. (2006). Our constraints on the empirical nuclear matter properties (same as those on SLy4) that are absent in these HFB mass fits might be the deeper reason for this conflict.

Large tensor-term coupling constants straighten the arches in the mass residuals in the heavy Sn and Pb isotopic chains, but the improvements are not completely satisfactory. Large, combined proton-neutron and like-particle coupling constants tend to transform the arch for the tin isotopic chain into a an s-shaped curve, which is not very realistic from the standpoint of expected corrections through collective effects. It can again be assumed that the deficiencies of the single-particle spectra pointed out in Fig. 17 are responsible, where the $\nu$ $1h_{11/2}$ and $\pi$ $1g_{9/2}$ are placed too high above the rest of the single-particle spectra in heavy Sn isotopes. For Pb isotopes, large values of the tensor terms tend to overbind the neutron-deficient isotopes. It is noteworthy that the tensor terms seem to not much affect the mass residuals of the heavy Pb isotopes above $N=126$, which are on the flank of a very deep ravine that becomes visible when going towards heavier elements, cf. the SLy4 results in Ref. Bender et al. (2006a).

It has been often noted that effective interactions that give a similar satisfying description of masses close to the valley of stability give diverging predictions when extrapolated to exotic nuclei. The standard example is the two-neutron separation energy $S_{2n}(N,Z)=E(N,Z-2)-E(N,Z)$ for the chain of Sn isotopes. Results obtained with a subset of our parameterizations are shown in Fig. 30. It is noteworthy that the differences for neutron-rich nuclei beyond $N=82$ are not larger than those for the isotopes closer to stability. Around the valley of stability, increasing the coupling constants of tensor terms, in particular the like-particle ones, tilts the curve, pushing it up for light isotopes and pulling it down it for heavy ones, which reflects of course the position of the $\nu$ $1h_{11/2}$ level that is pushed into the $N=82$ gap, see Fig. 17. For the neutron-rich isotopes, small differences appear around $N=90$, which reflects the change of level structure above the $\nu$ $2f_{7/2}$ level and at the drip line, but they are much smaller than the differences seen between parameterizations obtained with different fit protocols, see Fig. 5 of Ref. Bender et al. (2003a).

IV.3.2 Systematics

In the preceding section we showed how the ${\mathbf{J}}^{2}$ terms in the energy functional modify the trends of mass residuals along isotopic and isotonic chains, in particular the amplitude of the arches between doubly-magic nuclei. In this section, we want to examine how this translates into quality criteria for the overall performance of the parameterizations for masses.

Figure 31 displays the root-mean-square deviation of the mass residuals for all our 36 parameterizations, evaluated for a set of 134 nuclei predicted to have spherical mean-field ground states when calculated with the parameterizations SLy4 Bender et al. (2006a). One observes a clear minimum around T46, i.e. $(\alpha,\beta)=(240,120)$, with $(E_{\text{th}}-E_{\text{exp}})_{\text{r.m.s.}}=1.96$ MeV, compared with $3.44$ MeV for T22 ($\alpha=\beta=0$). We found even slightly better values with even more repulsive isoscalar and isovector coupling constants, but the single-particle spectra of these interactions turn out to be quite unrealistic, cf. Sect. IV.2.1. This already demonstrates that in the presence of the ${\mathbf{J}}^{2}$ terms a good fit of masses does not necessarily lead to satisfactory single-particle spectra.

Figure 32 demonstrates how the distribution of the mass residuals $E_{\text{th}}-E_{\text{exp}}$ affects the evolution of their r.m.s. value for a subset of 9 parameterizations. For T22 ($\alpha=\beta=0$), the distribution is centered at positive mass residuals, with only very few nuclei being overbound. Increasing $\beta$ to 120 MeV fm5 (T42) or even 240 MeV fm5 (T62) shifts the median of the distribution to smaller values, which yields more and more overbound nuclei. For large values of $\beta$, the distribution spreads out more, which diminishes the improvement from centering the distribution closer to zero. For given $\beta$, increasing $\alpha$ mainly shifts the median of the distribution without spreading out its overall shape, which is preferable to optimize the r.m.s. value.

These considerations, however, have to be taken with caution. As said above, we aim at a model where certain correlations beyond the mean-field are treated explicitly, which asks for a distribution of mean-field mass residuals with an asymmetric distribution towards positive mass residuals, and a width that is similar to the difference between the maximum and minimum correlation energies to be found.

IV.4 Radii

The evolution of nuclear charge radii along isotopic chains reflects how the mean field of the protons changes when neutrons are added in the system. In the simplistic liquid-drop model, it just follows the geometrical growth of the nucleus $\sim A^{1/3}$, but data show that there are many local deviations from this global trend. On the one hand, radii are of course subject to correlations beyond the mean field Reinhard and Drechsel (1979); Girod and Reinhard (1982); Bonche et al. (1991); Heenen et al. (1993); Bender et al. (2006a) On the other hand, they are also sensitive to the detailed shell structure, which, in turn, might be influenced by tensor terms. We will concentrate here on two anomalies of the evolution of charge radii, both of which are not much influenced by collective correlations beyond the mean-field (at least in calculations with the Skyrme interaction SLy4) Bender et al. (2006a): that the root-mean-square (r.m.s.) charge radius of 48Ca is almost the same as the one of the lighter 40Ca or possibly slightly smaller, and the kink in the isotopic shifts of mean-square (m.s.) charge radii in the Pb isotopes, where Pb isotopes above 208Pb are larger than what could be expected from liquid-drop systematics. In both cases it is plausible that shell effects are the determining factor, although alternative explanations that involve pairing effects have been put forward for the latter case as well Tajima et al. (1993); Fayans et al. (2000).

Charge radii have been calculated with the approximation used in Ref. Chabanat et al. (1997)333There is a typographical error in Eq. (4.2) in Ref. Chabanat et al. (1997), that was copied to Eq. (110) in Ref. Bender et al. (2003a): the $\hbar/mc$ factor should be squared, as is trivially found by dimensional analysis and confirmed by Ref. Bertozzi et al. (1972). and derived from Ref. Bertozzi et al. (1972)

[TABLE]

where the mean-square (m.s.) radius of the point-proton distribution $\langle r^{2}\rangle_{\mathrm{p}}$ is corrected by three terms: the first two estimate the effects of the intrinsic charge distribution of the free proton and neutron (with m.s. radii $r^{2}_{\mathrm{p}}$ and $r^{2}_{\mathrm{n}}$) and the third adds a correction from the magnetic moments of the nucleons. Since we will consider the shift of charge radii for different isotopes of the same series, the actual value of $r^{2}_{\mathrm{p}}$ cancels out. For the second correction term, which is independent from the interaction, we take $r^{2}_{\mathrm{n}}=-0.117$ fm2 Bender et al. (2003a). Finally, the magnetic correction can only depend weakly on the details of the interaction through the occupation factors $v^{2}_{i}$ when non-magic nuclei are considered. The same expressions had been used during the fit of our parameterizations.

We begin with the Ca isotopes. Most parameterizations of Skyrme’s interaction are not able to reproduce that the charge radius of 48Ca has about the same size as that of 40Ca, see Fig. 11 in Ref. Bender et al. (2003a). The middle panel of Fig. 33 shows the difference of the the m.s. radii of 48Ca and 40Ca in dependence of the tensor term coupling constants $\alpha$ and $\beta$. First, this difference is almost independent of $\alpha$, the strength of the like-particle tensor terms. Second, it is strongly correlated with $\beta$, the strength of the proton-neutron tensor term, with large positive values of $\beta$ bringing the difference of radii into the domain of experimentally acceptable values Otten (1989) or even below, with a best match obtained for $\beta=80$ MeV fm5. This effect can be explained by looking at the proton single-particle spectra of 40Ca (Fig. 19) and 48Ca (Fig. 20). Indeed, one observes that a positive neutron-proton tensor coupling constant decreases the strength of the proton spin-orbit field in 48Ca, which in turn lowers the $\pi$ $1d_{3/2}$ level in 48Ca (compare the parameterizations T$IJ$ in Fig. 20 with increasing $I$ for given $J$). As a consequence, the m.s. radius of this state decreases as it sinks deeper into the potential well of 48Ca. At the same time, this level is pushed up in 40Ca, which slightly increases the contribution of this state to the charge m.s. radius of this nucleus. This effect is demonstrated in the top panel of Fig. 33, which displays the degeneracy-weighted and normalized change of the m.s. radii of proton hole states between 40Ca and 48Ca as a function of the proton-neutron tensor term coupling constant $\beta$ for forces with a like-particle tensor term coupling constant $\alpha=120$ MeV fm5. Indeed, the decreasing contribution from the $\pi 1d_{3/2}$ state to the m.s. radius significantly decreases the isotopic shift between both Ca isotopes. It has to be noted that the m.s. value of the charge radii of 40Ca and 48Ca are almost independent of alpha and that their absolute values are not reproduced for any of our parameterizations.

The latter study demonstrates the correlation between the isotopic shift of m.s. charge radius between 40Ca and 48Ca and the absolute single-particle energy of the proton $1d_{3/2}$ state. This level can be moved around within the single-particle spectrum with the ${\mathbf{J}}^{2}$ terms. However, the agreement of the calculated single-particle energy of the proton $1d_{3/2}$ state in both nuclei with experiment is not necessarily improved for the parameterizations that reproduce the isotopic shift of the m.s. charge radius. Furthermore, a good reproduction of the isotopic shift does not guarantee that the absolute values of the charge radii are well reproduced, see the bottom panel in Fig. 33. In fact, they are predicted too large for all of our parameterizations, which again points to deficiencies of the central field. Altogether, this suggests that in spite of its sensitivity to the coupling constants of the ${\mathbf{J}}^{2}$ terms, the isotopic shift of m.s. charge radius between 40Ca and 48Ca should not be used to constrain them before one has gained sufficient control over the central interaction.

A few further words of caution are in place. The charge radii of all light nuclei are significantly increased by dynamical quadrupole correlations, see Fig. 23 of Ref. Bender et al. (2006a). Correlations beyond the static self-consistent mean field are also at the origin of the arch of the ms charge radii between 40Ca and 48Ca that is neither reproduced by any pure mean-field model, see again Fig. 11 in Ref. Bender et al. (2003a), nor by the beyond-mean-field calculations with SLy4 of Ref. Bender et al. (2006a), while the shell model allows for a satisfactory description Caurier et al. (2001).

Many explanations have been put forward to explain the kink in the isotopic shifts of Pb radii. As it qualitatively appears in relativistic mean-field models, but not in non-relativistic ones using the standard spin-orbit interaction (16), it has been used as a motivation to generalize the isospin mix of the standard spin-orbit energy density functional, Eq. (II.2.3), to simulate the isospin dependence of the relativistic Hartree models Sharma et al. (1995); Reinhard and Flocard (1995). The resulting parameterizations are not completely satisfactory, as the price for the improvement of the radii is a further deterioration of spin-orbit splittings Bender et al. (1999), while the relativistic mean field gives a satisfactory description of both. Some standard Skyrme interactions that take the tensor terms from the central Skyrme force into account also give a kink, but it is by far too small to reproduce the experimental values Chabanat et al. (1998).

Plotting the m.s. radii along the chain of Pb isotopes as a function of $N$, the slopes are nearly linear when looking separately at the isotopes below and above 208Pb. We will concentrate on the change in the slope at 208Pb that is brought by the tensor terms, which can be quantified through the second finite difference of the m.s. radii at 208Pb

[TABLE]

There are two conflicting values to be found in the literature, either $46.4\pm 1.4$ fm2 Otten (1989) and the significantly larger $59\pm 3$ fm2 Angeli (2004). Figure 34 shows the change of slope around 208Pb as defined through Eq. (IV.4) as a function of the like-particle tensor coupling constant $\alpha$ and for three different values of $\beta$. It is striking to see that this quantity is almost independent of the neutron-proton tensor coupling constant $\beta$, so the change is mainly induced by the tensor interaction between particles of the same kind. It has been noted before that the kink in the isotopic shift of the charge radii in Pb isotopes is correlated to the single-particle spectrum of neutrons above $N=126$, in particular the position of the $1i_{11/2}$ level. (This has to be contrasted with the Ca isotopic chain discussed above, where the difference of charge radii between 40Ca and 48Ca appears to be particularly sensitive to the single-particle spectrum of the protons.) The closer the $1i_{11/2}$ level is to the $2g_{9/2}$ level that is filled above $N=126$, the more the $1i_{11/2}$ becomes occupied through pairing correlations. Through the shape of its radial wave function, the partial filling of the nodeless $1i_{11/2}$ increases the neutron radius faster than filling only the $2g_{9/2}$, and in particular faster than for the isotopes below $N=126$. As the protons follow the density distribution of the neutrons, the charge radius grows rapidly beyond $N=126$. This offers an explanation why the kink increases with the like-particle tensor term coupling constant $\alpha$: for large values of the weight $\alpha$ of the neutron spin-orbit current in the neutron spin-orbit potential, Eq. (35), the spin-orbit splitting of the $\nu$ $1i$ levels is reduced such that the $1i_{11/2}$ approaches the $2g_{9/2}$ level in 208Pb, see Fig. 18.

While the kink is clearly sensitive to the tensor terms, they cannot be responsible for the entire effect, as even for extreme parameterizations that give unrealistic single-particle spectra the calculated kink hardly reaches about three quarters of its experimental value.

V Summary and conclusions

We have reported a systematic study of the effects of the ${\mathbf{J}}^{2}$ (tensor) terms in the Skyrme energy functional for spherical nuclei. The aim of the present study was not to obtain a unique best fit of the Skyrme energy functional with tensor terms, but to analyze the impact of the tensor terms on a large variety of observables in calculations at a pure mean-field level and to identify, if possible, observables that are particularly, even uniquely, sensitive to the ${\mathbf{J}}^{2}$ terms. To reach our goal, we have built a set of 36 parameterizations that cover the two-dimensional parameter space of the coupling constants of the ${\mathbf{J}}_{t}^{2}$ terms that does not give obviously unphysical predictions for a wide variety of observables we have looked at. The fits were performed using a protocol very similar to that of the SLy parameterizations Chabanat et al. (1997, 1998). The 36 actual sets of parameters can be found in the Physical Review archive EPA .

We use a formalism that explicitly relates the tensor terms in the energy functional to underlying effective density-dependent central, spin-orbit and tensor forces (or vertices) in the particle-hole channel. As has been known for long, a zero-range tensor force gives no qualitatively new terms for spherical mean-field states when combined with a central Skyrme force, but solely modifies the coupling constants of the ${\mathbf{J}}^{2}$ terms that are already present. The contribution from the central Skyrme force to the coupling constants of the ${\mathbf{J}}^{2}$ terms depends on the same parameters $t_{1}$, $x_{1}$, $t_{2}$ and $x_{2}$ that determine the effective mass and contribute to the surface terms. As the latter terms are much more important for the description of bulk properties than the ${\mathbf{J}}^{2}$ terms, the coupling constants of the ${\mathbf{J}}^{2}$ terms are confined to a very small region of the parameter space. From this point of view, adding a tensor force is necessary to explore it fully.

There is, however, the alternative interpretation of the Skyrme energy functional from the density matrix expansion, which in the absence of ab-initio realizations so far is used as a motivation to set up energy functionals with independent, and phenomenologically fitted, coupling constants of all terms not constrained by symmetries. In particular, this can be used to set unwanted or underconstrained terms to zero, as it is done for many existing parameterizations of the (central) Skyrme interaction. For the ground states of spherical nuclei, as discussed here, the frameworks cannot be distinguished. For deformed nuclei and, in particular, polarized nuclear matter, this choice will make a difference.

As a result of our study, we have obtained a long list of potential deficiencies of the Skyrme energy functional, most of which can be expected to be related to the properties of the central and spin-orbit interactions used. In fact, these deficiencies become more obvious the moment one adds a tensor force, as it appears that the presence of a tensor force unbalances a delicate compromise within various terms of the Skyrme interaction that permits to get the global trend of gross features of the shell structure right.

Our conclusions, however, have to be taken with a grain of salt. On the one hand, some might depend on the fit protocol; and on the other hand, we have to stress that (within the framework of our study – and all others available so far using mean-field methods) the comparison between calculated and empirical single-particle energies is not straightforward and without the risk of being misled. However, without even looking at single-particle spectra, we find that

The presence of the tensor terms leads to a strong rearrangement of the other coupling constants, most notably that of the spin-orbit force. In fact, we find that the variation of the spin-orbit strength $W_{0}$ provoked by the presence of tensor terms has a larger impact on the global systematics of single-particle spectra than the tensor terms themselves. The rearrangement of the parameters of the central and spin-orbit parts of the effective interaction suggests that perturbative studies of the tensor terms, in which they are added to an existing parameterization without readjustment, allow only very limited conclusions. 2. 2.

In the Skyrme energy functional, the combined coupling constants of the spin-orbit and tensor terms are nearly exclusively fixed by the mass differences between 40Ca, 48Ca and 56Ni. This correlation appears to be (at least partly) spurious, the rapidly varying spin-orbit and tensor terms being misused to simulate missing physics in the standard Skyrme functional. 3. 3.

The cost function $\chi^{2}$ used in our fit protocol prefers parameterizations with $\beta=0$, i.e. pure like-particle tensor terms $\sim({\mathbf{J}}^{2}_{\mathrm{n}}+{\mathbf{J}}^{2}_{\mathrm{p}})$, without giving a clear preference for a value of the corresponding coupling constant $\alpha$. By contrast, the mass residuals of 134 spherical even-even nuclei are minimized for interactions with large $\alpha$ and $\beta$. However, and as we will discuss in Bender et al. (2007), the deformation properties of many nuclei obtained with the latter parameterizations are unrealistic, which disfavors this region of the parameter space. 4. 4.

The difference of the charge radii of 40Ca and 48Ca turns out to be particularly sensitive to the absolute single-particle energy of the proton $1d_{3/2}$ level, which can be moved around by the ${\mathbf{J}}^{2}$ terms. As the parameterizations that give the best agreement for the absolute placement of this level do not necessarily give the best overall single-particle spectra for these two nuclei, this quantity should not be used to constrain the ${\mathbf{J}}^{2}$ terms.

Concerning the global properties of the spin-orbit current ${\mathbf{J}}$ and its contribution to the spin-orbit potential, we have shown that

The spin-orbit current ${\mathbf{J}}$ in non-spin-saturated doubly-magic nuclei as 56Ni, 100Sn, 132Sn or 208Pb is dominated by the nodeless intruder orbitals. Through the contribution of the tensor terms to the spin-orbit field, the feedback effect on their own spin-orbit splitting is maximized. 2. 2.

In light nuclei, ${\mathbf{J}}$ and consequently the contribution of the ${\mathbf{J}}^{2}$ terms to the binding energy and the spin-orbit potential, vary rapidly between near-zero and very large values when adding just a few nucleons to a given nucleus. In heavy spherical nuclei, the variation becomes much slower and smoother as on the one hand one does not encounter spin-saturated configurations anymore, and on the other hand there are more and more high-$\ell$ states with large degeneracy that require more nucleons to be filled. 3. 3.

The contribution from the zero-range spin-orbit force to the spin-orbit potential is peaked at the nuclear surface, as it is proportional to the gradient of the density. By contrast, the contribution from the zero-range tensor terms is peaked further inside of the nucleus, modifying the width of the spin-orbit potential with varying nucleon numbers. As shown in Ref. Brown et al. (2006), experimental data tend to dislike such a modification. 4. 4.

Large negative coupling constants of the tensor terms will lead to instabilities, where a nucleus gains energy separating the levels from many spin-orbit partners on both sides of the Fermi energy. This process leads to unphysical single-particle spectra and rules out a large part of the parameter space. In particular cases, one might even obtain a (probably spurious) coexistence of two spherical configurations with different shell structure in the same nucleus, which are separated by a barrier.

The main motivation to add ${\mathbf{J}}^{2}$ terms is of course to improve the single-particle spectra. All observations and conclusions concerning those have to be taken with care, as in this study we compare the eigenvalues of a spherical single-particle Hamiltonian with the separation energy to low-lying states in the odd-$A$ neighbors of doubly and semi-magic nuclei (as was done in all existing earlier studies). When looking at the single-particle spectra in doubly-magic nuclei (or semi-magic nuclei combined with a strong subshell closure of the other species) we find that

The relative error of the spin-orbit splittings depends strongly on the principal quantum number of the orbitals within a given shell, such that for parameterizations without the tensor terms the splittings of the intruder state (without nodes in the radial wave function) is tentatively too small, while it becomes too large with increasing number of nodes. Adding the tensor terms further increases the discrepancy. This problem can only be resolved by an improved control over the shape of the spin-orbit potential. Indeed, the size of the spin-orbit splittings is related to the overlap of the radial wave function of a given single-particle state with the spin-orbit potential. The tensor terms modify the width of the spin-orbit potential, but to cure this deficiency calls for a large negative like-particle tensor coupling constant $\alpha$, which is not consistent with the evolution of spin-orbit splittings along chains of semi-magic nuclei, and will lead to instabilities. 2. 2.

We also find that, in a given nucleus, the predicted spin-orbit splittings of neutron levels are larger than those of the protons when both are compared to experiment, which hints at an unresolved isospin trend in the spin-orbit interaction. 3. 3.

For spin-saturated doubly-magic nuclei as 16O and 40Ca, the spin-orbit splittings of the spin-saturated species of nucleons depends strongly on the coupling constants of the ${\mathbf{J}}^{2}$ terms, although they do not contribute to the spin-orbit field. This is a consequence of the strong correlation between the spin-orbit and tensor term coupling constants, which try to compensate each other in spin-unsaturated nuclei. For parameterizations with strong tensor-term coupling constants, the resulting spin-orbit force leads to unrealistic single-particle spectra of spin-saturated configurations. 4. 4.

The centroid of the spin-orbit partners that give the intruder state is tentatively too high compared to the major shell below.

The main effect of the tensor terms, that most of the recent studies concentrate on, is the evolution of spin-orbit splittings with $N$ and $Z$. Unfortunately, there are no data for the splittings themselves, such that one relies on data for the evolution of the distance of two levels with different $\ell$. The comparison is compromised by the global deficiencies of the single-particle spectra listed above. Still, a careful comparison of calculations and experiment suggests that

The evolution of the proton $1h_{11/2}$, $1g_{7/2}$ and $2d_{5/2}$ levels in the chain of Sn isotopes and that of the proton $1f_{5/2}$ and $2p_{3/2}$ levels in Ni isotopes call for a positive proton-neutron tensor coupling constant $\beta$ with a value around 120 MeV fm5, consistent with the findings of Refs. Brown et al. (2006); Colò et al. (2007); Brink and Stancu (2007). 2. 2.

The evolution of the neutron $1d_{3/2}$ and $2s_{1/2}$ levels between 40Ca and 48Ca calls for a like-particle tensor coupling constant $\alpha$ with a similar value around 120 MeV fm5. This it at variance to the findings of Refs. Brown et al. (2006); Colò et al. (2007); Brink and Stancu (2007), but in qualitative agreement with the parameterization skxta of Brown et al. Brown et al. (2006) for which the tensor terms were derived from a realistic interaction but disregarded thereafter because of its poor description of spin-orbit splittings. 3. 3.

Combined this leads to a dominantly isoscalar tensor term with a coupling constant $C^{J}_{0}$ around 120 MeV fm5, while the isovector coupling constant will have a small, near-zero, value.

Our study is obviously only a stepping stone towards improved parameterizations of the Skyrme energy density functional. There are a number of necessary further studies and future theoretical developments

The deformation properties of selected parameterizations T$IJ$ from this study will be discussed in a forthcoming paper Bender et al. (2007). 2. 2.

The influence of the terms depending on time-odd densities and currents in the complete energy functional (23) on nuclear matter and finite nuclei (rotational bands etc) is under investigation as well. The existing stability criteria of polarized matter have to be generalized as the tensor force introduces new unique terms, for example in the Landau parameters Haensel and Jerzak (1982). 3. 3.

It is well known that the strength of the spin-orbit force has to scale with the effective mass of an interaction, which in turn determines the average density of single-particle levels. All parameterizations discussed here have a similar effective mass close to $m^{*}_{0}/m=0.7$ that was already used for the SLy parameterizations. This value is somewhat smaller than the one obtained from ab-initio calculations. We have checked that increasing the effective isoscalar mass to the more realistic $m^{*}_{0}/m=0.8$ (which within our fit protocol requires to use two density dependent terms Lesinski et al. (2006)) does not significantly affect any of our conclusions. 4. 4.

It is evident that improvements of the central and spin-orbit parts of the energy density functional are necessary, which will require a generalization of its functional form. Other motivations were found recently to perform such a generalization Lesinski et al. (2006). 5. 5.

The only quantity that we found sufficiently sensitive to the tensor terms is the evolution of the distance between single-particle levels in isotopic or isotonic chains of semi-magic nuclei. The distance between the levels that can be used for such studies is so large, that it might be compromised by their coupling to collective excitations. Reliable calculations including pairing, polarization as well as particle-vibration coupling effects Bernard and Nguyen Van Giai (1980); Litvinova and Ring (2006) along isotopic and isotonic chains are needed to test the quality, reliability and limits of the simplistic identification of the eigenvalues of the spherical mean-field Hamiltonian in an even-even nucleus with the separation energy to or from low-lying states in the adjacent odd-$A$ nuclei.

Acknowledgments

We thank P. Bonche, H. Flocard, P.-H. Heenen and B. A. Brown for stimulating and encouraging discussions. Work by M. B. and K. B. was performed within the framework of the Espace de Structure Nucléaire Théorique (ESNT). T. L. acknowledges the hospitality of the SPhN and ESNT on many occasions during the realization of this work. This work was supported by the U.S. National Science Foundation under Grant No. PHY-0456903.

Appendix A Coupling constants of the Skyrme energy functional

The coupling constants of the central Skyrme energy density functional in terms of the parameters of the central Skyrme force are given by

[TABLE]

The coupling constants of the spin-orbit energy density functional in terms of the parameters of the spin-orbit force are given by

[TABLE]

The coupling constants of the tensor energy density functional in terms of the parameters of Skyrme’s tensor force are given by (Table I in Perlińska et al. (2004))

[TABLE]

Appendix B Phase transitions

The densities $\rho$ and $\tau$ entering the energy functional (28) vary smoothly with nucleon numbers as they follow the geometric growth of the nucleus. As a result, a functional depending only on $\rho$ and $\tau$ usually shows a unique minimum for given $N$, $Z$ and shape. The situation is quite different when the tensor terms are taken into account. Indeed, the amplitude of the spin-orbit current density ${\mathbf{J}}$ (27) depends on the number of spin-unsaturated single-particle states in the nucleus; it varies from (almost) zero in spin-saturated nuclei to large finite values as a consequence of shell and finite-size effects, see Fig. 5.

This behavior poses the risk of an instability, which was already reported in Beiner et al. (1975): multiplying ${\mathbf{J}}$ with a large coupling constant in the spin-orbit potential (35) might, for certain combinations of the signs of the coupling constant and the spin-orbit currents of protons and neutrons, increase the spin-orbit splittings. In some nuclei, this will cause two levels originating from different $\ell$ shells to approach the Fermi energy, one from above and the other from below, or even to cross. In that situation, their occupation numbers will change such that ${\mathbf{J}}$ increases further, which feeds back onto the spin-orbit potential and ultimately leads to a dramatic rearrangement of the single-particle spectrum.

We faced this problem when attempting to fit parameter sets with large negative $C_{0}^{J}$ and $C_{1}^{J}$. During the fit, some nuclei sometimes fell into the instability, depending on the values of the other coupling constants. As this is a highly nonlinear threshold effect that results in a very large energy gain from tiny modifications of the coupling constants, the corresponding fits did not, and could not, converge.

In special cases, one might even run into a situation with two coexisting minima, where as a function of a suitable coordinate the configuration with regular shell structure is separated from a configuration with unphysical large spin-orbit splittings by a barrier. In such a case, a calculation of the ground state might converge into one or the other minimum depending on the initial conditions chosen for the iterative solution of the HFB equations. In a calculation along an isotopic or isotonic chain, the coexistence will reveal itself through a large scattering of the mass residuals, which will fall on two distinct curves. We illustrate this phenomenon in Fig. 35 for 120Sn using a parameter set denoted “T$XX$” with $C_{0}^{J}=-157.57$ MeV fm5 and $C_{1}^{J}=-114.88$ MeV fm5, which is located outside the parameter space shown in Fig. 1, to its lower left. Among the various possible recipes for a constraint on the spin-orbit current density, we chose to minimize the following quantity:

[TABLE]

where $\mu$ is a Lagrange parameter and $C$ is a constant used to tune the constraint. The energy curve exhibits two minima denoted $(a)$ and $(b)$. The corresponding single-particle spectra are shown in Fig. 36 along with those obtained for SLy5. The minimum $(a)$ corresponds to an almost spin saturated neutron configuration where both spin partners are either occupied or empty,444 Note that the spin-saturation is a consequence of the wrong ordering of the $\nu$ $1h_{11/2}$, $2d_{3/2}$ and $1s_{1/2}$ levels compared to empirical data that is found for all our parameterizations, see Fig. 16, and in fact virtually all mean-field models Bender et al. (2003a).

which is very similar to what is found using SLy5. In the minimum $(b)$, which is deeper by more than 7 MeV, the single-particle spectrum is completely reorganized in order to maximize the spin-orbit current density and take advantage of its contribution in the functional. In this situation the neutron spin doublets $2d$, $1g$ and $1h$ split on both sides of the Fermi surface and generate a large spin-orbit current density.

This clearly shows that the parameter sets with large and negative coupling constants of the ${\mathbf{J}}^{2}$ terms must be discarded since for many nuclei they lead to ground states with unrealistic single-particle structure.

Note that this kind of instability does not appear for the spin-orbit term: although its contribution to the energy functional (28) also varies between small, sometimes near-zero, and very large values, see Figs. 28 and 29, it is only linear in ${\mathbf{J}}$. As a consequence, its contribution to the spin-orbit potential (35) lacks the feedback mechanism outlined above as it does not scale with ${\mathbf{J}}$. Still, its contribution to the total energy is usually much larger than that of the ${\mathbf{J}}^{2}$ terms, so it plays a decisive role for the absolute energy gained when varying ${\mathbf{J}}$.