EQUATION OF STATE FOR NUCLEONIC AND HYPERONIC NEUTRON STARS WITH MASS AND RADIUS CONSTRAINTS

Neutron stars are the most compact known objects without event horizons. They are formed in the aftermath of core-collapse supernovae and are usually observed as pulsars. Their features, such as the mass and radius, strongly depend on the properties of their dense interior. Thus, neutron stars serve as a unique laboratory for dense matter physics.

With more than 2000 pulsars known up to date, one of the best determined pulsar masses is that of the Hulse–Taylor pulsar of 1.4 M_⊙ (Hulse & Taylor 1975). Until very recently, the most precise measurements of neutron star masses clustered around this canonical value. Higher masses in neutron star binary systems have been measured in recent years with very high precision, using post-Keplerian parameters. This is the case of the binary milisecond pulsar PSR J1614-2230 of M = 1.97 ± 0.04 M_⊙ (Demorest et al. 2010) and the PSR J0348+0432 of M = 2.01 ± 0.04 M_⊙ (Antoniadis et al. 2013).

While the measurement of neutron star masses is accurate, the observational determination of their radii is more difficult and, as a consequence, comparably accurate values of radii do not yet exist. The radius of a neutron star can be extracted from the analysis of X-ray spectra emitted by the neutron star atmosphere. The modeling of the X-ray emission strongly depends on the distance to the source, its magnetic field and the composition of its atmosphere, thus making the determination of the radius a difficult task. As a result, different values for the stellar radii have been derived (Verbiest et al. 2008; Ozel et al. 2010, 2016; Suleimanov et al. 2011; Bogdanov 2013; Guillot et al. 2013; Guver & Ozel 2013; Lattimer & Lim 2013; Steiner et al. 2013; Guillot & Rutledge 2014; Heinke et al. 2014; Lattimer & Steiner 2014; Poutanen et al. 2014; Ozel & Psaltis 2015; Lattimer & Prakash 2016; Ozel & Freire 2016). In general, the extractions based on the spectral analysis of X-ray emission from quiescent X-ray transients in low-mass binaries (QLMXBs) favor small stellar radii in the 9–12 km range, whereas the determinations from neutron stars with recurring powerful bursts may lead to larger radii, of up to 16 km, although they are subject to larger uncertainties and controversy (see the discussion in the analysis of Fortin et al. 2015). The very recent work of Lattimer & Prakash (2016) indicates that the realistic range for radii of canonical neutron stars should be 10.7 km to 13.1 km. This analysis is based on observations of pulsar masses and estimates of symmetry properties derived from neutron matter studies and nuclear experiments. It is expected that robust observational upper bounds on stellar radii will be within reach in the near future. With space missions such as NICER (Neutron star Interior Composition ExploreR; Arzoumanian et al. 2014), high-precision X-ray astronomy will be able to offer precise measurements of masses and radii (Watts et al. 2016), while gravitational-wave signals from neutron star mergers hold promise to determine neutron star radii with a precision of 1 km (Bauswein & Janka 2012; Lackey & Wade 2015).

In anticipation that these upcoming astrophysical determinations could confirm small neutron star sizes, it is important and timely to explore the smallest radii that can be delivered by the theoretical models of compressed matter that are able to fulfil the 2 M_⊙ maximum mass constraint, while reproducing at the same time the phenomenology of atomic nuclei. The masses and radii of neutron stars are linked to the physics of their interior, that is, the equation of state (EoS) of dense matter (Lattimer & Prakash 2004, 2007; Oertel et al. 2016). Many of the current nuclear models for the EoS are able to satisfy the 2 M_⊙ constraint required by the discovery of massive neutron stars (Demorest et al. 2010; Antoniadis et al. 2013). However, the possible existence of neutron stars with small radii suggested by recent astrophysical analyses (Guillot et al. 2013; Guver & Ozel 2013; Guillot & Rutledge 2014; Heinke et al. 2014; Lattimer & Steiner 2014; Lattimer & Prakash 2016; Ozel & Freire 2016; Ozel et al. 2016) poses a difficult challenge to most of the nuclear models (Chen & Piekarewicz 2015a; Dexheimer et al. 2015; Jiang et al. 2015; Ozel & Freire 2016).

A small neutron star radius for a canonical neutron star requires a certain softening of the pressure of neutron matter, and hence of the nuclear symmetry energy, around 1–2 times saturation density n₀ (n₀ ≈ 0.16 fm⁻³) (Lattimer & Prakash 2007; Tsang et al. 2012; Ozel & Freire 2016). The star radius could also be reduced by decreasing the pressure of the isospin-symmetric part of the EoS in the intermediate-density region, but this is only possible with severe limitations due to the saturation properties of nuclear matter and the constraints on the EoS of dense nuclear matter extracted from nuclear collective flow (Danielewicz et al. 2002) and kaon production (Fuchs et al. 2001; Lynch et al. 2009) in high-energy heavy-ion collisions (HICs). Moreover, the requirement of maximum masses of 2 M_⊙ does not allow a significant reduction of the total pressure. Indeed, very few models seem to exist that can meet both constraints (small radius and large mass) simultaneously, and fewer such models can in addition render accurate descriptions of the finite nuclei properties (Horowitz & Piekarewicz 2001a, 2001b; Chen & Piekarewicz 2015a; Jiang et al. 2015; Sharma et al. 2015).

It has also been long known that the transition from nuclear to hyperonic matter is energetically favored as the density increases inside neutron stars (Ambartsumyan & Saakyan 1960). The opening of hyperon degrees of freedom leads to a considerable softening of the EoS (Glendenning 1982). As a consequence, the maximum neutron star masses obtained are usually smaller than the 2 M_⊙ observations (Demorest et al. 2010; Antoniadis et al. 2013). The solution of this so-called "hyperon puzzle" is not easy, and requires a mechanism that could provide additional repulsion to make the EoS stiffer. Possible mechanisms could be: (1) stiffer hyperon–nucleon and/or hyperon–hyperon interactions (see the recent works by Bednarek et al. 2012; Weissenborn et al. 2012; Maslov et al. 2015; Oertel et al. 2015); (2) the inclusion of three-body forces with one or more hyperons (see Vidana et al. 2011; Yamamoto et al. 2014; Lonardoni et al. 2015 for recent studies); (3) the appearance of other hadronic degrees of freedom such as the Δ isobar (Drago et al. 2014) or meson condensates that push the onset of hyperons to higher densities; and (4) the appearance of a phase transition to deconfined quark matter at densities below the hyperon threshold (Alford et al. 2007; Klahn et al. 2013; Zdunik & Haensel 2013). For a detailed review on the "hyperon puzzle," we refer the reader to Chatterjee & Vidana (2016) and references therein.

The presence of strong magnetic fields inside neutron stars is another possible source for a stiffer EoS that could sustain masses of 2 M_⊙. Anomalous X-ray pulsars and soft γ-ray repeaters are identified with highly magnetized neutron stars with a surface magnetic field of ∼10¹⁴–10¹⁵ G (Vasisht & Gotthelf 1997; Kouveliotou et al. 1998; Woods et al. 1999). This class of compact objects has been named "magnetars," i.e., neutron stars with magnetic fields several orders of magnitude larger than the canonical surface dipole magnetic fields B ∼ 10¹²–10¹³ G of the bulk of the pulsar population (Mereghetti 2008; Rea & Esposito 2011; Turolla et al. 2015). It has been shown that magnetic fields larger than $B/{B}_{c}^{e}={10}^{5}$, with ${B}_{c}^{e}=4.414\times {10}^{13}$ G being the critical magnetic field at which the electron cyclotron energy is equal to the electron mass, will affect the EoS of dense nuclear matter (Chakrabarty et al. 1997; Bandyopadhyay et al. 1998; Broderick et al. 2000; Suh & Mathews 2001; Harding & Lai 2006; Chen et al. 2007; Rabhi et al. 2008; Dexheimer et al. 2012; Strickland et al. 2012). The study of the effects upon the EoS of hyperonic matter of very strong magnetic fields (B ∼ 10¹⁸–10¹⁹ G in the star center) was initiated by Broderick et al. (2002) and has been recently addressed (Rabhi & Providencia 2010; Lopes & Menezes 2012; Sinha et al. 2013; Gomes et al. 2014).

In the present paper we reconcile the 2 M_⊙ mass observations with the recent analyses of radii below 13 km for neutron stars, while fulfilling the constraints from the properties of nuclear matter, nuclei, and HICs at high energy. This is accomplished for neutron stars with nucleonic and hyperonic cores. The formalism is based on the covariant field-theoretical approach to hadronic matter (see for example Serot & Walecka 1986, 1997, and chapter 4 of Glendenning 2000, and references therein). The nucleonic EoS is obtained as a new parameterization of the nonlinear realization of the relativistic mean-field (RMF) model (Serot & Walecka 1986, 1997; Glendenning 2000; Chen & Piekarewicz 2014). Starting from the recent RMF parameter set FSU2 (Chen & Piekarewicz 2014), we find that by softening the pressure of neutron star matter in the neighborhood of saturation one can accommodate smaller stellar radii, while the properties of nuclear matter and finite nuclei are still fulfilled. Moreover, we are able to keep the pressure at high densities in agreement with HIC data and sufficiently stiff such that it can sustain neutron stars of ≈2 M_⊙. We denote the new parameterization by FSU2R. Next we introduce hyperons in our calculation and fit the hyperon couplings to the value of the hyperon–nucleon and hyperon–hyperon optical potentials extracted from the available data on hypernuclei. Whereas the radius of the neutron stars is insensitive to the appearance of the hyperons, we find a reduction of the maximum mass of the neutron stars due to the expected softening of the EoS. However, we find that the 2 M_⊙ constraint is still fulfilled when hyperons are considered by means of a slight modification of the parameters of the model, denoted as FSU2H, compatible with astrophysical observations and empirical data. We also analyze the effect of strong magnetic fields in the mass and radius of neutron stars. The origin of the intense magnetic fields in magnetars is still open to debate and the strength of the inner values is still unknown (Thompson & Duncan 1993; Ardeljan et al. 2005; Vink & Kuiper 2006). Nevertheless, it is worth exploring the modification on the EoS and on the neutron star properties induced by magnetic fields that are as large as the upper limit imposed by the scalar virial theorem (Chandrasekhar & Fermi 1953; Shapiro & Teukolsky 1983), which is of the order of $B\sim 2\times {10}^{8}{(M/{M}_{\odot })(R/{R}_{\odot })}^{-2}$. For a star of R ∼ 10 km and M ∼ 2 M_⊙ the magnetic field could then reach around 2 × 10¹⁸ G. In our study we have magnetic fields close to this value only at the very center of the star and assume a magnetic field profile toward a value of 10¹⁵ G at the surface, hence fulfilling the stability constraint. From the calculations with the proposed EoS we conclude that nucleonic and hyperonic magnetars with a surface magnetic field of ∼10¹⁵ G and with magnetic fields values of ∼10¹⁸ G in the interior can reach maximum masses of 2 M_⊙ with radii in the 12–13 km range.

The paper is organized as follows. In Section 2 we present the RMF model and the inclusion of magnetic fields for the determination of the EoS in beta-equilibrated matter. In Section 3 we show how we calibrate the nucleonic model FSU2R by fulfilling the constraints of 2 M_⊙ mass observations and small neutron star radii, as well as the properties of nuclear matter, nuclei, and HICs at high energy. Then, in Section 4 we introduce hyperons and magnetic fields and provide a slightly changed parameterization, FSU2H, that also fulfils the observational and experimental requirements while allowing for maximum masses of 2 M_⊙. We finally summarize our results in Section 5.

Our starting point is the RMF model of matter, where baryons interact through the exchange of mesons and which provides a covariant description of the EoS and nuclear systems. The Lagrangian density of the theory can be written as (Serot & Walecka 1986, 1997; Glendenning 2000; Chen & Piekarewicz 2014)

$$\begin{eqnarray}&&{ \mathcal L }=\displaystyle \sum _{b}{{ \mathcal L }}_{b}+{{ \mathcal L }}_{m}+\displaystyle \sum _{l=e,\mu }{{ \mathcal L }}_{l},\end{eqnarray} \tag{ 1 }$$

with the baryon (b), lepton (l = e, μ), and meson (m = σ, ω, ρ, and ϕ) Lagrangians given by

$$\begin{eqnarray}{{ \mathcal L }}_{b} & = & {\bar{{\rm{\Psi }}}}_{b}(i{\gamma }_{\mu }{\partial }^{\mu }-{q}_{b}{\gamma }_{\mu }{A}^{\mu }-{m}_{b}\\ & & +{g}_{\sigma b}\sigma -{g}_{\omega b}{\gamma }_{\mu }{\omega }^{\mu }-{g}_{\phi b}{\gamma }_{\mu }{\phi }^{\mu }-{g}_{\rho b}{\gamma }_{\mu }{{\boldsymbol{I}}}_{b}\,{{\boldsymbol{\rho }}}^{\mu }){{\rm{\Psi }}}_{b},\\ {{ \mathcal L }}_{l} & = & {\bar{\psi }}_{l}(i{\gamma }_{\mu }{\partial }^{\mu }-{q}_{l}{\gamma }_{\mu }{A}^{\mu }-{m}_{l}){\psi }_{l},\\ {{ \mathcal L }}_{m} & = & \displaystyle \frac{1}{2}{\partial }_{\mu }\sigma {\partial }^{\mu }\sigma -\displaystyle \frac{1}{2}{m}_{\sigma }^{2}{\sigma }^{2}-\displaystyle \frac{\kappa }{3!}{({g}_{\sigma N}\sigma )}^{3}-\displaystyle \frac{\lambda }{4!}{({g}_{\sigma N}\sigma )}^{4}\\ & & -\displaystyle \frac{1}{4}{{\rm{\Omega }}}^{\mu \nu }{{\rm{\Omega }}}_{\mu \nu }+\displaystyle \frac{1}{2}{m}_{\omega }^{2}{\omega }_{\mu }{\omega }^{\mu }+\displaystyle \frac{\zeta }{4!}{({g}_{\omega N}{\omega }_{\mu }{\omega }^{\mu })}^{4}\\ & & -\displaystyle \frac{1}{4}{{\boldsymbol{R}}}^{\mu \nu }{{\boldsymbol{R}}}_{\mu \nu }+\displaystyle \frac{1}{2}{m}_{\rho }^{2}{{\boldsymbol{\rho }}}_{\mu }{{\boldsymbol{\rho }}}^{\mu }+{{\rm{\Lambda }}}_{\omega }{g}_{\rho N}^{2}{{\boldsymbol{\rho }}}_{\mu }{{\boldsymbol{\rho }}}^{\mu }{\,g}_{\omega N}^{2}{\omega }_{\mu }{\omega }^{\mu }\\ & & -\displaystyle \frac{1}{4}{P}^{\mu \nu }{P}_{\mu \nu }+\displaystyle \frac{1}{2}{m}_{\phi }^{2}{\phi }_{\mu }{\phi }^{\mu }-\displaystyle \frac{1}{4}{F}^{\mu \nu }{F}_{\mu \nu },\end{eqnarray} \tag{ 2 }$$

where Ψ_b and ψ_l are the baryon and lepton Dirac fields, respectively. The mesonic and electromagnetic field strength tensors are ${{\rm{\Omega }}}_{\mu \nu }={\partial }_{\mu }{\omega }_{\nu }-{\partial }_{\nu }{\omega }_{\mu }$, ${{\boldsymbol{R}}}_{\mu \nu }={\partial }_{\mu }{{\boldsymbol{\rho }}}_{\nu }-{\partial }_{\nu }{{\boldsymbol{\rho }}}_{\mu }$, ${P}_{\mu \nu }={\partial }_{\mu }{\phi }_{\nu }-{\partial }_{\nu }{\phi }_{\mu }$, and ${F}_{\mu \nu }={\partial }_{\mu }{A}_{\nu }-{\partial }_{\nu }{A}_{\mu }$. The electromagnetic field is assumed to be externally generated and, as we will discuss below, we do not consider the coupling of the particles to the electromagnetic field tensor via the baryon anomalous magnetic moments. The strong interaction couplings of a meson to a certain baryon are denoted by g (with N indicating nucleon), the electromagnetic couplings by q and the baryon, meson, and lepton masses by m. The vector ${{\boldsymbol{I}}}_{b}$ stands for the isospin operator.

The Lagrangian density (2) incorporates scalar and vector meson self-interactions as well as a mixed quartic vector meson interaction. The nonlinear meson interactions are important for a quantitative description of nuclear matter and finite nuclei, as they lead to additional density dependence that represents in an effective way the medium dependence induced by many-body correlations. The scalar self-interactions with coupling constants κ and λ, introduced by Boguta & Bodmer (1977), are responsible for softening the EoS of symmetric nuclear matter (SNM) around saturation density and allow one to obtain a realistic value for the compression modulus of nuclear matter (Boguta & Bodmer 1977; Boguta & Stoecker 1983). The quartic isoscalar–vector self-interaction (with coupling ζ) softens the EoS at high densities (Mueller & Serot 1996), while the mixed quartic isovector–vector interaction (with coupling Λ_ω) is introduced (Horowitz & Piekarewicz 2001a, 2001b) to modify the density dependence of the nuclear symmetry energy, which measures the energy cost involved in changing the protons into neutrons in nuclear matter.

The Dirac equations for baryons and leptons are given by

$$\begin{eqnarray}&&(i{\gamma }_{\mu }\,{\partial }^{\mu }-{q}_{b}\,{\gamma }_{\mu }\,{A}^{\mu }-{m}_{b}^{* }\\ &&\quad -{g}_{\omega b}\,{\gamma }_{0}\,{\omega }^{0}-{g}_{\phi b}\,{\gamma }_{0}\,{\phi }^{0}-{g}_{\rho b}\,{I}_{3b}\,{\gamma }_{0}\,{\rho }_{3}^{0}){{\rm{\Psi }}}_{b}=0,\end{eqnarray} \tag{ 3 }$$

$$\begin{eqnarray}&&(i{\gamma }_{\mu }\,{\partial }^{\mu }-{q}_{l}\,{\gamma }_{\mu }\,{A}^{\mu }-{m}_{l}){\psi }_{l}=0,\end{eqnarray} \tag{ 4 }$$

where the effective baryon masses are defined as

$$\begin{eqnarray}&&{m}_{b}^{* }={m}_{b}-{g}_{\sigma b}\sigma .\end{eqnarray} \tag{ 5 }$$

The field equations of motion follow from the Euler–Lagrange equations. In the mean-field approximation, the meson fields are replaced by their respective mean-field expectation values, which are given in uniform matter as $\bar{\sigma }=\langle \sigma \rangle$, $\bar{\rho }=\langle {\rho }_{3}^{0}\rangle$, $\bar{\omega }=\langle {\omega }^{0}\rangle$, and $\bar{\phi }=\langle {\phi }^{0}\rangle$. Thus, the equations of motion for the meson fields in the mean-field approximation for the uniform medium are

$$\begin{eqnarray}&&{m}_{\sigma }^{2}\,\bar{\sigma }+\displaystyle \frac{\kappa }{2}{g}_{\sigma N}^{3}{\bar{\sigma }}^{2}+\displaystyle \frac{\lambda }{3!}{g}_{\sigma N}^{4}{\bar{\sigma }}^{3}=\displaystyle \sum _{b}{g}_{\sigma b}{n}_{b}^{s},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&{m}_{\omega }^{2}\,\bar{\omega }+\displaystyle \frac{\zeta }{3!}{g}_{\omega N}^{4}{\bar{\omega }}^{3}+2{{\rm{\Lambda }}}_{\omega }{g}_{\rho N\,}^{2}{g}_{\omega N}^{2}{\bar{\rho }}^{2}\bar{\omega }=\displaystyle \sum _{b}{g}_{\omega b}{n}_{b},\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{m}_{\rho }^{2}\,\bar{\rho }+2{{\rm{\Lambda }}}_{\omega }{g}_{\rho N}^{2}{g}_{\omega N}^{2}{\bar{\omega }}^{2}\bar{\rho }=\displaystyle \sum _{b}{g}_{\rho b}{I}_{3b}{n}_{b},\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&{m}_{\phi }^{2}\bar{\phi }=\displaystyle \sum _{b}{g}_{\phi b}{n}_{b},\end{eqnarray} \tag{ 9 }$$

where ${I}_{3b}$ represents the third component of isospin of baryon b, with the convention ${I}_{3p}=1/2$. The quantities ${n}_{b}^{s}=\langle {\bar{{\rm{\Psi }}}}_{b}{{\rm{\Psi }}}_{b}\rangle$ and ${n}_{b}=\langle {\bar{{\rm{\Psi }}}}_{b}{\gamma }^{0}{{\rm{\Psi }}}_{b}\rangle$ are the scalar and baryon density for a given baryon, respectively.

In the presence of a magnetic field, the single-particle energy of the charged baryons and leptons is quantized in the perpendicular direction to the magnetic field. Taking the magnetic field in the z-direction, ${\boldsymbol{B}}=B\hat{z}$, the single particle energies of all baryons and leptons are given by (Broderick et al. 2000)

$$\begin{eqnarray}{E}_{\nu }^{{cb}} & = & \sqrt{{k}_{z}^{2}+{m}_{{cb}}^{* 2}+2\nu | {q}_{{cb}}| B}\\ & & +{g}_{\omega {cb}}\,\bar{\omega }+{g}_{\rho {cb}}\,{I}_{3b}\,\bar{\rho }+{g}_{\phi {cb}}\,\bar{\phi },\end{eqnarray} \tag{ 10 }$$

$$\begin{eqnarray}{E}^{{ub}} & = & \sqrt{{k}^{2}+{m}_{{ub}}^{* 2}}\\ & & +{g}_{\omega {ub}}\,\bar{\omega }+{g}_{\rho {ub}}\,{I}_{3b}\,\bar{\rho }+{g}_{\phi {ub}}\,\bar{\phi },\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{E}_{\nu }^{l}=\sqrt{{k}_{z}^{2}+{m}_{l}^{2}+2\nu | {q}_{l}| B},\end{eqnarray} \tag{ 12 }$$

with cb denoting charged baryons and ub uncharged baryons. The quantity $\nu =\left(n+\tfrac{1}{2}-\tfrac{1}{2}\tfrac{q}{| q| }{\sigma }_{z}\right)=0,1,2,\ldots$, with n being the principal quantum number and σ_z the Pauli matrix, indicates the Landau levels of the fermions with electric charge q.

As mentioned above, we have omitted the coupling of the baryons to the electromagnetic field tensor via their anomalous magnetic moments. The interaction of the baryon anomalous magnetic moments with the field strength has been found to partly compensate for the effects on the EoS associated with Landau quantization (Broderick et al. 2000). However, to see some appreciable changes in the EoS and the neutron star composition, intense fields of the order of 5 × 10¹⁸ G are needed. Moreover, those effects are mostly concentrated at low densities (≲2n₀) for such a field strength (Broderick et al. 2000; Rabhi et al. 2008). Therefore, neglecting the effects associated to the anomalous magnetic moments is a reasonable approximation in the present work since we consider neutron stars with magnetic fields at the core of at most 2 × 10¹⁸ G and magnetic field profiles that do not reach 5 × 10¹⁷ G in the region ≲2n₀.

The Fermi momenta of the charged baryons, ${k}_{F,\nu }^{{cb}}$, uncharged baryons, ${k}_{F}^{{ub}}$, and leptons, ${k}_{F,\nu }^{l}$, are related to the Fermi energies ${E}_{F}^{{cb}}$, ${E}_{F}^{{ub}}$ and ${E}_{F}^{l}$ as

$$\begin{eqnarray}{k}_{F,\nu }^{{cb}} & = & \sqrt{{({E}_{F}^{{cb}})}^{2}-({m}_{{cb}}^{* 2}+2\nu | {q}_{{cb}}| B)},\\ {k}_{F}^{{ub}} & = & \sqrt{{({E}_{F}^{{ub}})}^{2}-{m}_{{ub}}^{* 2}},\\ {k}_{F,\nu }^{l} & = & \sqrt{{({E}_{F}^{l})}^{2}-({m}_{l}^{2}+2\nu | {q}_{l}| B)},\end{eqnarray} \tag{ 13 }$$

while the chemical potentials of baryons and leptons are defined as

$$\begin{eqnarray}&&{\mu }_{b}={E}_{F}^{b}+{g}_{\omega b}\,\bar{\omega }+{g}_{\rho b}\,{I}_{3b}\,\bar{\rho }+{g}_{\phi b}\,\bar{\phi },\end{eqnarray} \tag{ 14 }$$

$$\begin{eqnarray}&&{\mu }_{l}={E}_{F}^{l}.\end{eqnarray} \tag{ 15 }$$

The largest value of ν is obtained by imposing that the square of the Fermi momentum of the particle is still positive, i.e., by taking the closest integer from below defined by the ratio

$$\begin{eqnarray*}&&{\nu }_{\max }=\left[\displaystyle \frac{{({E}_{F}^{l})}^{2}-{m}_{l}^{2}}{2| {q}_{l}| \,B}\right],\quad {\rm{leptons}}\end{eqnarray*}$$

$$\begin{eqnarray*}&&{\nu }_{\max }=\left[\displaystyle \frac{{({E}_{F}^{{cb}})}^{2}-{{m}_{{cb}}^{* }}^{2}}{2| {q}_{{cb}}| \,B}\right],\quad {\rm{charged}}\,{\rm{baryons}}.\end{eqnarray*}$$

With all these ingredients, the scalar and vector densities for baryons and leptons are given by (Broderick et al. 2000)

$$\begin{eqnarray}{n}_{{cb}}^{s} & = & \displaystyle \frac{| {q}_{{cb}}| {{Bm}}_{{cb}}^{* }}{2{\pi }^{2}}\displaystyle \sum _{\nu =0}^{{\nu }_{\max }}{r}_{\nu }\ \mathrm{ln}\displaystyle \frac{{k}_{F,\nu }^{{cb}}+{E}_{F}^{{cb}}}{\sqrt{{m}_{{cb}}^{* 2}+2\nu | {q}_{{cb}}| B}},\\ {n}_{{ub}}^{s} & = & \displaystyle \frac{{m}_{{ub}}^{* }}{2{\pi }^{2}}\left[{E}_{F}^{{ub}}{k}_{F}^{{ub}}-{m}_{{ub}}^{* 2}\mathrm{ln}\displaystyle \frac{{k}_{F}^{{ub}}+{E}_{F}^{{ub}}}{{m}_{{ub}}^{* }}\right],\\ {n}_{{cb}} & = & \displaystyle \frac{| {q}_{{cb}}| B}{2{\pi }^{2}}\displaystyle \sum _{\nu =0}^{{\nu }_{\max }}{r}_{\nu }\ {k}_{F,\nu }^{{cb}},\\ {n}_{{ub}} & = & \displaystyle \frac{{({k}_{F}^{{ub}})}^{3}}{3{\pi }^{2}},\\ {n}_{l} & = & \displaystyle \frac{| {q}_{l}| B}{2{\pi }^{2}}\displaystyle \sum _{\nu =0}^{{\nu }_{\max }}{r}_{\nu }\ {k}_{F,\nu }^{l},\end{eqnarray} \tag{ 16 }$$

where r_ν is the degeneracy of the ν Landau level, which is 1 for ν = 0 and 2 for $\nu \ne 0$.

We can now obtain the energy density ε and pressure P of the system. The energy density of matter, ε_matt, is given by

$$\begin{eqnarray}{\varepsilon }_{\mathrm{matt}} & = & \displaystyle \sum _{b}{\varepsilon }_{b}+\displaystyle \sum _{l}{\varepsilon }_{l}\\ & & +\displaystyle \frac{1}{2}{m}_{\sigma }^{2}{\bar{\sigma }}^{2}+\displaystyle \frac{1}{2}{m}_{\omega }^{2}{\bar{\omega }}^{2}+\displaystyle \frac{1}{2}{m}_{\rho }^{2}{\bar{\rho }}^{2}+\displaystyle \frac{1}{2}{m}_{\phi }^{2}{\bar{\phi }}^{2}\\ & & +\displaystyle \frac{\kappa }{3!}{({g}_{\sigma }\bar{\sigma })}^{3}+\displaystyle \frac{\lambda }{4!}{({g}_{\sigma }\bar{\sigma })}^{4}\\ & & +\displaystyle \frac{\zeta }{8}{({g}_{\omega }\bar{\omega })}^{4}+3{{\rm{\Lambda }}}_{\omega }{({g}_{\rho }{g}_{\omega }\bar{\rho }\bar{\omega })}^{2},\end{eqnarray} \tag{ 17 }$$

where the energy densities of baryons and leptons have the following expressions

$$\begin{eqnarray}&&{\varepsilon }_{{cb}}=\displaystyle \frac{| {q}_{{cb}}| B}{4{\pi }^{2}}\displaystyle \sum _{\nu =0}^{{\nu }_{\max }}{r}_{\nu }\\ &&\quad \times \left[{k}_{F,\nu }^{{cb}}{E}_{F}^{{cb}}+({m}_{{cb}}^{* 2}+2\nu | {q}_{{cb}}| B)\mathrm{ln}\displaystyle \frac{{k}_{F,\nu }^{{cb}}+{E}_{F}^{{cb}}}{\sqrt{{m}_{{cb}}^{* 2}+2\nu | {q}_{{cb}}| B}}\right],\\ &&{\varepsilon }_{{ub}}=\displaystyle \frac{1}{8{\pi }^{2}}\left[{k}_{F}^{{ub}}{({E}_{F}^{{ub}})}^{3}+{({k}_{F}^{{ub}})}^{3}{E}_{F}^{{ub}}-{m}_{{ub}}^{* 4}\mathrm{ln}\displaystyle \frac{{k}_{F}^{{ub}}+{E}_{F}^{{ub}}}{{m}_{{ub}}^{* }}\right],\\ &&{\varepsilon }_{l}=\displaystyle \frac{| {q}_{l}| B}{4{\pi }^{2}}\displaystyle \sum _{\nu =0}^{{\nu }_{\max }}{r}_{\nu }\\ &&\quad \times \left[{k}_{F,\nu }^{l}{E}_{F}^{l}+({m}_{l}^{2}+2\nu | {q}_{l}| B)\mathrm{ln}\displaystyle \frac{{k}_{F,\nu }^{l}+{E}_{F}^{l}}{\sqrt{{m}_{l}^{2}+2\nu | {q}_{l}| B}}\right].\end{eqnarray} \tag{ 18 }$$

The pressure of matter, P_matt, is obtained using the thermodynamic relation

$$\begin{eqnarray}&&{P}_{\mathrm{matt}}=\displaystyle \sum _{i}{\mu }_{i}{n}_{i}-{\varepsilon }_{\mathrm{matt}}.\end{eqnarray} \tag{ 19 }$$

While the contribution from the electromagnetic field to the energy density is ${B}^{2}/8\pi$, we use the so-called "chaotic field" prescription for the calculation of the pressure of the system (Menezes & Lopes 2016), so that we have

$$\begin{eqnarray}&&\varepsilon ={\varepsilon }_{\mathrm{matt}}+\displaystyle \frac{{B}^{2}}{8\pi },\end{eqnarray} \tag{ 20 }$$

$$\begin{eqnarray}&&P={P}_{\mathrm{matt}}+\displaystyle \frac{{B}^{2}}{24\pi }.\end{eqnarray} \tag{ 21 }$$

2.1. Neutron Star Matter in β-equilibrium

In order to determine the structure of neutron stars one needs to obtain the EoS over a wide range of densities. For the inner and outer crust of the star we employ the EoS of Sharma et al. (2015), which has been obtained from microscopic calculations. In the core of neutron stars, we find β-equilibrated globally neutral, charged matter. Consequently, the chemical potentials, μ_i, and particle densities, n_i, satisfy the conditions

$$\begin{eqnarray}{\mu }_{i} & = & {b}_{i}{\mu }_{n}-{q}_{i}{\mu }_{e},\\ 0 & = & \displaystyle \sum _{{cb},l}{q}_{i}\,{n}_{i},\\ n & = & \displaystyle \sum _{{cb},{ub}}{n}_{i},\end{eqnarray} \tag{ 22 }$$

with b_i the baryon number and q_i the charge of the particle i. These relations together with Equations (3), (4) and the field Equations (6)–(9) for σ, ω, ρ and ϕ have to be solved self-consistently for total baryon density n in the presence of a magnetic field. In this way, we obtain the chemical potential and the corresponding density of each species for a given n, so that we can determine the energy density and pressure of the neutron star matter at each density.

Once the EoS is known, the mass M and the corresponding radius R of the neutron star are obtained from solving the Tolman–Oppenheimer–Volkoff (TOV) equations (Oppenheimer & Volkoff 1939)

$$\begin{eqnarray}\displaystyle \frac{{dP}(r)}{{dr}} & = & -\,\displaystyle \frac{G}{{r}^{2}}[\varepsilon (r)+P(r)]\\ & & \times [M(r)+4\pi {r}^{3}P(r)]{\left[1-\displaystyle \frac{2{GM}(r)}{r}\right]}^{-1},\\ \displaystyle \frac{{dM}(r)}{{dr}} & = & 4\pi \varepsilon (r){r}^{2},\end{eqnarray} \tag{ 23 }$$

where r is the radial coordinate. To solve these equations one needs to specify the initial conditions, namely the enclosed mass and the pressure at the center of the star, M(r = 0) = 0 and P(r = 0) = P_c, while the energy density is taken from the assumed EoS. The integration of the TOV equations over the radial coordinate ends when P(r = R) = 0.

3.1. Equation of State, Stellar Masses, and Stellar Radii

We start our analysis by defining the baseline model for nuclear matter to compute masses and radii of neutron stars. Nuclear models that perform similarly well in the description of finite nuclei often extrapolate very differently at high densities, as usually no information on the high-density sector of the EoS has been incorporated in the fitting of the model. In this work we are interested in a model that gives neutron star radii as small as possible and massive enough neutron stars, in order to reconcile in a unified formalism the new astrophysical indications of small stellar radii and the existence of stars of 2 M_⊙ masses, while still meeting the constraints from the nuclear data of terrestrial laboratories.

For the nuclear model we start from the Lagrangian density of Equations (1), (2) by only considering nucleons and mesons. As mentioned in Section 2, the ζ self-coupling of the ω meson (see Equation (2)) is efficient in softening the EoS at supranormal densities while the Λ_ω cross-coupling of the ω and ρ mesons (Equation (2)) regulates the density dependence of the symmetry energy. In order to show the effect of these nonlinear contributions to the EoS, in Figure 1 we plot for selected interactions the pressures of SNM in the upper panel and of pure neutron matter (PNM) in the lower panel. The two shaded areas in the SNM panel depict the regions that are compatible with the data on collective flow (Danielewicz et al. 2002) (gray area) and on kaon production (Fuchs et al. 2001; Lynch et al. 2009) (turquoise region) according to the modeling of energetic HICs. The shaded areas in the PNM panel correspond to the constraints from the flow data supplemented by a symmetry energy with weak (gray area) or strong (brown area) density dependence (Danielewicz et al. 2002).

Zoom In Zoom Out Reset image size

We first consider the well-known parameter sets NL3 (Lalazissis et al. 1997) and FSU (also called FSUGold) (Todd-Rutel & Piekarewicz 2005). NL3 has ζ = Λ_ω = 0 while FSU has ζ = 0.06 and Λ_ω = 0.03 (the full set of parameters of the models can be found in Table 1). Both NL3 and FSU reproduce quite well a variety of properties of atomic nuclei. However, they render two EoSs in SNM with different behavior at supranormal densities due to the different ζ value (we recall that the Λ_ω coupling does not contribute in SNM). We can see in Figure 1(upper panel) that above density n ∼ 1.5–2n₀ the FSU model with ζ = 0.06 (dot-dashed blue line) yields a much softer SNM pressure than the NL3 model with ζ = 0 (dotted magenta line). In PNM, the isovector coupling Λ_ω tunes the change with density of the EoS, as it softens the symmetry energy. Indeed, if we compare the same models FSU and NL3 in PNM (see Figure 1(lower panel)), FSU (Λ_ω = 0.03) has its pressure strongly further reduced with respect to NL3 (Λ_ω = 0) in the density window from around saturation density n₀ up to n ∼ 1.5–2n₀. For densities above 2n₀ the softening effect from Λ_ω is less prominent and the PNM pressures of FSU and NL3 show comparable differences with the case of SNM. We therefore note, consistently with the systematics in earlier works (Horowitz & Piekarewicz 2001a, 2001b; Carriere et al. 2003; Chen & Piekarewicz 2014), that the Λ_ω and ζ couplings have a complementary impact on the EoS by each one influencing almost separate density sectors. This will be important for our goals for stellar radii and masses as we shall see below.

Table 1.  Parameters for the Models NL3 (Lalazissis et al. 1997), FSU (Todd-Rutel & Piekarewicz 2005), FSU2 (Chen & Piekarewicz 2014), FSU2R (This Work), and FSU2H (This Work, Section 4)

Models	NL3	FSU	FSU2	FSU2R	FSU2H
mσ [MeV]	508.194	491.500	497.479	497.479	497.479
mω [MeV]	782.501	782.500	782.500	782.500	782.500
mρ [MeV]	763.000	763.000	763.000	763.000	763.000
${g}_{\sigma N}^{2}$	104.3871	112.1996	108.0943	107.5751	102.7200
${g}_{\omega N}^{2}$	165.5854	204.5469	183.7893	182.3949	169.5315
${g}_{\rho N}^{2}$	79.6000	138.4701	80.4656	247.3409	247.3409
κ	3.8599	1.4203	3.0029	3.0911	4.0014
λ	−0.015905	0.023762	−0.000533	−0.001680	−0.013298
ζ	0.00	0.06	0.0256	0.024	0.008
Λω	0.00	0.03	0.000823	0.05	0.05
n0 [fm−3]	0.1481	0.1484	0.1505	0.1505	0.1505
E/A [MeV]	−16.24	−16.30	−16.28	−16.28	−16.28
K [MeV]	271.5	230.0	238.0	238.0	238.0
${m}_{N}^{* }/{m}_{N}$	0.595	0.610	0.593	0.593	0.593
${E}_{\mathrm{sym}}({n}_{0})$ [MeV]	37.3	32.6	37.6	30.2	30.2
L [MeV]	118.2	60.5	112.8	44.3	41.0
$P({n}_{0})[\mathrm{MeV}\,{\mathrm{fm}}^{-3}]$	5.99	3.18	5.81	2.27	2.06

Note.  The mass of the nucleon is set to 939 MeV. Also shown are the corresponding energy per particle (E/A), compression modulus (K), and effective nucleon mass ${m}_{N}^{* }/{m}_{N}$ at saturation density n₀, as well as the symmetry energy (E_sym), slope of the symmetry energy (L) and PNM pressure (P) at n₀.

Download table as:  ASCIITypeset image

Next, we obtain the mass–radius (M–R) relation of neutron stars for a given EoS by solving the TOV equations (Oppenheimer & Volkoff 1939). As mentioned in Section 2, for the crust region of the star we have employed the EoS recently derived in Sharma et al. (2015).⁴ In this section we focus on neutron stars with cores of purely nucleonic matter; hence, we compute the EoS of the core assuming a β-equilibrated and charge neutral uniform liquid of neutrons, protons, and leptons (electrons and muons). As expected from its stiff EoS, the NL3 set predicts a large maximum mass (M_max ≈ 2.8 M_⊙) and large stellar radii (≈13 km for M_max and ≈15 km for a typical neutron star of 1.5 M_⊙); see the M–R relations plotted in Figure 2 and the values in Table 2. In comparison, the soft EoS of the FSU model brings in a dramatic reduction of the stellar masses and radii. The two shaded bands in Figure 1 portray the observed masses of the heaviest neutron stars known, i.e., M = 1.97 ± 0.04 M_⊙ in the pulsar PSR J1614–2230 (Demorest et al. 2010) and M = 2.01 ± 0.04 M_⊙ in the pulsar PSR J0348+0432 (Antoniadis et al. 2013). These two astrophysical measurements are arguably the most accurate constraints available so far to validate or defeat the model predictions for the high-density EoS. The recently formulated relativistic parameter set FSU2 (Chen & Piekarewicz 2014)—based on the same Lagrangian we are discussing—is one of the first best-fit models to take into account the condition of a limiting stellar mass of 2 M_⊙ in the calibration of the parameters (also see Erler et al. 2013; Chen & Piekarewicz 2015b). The FSU2 model has been optimized to accurately reproduce the experimental data on a pool of properties of finite nuclei with the maximum neutron star mass observable included in the fit (Chen & Piekarewicz 2014). The resulting FSU2 set has ζ = 0.0256 and Λ_ω = 0.0008, see Table 1. In consonance with these values, we can appreciate in Figure 1 that FSU2 predicts an intermediate EoS between the stiff EoS of the NL3 set (ζ = Λ_ω = 0) and the soft EoS of the FSU set (ζ = 0.06, Λ_ω = 0.03). Accordingly, FSU2 produces a neutron star M–R relation located in between the curves of NL3 and FSU in Figure 2. FSU2 yields a heaviest stellar mass ${M}_{\max }=2.07\,{M}_{\odot }$ with a radius of 12.1 km, and predicts 1.5 M_⊙ stars with a radius of 14 km, see Table 2.

Zoom In Zoom Out Reset image size

Table 2.  Neutron Star Properties Obtained for the Different Nuclear Models Discussed in This Work

Composition	Models	${M}_{\max }/{M}_{\odot }$	R(Mmax)	nc(Mmax)/n0	R(1.5 M⊙)	Y onset
			(km)		(km)	(n/n0)
	NL3	2.77	13.3	4.5	14.8	⋯
	FSU	1.72	10.8	7.8	12.2	⋯
pneμ	FSU2	2.07	12.1	5.9	14.0	⋯
	FSU2R	2.05	11.6	6.3	12.8	⋯
	FSU2H	2.38	12.3	5.3	13.2	⋯
	NL3	2.27	12.9	5.3	14.8	1.9
pnYeμ	FSU2	1.76	12.1	6.3	13.9	2.1
	FSU2R	1.77	11.6	6.5	12.8	2.4
	FSU2H	2.03	12.0	5.8	13.2	2.2
pneμ	FSU2R	2.11	11.6	6.1	12.8
(Bc = 2 × 1018 G)	FSU2H	2.42	12.3	5.2	13.2
pnYeμ	FSU2R	1.88	11.6	6.3	12.8	2.4
(Bc = 2 × 1018 G)	FSU2H	2.15	12.3	5.3	13.2	2.2

Note.  Results are shown for nucleonic-only (pneμ) or hyperonic (pnYeμ) stars, and including or not a magnetic field having the profile of the solid line in Figure 6 (Section 4). The quantity n_c(M_max)/n₀ denotes the central baryonic density at the maximum mass, M_max, normalized to the corresponding saturation density, n₀, whereas Y onset is the onset of appearance of hyperons normalized to n₀.

Download table as:  ASCIITypeset image

While the limiting stellar mass is governed by the stiffness of the EoS above several times the saturation density n₀ (see column ${n}_{c}({M}_{\max })/{n}_{0}$ in Table 2), the radius of a canonical neutron star is dominated by the density dependence of the EoS of PNM at 1–2 times n₀ (Lattimer & Prakash 2007; Ozel & Freire 2016). Thus, observational information on masses and radii of neutron stars has the potential to uniquely pin down the nuclear EoS in a vast density region. As mentioned in the Introduction, several of the recent astrophysical analyses for radii (Guillot et al. 2013; Guver & Ozel 2013; Guillot & Rutledge 2014; Heinke et al. 2014; Lattimer & Steiner 2014; Ozel & Freire 2016; Ozel et al. 2016) are converging in the 9–12 km range (also see Fortin et al. 2015 for a detailed discussion). The review study of Lattimer & Prakash (2016) indicates a similar range around 11–13 km for the radii of canonical neutron stars. The possibility that neutron stars have these small radii is as exciting as it is deeply challenging for nuclear theory. Note that small radii demand a sufficiently soft EoS below twice the saturation density, while the observed large masses require that the same EoS must be able to evolve into a stiff EoS at high densities. It is therefore timely to explore whether such small radii can be obtained by the EoS of the covariant field-theoretical Lagrangian (1), (2), while fulfilling at the same time the maximum mass constraint of 2 M_⊙ and the phenomenology of the atomic nucleus.

To construct the new EoS we start from the FSU2 model and increase the Λ_ω coupling. This softens the PNM pressure especially up to densities of 1.5–2n₀. For a given stellar mass there is less pressure to balance gravity, thereby leading to a more compact object of smaller radius. The increase of Λ_ω also produces a certain reduction of the PNM pressure in the high-density sector. This may spoil the 2 M_⊙ maximum mass but can be counteracted by a decrease of the strength of the ζ coupling. During the change of the (Λ_ω, ζ) couplings, we refit the remaining couplings ${g}_{\sigma N},{g}_{\omega N},{g}_{\rho N},\kappa$, and λ of the nucleon–meson Lagrangian (1), (2) by invoking the same saturation properties of FSU2 in SNM (i.e., same saturation density n₀, energy per particle E/A, compression modulus K, and effective nucleon mass ${m}_{N}^{* }$) and a symmetry energy E_sym = 25.7 MeV at subsaturation density n = 0.10 fm⁻³. The last condition arises from the fact that the binding energies of atomic nuclei constrain the symmetry energy at an average density of nuclei of ∼0.10 fm⁻³ better than the symmetry energy at normal density n₀ (Horowitz & Piekarewicz 2001a; Centelles et al. 2009). We found that under this protocol a noteworthy decrease of neutron star radii is achieved with Λ_ω = 0.05 and ζ = 0.024. We refer to the resulting model as FSU2R. The coupling constants and several bulk properties of FSU2R are collected in Table 1.

We observe in Figure 1 that the EoS of the new FSU2R model is within the boundaries deduced in the studies of energetic HICs (Fuchs et al. 2001; Danielewicz et al. 2002; Lynch et al. 2009). It is worth noting that FSU2R features a soft PNM EoS at $n\lesssim 1.5\mbox{--}2{n}_{0}$ and a stiff PNM EoS at n ≳ 2n₀—apparently a necessary condition to satisfy small radii and heavy limiting neutron star masses. The reduction of the stellar radii in FSU2R compared with the other parameterizations of the Lagrangian (1), (2) is very clear from Figure 2, also see Table 2. The maximum mass of 2.05 M_⊙ calculated with FSU2R is compatible with the heaviest neutron stars (Demorest et al. 2010; Antoniadis et al. 2013) and is characterized by a radius of 11.6 km. For canonical neutron stars of 1.4–1.5 solar masses, FSU2R predicts radii of ≈12.8 km, which are more compact than in the other EoSs, see Table 2. Hence, the smaller radii reproduced by the new model point toward the reconciliation between the nuclear EoS, the largest neutron star masses (Demorest et al. 2010; Antoniadis et al. 2013), and the recent extractions of small neutron star sizes from the astrophysical observations of quiescent low-mass X-ray binaries (Guillot & Rutledge 2014) and X-ray bursters (Guver & Ozel 2013) (also see Guillot et al. 2013; Heinke et al. 2014; Lattimer & Steiner 2014; Lattimer & Prakash 2016; Ozel & Freire 2016; Ozel et al. 2016). We are only aware of similar models RMF012 and RMF016 (also called FSUGarnet) introduced in a recent work (Chen & Piekarewicz 2015b). The RMF012 and RMF016 models were fitted with the same procedure of the FSU2 model of Chen & Piekarewicz (2014) but requiring values for the neutron skin thickness of the ²⁰⁸Pb nucleus of, respectively, 0.12 fm and 0.16 fm. As reported in Chen & Piekarewicz (2015a, 2015b), the RMF016 model supports 2 M_⊙ neutron stars and leads to a radius of 13 km for a 1.4 M_⊙ star, similarly to the predictions we obtain with our FSU2R model.

3.2. Implications for Finite Nuclei: Symmetry Energy, Slope of the Symmetry Energy, and Neutron Skin Thickness

Once the new EoS has been calibrated for neutron stars, it is important to review its implications for the physics regime of atomic nuclei since this regime is accessible in laboratory experiments. We first verify that the new model FSU2R is able to provide a satisfactory description of the best known properties of nuclei, i.e., nuclear ground-state energies and sizes of the nuclear charge distributions. We display in Figure 3 our results for the energies and charge radii of a set of representative nuclei ranging from the light ¹⁶O to the heavy ²⁰⁸Pb. The experimental data of these same nuclei were used in the fit of the FSU2 model in Chen & Piekarewicz (2014). In Figure 3, we show the predictions of our FSU2R model alongside the experimental values and the results from the parameter sets NL3, FSU, and FSU2. It can be seen that the four models successfully reproduce the energies and charge radii across the mass table. The agreement of FSU2R with experiment is overall comparable to the other models. We find that the differences between FSU2R and the experimental energies and radii are at the level of 1% or smaller. We mention that we have not drawn error bars of the experimental data in Figure 3, because the nuclear masses and charge radii are measured so precisely (Wang et al. 2012; Angeli & Marinova 2013) that the experimental uncertainties cannot be resolved in the plot.

Zoom In Zoom Out Reset image size

For our purposes, of special relevance is the fact that the neutron density distributions and other isospin-sensitive observables of atomic nuclei are closely related to the density dependence of the symmetry energy, which in FSU2R has been tailored to supply small stellar radii. The stiffness of the symmetry energy with density is conventionally characterized by its density slope L at the saturation point: $L=3{n}_{0}{\left(\tfrac{\partial {E}_{\mathrm{sym}}(n)}{\partial n}\right)}_{{n}_{0}}$. The L parameter and the pressure P(n₀) of PNM at saturation density are related as $P({n}_{0})=\tfrac{1}{3}{n}_{0}L$ (Lynch et al. 2009; Lattimer & Prakash 2016). The new FSU2R EoS yields ${E}_{\mathrm{sym}}({n}_{0})=30.2\,\mathrm{MeV}$ for the symmetry energy at saturation and a slope parameter L = 44.3 MeV, which corresponds to a mildly soft nuclear symmetry energy. The PNM pressure at saturation is P(n₀) = 2.27 MeV fm⁻³. We have collected these values in Table 1 along with the results for E_sym(n₀), L, and P(n₀) from the other discussed EoSs; now, large differences can be appreciated among the models.

Despite the fact that a precise knowledge of the density dependence of the symmetry energy remains elusive, the windows of values for E_sym(n₀) and the slope parameter L have been continuously narrowed down as the empirical and theoretical constraints have improved over recent years (see, e.g., Li et al. 2014 for a topical review). Remarkably, the values of 30.2 MeV for E_sym(n₀) and 44.3 MeV for L that we find after constraining the EoS to reflect small neutron star radii turn out to be very consistent with the newest determinations of the symmetry energy and its slope at saturation, see Figure 4. Indeed, the quoted FSU2R values overlap with the ranges 30 ≲ E_sym(n₀) ≲ 35 MeV and 20 ≲ L ≲ 66 MeV extracted in Roca-Maza et al. (2015) from the recent high-resolution measurements at RCNP and GSI of the electric dipole polarizability α_D in the nuclei ²⁰⁸Pb (Tamii et al. 2011), ¹²⁰Sn (Hashimoto et al. 2015), and ⁶⁸Ni (Rossi et al. 2013). We note that the dipole polarizability α_D, related to the response of nuclei to an external electric field, has been identified as one of the strongest isovector indicators (Reinhard & Nazarewicz 2010). Also note that, compared to hadronic experiments used to probe the symmetry energy, the electromagnetic reactions involved in the measurements of the α_D observable (Tamii et al. 2011; Rossi et al. 2013; Hashimoto et al. 2015) are particularly suited because they are not hindered by large or uncontrolled uncertainties. The FSU2R predictions for E_sym(n₀) and L also fit within the windows 29 ≲ E_sym(n₀) ≲ 33 MeV and 40 ≲ L ≲ 62 MeV obtained in Lattimer & Lim (2013) from the combined analysis of a variety of empirical nuclear constraints and astrophysical information, which are in line with similar windows obtained in other recent studies (Tsang et al. 2012; Li et al. 2014; Lattimer & Prakash 2016). It also deserves to be mentioned that the E_sym(n₀) and L values of the FSU2R EoS are quite compatible with the theoretical ranges 25.2 ≲ E_sym(n₀) ≲ 30.4 MeV and 37.8 ≲ L ≲ 47.7 MeV that have been derived from the latest progress in ab initio calculations of nuclear systems with chiral forces (Hagen et al. 2015).

Zoom In Zoom Out Reset image size

The neutron skin thickness ${\rm{\Delta }}{r}_{{np}}={r}_{n}-{r}_{p}$ (difference between the neutron and proton matter radii) of a heavy nucleus such as ²⁰⁸Pb, also provides strong sensitivity to the symmetry energy and the pressure of neutron-rich matter near saturation (Alex Brown 2000; Horowitz & Piekarewicz 2001a, 2001b; Centelles et al. 2009). Basically, the same nuclear pressure that is responsible for determining the radius of a canonical neutron star determines how far neutrons extend out further than protons in a nucleus. By the same token, models that produce smaller stellar radii are expected to predict thinner neutron skins. We find that the FSU2R model, constrained to small neutron star radii, predicts Δr_np = 0.133 fm in ²⁰⁸Pb. Unfortunately, neutron skins are difficult to extract from experiments in a model-independent fashion. The new experiments to measure neutron skins are being designed with electroweak and electromagnetic probes where, unlike hadronic experiments, the interactions with the nucleus (Abrahamyan et al. 2012), or at least the initial state interactions (Tarbert et al. 2014), are not complicated by the strong force. The challenging, purely electroweak (nearly model-independent) measurement of the neutron skin of ²⁰⁸Pb by parity violating electron scattering at JLab (Abrahamyan et al. 2012; Horowitz et al. 2012) has been able to provide Δr_np = 0.302 ± 0.177 fm for this isotope (Horowitz et al. 2012), although the data are not conclusive due to the large error bars (a follow-up measurement at JLab with better statistics has been proposed). The recent measurement of the neutron skin of ²⁰⁸Pb at the MAMI facility from coherent pion production by photons (Tarbert et al. 2014) has obtained Δr_np = 0.15 ± 0.03 fm. A similar range 0.13 ≲ Δr_np ≲ 0.19 fm for ²⁰⁸Pb is extracted (Roca-Maza et al. 2015) by comparing theory with the accurately measured electric dipole polarizability in ²⁰⁸Pb (Tamii et al. 2011), ¹²⁰Sn (Hashimoto et al. 2015), and ⁶⁸Ni (Rossi et al. 2013). Thus, the FSU2R prediction of a neutron skin of 0.133 fm in ²⁰⁸Pb turns out to be fairly compatible within error bars with the recent determinations of this isospin-sensitive observable.

In summary, when the nuclear EoS has been constrained to encode the recent astrophysical indications of small neutron star radii, yet without compromising massive stars, a high degree of consistency has emerged between the predictions of the model and the latest terrestrial information on the symmetry energy, its density dependence, and neutron skins, as well as with the constraints inferred from state-of-the-art ab initio microscopic calculations (Hagen et al. 2015). All in all, we believe that the present findings make a compelling case in favor of the prospect that neutron stars may have small, or moderate-to-small, radii.

Having calibrated the nuclear model to produce small neutron star radii and fulfil maximum masses of 2 M_⊙, while at the same time reproducing the phenomenology of atomic nuclei and the empirical constraints from collective flow and kaon production in HICs, we explore in this section the effect on the EoS and neutron stars of including hyperons and magnetic fields.

We should first determine the value of the hyperon couplings in our RMF model. Those couplings are calculated by fitting the experimental data available for hypernuclei, in particular, the value of the optical potential of hyperons extracted from these data. In our model, the contribution to the potential of a hyperon i in j-particle matter is given by

$$\begin{eqnarray}&&{U}_{i}^{(j)}({n}_{j})=-{g}_{\sigma i}\,{\bar{\sigma }}^{(j)}+{g}_{\omega i}\,{\bar{\omega }}^{(j)}+{g}_{\rho i}\,{I}_{3i}\,{\bar{\rho }}^{(j)}+{g}_{\phi i}\,{\bar{\phi }}^{(j)},\end{eqnarray} \tag{ 24 }$$

where ${\bar{\sigma }}^{(j)}$, ${\bar{\omega }}^{(j)}$, ${\bar{\rho }}^{(j)}$ and ${\bar{\phi }}^{(j)}$ are the values of the meson fields in the j-particle matter and I_3i stands for the third component of the isospin operator.

The couplings of the hyperons to the vector mesons are related to the nucleon couplings, ${g}_{\omega N}$ and ${g}_{\rho N}$, by assuming SU(3)-flavor symmetry, the vector dominance model and ideal mixing for the physical ω and ϕ mesons, as e.g., employed in many recent works (Schaffner & Mishustin 1996; Weissenborn et al. 2012; Colucci & Sedrakian 2013; Miyatsu et al. 2013; Banik et al. 2014). This amounts to assuming the following relative coupling strengths:

$$\begin{eqnarray}{g}_{\omega {\rm{\Lambda }}}\,:{g}_{\omega {\rm{\Sigma }}}\,:{g}_{\omega {\rm{\Xi }}}\,:{g}_{\omega N} & = & \displaystyle \frac{2}{3}\,:\displaystyle \frac{2}{3}\,:\displaystyle \frac{1}{3}\,:1\\ {g}_{\rho {\rm{\Lambda }}}\,:{g}_{\rho {\rm{\Sigma }}}\,:{g}_{\rho {\rm{\Xi }}}\,:{g}_{\rho N} & = & 0\,:1\,:1\,:1\\ {g}_{\phi {\rm{\Lambda }}}\,:{g}_{\phi {\rm{\Sigma }}}\,:{g}_{\phi {\rm{\Xi }}}\,:{g}_{\omega N} & = & -\displaystyle \frac{\sqrt{2}}{3}\,:-\displaystyle \frac{\sqrt{2}}{3}\,:-\displaystyle \frac{2\sqrt{2}}{3}\,:1,\end{eqnarray} \tag{ 25 }$$

and ${g}_{\phi N}=0$. Note that the isospin operator ${I}_{3i}$ appearing in the definition of the potentials in Equation (24) implements the relative factor of 2 missing in the 1:1 relation between ${g}_{\rho {\rm{\Sigma }}}$ and ${g}_{\rho N}$ displayed in Equation (25), so that the effective coupling of the ρ meson to the Σ hyperon (I₃ = −1, 0, +1) is twice that to the nucleon (I₃ = −1/2, +1/2), as required by the symmetries assumed.

The coupling of each hyperon to the σ field is adjusted to reproduce the hyperon potential in SNM derived from hypernuclear observables (see, e.g., Hashimoto & Tamura 2006; Gal et al. 2016). The Λ binding energy of Λ-hypernuclei is well reproduced by an attractive Woods–Saxon potential of depth ${U}_{{\rm{\Lambda }}}^{(N)}({n}_{0})\sim -28\,\mathrm{MeV}$ (Millener et al. 1988). The analyses of the (π⁻, K⁺) reaction data on medium to heavy nuclei (Noumi et al. 2002) performed in Harada & Hirabayashi (2006) and Kohno et al. (2006) revealed a moderately repulsive Σ-nuclear potential in the nuclear interior of around 10–40 MeV, while the fits to Σ⁻ atomic data (Friedman & Gal 2007) indicate a clear transition from an attractive Σ potential in the surface, to a repulsive one in the interior, although the size of the repulsion cannot be precisely determined. As for the strangeness −2 systems, the Nagara event (Takahashi et al. 2001) and other experiments providing consistency checks established the size of the ΛΛ interaction to be mildly attractive, ${\rm{\Delta }}{B}_{{\rm{\Lambda }}{\rm{\Lambda }}}{(}_{{\rm{\Lambda }}{\rm{\Lambda }}}^{6}\mathrm{He})=0.67\pm 0.17\,\mathrm{MeV}$ (Ahn et al. 2013), while the knowledge obtained for the Ξ–N interaction is more uncertain. Analyses of old emulsion data indicate a sizable attractive Ξ-nucleus potential of ${U}_{{\rm{\Xi }}}^{(N)}({n}_{0})=-24\pm 4\,\mathrm{MeV}$ (Dover & Gal 1983), while the missing-mass spectra of the $({K}^{-},{K}^{+})$ reaction on a ¹²C target suggest a milder attraction of ≳−20 MeV (Fukuda et al. 1998) or ∼−14 ± 2 MeV (Khaustov et al. 2000). These values are compatible with the recent analysis of the nuclear emulsion event KISO, claiming to have observed a nuclear bound state of the Ξ⁻–¹⁴N system with a binding energy of 4.38 ± 0.25 MeV (Nakazawa et al. 2015). From the above considerations, we fix the hyperon potentials in SNM to the following values:

$$\begin{eqnarray}{U}_{{\rm{\Lambda }}}^{(N)}({n}_{0}) & = & -\,28\,\mathrm{MeV}\\ {U}_{{\rm{\Sigma }}}^{(N)}({n}_{0}) & = & +\,30\,\mathrm{MeV}\\ {U}_{{\rm{\Xi }}}^{(N)}({n}_{0}) & = & -\,18\,\mathrm{MeV},\end{eqnarray} \tag{ 26 }$$

which allow us to determine the couplings ${g}_{\sigma {\rm{\Lambda }}}$, ${g}_{\sigma {\rm{\Sigma }}}$ and ${g}_{\sigma {\rm{\Xi }}}$, from Equation (24). We finally note that the coupling of the ϕ meson to the Λ baryon is reduced by 20% with respect to its SU(3) value in order to obtain a ΛΛ bond energy in Λ matter at a density n_Λ ≃ n₀/5 of ${\rm{\Delta }}{B}_{{\rm{\Lambda }}{\rm{\Lambda }}}({n}_{0}/5)=0.67\,\mathrm{MeV}$, thereby reproducing the Nagara event (Ahn et al. 2013).

Let us comment on the fact that the presence of hyperons in the neutron star interior and their influence on the EoS suffer from uncertainties tied to our lack of knowledge of the hyperon–nucleon and hyperon–hyperon interactions around the hyperon onset density of ∼2n₀ and beyond. This freedom has been exploited by different groups to build up RMF models that ensure the existence of neutron stars with masses larger than 2 M_⊙ even with the presence of hyperons (see for example Bednarek et al. 2012; Weissenborn et al. 2012; van Dalen et al. 2014; Oertel et al. 2015). While the radii of neutron stars are essentially determined by the nucleonic part of the EoS, the uncertainties in the hyperon interactions reflect on maximum masses that are scattered within a 0.3 M_⊙ band, as can be seen from the thorough analysis of various models done by Fortin et al. (2015), which is consistent with admitting a deviation by at most 30% of the symmetries assumed to determine the hyperon coupling constants (Weissenborn et al. 2012). These results provide an estimate of the uncertainties that one must admit in the hyperonic sector until data on the hyperon interactions at higher densities, coming for instance from HIC experiments (Morita et al. 2015), become available. As explained in the preceding paragraph, in this work we have simply made a minimal adjustment of the hyperon parameters away from the symmetry constraints imposed by Equation (25) in order to reproduce the known hypernuclear properties.

In Figure 5 we show how the presence of hyperons affects the M–R relationship for some representative nuclear EoSs selected from the previous section: the highly stiff EoS of the NL3 model and the FSU2 and FSU2R EoSs. These models differ essentially on the lower density ($n\lesssim 2{n}_{0}$) and/or the higher density ($n\gtrsim 2{n}_{0}$) stiffness of the EoS. As already noted (see Horowitz & Piekarewicz 2001b; Chen & Piekarewicz 2014, 2015a) and also discussed in the previous section, models with a larger value of the Λ_ω coupling produce a softer symmetry energy and, in consequence, become more compressible leading to stars with higher central densities and smaller radii. The presence of hyperons softens the EoS by essentially releasing Fermi pressure. Thereby, the stars get further compressed than their nucleonic counterparts, and the maximum masses get reduced by about 15%, as seen by the thick lines in Figure 5. It is also seen from this figure that the occurrence of hyperons leaves the stellar radii almost unaffected.

Zoom In Zoom Out Reset image size

Except for the NL3 model, which has shown to be exceedingly stiff at supranormal densities, the maximum masses of hyperonic stars attained by the other two models are too low, of about M_max = 1.8 M_⊙, to reproduce the 2 M_⊙ limit. The specific values of the maximum masses of hyperonic stars for all these models can be seen in Table 2 in the "pneYμ" section. We observe a weak sensitivity to the slope of the symmetry energy, as FSU2 and FSU2R produce essentially the same maximum masses. This was already noted in the context of nucleonic-only EoSs in Horowitz & Piekarewicz (2001b). As the symmetry energy softens, the star simply becomes more compressed and attains a larger central density, but reaches a similar maximum mass value. This phenomenology remains when hyperons are present, as was also found in Providencia & Rabhi (2013), the only difference being that the hyperonic stars attain a lower maximum mass and have a higher central density than their nucleonic-only counterparts, as expected for a softer EoS. This can be seen upon comparing the values shown in the "pneμ" and "pneYμ" sections of Table 2.

Since the hyperonic EoS based on the FSU2R model does not produce M_max > 2 M_⊙, we tense the parameters of this nuclear model a little further so as to make it stiffer. We essentially reduce the value of ζ from 0.024 to 0.008, which stiffens the EoS at densities larger than twice the saturation density, i.e., around the region where hyperons start appearing (see the hyperon onset density for the different models in Table 2). The remaining parameters of the model are refitted so as to reproduce the SNM saturation properties of the FSU2 model and a symmetry energy E_sym = 26.2 MeV at density n = 0.10 fm⁻³. The values of the parameters of this new interaction, named FSU2H, are listed in Table 1, together with the predicted E_sym value at saturation density and its slope L, which fall comfortably within the newest empirical and theoretical constraints of these quantities, as can be seen in Figure 4. The couplings of the hyperons to the different vector mesons can be readily obtained from Equation (25), and those to the σ meson, determined from fixing the hyperon potentials in SNM, are ${g}_{\sigma {\rm{\Lambda }}}=0.6113$, ${g}_{\sigma {\rm{\Sigma }}}=0.4665$ and ${g}_{\sigma {\rm{\Xi }}}=0.3157$.

We note that the FSU2H interaction produces a certain overpressure in SNM at n ≳ 2n₀, since the pressure falls above the allowed region obtained from the modeling of collective flow in HICs, as seen by the long dashed line in the upper part of Figure 1. Nevertheless, the EoS for PNM, seen in the lower panel of this figure, falls within the PNM extrapolated band compatible with collective flow. Since neutron-star matter in beta equilibrium is highly asymmetric we consider this model to be sufficiently realistic to describe neutron stars, whose properties are presented in Table 2. We observe that the maximum mass of 2.38 M_⊙ obtained for a "pneμ" neutron star with the FSU2H model gets reduced to 2.03 M_⊙, with a radius of 12 km, when hyperons are present. We also observe that the radius of a canonical star of ∼1.5 M_⊙ gets slightly enhanced from 12.8 km for FSU2R to 13.2 km for FSU2H, which is the price one pays for having stiffened the EoS.

On comparing the "pneμ" with the "pneYμ" parts of Table 2 we essentially see, as in Figure 5, a reduction of about 15% on the maximum mass when hyperons are allowed to appear in the neutron star cores. Since the hyperonic EoSs become more compressible, the "pnYeμ" stars attain higher central densities, but the radii of the maximum-mass stars stay rather similar to their nucleonic-only counterparts.

We note that the FSU2H parameterization, which produces M_max > 2 M_⊙ even in the presence of hyperons, fulfils the pressure constraint in neutron star matter at saturation density n₀:

$$\begin{eqnarray}&&1.7\,{{\rm{MeV}}{\rm{fm}}}^{-3}\lt P({n}_{0})\lt 2.8\,{{\rm{MeV}}{\rm{fm}}}^{-3},\end{eqnarray} \tag{ 27 }$$

estimated in Fortin et al. (2015) from the results shown in Hebeler et al. (2013), which were obtained from microscopic calculations of PNM based on chiral two-nucleon and three-nucleon interactions, and which are in remarkable agreement with the quantum Monte Carlo results of Gandolfi et al. (2012), obtained from the phenomenological Argonne v18 NN potential plus three-nucleon forces. It is argued in Fortin et al. (2015) that nearly all hyperonic EoS models that are able to sustain M_max > 2 M_⊙ produce large PNM pressures of about 5 MeV fm⁻³ at saturation density, leading to an overpressure of the nucleonic (pre-hyperon) segment and resulting in large radii of around 14 km or more for neutron stars in the range 1 < M/M_⊙ < 1.6. Our FSU2H model does not encounter this problem, since it gives a PNM pressure of ∼2 MeV fm⁻³ at n₀ (see Table 1), well within the constraint of Equation (27), and as a consequence is able to reach a smaller radius of 13 km. We note that the symmetry energy slope parameter L of the hyperonic models analyzed in Fortin et al. (2015) lies in the range of values 67–118 MeV, which deviate considerably from the recent constraints displayed in Figure 4.

We now discuss the effect of including a magnetic field in our nucleonic and hyperonic stars. We consider a density-dependent magnetic field with a profile of the type

$$\begin{eqnarray}&&B(n)={B}_{s}+{B}_{c}\{1-\exp [-\beta {(n/{n}_{0})}^{\gamma }]\},\end{eqnarray} \tag{ 28 }$$

introduced in Chakrabarty et al. (1997) and employed in several other works (Rabhi & Providencia 2010; Lopes & Menezes 2012; Sinha et al. 2013). We take a surface magnetic field value of B_s = 10¹⁵ G, consistent with the surface magnetic fields of observed magnetars (Vasisht & Gotthelf 1997; Kouveliotou et al. 1998; Woods et al. 1999) and a core magnetic field value of B_c = 2 × 10¹⁸ G, which is sufficiently strong to produce distinguishable effects on the properties of neutron stars. The parameters β and γ control the density where the magnetic field saturates and the steepness of the transition from the surface to the core field, respectively. We take β = 0.0065 and γ = 3.5 which ensure that the magnetic field has practically saturated to its maximum value at around $5\mbox{--}6{n}_{0}$, a range that covers the typical central densities of the maximum mass neutron stars explored in this work. Moreover, the indicated β and γ values produce moderate field values below saturation density, as can be seen by the solid line in Figure 6. We note that this field profile does not incur on instabilities of the parallel component of the pressure ${P}_{| | }$ associated to rapidly rising magnetic field toward relatively strong central values (Sinha et al. 2013).

Zoom In Zoom Out Reset image size

The effect of this magnetic field on the M–R relationship is displayed in Figure 7. On the left (right) panel we show the results obtained for the EoS employing the FSU2R (FSU2H) model. The solid lines correspond to vanishing magnetic field, while the dashed lines include the effects of the magnetic field with the density profile discussed above. The thin black lines show the results for nucleonic neutron stars and the thick red lines correspond to the hyperonic stars. As observed in earlier works (Lopes & Menezes 2012), including the magnetic field produces stars with larger maximum masses. This is essentially a consequence of the increase in the total pressure which, apart from the matter pressure P_matt, also includes the extra average field pressure component, as seen in Equation (21). The size of this enhancement is larger for the hyperonic than for the nucleonic stars, which is essentially due, as we will show below, to the additional effect of de-hyperonization that takes place in the presence of a magnetic field. The reduction of hyperons is responsible for enhancing the value of the matter pressure P_matt, since the Fermi contributions of the other species are larger than in the B = 0 case. Nevertheless, the increase in the maximum mass induced by magnetic field effects is not enough to produce hyperonic star masses of the order of 2 M_⊙ in the case of the FSU2R model, as the dashed red line on the left panel does not reach the observational bands. The effects of the magnetic field on the M–R relationships obtained with the FSU2H EoS (right panel) are similar to those for the FSU2R EoS, the only difference being that the constraint M_max > 2 M_⊙ is now amply fulfilled, since the FSU2H model served this purpose even in the absence of a magnetic field.

Zoom In Zoom Out Reset image size

We now explore the effect of employing different magnetic field profiles having the same surface and central values, B_s = 10¹⁵ G and B_c = 2 × 10¹⁸ G, but different β and γ parameters. To this end, we consider, in addition to the profile obtained with the parameters (β, γ) = (0.0065, 3.5) chosen in this work, the profiles with (β, γ) = (0.05, 2) (Rabhi & Providencia 2010) and (β, γ) = (0.1, 1) (Sinha et al. 2013), which are represented, respectively, by the dashed–dotted and dashed lines in Figure 6. We observe that our parameterization produces a substantially lower magnetic field in the n < 2n₀ region and reaches 90% of the saturation value around 5n₀, while the dashed-dotted parameterization does it right after 6n₀. The parameterization of the dashed line does not even reach the value B = 10¹⁸ G within the densities of interest (n ≲ 6n₀).

In Figure 8 we display the M–R relationships obtained with these profiles, together with the zero magnetic field case, represented by a thin solid line. A noticeable dependence of the M–R relationship on the magnetic field profile is observed. The results for the (β, γ) = (0.05, 2) case (thick dashed–dotted line) are similar to those obtained with our (β, γ) = (0.0065, 3.5) parameterization (thick solid line), but the stars are produced with a somewhat larger radius since the magnetic field and, hence, the total pressure are larger in the pre-hyperon region. This is even more evident for the M–R relationship obtained with the (β, γ) = (0.1, 1) profile (thick dashed line), which produces stars that are ∼0.5 km wider than the other two cases and deviates from the 13 km maximum radius constraint. The reason is that this profile clearly gives larger magnetic fields in the n ≲ n₀ region, hence producing a larger total pressure and making the star less compressible.

Zoom In Zoom Out Reset image size

The particle fractions for beta-stable neutron star matter obtained using the FSU2H EoS are shown in Figure 9 as functions of the baryonic density. The upper panel displays the fractions in the absence of magnetic field, while the other panels implement the magnetic field with the three different profiles shown in Figure 6. Landau oscillations are seen in the charged particle fractions when a magnetic field is applied, reflecting the successive filling of the Landau levels as the quantity $({E}_{F}^{2}-{{m}^{* }}^{2})/2| q| B$ reaches integer values. For a fixed density, smaller magnetic fields accommodate more Landau levels and, correspondingly, more oscillations are observed, as seen for instance when comparing the three $B\ne 0$ panels in the n < 2n₀ density region, where the smallest field corresponds to the (β, γ) = (0.0065, 3.5) case. As density increases, so does the magnetic field in all the considered profiles, eventually needing only one Landau level to accommodate the population of the charged particles. The oscillations then tend to smooth out and disappear with increasing density. As is evident, the magnetic field mostly affects the charged particles, which in general increases their population with respect the B = 0 case. At low and intermediate densities up to n ∼ 4n₀, we clearly observe an increase in the occupation of negatively charged electrons and muons. This delays the appearance of the negatively charged hyperons, an effect that is especially visible for the Σ⁻ baryon, whose onset moves to n ≳ 4n₀ for all the considered magnetic field profiles.

Zoom In Zoom Out Reset image size

According to the results shown in the second panel of Figure 9, in the case of a magnetized hyperonic star having a mass of about 2 M_⊙ with a maximum density of about 5n₀ (see Table 2), the baryon fractions at the center would be 38% for n, 28% for Λ, 26% for p, 6% for Ξ⁻ and 2% for Σ⁻. In the B = 0 case (upper panel), these fractions would be 45% for n, 31% for Λ, 13% for p, 6% for Ξ⁻ and 5% for Σ⁻. We then see that the proton abundance can be twice as large in a magnetar as it is in a field-free star. Our results are qualitatively consistent with those obtained by other works in the literature studying the effect of magnetic fields in hyperonic stars (Broderick et al. 2002; Yue et al. 2009; Rabhi & Providencia 2010; Lopes & Menezes 2012; Sinha et al. 2013). We can conclude that, in general, hyperonic magnetars re-leptonize and de-hyperonize with respect to zero-field stars, while the proton abundance increases substantially. This might facilitate direct Urca processes, drastically altering the cooling evolution of the star.

We have obtained a new EoS for the nucleonic inner core of neutron stars that fulfils the constraints coming from recent astrophysical observations of maximum masses and determinations of radii, as well as the requirements from experimental nuclear data known from terrestrial laboratories. This EoS results from a new parameterization of the FSU2 force (Chen & Piekarewicz 2014), the so-called FSU2R model, that reproduces: (i) the 2 M_⊙ observations (Demorest et al. 2010; Antoniadis et al. 2013), (ii) the recent determinations of radii below 13 km region (Guillot et al. 2013; Guillot & Rutledge 2014; Heinke et al. 2014; Lattimer & Steiner 2014; Lattimer & Prakash 2016; Ozel et al. 2016), (iii) the saturation properties of nuclear matter and finite nuclei (Tsang et al. 2012; Chen & Piekarewicz 2014) and (iv) the constraints extracted from nuclear collective flow (Danielewicz et al. 2002) and kaon production (Fuchs et al. 2001; Lynch et al. 2009) in HICs.

The FSU2R model is obtained by modifying the Λ_ω and ζ couplings of the Lagrangian simultaneously, while recalculating the couplings ${g}_{\sigma N},{g}_{\omega N},{g}_{\rho N},\kappa$, and λ to grant the same saturation properties of FSU2 in SNM and a symmetry energy of 25.7 MeV at n = 0.10 fm⁻³. On the one hand, radii of 12–13 km are obtained, owing to the fact that we softened the symmetry energy and, consequently, the pressure of PNM at densities ∼1.5–2n₀, while reproducing the properties of nuclear matter and nuclei. Indeed, we obtain E_sym = 30 MeV and L = 44 MeV, which lie within the limits of recent determinations (see Lattimer & Lim 2013; Hagen et al. 2015; Roca-Maza et al. 2015). Moreover, the FSU2R model predicts a neutron skin thickness of 0.133 fm for the ²⁰⁸Pb nucleus, which is compatible with recent experimental results (Abrahamyan et al. 2012; Horowitz et al. 2012; Tarbert et al. 2014; Roca-Maza et al. 2015). On the other hand, we have stiffened the EoS above twice the saturation density, which satisfies the constraints of HICs (Fuchs et al. 2001; Danielewicz et al. 2002; Lynch et al. 2009) and allows for maximum masses of 2 M_⊙ (Demorest et al. 2010; Antoniadis et al. 2013). All in all, the FSU2R parameterization allows for a compromise between small stellar sizes and large masses, a task that seemed difficult to achieve in up-to-date RMF models.

We also analyze the consequences of the appearance of hyperons inside the core of neutron stars. The values of the hyperon couplings are determined from the available experimental information on hypernuclei, in particular by fitting to the optical potential of hyperons extracted from the data. On the one hand, we find that the radii of the neutron stars are rather insensitive to the appearance of the hyperons and, thus, still respect the observations of radii <13 km. On the other hand, we obtain a reduction of the maximum mass below 2 M_⊙ once hyperons appear due to the expected softening of the EoS. However, by refitting the parameters of the FSU2R model slightly, the new parameterization FSU2H fulfils the 2 M_⊙ limit while still reproducing the properties of nuclear matter and nuclei. The price to pay is a stiffer EoS in SNM as compared to the constraint derived from the modeling of HICs. Nonetheless, the HICs estimate in PNM is still satisfied by the FSU2H parameterization (Danielewicz et al. 2002).

We finally study the effect of high magnetic fields on the nucleonic and hyperonic EoSs. This is of particular interest for understanding the behavior of highly magnetized neutron stars, the so-called magnetars. Employing magnetic fields with crustal and interior values of ∼10¹⁵ G and ∼10¹⁸ G, respectively, we find EoSs that are stiffer and produce larger maximum-mass stars, while keeping radii in the 12–13 km range, both for nucleonic and hyperonic magnetars, as long as the magnetic field does not reach values larger than about 10¹⁷ G at saturation density. The particle fractions in the interior of the stars depend on the specific profile of the magnetic field, but the general trend with respect to zero-field stars is that hyperonic magnetars re-leptonize and de-hyperonize, while the amount of protons may double, a fact that may trigger direct Urca processes affecting the cooling and other transport properties of the star.

We are most grateful to J. Piekarewicz for a careful reading of the manuscript and for valuable comments. L.T. acknowledges support from the Ramón y Cajal research programme, FPA2013-43425-P Grant from Ministerio de Economia y Competitividad (MINECO) and NewCompstar COST Action MP1304. M.C. and A.R. acknowledge support from Grant No. FIS2014-54672-P from MINECO, Grant No. 2014SGR-401 from Generalitat de Catalunya, and the project MDM-2014-0369 of ICCUB (Unidad de Excelencia María de Maeztu) from MINECO. L.T. and A.R. acknowledge support from the Spanish Excellence Network on Hadronic Physics FIS2014-57026-REDT from MINECO.