Studying the parameters of the extended $\sigma$-$\omega$ model for neutron star matter

In this work we study the parameters of the extended $\sigma$-$\omega$ model for neutron star matter by a Bayesian analysis on state-of-the-art multi-messenger astronomy observations, namely mass, radius and tidal deformabilities. We have considered three parameters of the model, the Landau mass $m_L$, the nuclear compressibility $K_0$, and the value of the symmetry energy $S_0$, all at saturation density $n_0$. As a result, we are able to estimate the best values of the Landau mass of $m_L \approx 0.73$ GeV, whereas the values of $K_0$ and $S_0$ fall within already known empirical values. Furthermore, for neutron stars we find the most probable value of 13 km $<R_{1.4}<$ 13.5 km and the upper mass limit of $M_{max} \approx 2.2$ M$_{\odot}$.


Introduction
During the recent years the uncertainty on the determination of the equation of state (EoS) of neutron stars has been reduced by both laboratory experiments and by multi-messenger astronomy observations. In the neutron star interior matter is as dense as several times saturation density n 0 , the mean density of heavy atomic nuclei. Studying this broad range of densities represents a great challenge for both experiments and theories. The properties of nuclear matter around saturation density have been measured in different experiments like is the case of relativistic heavy ion collisions, giant dipole resonances of nuclei, parity violating weak neutral interaction scattering on heavy nuclei, among others [1,2]. Higher densities are hard to study in the laboratory particularly in the case of cold and neutron rich nuclear matter inside compact stars. On the contrary, atomic nuclei are mostly symmetric in their proton and neutron content. Observations of neutron stars may provide estimates of their masses and radii as well as other properties like their deformation, which is intimately related to gravitational wave emissions, see [3] for an overview. From the theoretical front, nuclear matter at high densities harbors the possibility that neutron Send offprint requests to: arXiv:2006.03676v1 [astro-ph.HE] 5 Jun 2020 star interiors bear exotic content beyond protons and neutrons, like hyperons, meson condensates or deconfined quark matter, resulting in different predictions for mass and radius relations. It is in this high density part where lies the discrepancy between the many EoS models, with many extrapolating from their low density results. In our study, we consider canonical neutron stars composed of protons and neutrons as well as leptons that account for charge neutrality. As a theoretical tool, the extended σ-ω relativistic model [4] is solved in its mean field version. We rely on the simplicity of this EoS model in order to obtain a robust framework where parameter changes can control the high density EoS extrapolation thus able to cover regions of interest in the mass-radius diagram for neutron stars. We choose three parameters to characterise the neutron star EoS which we estimate through a Bayesian analysis based on astronomical observations: the mass measurements of the PSR J0740+6620 [5] and PSR J0348+0432 [6], the tidal deformability estimates from the gravitational wave event GW170817 [7,8] and maximum mass boundary from its electromagnetic signal counterpart [9], as well as the X-ray measurement of the PSR J0030+0451 by the NICER detector [10,11].
The EoS parameters we consider here are the compressibility of nuclear matter K 0 , and the nuclear symmetry S 0 , and the Landau mass m L , all evaluated at saturation density. These first two quantities have already been measured in the laboratory with resulting values of about K 0 ≈ 240 ± 20 MeV [12,13] and S 0 = 31.7 ± 3.2 MeV [14,15]. We shall compare these empirical values with the ones resulting from our Bayesian analysis based on neutron star observations. Moreover, the Landau mass is an effective mass resulting from the interaction of nucleons in the dense medium. In a primer study [16] we have estimated its value to be of m L ≈ 750±15 MeV where we varied only one parameter keeping the others fixed with values of n 0 = 0.156 fm −3 , K 0 = 240 MeV and S 0 = 32.5 MeV was performed. Consequently, a generalization of that method to three parameters is presented in this work.
The compressibility of nuclear matter K has a dominant effect on the stiffness of the EoS which is translated into the maximum neutron star mass M max . While the symmetry energy S strongly affects the size of neutron stars, it also fully determines the proton fraction in their interiors. Cooling of neutron stars may undergo the fast DUrca cooling [17] if the proton fraction goes above some threshold value. From population observations, this fast cooling is not favoured [18]. Thus, the link between the symmetry energy and temperature observations of compact stars, see [19,20,21] and references therein. Most importantly, a correlation between the neutron star radius and the tidal deformability has been found in [22,23].
Bayesian methods have been widely used for estimation of model parameters as well as for statistical inference. In the case of the neutron star EoS studies, we can briefly summarise a few works before and after the multi-messenger era, i.e. before and after the GW170817 detection when many electromagnetic counterparts complemented the gravitational wave observation. As for earlier works it is worth to mention the seminal work which performed an X-ray bursters Bayesian analysis [24,25]. Raithel et al. [25] introduced the popular multi-polytrope EoS parameterization (MPP) [26,27,28] where as [29] exploited realistic EoS instead. Moreover, Lackey and Wade [30] considered future detections of gravitation radiation from binary neutron stars populations with a four parameter polytrope formulation EoS concluding that the high density EoS can be more accurately determined. With the advent of the multi-messenger astronomy era many other Bayesian studies appeared, leading to compatible results. For instance [9], presents an exhaustive study using the MPP approach to cover the phase space of EoS parameters by imposing constraints on the lower bound of the maximum neutron star mass and tidal deformability values. Later on, the mass-radius measurement of PSR J0030+0451 by NICER would help to provide tighter constraints on the resulting EoS regions, see [31,32,33,34,35,36,37,38,39,40] where different considerations and constraints for example different Bayesian priors, have been taken into account. Most importantly, the EoS study parameters considered here have been studied in a few other Bayesian works that consider different EoS models, as is the case of Miller et al. [32] that use the MPP EoS parameterization to study the impact of laboratory measured the symmetry energy on the compact star EoS, below nuclear saturation densities.
The Bayesian method used in this work is a continuation of previous works [29,31] where the EoS described hybrid stars with hadron matter undergoing a phase transition into deconfined quark matter described both by a Maxwell construction or a mixed phase. A very important aspect shared in these Bayesian studies is the implementation of an unrestricted lower bound for the maximum neutron star mass, a quantity that keeps changing due to the observational updates, therefore influencing the posterior probabilities results [29]. The accompanying analysis [41] to our study here a) focuses on EoS models that show a third branch in the M -R diagram [42,43,44] and, b) uses the same observational constraints and Bayesian methodology used here. In addition, a more detailed overview of the state-of-the-art neutron star Bayesian studies can be found there. Alternatively, an EoS study using deep neural networks has been presented in [45] bearing similar results to most Bayesian analyses.
This paper is organized as follows. In section 2 we introduce the extended σ-ω model for the neutron star EoS. In section 3 we briefly describe the methods used to compute neutron star parameters followed by section 4 where the Bayesian methods together with their corresponding observational inputs are presented. We present our results and conclude with a summary and outlook in the last two sections.

The σ-ω model equation of state for neutron star matter
In this study we employed an extended version of the σ-ω model to describe the interior of the neutron star [4]. The Lagrangian of the model is given by, where Ψ = (Ψ n , Ψ p ) is the vector of proton and neutron fields, m N m σ m ω are the σ and ω meson masses and g σ , g ω , and g ρ are the Yukawa couplings corresponding to the σ-nucleon, ω-nucleon and ρ-nucleon interactions respectively. The kinetic terms corresponding to the ω and ρ meson can be written as: In eq. (1) U i (σ) is a self interaction term for the σ-meson and it has the following form: The nuclear interaction is approximated by the meson fields introduced in (1). The ω meson describes the repulsive nature of the interaction while the σ meson is responsible for the attraction between the nucleons. The ρ meson takes into account the difference between the proton and the neutron. To describe the neutron rich matter present inside compact stars the protons, electrons and neutrons are considered to be in β-equilibrium: This assumption provides a connection between chemical potential of the nucleons and the electrons: In this study, we consider the σ-ω model in the mean field approximation. This is a reasonable assumption because quantum fluctuations of nuclear matter do not influence significantly the observable properties of compact stars considering the sensitivity of the currently available measurement techniques [46,47]. Since the temperature of neutron stars is negligible compared to the energy scale corresponding to nuclear matter it is reasonable to assume while calculating the thermodynamics of the system that the nuclear matter is at zero temperature [48,49]. Using these assumptions eq. (1) simplifies considerably. All of the kinetic terms become zero and only the following components of the mesons have non-zero values: ω 0 = ω and ρ 3 0 = ρ. With these assumptions the free energy of the model can be calculated as it is described in for example in Ref. [50]: where µ p , µ n and µ e are the proton, neutron, and electron chemical potential respectively. The f F function gives the free energy contribution corresponding to one fermionic degree of freedom it is given as: where E 2 k = k 2 + m 2 . As described above in the case of cold nuclear matter of neutron stars one has to take T → 0 in eq. (7). This means that that the fermionic free energy has only two variables f F (m, µ).
The couplings in eq. (7) are determined by nuclear saturation data [48,51]. The values of the nuclear parameters can be found in Table 1. The definition of the Landau where k = k F is the Fermi-momentum and E k is the dispersion relation corresponding to the nucleons. The Landau mass is not independent of the effective mass in the mean field approximation. Substituting the energy of nucleons from eq. (7) into the definition of the Landau mass yields the above relation: This connection makes it impossible to fit the value of the Landau mass and the effective mass simultaneously [51]. One of them has to be chosen as a free parameter and only after fixing its value the other one becomes completely determined. The compression modulus is defined as usual for example in Refs. [49,48]: The asymmetry energy is defined as it is described for example in [48]: Here t = nn−np n B is the relative difference between the number of neutrons and protons. To describe the crust of the neutron star we complemented the extended σ-ω model with a low density EoS. The two EoS are joined at the point where both models hold the same pressure. For the crust EoS we used the widely known BPS equation of state [52].
In this study our aim is to use astrophysical data to estimate the most relevant parameters of the model, some of them already measured in the laboratory. Earlier works suggest that from the nuclear parameters listed in Table 1 whereas Table 2 shows parameter values for other EoS models for the sake of comparison and shows the parameter similarity of our model with other relativistic field models like DD. The Landau mass and the compression modulus has the largest influence on neutron star observables [4,53]. To provide data for the Bayesian-analysis we fitted the extended σ-ω model using different values for the above mentioned two parameters at saturation density, K 0 ≡ K(n 0 ) and S 0 ≡ S(n 0 ) ≡ a sym (n 0 ). Using their values in table 1 as reference points we varied the given parameter between 80 % and 110 % of these values. Thus, we have chosen 10 points in the above mentioned range for each parameter and ended up with 1000 EoS models after calculating the model for every combination.

Neutron Star Sequences
Predictions for the properties of static neutron stars from the extended σ-ω EoS sets are obtained by solving the Tolman-Oppenheimer-Volkoff equations [54,55]: dm(r) dr = 4πr 2 ε(r).
from which the mass is derived. The radial functions ε(r) and p(r) represent the energy density and pressure in the interior of the star. The boundary conditions for the system are zero pressure at the star surface p(r = R) = 0 as well as null mass at the origin m(r = 0) = 0. The total mass M is defined as M = m(r = R). To determine the mass and radius of a single star it is necessary to specify its central density ε c at the origin. To obtain the entire sequence of neutron stars displayed in the corresponding M -R diagram, ε c is increased starting from a low density value around nuclear saturation. In addition, the dimensionless tidal deformability of neutron stars Λ is computed following the approach introduced in [56]. It is related to the Love number k 2 : whose computation involves a perturbation treatment of the spherical metric of the neutron star, see [57,58,59,60] for details.

Bayesian inference Formalism
Let's introduce a notation for a vector of parameters, so that each vector represents one EoS from considered EoS model class: where N 2 and N 3 denote the cardinalities of te ordered sets of values of the parameters m L , K 0 and S 0 respectively. The full likelihood for the given − → π q can be calculated as a product of all likelihoods, since the considered constraints are independent of each other where w is an index of the constraints. The posterior distribution is given by Bayes theorem where P ( − → π q ) is a prior distribution of a models taken to be uniform: P ( − → π q ) = 1/N .

Likelihood of a model for the GW170817 (Λ 1 -Λ 2 constraint)
In order to implement the tidal deformability constraint on the compact star EoS, reflected on the Λ 1 -Λ 2 diagram that includes probability regions from GW170817 event [7,8], we employ the following formula for the likelihood: where l are the length of the line at Λ 1 -Λ 2 diagram produced by − → π q , and n c is the central density of a star. β(Λ 1 , Λ 2 ) is the Probability Distribution Function (PDF) that has been reconstructed (as previously in [31]) by the method of Gaussian kernel density estimation with Λ 1 -Λ 2 data given at the LIGO web-page [61].

Likelihood of a model for the maximum mass (maximum mass constraint)
The likelihood of the model is a conditional probability of the maximum mass for the given model under the considered maximum mass constraints. Here we have considered the recent estimation of the mass of the heaviest known pulsar PSR J0740+6620 given by the measurements presented in [5] 2.14 +0.10 −0.09 M . Since it has relatively large error deviation, the mass estimation of PSR J0348+0432 2.01 +0.04 −0.04 M [6] has been considered as well. The last one is lower, but it has a tiny standard deviation of the error, so that it helps to constraints the maximum mass from below starting from the 3σ region. Additionally, the constraint for the upper limit on the maximum mass 2.16 +0. 17 −0.15 M [9] has been added to the analysis. Note, that all those numbers are given at the 68.3% confidence level. The likelihood for the maximum mass constraint has been introduced in the following form: here M q is maximum mass of the given by π q , and Φ(M, µ, σ) and N (M, µ, σ) are Comulative Distribution Function (CDF) and Probability Distribution Function (PDF) for normal distribution respectively. The mass measurement for PSR J0740+6620 has been approximated with (µ C = 2.14, σ C = 0.095). For PSR J0348+0432 the original mean µ A = 2.01 and deviation σ A = 0.04 have been used. The upper limit constraint has been approximated by (µ U = 2.16, σ U = 0.16).

Likelihood of a model for mass and radius (M -R constraint)
The results of the Neutron Star Interior Composition Explorer (NICER) observation of PSR J0030+0451 have been recently reported in a collection of publications, see for instance Ref. [10,11]. There were two estimates of the mass and equatorial radius based on mutually exclusive assumptions about the uniform-temperature emitting spots. The one radius and mass estimates are M 1 = 1.44 +0.15 −0.14 M and R 1 = 13.02 +1.24 −1.06 km [10] and the second estimates are M 2 = 1.34 +0.15 −0.16 M and R 2 = 12.71 +1.14 −1.19 km [11]. Here we introduced approximations of those estimates by a bivariate normal distribution PDF exactly as in [41]. Assuming, those two estimates are equiprobable, since they are mutually exclusive, a likelihood for M -R constraint has been introducing as follows where (µ M , µ R ) is a mathematical expectation at M -R plot, σ M and σ R are standard deviation for mass and radius respectively, and ρ is an angle of an ellipse rotation Fig. 1, which corresponds to correlation between mass and radius, here l is a length of the curve on M -R plot for particular EoS and n c , as mentioned above, is a central density of a star. The approximations of [10] and [11] by bivariate normal distribution are

Marginalization
For the convenience of analysis of the Bayesian inference results, the marginalisation procedure has been used. In order to show posterior distribution of one parameter marginalization over other parameters is required: For a 2D posterior distribution is obtained by marginalization over one parameter: Fig. 1 shows the input constraints for our Bayesian Analysis in the M − R diagram of neutron stars where also one of the most probable EoS sequence is displayed for completeness.

Results
In  Fig. 2, where several plots that relate mass, radius, central density and tidal deformability for the entire set of EoS models are included. As it can be seen in the M -R diagram of Fig. 2, we find that the resulting neutron star sequences cover radius values of between 11 km and 14.5 km, for masses above half a solar mass and a maximum mass value of about 2.7 M for the stiffest EoS. It can also be seen that matter can be compressed up to 6 times saturation density for the softer EoS which are unable to support the observed 2 M neutron star. In all the plots but the Λ 1 − Λ 2 diagram, sequences of neutron stars appear to be grouped in sets being separated by gap regions. The sequences in these groups share the same values of m L , an observation that points at the influence of the Landau mass over the other two parameters, especially for the highest masses. The black line sequences represent the top 10 EoS models with the highest probabilities with respect to the rest which has been greyed out. Fig. 3 shows the posterior probabilities resulting from our Bayesian Mass-radius diagram for neutron stars. Shaded regions correspond to the different measurements used as inputs for the Bayesian analysis. The two horizontal bands correspond to the pulsar masses M=2.14 +0.10 −0.09 M of PSR J0740+6620 [5], and M=2.01 +0.04 −0.04 M of PSR J0348+0432 [6]. The grey and brown regions above and below the mass value of M=1.365 M correspond to the mass-radius estimates for the two compact stars of the merger in GW170817 [7,8]. The elliptical regions denoted by dashes lines correspond to the mass and radius measurement of the pulsar PSR J0030+0451: M=1.44 +0. 15 −0.14 M with R=13.02 +1.24 −1.06 km [10] or M=1.34 +0.15 −0.16 M and R=12.71 +1.14 −1.19 km [11]. The forbidden upper region above M > 2.33M has been derived from GW170817 together with the associated kilonova AT2017gfo in Ref. [9].
analysis. The upper figures show the marginalised probabilities for each of our EoS parameters under study (see Eq. 22) whereas the lower LEGO plots represent the 2D posterior distributions for each pair of parameters that have been marginalized over a third one. These probability regions correspond to the two parameters K 0 and S 0 with marginalized m L . The effects of the EoS parameters on the maximum mass have been found to follow a variation within one order as follows [4,62] ∆M max (δm L ) Most importantly, in order to corroborate the major effect of m L over the remaining parameters on neutron stars we have reduced our Bayesian study to two EoS parameters, keeping the symmetry energy fixed. Therefore, we vary the nuclear compressibility and the Landau mass since these quantities significantly contribute to the stiffness of the EoS whereas the symmetry energy has a major influence on the stellar radius. We present the setup and results in Fig. 4 where the M -R diagram clearly shows the subsets of sequences with the same m L value. The posterior probabilities are displayed in the lower plots which indicate the same value for the most probable m L ≈ 0.73 GeV which coincides with the previous result from the three EoS parameter study as well as being compatible with the one parameter study in Ref. [16].

Summary and conclusions
In this work we have performed a two and three dimensional Bayesian study for parameter estimation of the extended σ-ω EoS model of neutron stars. The chosen EoS parameters under study are the compressibility of nuclear matter K 0 , the symmetry energy S 0 and the Landau mass m L , properties of saturation density dense matter. These parameters where varied within an acceptable range of empirical values. It is well known that compressibility of nuclear matter dominates the stiffness of neutron star matter therefore significantly contributes to the determination of the maximum neutron star mass M max . On the contrary, the symmetry energy has a greater influence on the determination of the neutron star radius R [63,64,65] rather than on M max . For our analysis we have used as observational inputs for our Bayesian study the multi-messenger astronomy measurements: mass measurements of PSR J0740+6620 [5] and PSR J0348+0432 [6], tidal deformability data from GW170817 [7,8], M max boundaries [9], and combined mass-radius measurements for PSR J0030+0451 [10,11]. An important element in our Bayesian analysis is that the choice of the corresponding prior distributions for these astrophysical measurements have been taken as uniform. We have observed that pulsar maximum mass measurements favour any stiff EoS models with a high maximum mass value above 2 M . On the other hand, gravitational wave data both in the form of neutron star tidal deformabilities as well as the excluded high mass region do not favour the stiffest EoS models in our sample. Thus, the interplay existence between these constraints.
The Landau mass m L is the best determined parameter whereas the remaining two have probabilities distributions that peak at the lowest values we have considered. Therefore, the values of K 0 and S 0 found by our analysis have been compared with empirical values found in the laboratory. We find that probability distributions do not exactly peak at their measured value, however they lie within the one sigma confidence region of their posterior distribution. Moreover, their influence on the neutron star properties results to be negligible with respect to m L . Our one and two parameter study corroborates that it is reliable to set K 0 and S 0 within their empirical values, therefore entrusting the impact on neutron stars to the Landau mass m L whose values are inversely proportional to M max . Moreover, the Bayesian study of [66] considers an inversion of the TOV equations in order to study a set of nuclear parameters derived from an expansion of the EoS functional under three different scenarios for DUrca cooling activation. For the first set DUrca is activated for masses below 2 M , for their intermediate set DUrca holds within 1.8 < M/M < 2, and within 1.8 < M/M < 2 for the last one. They corroborate the hypothesis of a universal contribution from the symmetry energy under the DUrca constraint [67] and find the corresponding results for these three categories: a) K 0 = 232. All in all, we find that for the most probable neutron star sequences the values of 13 km < R 1.4 < 13.5 km and M max ≈ 2.2 M for the stiffest EoS models. The latter quantity is in agreement with estimations on the bound for the maximum neutron