Maximum gravitational mass $M_{\mathrm{TOV}}=2.25_{-0.07}^{+0.08}{M}_{\odot }$ inferred at about 3% precision with multimessenger data of neutron stars

The maximal gravitational mass of nonrotating neutron stars ($M_{\mathrm{TOV}}$) is one of the key parameters of compact objects and only loose bounds can be set based on the first principle. With reliable measurements of the masses and/or radii of the neutron stars, $M_{\mathrm{TOV}}$ can be robustly inferred from either the mass distribution of these objects or the reconstruction of the equation of state of the very dense matter. For the first time we take the advantages of both two approaches to have a precise inference of $M_{\mathrm{TOV}}=2.25_{-0.07}^{+0.08}{M}_{\odot }$ (68.3% credibility), with the updated neutron star mass measurement sample, the mass-tidal deformability data of GW170817, the mass-radius data of PSR $\mathrm{J}0030+0451$ and PSR $\mathrm{J}0740+6620$, as well as the theoretical information from the chiral effective theory and perturbative quantum chromodynamics (pQCD) at low- and very high-energy densities, respectively. This narrow credible range is benefited from the suppression of the high $M_{\mathrm{TOV}}$ by the pQCD constraint and the exclusion of the low $M_{\mathrm{TOV}}$ by the mass function. Three different EoS reconstruction methods are adopted separately, and the resulting $M_{\mathrm{TOV}}$ and $R_{\mathrm{TOV}}$ are found to be almost identical, where $R_{\mathrm{TOV}}=11.90_{-0.60}^{+0.63}\mathrm{km}$ is the radius of the most massive nonrotating neutron star. This precisely evaluated $M_{\mathrm{TOV}}$ suggests that the EoS of neutron star matter is just moderately stiff and the $\sim 2.5–3M_{\odot }$ compact objects detected by the second-generation gravitational wave detectors are most likely the lightest black holes.

Figure 1

The red line represents the best-fit mass distribution, i.e., a two-component Gaussian mixture with a sharp cutoff of $M_{\max }=2.28M_{\odot }$, of the 136 neutron stars with gravitational mass measurements. Here we take 1000 independent posterior samples (the gray lines) to give a visual guide for the uncertainties. The inset shows $P(M_{\max })$, the posterior distribution of $M_{\max }$.

Figure 2

The EoS reconstruction results for the single-layer FFNN model of Han et al. [26] (the blue dash-dotted curve), the Gaussian process method mainly following Gorda et al. [25] (the orange dash-dotted line), and the flexible piecewise method of Jiang et al. [35] (the purple dashed curve) in the case of $n_{\mathrm{L}}=10n_{\mathrm{s}}$. In the plot the regions are for 90% confidence level (the vertical lines mark the regions of inferred $n_{\mathrm{c},\mathrm{TOV}}$).

Figure 3

Almost identical posterior distributions of $M_{\mathrm{TOV}}$ and $R_{\mathrm{TOV}}$ found in three independent EoS reconstruction approaches for $n_{\mathrm{L}}=10n_{\mathrm{s}}$, where we have taken into account the mass/radius data of four NSs (including PSR $\mathrm{J}0030+0451$, PSR $\mathrm{J}0740+6620$ and the two involved in GW170817), the posterior distribution of the maximum cutoff of the neutron star mass function, as well as the theoretical $\chi \mathrm{EFT}$ and pQCD constraints. The orange, purple and blue lines represent the results for the GP, PWLS, and FFNN models, while the black gives their average.

Figure 4

$M_{\mathrm{TOV}}$ obtained with different data sets and/or conditions. The gray dash-dotted line represents $P(M_{\max })$. The black dotted line represents that deduced from the $\chi \mathrm{EFT}$ constraint as well as the four NSs with both mass and radius measurements. The blue dashed curve employs also the pQCD constraint. The red solid line further takes into account the information of the mass distribution of NSs. The last three lines are calculated with the FFNN approach.

Figure 5

Top panel: The posterior distributions of the $M_{\mathrm{TOV}}$ found in the three EoS reconstruction approaches in the case of $n_{\mathrm{L}}=n_{\mathrm{c},\mathrm{TOV}}$, which are moderately weaker than the case of $n_{\mathrm{L}}=10n_{\mathrm{s}}$ (i.e., the red line). Middle panel: The posterior distribution of $M_{\mathrm{TOV}}$ obtained in setting the pQCD constraint at a density uniformly distributed between ${\rho }_{\mathrm{TOV}}$ and $10n_{\mathrm{s}}$ (i.e., the random check). Bottom panel: The posterior distribution of $M_{\mathrm{TOV}}$ obtained in the GP model with different $\overline{l}$.

Figure 6

The top panel is the reconstructed $c_{\mathrm{s}}^2$ as a function of $n/n_{\mathrm{c},\mathrm{TOV}}$ for $n_{\mathrm{L}}=n_{\mathrm{c},\mathrm{TOV}}$, where the part of $n>n_{\mathrm{c},\mathrm{TOV}}$ is a simple extrapolation of the low density EoS. The left part of the middle panel is the same as that of the top panel, while the right part is derived by requiring the extrapolated EoS at higher densities (i.e., $\ge n_{\mathrm{c},\mathrm{TOV}}$) to also satisfy the pQCD likelihood at $10n_{\mathrm{s}}$ with three constant-sound-speed segment approximation (i.e., a specific PW approach). The bottom panel shows the constrained EoS with the pQCD likelihood at $n_{\mathrm{L}}=10n_{\mathrm{s}}$ for a comparison. The red dash-dotted curves and orange dashed curves denote the 50% and 68.3% intervals, respectively.

Figure 7

Similar to Fig. 1. The top row presents the outcomes derived from the two-component Gaussian mixture model, whereas the bottom row corresponds to the results utilizing the Gaussian+Cauchy-Lorentz model. The left column excludes the data from GW170817, PSR $\mathrm{J}0030+0451$, and PSR $\mathrm{J}0740+6620$, whereas the right column presents the results obtained after omitting some high-mass NSs.

Figure 8

Posterior distribution of mass distribution parameters derived from various tests. The parameters ${\mu }_1$ and ${\sigma }_1$ represent the location and scale, respectively, for the first component, whereas ${\mu }_2$ and ${\sigma }_2$ correspond to those for the second component. The variable $r$ denotes the weight attributed to the first component.

Figure 9

Posterior distributions of $M_{\mathrm{TOV}}$ derived from different $P(M_{\max })$ results, based on the GP EoS reconstruction approach in the case of $n_{\mathrm{L}}=10n_s$.