THE NEUTRON STAR MASS–RADIUS RELATION AND THE EQUATION OF STATE OF DENSE MATTER

The masses of several neutron stars have been precisely measured using pulsar timing (Lorimer 2008); simultaneous mass and radius measurements, however, are considerably less certain. The leading candidates for such measurements are bursting neutron stars that show photospheric radius expansion (PRE; van Paradijs 1979; for a review, see Lewin et al. 1993) and transiently accreting neutron stars in quiescence (Rutledge et al. 1999).

Observations have already begun to determine the M–R relation for neutron stars and to place strong constraints on the equation of state (EOS) of dense matter. Previously derived constraints (Read et al. 2009; Özel et al. 2010; Steiner et al. 2010; Steiner & Gandolfi 2012) have several limitations, including the use of fixed parameterizations for the EOS. Those works neither showed their results to be independent of their parameterizations, including the possibility of deconfined quark matter, nor did they address the full set of sources analyzed here. Özel et al. (2010) considered only PRE sources, but these may be subject to considerable systematic errors (Steiner et al. 2010; Boutloukos et al. 2010; Suleimanov et al. 2011). Steiner et al. (2010) considered both types of sources but used a smaller data set. Hebeler et al. (2010) and Steiner & Gandolfi (2012) took advantage of modern predictions for pure neutron matter near the saturation density to constrain the M–R relation, but ignored systematic uncertainties associated with the observations. These systematic uncertainties could invalidate M and R results for all the PRE burst sources and/or some quiescent low-mass X-ray binary (qLMXB) sources, in which case the constraints on the EOS obtained by these works might be significantly altered. In this Letter, we address these limitations and attempt to make a more complete accounting of the systematic uncertainties due to modifications in either the EOS model or assumptions about the data set. Our constraints on, for example, the radius for a 1.4 solar mass neutron star are in some cases less stringent than those obtained previously, but are more robust because of the more careful accounting of the uncertainties described above.

At the lowest energy densities (⩽15 MeV fm⁻³ or mass densities ≲ 3 × 10¹³ g cm⁻³) the pressure–density relation is well understood (Haensel et al. 2007). Between 15 and 200–300 MeV fm⁻³, the EOS is well described by four parameters: the incompressibility, the skewness, the magnitude of the symmetry energy (S_v), and the parameter describing the density derivative of the symmetry energy (L), all evaluated at the nuclear saturation density (approximately 150 MeV fm⁻³). These parameters are constrained to varying degrees by experimental data (Lattimer & Lim 2012), including nuclear masses (Kortelainen et al. 2010), neutron skin thicknesses (Chen et al. 2010), giant dipole resonances and dipole polarizabilities (Trippa et al. 2008; Tamii et al. 2011; Piekarewicz et al. 2012), and heavy-ion collisions (Tsang et al. 2009). High-density matter is constrained by (1) causality (the speed of sound must not exceed the speed of light), (2) hydrodynamical stability, and (3) having a sufficient maximum mass (it must be greater than the largest well-determined neutron star mass, 1.97 ± 0.04 M_☉ for PSR J1614−2230; Demorest et al. 2010). In addition, we impose the constraint that implied neutron star masses cannot be less than what is thought to be achievable in supernova, about 0.8 M_☉.

The low-density part of all the EOS parameterizations (except for strange quark stars) is described as in Steiner et al. (2010). All of the parameters in the low-density EOS, as well as the parameters in the high-density EOS models, are described with uniform prior distributions. Our fiducial EOS model (A) parameterizes the high-density EOS as a set of two piecewise continuous power laws defining pressure P = ε^{1 + 1/n} as a function of energy density ε. Model A has four high-density parameters: the transition energy density between the low-density EOS and the first polytrope, the transition energy density between the first and second polytrope, and the two polytropic indices, n₁ and n₂. We also employ an EOS model (B) which is similar to model A, except that the exponents in the two polytropes are parameterized with uniform priors in Γ_i ≡ 1 + 1/n_i instead of n_i. The upper and lower limits of Γ_i correspond to the upper and lower limits of n_i from model A and Steiner et al. (2010). A third EOS model (C) is piecewise linear EOS: the low-density EOS is used up to energy densities of 200 MeV fm⁻³, and then consists of line segments that begin and end at four fixed energy densities: 400, 600, 1000, and 1400 MeV fm⁻³. The linear relation between energy densities of 1000 and 1400 MeV fm⁻³ is extrapolated to higher energy densities when necessary. The pressure on the energy density grid is allowed to be as large as 600 MeV.

A fourth EOS model (D) assumes that matter is a polytrope at intermediate densities and quark matter at high densities. The pressure of quark matter is described by Alford et al. (2005),

$$\begin{equation} P = \frac{3 a_4}{4 \pi ^2} \mu ^4 - \frac{3 a_2}{4 \pi ^2} \mu ^2 -B \,, \end{equation} \tag{ 1 }$$

where μ is the quark chemical potential. The quantity B is the bag constant and simulates confinement. The parameter a₄ describes corrections to the leading coefficient from a non-interacting Fermi gas (for which a₄ = 1). Corrections from perturbative quantum chromodynamics at high density suggest 0.6 < a₄ < 1. The parameter a₂ approximately describes corrections from the finite strange quark mass m_s and the quark superfluid gap Δ and is given by a₂ = m²_s − 4Δ² (see also Steiner & Gandolfi 2012). In order to include the effects of a possible mixed phase in a model-independent way, a polytrope is added in between the low-density EOS and quark matter. In this model, there are five parameters in total: the transition energy density between low densities and the polytrope, the polytropic index, the transition energy density between the polytrope and quark matter (used to fix the value of B), a₂, and a₄. To describe bare strange quark stars, model E applies Equation (1) at all densities, with neither a low-density EOS nor an intermediate polytrope. Model E thus has three parameters: a₂, a₄, and B. All parameter ranges (i.e., the boundaries of the prior distributions) are large enough to ensure that the results do not change significantly when the range is increased. Our results are sensitive to the shape of these prior distributions, and our variation in EOS models effectively probes variations of these shapes.

Our baseline data include eight neutron stars: four which produced PRE X-ray bursts (4U 1608–522 (Güver et al. 2010a), KS 1731–260 (Özel et al. 2011), EXO 1745–248 (Özel et al. 2009), and 4U 1820–30 (Güver et al. 2010b)) and four qLMXBs in the globular clusters M13 (Webb & Barret 2007), ω Cen (Webb & Barret 2007), 47 Tuc (Heinke et al. 2006), and NGC 6397 (Guillot et al. 2011). There are a handful of other neutron star mass and radius constraints, such as those for RXJ 1856 (Pons et al. 2002), but all of these may be subject to even larger uncertainties than the observations we use here and would thus not strongly modify our final results. Our baseline model interprets the PRE burst sources with an extended photosphere, and assumes the same distribution of distances as in Steiner et al. (2010).

We do find that neutron star radii depend (albeit weakly) on the EOS model used to analyze the data, causing the lower and upper 95% confidence limits to vary by about 0.8 km. This effectively tests the robustness of our constraints under variations of the shape of the prior distributions. We also consider several modifications to the baseline data set. Suleimanov et al. (2011) have suggested that the X-ray spectra for PRE sources are affected by accretion and the eclipse of the neutron star by the disk. This affects the normalization at late times and we take this into account by increasing the color correction factor f_C, taking 1.45 < f_C < 1.8 (modification I). Alternatively, Boutloukos et al. (2010) have suggested that the color correction factors are instead smaller as a result of magnetic confinement of the X-ray burst; this is considered in modification II in which we assume 1 < f_C < 1.35. Some previous works assumed that the photosphere of PRE bursts is coincident with the neutron star surface (Özel et al. 2010), and we also consider this scenario (III). Finally, we test the sensitivity of our results to the removal of any one source or class of sources by removing X7 (IV) or M13 (V), by removing all PRE sources (VI), and by removing all the qLMXBs (VII). To summarize, the standard models we examine include the baseline case (A), three variations of the EOS (B–D), and six other modifications of the baseline model varying the included data or its interpretation (A I–A VI). In addition, we also examine several more extreme scenarios which are described below.

We use the Bayesian method of Steiner et al. (2010), using marginal estimation to determine the posterior probability densities of quantities of interest. The marginal estimation integrals are performed using Markov Chain Monte Carlo.

In our baseline analysis, we find strong constraints on the M–R curve and on the dense matter EOS: the radius of a 1.4 M_☉ neutron star is between 11.2 and 12.3 km (95% confidence). The permissible radius range encompassing all variations of the EOS and interpretations of the astrophysical data, but not including the more extreme scenarios, is 10.4–12.9 km (95% confidence), only moderately larger than the baseline result. The 68% and 95% confidence ranges are displayed in the upper and middle portions of Table 1. We determine the M–R relation for a range of neutron star masses (see Figure 1). We also determine the 68% and 95% confidence intervals of the EOS of dense matter (Figure 2). The estimated uncertainty of the pressure is approximately 30%–50% at all densities achievable in neutron star interiors. In addition, the posterior distributions of the central energy density of the maximum mass star imply that the highest central density is ε ∼ 1200 MeV fm⁻³ (Lattimer & Prakash 2005).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Producing significantly different neutron star radii requires extreme assumptions regarding the EOS and the data. We now consider more extreme scenarios, which are presented in the bottom portion of Table 1 (see also Figure 3). To achieve significantly smaller radii, we must assume both that the color correction factor is anomalously small (⩽1.3) for all of the PRE sources and that the EOS has strong phase transitions (model C). In this case, we get radii as small as 9 km. Increasing the maximum mass constraint (modification VIII), as would be the case if the estimated most-likely mass of the pulsar B1957+20 is 2.4 M_☉ (van Kerkwijk et al. 2011), slightly increases radii relative to the baseline case.

Zoom In Zoom Out Reset image size

We obtain even larger radii, up to almost 14 km, if we add the long PRE burst source 4U 1724−307 and further assume, as suggested in Suleimanov et al. (2011; modification IX), that the short PRE bursts and the qLMXBs M13 and ω Cen not be considered because of modifications to their spectra due to accretion (Suleimanov et al. 2011). On the other hand, Güver et al. (2012) find that the long PRE burst of 4U 1724 does not fit modern atmosphere models as well as short bursts from the same source. A full resolution of this discrepancy is outside the scope of this work and may require more observational data to fully understand PRE bursts. Nevertheless, we have attempted to cover the most likely scenarios.

While we are able to significantly constrain the P–ε relation, determination of the composition of neutron star cores is not yet possible. To probe the core composition, we consider EOS model E, which describes the entire star by the high-density quark matter EOS used in model D, i.e., a self-bound strange quark star. In the mass range 1.4–2 solar masses, the radii are not significantly different from our baseline model so that there is no strong preference for either strange quark or hadronic stars; however, model E predicts radii significantly less than 10 km for low masses (⩽1.2 M_☉).

Our neglect of rotation is unlikely to affect our conclusions. Rotation increases the radius at the equator and decreases the radius at the poles, and this could be relevant for the interpretation of some PRE X-ray bursts: the rotation rate of 4U 1608−522 is 619 Hz. However, for EOSs that are likely to reproduce the observational data, this rotation rate increases the radius by less than 10% (Weber 1999). This introduces an uncertainty smaller than that due to variations in f_C, which we have already taken into account. The rotation rates for the qLMXBs in our sample are unknown. Assuming that they are similar to other qLMXBs, however, means that the effect of rotation is smaller than that of their distance uncertainties.

The relationship between pressure and energy density (Figure 2) that we determine from our baseline analysis from observations is consistent with effective field theory (Hebeler et al. 2010) and quantum Monte Carlo (Gandolfi et al. 2012; Steiner & Gandolfi 2012) calculations of low-density neutron matter. Note that these neutron matter results are incompatible with the Suleimanov et al. (2010) interpretation of 4U 1724−307 (Suleimanov et al. 2011) which suggested exclusion of short PRE bursts and qLMXBs M13 and ω Cen, also pointed out by Hebeler et al. (2010). Our results are also consistent with the high-density constraints on neutron matter from heavy-ion collisions (Danielewicz et al. 2002). In order to infer the constraints on neutron star matter from the neutron matter constraints in Danielewicz et al. (2002), we performed a small phenomenological correction for the small proton content using the method in Steiner & Gandolfi (2012). Also, we should note that the neutron matter constraints in Danielewicz et al. (2002) are not model-independent, and depend on assumptions about the high-density behavior of the nuclear symmetry energy.

Our results imply that over one-third of the modern Skyrme models studied in Stone et al. (2003) are inconsistent with observations. Covariant field-theoretical models that have symmetry energies that increase nearly linearly with density, such as the model NL3 (Lalazissis et al. 1997), are also inconsistent with our results, although they may still adequately describe isospin-symmetric matter in nuclei.

Our models do not place effective constraints on the symmetry parameter S_v, but do place significant constraints on the symmetry energy parameter L; these are summarized in Figure 4. The probability distribution for each model is renormalized to fix the maximum probability at unity and is then shifted upward by an arbitrary amount. The range that encloses all of the models and modifications to the data is 43.3–66.5 MeV to 68% confidence and 41.1–83.4 MeV to 95% confidence. The allowed values of L are substantially larger for model C because this parameterization more effectively decouples the low- and high-density behaviors of the EOS.

Zoom In Zoom Out Reset image size

Our preferred range for L is similar to that obtained from other experimental and observational studies (Tamii et al. 2011; Tsang et al. 2012; Steiner & Gandolfi 2012; Lattimer & Lim 2012) and experimental studies (e.g., Tsang et al. 2012; Tamii et al. 2011). Our results suggest that the neutron skin thickness of ²⁰⁸Pb (Typel & Brown 2000; Steiner et al. 2005) is less than about 0.20 fm. This result is independent of the EOS models (which include possible phase transitions) and data modifications described above. It is compatible with experiment (Horowitz et al. 2001) and also with measurements of the dipole polarizability of ²⁰⁸Pb (Reinhard & Nazarewicz 2010).

While we have endeavored to take into account some systematic uncertainties in our analysis, we cannot rule out corrections due to the small number of sources and to possible drastic modifications of the current understanding of low-mass X-ray binaries. Nevertheless, it is encouraging that these astrophysical considerations agree not only with nuclear physics experiments but also with theoretical studies of neutron matter at low densities and heavy-ion experiments at higher densities.

We thank G. Bertsch, J. Linnemann, and S. Reddy for useful discussions. This work is supported by Chandra grant TM1-12003X (A.W.S.), by the Joint Institute for Nuclear Astrophysics at MSU under NSF PHY grant 08-22648 (A.W.S. and E.F.B.), by NASA ATFP grant NNX08AG76G (A.W.S. and E.F.B.), and by DOE grants DE-FG02-00ER41132 (A.W.S.) and DE-AC02-87ER40317 (J.M.L.). J.M.L. and E.F.B thank the Institute for Nuclear Theory at the University of Washington for partial support during this work. E.F.B. is a member of an International Team in Space Science on type I X-ray bursts sponsored by the International Space Science Institute (ISSI) in Bern, and thanks the ISSI for hospitality during part of this work. Some of the computer code for this project is available at http://www.bamr.sourceforge.net/.