Approximate Universal Relations among Tidal Parameters for Neutron Star Binaries
Abstract
One of largest uncertainties in nuclear physics is the relation between the pressure and density of supranuclear matter: the equation of state. Some of this uncertainty may be removed through future gravitational wave observations of neutron star binaries by extracting the tidal deformabilities (or Love numbers) of neutron stars, a novel way to probe nuclear physics in the highdensity regime. Previous studies have shown that only a certain combination of the individual (quadrupolar) deformabilities of each body (the socalled chirp tidal deformability) can be measured with secondgeneration, gravitational wave interferometers, such as Adv. LIGO, due to correlations between the individual deformabilities. To overcome this, we search for approximately universal (i.e. approximately equationofstate independent) relations between two combinations of the individual tidal deformabilities, such that once one of them has been measured, the other can be automatically obtained and the individual ones decoupled through these relations. We find an approximately universal relation between the symmetric and the antisymmetric combination of the individual tidal deformabilities that is equationofstateinsensitive to for binaries with masses less than . We show that these relations can be used to eliminate a combination of the tidal parameters from the list of model parameters, thus breaking degeneracies and improving the accuracy in parameter estimation. A simple (Fisher) study shows that the universal binary Love relations can improve the accuracy in the extraction of the symmetric combination of tidal parameters by as much as an order of magnitude, making the overall accuracy in the extraction of this parameter slightly better than that of the chirp tidal deformability. These new universal relations and the improved measurement accuracy on tidal parameters not only are important to astrophysics and nuclear physics, but also impact our ability to probe extreme gravity with gravitational waves and cosmology.
1 Introduction
The equation of state (EoS), the thermodynamic relation between state variables, is key in the description of fluids and solids, but in particular, it is critical in the description of stars. The barotropic equation of state of matter, the relation between pressure and density, at nuclear saturation density () has been wellconstrained by terrestrial experiments [1]. For example, heavyion collisions [2, 3, 4] and measurements of the neutron skins of nuclei [5] have constrained the linear density dependence of the nuclear symmetry energy. However, terrestrial experiments cannot constrain the EoS beyond saturation density, rendering the supranuclear EoS one of the largest uncertainties in nuclear physics.
Neutron stars (NSs) are a perfect testbed to probe nuclear physics in this highdensity regime. For example, one can use independent measurements of the NS mass and radius to constrain the EoS [6, 7, 8, 9]. Current observations of Xray bursters and quiescent lowmass Xray binaries have already placed some constraints on the EoS [10, 11, 12, 13, 14, 15, 16], though these may suffer from large systematic errors due to uncertainties in the astrophysical modeling of the NS sources. Future Xray pulse profile observations from a hot spot on the NS surface using NICER [17] and LOFT [18, 19] may be able to place stronger constraints with less systematics [20, 21, 22, 23, 24]. Future radio observations of the double binary pulsar J07373039 [25, 26, 27] are expected to measure the moment of inertia of the primary pulsar [28, 29], which would also allow for constraints on the EoS.
Given the recent direct detection of gravitational waves (GWs) from black hole binaries [30, 31], one expects GWs from NS binaries to also be detected soon. The latter can be a novel probe of nuclear physics because they encode information about the EoS through finitesize effects. During the early inspiral, when the orbital separation is large, the NSs barely feel the gravitational field of their companion, and thus, are not really deformed. As the inspiral proceeds and the orbital separation decreases, the gravitational tidal field of one star at the location of its companion increases in magnitude, creating a deformation on the latter (and vice versa). Such a deformation forces the gravitational field of the deformed stars to not be spherically symmetric any longer, which affects their orbital trajectory. Moreover, since the quadrupolar tidal deformation is time dependent, it also modifies the amount of GW energy carried away from the binary, and thus, the orbital decay rate. Such changes in the trajectory, in turn, imprint directly onto the gravitational waveform. The largest finitesize effect in the waveform is characterized by the electrictype, quadrupolar tidal deformability [32, 33, 34], which we here refer to as the tidal deformability , and which is related to the tidal Love number [35]. One can think of this tidal deformability as a parameter that quantifies the linear response of a NS’s quadrupole moment to the gravitational tidal field of its companion.
Measurability of the NS tidal effects and constraints on the NS EoS with GW observations have been studied in [36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49]. For example, Ref. [39] carried out an approximate Fisher analysis and found that a certain combination of the individual tidal deformabilities of the binary components can be measured to better than with secondgeneration GW interferometers, such as Adv. LIGO. Bayesian studies that followed [46] showed that the certain combination of the tidal deformabilities may be measurable to the accuracy mentioned above, for sufficiently high signaltonoise ratio events, although systematic error may also add to the error budget (see also [44, 45] for studies on systematic errors). References [42, 48] also showed that the detection of multiple NS binary signals could enhance even further the accuracy to which certain combinations of the tidal deformabilities can be measured.
The measurement of the tidal deformabilities can also be useful in probing cosmology [50, 51, 52] by using GW sources as standard sirens [53, 54, 55, 56, 57, 58, 59, 60, 61, 62]. The main idea here is to measure the luminosity distance (from the GW amplitude) and the redshift of the source independently from a set of GW and electromagneticwave observations; since these two quantities are related by cosmological evolution equations, they can be used to estimate cosmological parameters. The original idea was to obtain the redshift from host galaxy identification, provided the GW skylocalization is accurate enough. Typically, however, this is not expected to be the case for a large number of galaxies with secondgeneration detectors. Another idea is to obtain the redshift from precise measurements of the tidal deformabilities [50, 51, 52]. These measurements allow us to construct the intrinsic mass of the source, provided one knows the correct EoS a priori; comparing this mass to the observed (redshifted) mass, one can then estimate the redshift [50, 51, 52]^{1}^{1}1See e.g. [63, 64, 65, 66, 67, 68, 69, 70] for other possibilities of probing cosmology with GW observations alone..
Can one independently extract the individual dimensionless tidal deformabilities of each binary component with an Adv. LIGO observation? Previous literature claims this is not possible with the sources expected to be detected with Adv. LIGO because and are strongly correlated. One can, however, reparameterize the templates with two new tidal parameters, constructed from independent combinations of and , in such a way as to diminish the correlations. For example, the socalled and parameters partially break the degeneracies between and because they enter at different postNewtonian (PN) orders in the waveform^{2}^{2}2The PN expansion is a series in powers of the ratio of the orbital velocity of the binary to the speed of light . A term of PN order in a PN series corresponds to one that is of , relative to the leading order term in the series. [37, 38, 44, 46, 47]. Nevertheless, one can only measure in this parametrization, because enters at too high a PN order to affect the template’s phase sufficiently [46]. Since is the dominant tidal parameter entering the waveform, we shall call it the dimensionless chirp tidal deformability, in analogy with the chirp mass being the dominant mass parameter in the waveform.
In this paper, we search for a way to overcome this problem, i.e. a way to extract both tidal deformabilities of the binary from an Adv. LIGO GW observation. We accomplish this by finding approximately universal relations, i.e. relations that are approximately EOSinsensitive, among independent tidal parameters. In general, these relations serve two important purposes:

To Improve Parameter Estimation. The universal relations can be used to eliminate some of the tidal parameters from the template parameter vector in a data analysis study. This breaks degeneracies between the parameters removed and those left in the template and improves the latter’s parameter estimation accuracy.

To Extract Both Tidal Deformabilities. Given the measurement of a given combination of tidal parameters (e.g. ), these relations allow us to automatically obtain the other combination (e.g. ) to the accuracy of the approximate universality. These two independent combinations can then be easily decoupled to find and .
Universal relations in NSs are not really new (see the recent review [71] and references therein for various NS universal relations). In particular, universal relations among the NS moment of inertia, the tidal deformability and the quadrupole moment (the socalled ILoveQ relations) [72, 73] can be used to eliminate the individual quadrupole moments from the template parameter vector. This breaks degeneracies between the quadrupole moment and the individual NS spins, allowing us to extract the latter more accurately. Such relations also allow us to probe extreme gravity without having to know what the correct EoS is. For example, by combining future measurements of the moment of inertia with the double binary pulsar and the tidal deformability with GW observations, one can use the ILove relation to place constraints on the parityviolating sector of gravity that are six orders of magnitude stronger than the current bound [72, 73]. Universal relations among the various tidal deformabilities of an isolated NS (multipole Love relations) were found in [74], while relations among these deformabilities and certain binary parameters were studied in [75].
1.1 Executive Summary
In this paper, we find three different, approximately universal, binary Love relations, among the following quantities:

the symmetric tidal parameter (the average sum of and ) and the antisymmetric tidal parameter (the average difference of and ),

The and tidal parameters mentioned earlier,

The coefficients in the Taylor expansion of the dimensionless tidal deformabilities around a fiducial mass [42].
One can implement the first and second parameterizations for any NS binary systems, while the third parameterization can only be applied, in practice, to systems whose NS masses are close to . Otherwise, the systematic error on the leading tidal coefficient due to mismodeling the tidal deformability can dominate the statistical one. If one insists on retaining many terms in the Taylor expansion to reduce such systematic error, then correlations among the many parameters may increase the statistical error. If enough systems with masses that cluster together are observed, however, the third parameterization has an advantage that it allows one to combine such observations to increase the measurement accuracy of that is common to all systems [42, 48].
We find these new universal relations as follows. We first calculate the tidal deformability of an isolated NS with various realistic EoSs for a set of different NS masses. In doing this, we follow [32] and extract the tidal deformability from the asymptotic behavior of the gravitational potential of tidallydeformed NS solutions, treating the tidal deformations perturbatively. When considering the first and second relations, we then calculate the two independent tidal parameters in each relation, such as and , or and . Given these, we finally rewrite one of the tidal parameters in terms of the other and in terms of the mass ratio (with ), e.g. or , and study its EoS variation. When considering the third relation, we calculate the th coefficients in the Taylor expansion of about for the dimensionless tidal deformability of an isolated NS . We then eliminate using the zerothcoefficient in the expansion, , so as to find , and then we study its EoS variation.
Figure 1 shows the – relation, namely as a function of and , for various EoSs. The plane in the figure is simply a fit to all the numerical data. Observe that all the data points lie approximately on this approximately universal plane irrespective of the EoSs. The relative fractional difference between any data point and the plane never exceeds for NS masses less than . Similar figures can be constructed for the – and the – relations, although the universality for the former is weaker compared to the – universality while the latter does not depend on . The – relation is particularly useful when studying the possibility of using tidal deformabilities to probe cosmology with GWs, as the knowledge of the correct EoS, or the mass dependence of the tidal deformability, is crucial in such an analysis. These results show that once one measures one of the tidal parameters in the relations, one can automatically obtain the other. This in turn allows one to determine the tidal deformability of each body independently, allowing us to extract more astrophysical information beyond just the masses and spins of the binary.
The universality seems to become stronger as one decreases the mass ratio for fixed total mass, but it deteriorates as one increases the total mass for a fixed mass ratio. The former can be easily understood by considering the – relation in the limit, since then and this is exactly EoS universal. The latter, however, is somewhat more difficult to explain, because one intuitively expects the universal relations to improve as one approaches the black hole (BH) limit, i.e. as the NS mass, and thus the NS compactness, increase. Although it is true that and approach their BH values in this limit, their averaged difference is sensitive to how this limit is approached, i.e. to the slope of the – relation. This slope decreases with increasing mass, essentially because decreases very fast with stellar compactness, which is simply because massive stars deform less than light stars. Recalling that the relative fractional difference is a function of the difference in with different EoSs divided by for a reference EoS, the fact that the universality deteriorates in the BH limit is then simply a consequence of the slope of the – relation decreasing faster than the EoS variability in .
Figure 2 shows a measure of the improvement in parameter accuracy of when one uses the – relation. More precisely, we carried out two Fisher analyses: one with a template family that includes both tidal parameters (red dashed curve), and one with a template family that uses the approximately universal relations to eliminate one of the tidal parameters from the parameter vector (red solid curve), in both cases assuming a fiducial realistic NS EoS. Observe that parameter estimation accuracy improves by roughly an order of magnitude when using the universal relations. This is because the relation allows us to eliminate from the search parameters, which breaks the degeneracy between and . We have checked that the fractional systematic error on , due to the EoS variation in the – relation, is smaller than , and thus negligible in this figure. This is mainly because the term in the GW phase that depends on is proportional to the difference in the masses of the two NSs, and hence, the systematic error is suppressed in the near equalmass regime, precisely where the binary Love relations are less EoSuniversal. On the other hand, we find that the – relation does not really affect the measurement accuracy of the chirp tidal deformability (blue dotdashed curve). This is because the correlation between and is weaker than that between and , since and enter at different PN order in the waveform phase. We also found that the measurement accuracy of improves by a few times when one uses the – relations.
These findings have an impact on various branches of physics. On the astrophysics front, the individual tidal deformability of each NS brings important additional information, in addition to the mass and spin of NSs. On the nuclear physics front, the improved measurement accuracy of tidal parameters increases our ability to constrain the EoS with future GW observations. On the experimental relativity front, such an improvement in the measurement accuracy of the tidal deformability strengthens projected constraints on alternative theories of gravity through universal relations between e.g. the tidal deformability and the moment of inertia. On the cosmology front, the – relation improves our knowledge of how the tidal deformability depends on the intrinsic NS mass, which helps breaking the degeneracy between the mass and redshift, allowing us to probe cosmology with GW observations alone.
Some of the results discussed above were already presented in the letter [78], which we here explain in much more detail and extend by finding new results. For example, the binary Love relation between and was not shown in [78]. We also carry out analytic calculations in the Newtonian limit to obtain a better understanding of the universality and to create a base function for fits. Furthermore, we estimate not only statistical errors but also systematic errors on tidal parameters due to (i) the EoS variation in the binary Love relations and (ii) not including one of the tidal parameters in the search parameter set that can be eliminated with the binary Love relations. We also extend the binary Love relations for NSs to quark stars (QSs) and to ILoveQ relations for a binary system (binary ILoveQ relations).
1.2 Organization
The remainder of this paper presents the details of the results discussed above and it is organized as follows. In Sec. 2, we derive universal relations among two (or more) independent tidal parameters in the gravitational waveform of NS binaries. In Sec. 3, we carry out a parameter estimation study using the Fisher matrix and estimate the measurement accuracy of tidal parameters with Adv. LIGO. In Sec. 4, we discuss the possibility of applying the universal relations to probe astrophysics, nuclear physics, experimental relativity and cosmology. We end in Sec. 5 by presenting possible avenues for future work. We use geometric unit throughout.
2 Universal Binary Love Relations
In this section, we focus on relations among various tidal parameters that enter the gravitational waveform of a nonspinning NS/NS binary with masses , where the subscript denotes the th body. We assume throughout the paper. We adopt eleven realistic EoSs that can support a NS with a mass above : AP3 [76], AP4 [76], SLy [79], WFF1 [80], WFF2 [80], ENG [81], MPA1 [82], MS1 [83], MS1b [83], LS [84] with nuclear incompressibility of 220MeV (LS220) and Shen [85, 86], where for the latter two, we assume a temperature of 0.1MeV with a neutrinoless and betaequilibrium condition.
The tidal parameter that impacts the waveform the most is the electrictype, tidal deformability (or just tidal deformability for short), which we denote as ^{3}^{3}3For rotating NSs, tidal effects enter first through their quadrupole moment, which appear at lower PN order than in the gravitational waveform. The effect of the quadrupole, however, is suppressed for slowlyrotating NSs because it is proportional to the NS spin squared.. Let us consider a NS of mass that is immersed in the field of an external companion. The tidal deformability is then the linear response of the tidally induced quadrupole moment tensor of the NS due to the external quadrupolar tidal tensor induced by its companion [32, 33, 34]:
(1) 
These tensors are extracted from the asymptotic behavior of the gravitational potential, or the component of the metric of the NS, via
(2) 
where the NS is located at the origin and is a unit vector pointing from the NS to a field point.
One can calculate , and thus extract , by constructing tidally deformed NS solutions. One first constructs a background solution for a nonspinning, spherically symmetric NS and then adds a tidal deformation as a perturbation. One solves the background and perturbed Einstein equations in the interior region with a given EoS, imposing regularity at the stellar center. One then matches the interior solution to an analytic, exterior solution at the stellar surface to determine integration constants (see e.g. [32, 33, 73] for a detailed explanation of this procedure), which then determines modulo an overall constant. With this metric component at hand, one can then extract and by expanding far from the NS and reading off the coefficient of for the quadrupole moment and for the external tidal field, where stands for the symmetric and trace free operator. From these, one can easily extract from Eq. (1).
Finitesize effects encoded in the GWs emitted by NS binaries are then dominantly controlled by two tidal parameters, and . Two parameters enter the waveform because depends on the NS mass and NSs in binaries typically have different masses. For later convenience, we define the dimensionless tidal deformability [72, 73]. Higher (multipole) order tidal deformabilities also enter the waveform at higher PN order but they can be expressed in terms of through the multipole Love relations [74].
When constructing a GW template family for NS binaries, one must choose two independent tidal parameters to represent finitesize effects, but this choice is not unique as any function of would in principle be acceptable. We first consider two different sets [(,) [72, 73] in Sec. 2.1 and (,) [37, 38, 44, 46, 47] in Sec. 2.2] and look for universal relations, which as argued in the Introduction, will also depend on the mass ratio . We then consider a Taylor expansion of around a fiducial mass [50, 39, 42, 48] and look for EoSinsensitive relations among the coefficients in the expansion (Sec. 2.3).
2.1 – Relation
Let us first consider relations between two new dimensionless tidal parameters, and , defined by [72, 73]
(3) 
Notice that the symmetric tidal parameter corresponds to the mean dimensionless tidal deformability, while the antisymmetric tidal parameter is the mean difference between the two tidal deformabilities. For an equalmass system, they reduce to and respectively.
2.1.1 Newtonian Limit
Before studying the relativistic relation between and , let us investigate the relations in the socalled Newtonian limit. The latter, also known as the nonrelativistic limit, is nothing but the leadingorder expansion in stellar compactness , where and are the NS mass and radius respectively of the th star. One can carry out calculations analytically in such a limit for a polytropic EoS,
(4) 
where and correspond to the pressure and energy density respectively, is a constant and is the polytropic index. For such a polytrope, and are related to the compactness via [73]
(5) 
respectively, where the coefficient depends on the polytropic index . For an and polytrope, one finds and [73]. From Eq. (5) and the definition of and , one finds
(6) 
where the and correspond to and respectively. Since in the Newtonian limit , Eqs. (5) and (6) tell us that large values of or large values of also correspond to this limit. Solving for in terms of and then inserting this in Eq. (6) for , one finds the – relation in the Newtonian limit
(7) 
Observe that as vanishes in this case, but as then .
AP3  AP4  SLy  WFF1  WFF2  ENG  MPA1  MS1  MS1b  LS220  Shen  mean  
0.795  0.875  0.832  0.977  0.869  0.873  0.714  0.411  0.581  0.613  0.632  0.743 
Let us now investigate for values of that approximate more realistic EoSs in the Newtonian limit. The latter can be estimated by first constructing a sequence of NS solutions using a set of realistic EoSs with , and then computing the realistic – relation and fitting Eq. (7) to it. We show the best fit value of in Table 1, where notice that the mean value is . The top panel of Fig. 3 shows against with ’s that approximate the WFF1 (the softest EoS considered) and MS1 EoSs (the stiffest EoS considered), as well as the mean . Observe that is close to unity for , while it drops quickly to zero in the limit as explained earlier. The bottom panel shows the relative fractional difference between the WFF1 or the MS1 relations and the mean relation. Notice that the dependence is smaller in the region compared to that in the region. But even in the latter case, the Newtonian relation is approximately EoSinsensitive to an accuracy of .
One may think that the EoSuniversality should be exact in the limit, since then , but this contradicts the bottom panel of Fig. 3. Such a naive argument does not hold because Fig. 3 considers the relative fractional difference, (where is the mean of Table 1), and not the absolute difference. Although the numerator of this fraction vanishes in the limit, the denominator does as well. Thus, the fact that the fractional difference in does not vanish in the limit is a reflection of the ratio being finite in this limit. We can see this explicitly by expanding the relative fractional difference around :
(8) 
where we evaluated the second line with the data of Table 1 for the MS1 EoS. Notice that the correction is of and not of . Equation (2.1.1) agrees with the variation in the – relation in the Newtonian, limit shown in the bottom panel of Fig. 3.
0.743  0.743  1  1  

1  1  1/5  1/5  
0.07550  0.07319  0.1697  
0.8474  
10.45  
13.61  
7.941  
0.5658  
73.60  
0.9996  0.9939  0.9998  0.9971 
2.1.2 Relativistic Relations
We now turn our attention to the relativistic relation. We construct this as explained at the beginning of Sec. 2, without expanding in powers of compactness. Figure 1 presents against and from two different viewpoints for eleven different realistic EoSs with . Inspired by the Newtonian relation in Eq. (7), we created a fit for the relation given by
(9) 
with , and . The fitted coefficients , , and are given in the second column of Table 2. Using the Newtonian relation of Eq. (7) as a controlling factor in the fit has one big advantage: by construction, the relation reduces to (i) this equation in the Newtonian limit when , (ii) with and (iii) with (as ). The fit corresponds to a double (or bivariate) expansion in and , where we expand in this power of because in the Newtonian limit. Namely, we simultaneously expand asymptotically in small mass ratios and in the Newtonian region of the 3D plot in Fig. 1. Although the fit becomes less accurate as and , the fit is reasonably good as and . This fit is depicted as the blue approximately invariant plane of Fig. 1. Observe that the binary Love relation for each EoS lies nicely on this plane, which shows that the relation is quite EoSinsensitive. Let us quantify this statement further.
The top left panel of Fig. 4 shows the – relation for NSs with , 0.75 and 0.9 (constant hypersurfaces of Fig. 1), where correspond to the mass ratio of the currentlyknown double neutron star binary J0453+1559 [87], while roughly corresponds to that of J07373039 (the double binary pulsar) [25, 26, 27] and J17562251 [88]. The relation is given by a single curve for a fixed and a given EoS (see Fig. 1 of [78]), with the single parameter along the curve being the mass, or equivalently, the compactness. The red shaded region for each corresponds to the parameter space spanned by the relation with different EoSs. A smaller area corresponds to a stronger universality. We also show the fit of Eq. (9), as well as the Newtonian relation of Eq. (7) with the mean . The bottom panel of this figure shows the relative fractional difference between the numerical data and the fit. We selected four representative realistic EoSs: WFF1, AP4, LS220 and MS1, which lead to stars with a radius of 10.4km, 11.4km, 12.6km and 14.9km respectively when the NS mass is 1.4. We only considered cases in which the smallest NS has a mass larger than and the largest NS has a mass smaller than the maximum allowed to have a stable configuration given the EoS. For example, the LS220 EoS predicts a maximum NS mass of , and thus, when there is only a very small subset of binary masses that satisfy and with which one can compute the – relation.
The left panel of Fig. 4 contains several features that ought to be discussed further. First, observe that the numerical data approaches the analytic Newtonian relation in the large limit. This makes sense since , and thus, as then and one approaches the Newtonian limit. Second, observe that when , the difference between the Newtonian and numerical relation is small, while it becomes larger as one increases . This is because the relation reduces to in the limit both in the Newtonian and relativistic regimes, while the relativistic correction to the Newtonian relation becomes more important as one increases . This is one of the reasons that we chose the Newtonian relations as the controlling factor of the fit in Eq. (9).
Let us now consider how universal the relation is when exploring a large set of EoSs, a fact that cannot be inferred by considering only the four EoSs of the bottom left panel of Fig. 4. The contour plot in the right panel shows the maximum fractional difference between the fit and the numerical data obtained using all eleven realistic EoSs we considered in this paper. Observe that the approximate universality holds to if and for all . For example, the maximum fractional difference for is roughly . On the other hand, observe that the maximum fractional difference can reach a maximum of when . In fact, notice that the maximum fractional difference increases from to as one considers stars on the line of ever larger mass.
Why does the universality become worse as one increases the NS mass? One might find this puzzling since the stellar sequence approaches the BH limit as one increases the mass (or equivalently, as one increases the stellar compactness), and in this limit, the universality should become exact (because of the nohair theorems of General Relativity (GR)). Indeed, and do approach their BH values in the BH limit, but their averaged difference, , is sensitive to how this limit is approached, or in other words, to the slope of the – relation. Figure 5 shows (top) and (bottom) against for four different EoSs. Observe that becomes less sensitive to as one increases , since decreases. Physically, the slope of the – relation decreases with increasing mass because decreases very fast with stellar compactness, which, in turn, is simply due to the fact that more dense stars deform less than less dense stars [see Eq. (5)]. Therefore, when one calculates the relative fractional difference in , both the numerator and the denominator of this fraction decrease as increases. That the universality deteriorates in the BH limit is then simply a consequence of the slope of the – relation decreasing faster than the EoS variability in .
Let us end this section by studying the binary Love relations for QSs. The top panel of Fig. 6 presents the – relations using three different QS EoSs, SQM1, SQM2 and SQM3 [89]. We also show the fit created for the NS relations and the Newtonian relation with . The latter choice is made because the QS EoSs are similar to a constant density EoS (i.e. an polytropic EoS) in the lowpressure regime [6]. One sees that the QS relations are slightly different from the NS ones, in particular when is close to 1. The bottom panel of Fig. 6 shows the absolute fractional difference between the QS relations and the NS fit. When or 0.9, the difference between the two relations can be as large as in the large regime, which is comparable to the intrinsic EoSvariation in the NS relations. This shows that one can still apply the relations for NSs to QSs. That is, a separate set of QS universal relations is not needed at this stage.
2.2 – Relation
Another way to choose two independent tidal parameters in the gravitational waveform is through the coefficients in the PN expansion of the GW phase. For nonspinning compact binaries, finitesize effects enter the GW phase first at PN order relative to the leadingorder (Newtonian) term [37], with 1PN corrections to this entering at 6PN order. We can then choose the two independent tidal parameters to be the combination of tidal deformabilities that enter first as a coefficient of the 5PN (the chirp tidal deformability ) and the 6PN () terms in the phase, namely
(10)  
(11) 
where the symmetric mass ratio is given by
(12) 
Such a choice was made e.g. in [37, 38, 44, 46, 47]^{4}^{4}4 and in Eqs. (10) and (11) are equivalent to and in Eqs. (5) and (6) of [46]. The signs in front of terms proportional to are different because we assume (so that is positive) in this paper, while Ref. [46] assumed .. Observe that for an equal mass binary, these parameters reduce to and . Although such parametrization is less intuitive than , the former has the data analysis advantage that the correlation between and is smaller relative to that between and .
2.2.1 Newtonian Limit
As in Sec. 2.1, let us first investigate the relation in the Newtonian limit, i.e. to leading order in an expansion about small compactness. Such a relation can easily be obtained by substituting the Newtonian – relation of Eq. (7) into the definition of and in Eqs. (10) and (11) so as to obtain and . Next, we solve for to obtain . Finally, we substitute this into to obtain , which leads to the Newtonian relation
(13) 
where
(14) 
with
(15) 
In the small massratio limit, becomes , and hence the relation is exactly EoSuniversal. On the other hand, in the equalmass limit, the fractional difference between the relation with a polytropic index of and becomes
(16) 
where, in the last equality, we set and (see Table 1). Observe that the fractional difference is larger than the one for the Newtonian – relation by for . We thus expect the relativistic – relations to approach the Newtonian relations more slowly than the – ones do as the compactness increases.
2.2.2 Relativistic Relations
Let us now return to the relativistic relations. Once relativistic NSs have been constructed with a given EoS, one can easily obtain the relativistic – relations, just as we did in Sec. 2.1.2 for the – relation. The top left panel of Fig. 7 is the same as that of Fig. 4 but for the – relation. As in the latter figure, we focus on fixed , 0.75 and 0.9, and also plot the Newtonian relation with and a fit to the numerical data. The latter is created from the fitting function of Eq. (9), with , , , and the fitting coefficients of Table 2. As in the – case of Fig. 4, the numerical results approach the Newtonian relation in the large limit. The bottom left panel of Fig. 7 shows the absolute fractional difference between the numerical data and the fit. As in the – relation, the difference becomes larger as one increases .
Let us now study the maximum absolute fractional difference of the – relation. The right panel of Fig. 7 shows an contour plot of the maximum absolute fractional difference between the fit and all eleven realistic EoSs considered in this paper. Comparing this to the right panel of Fig. 4, one sees that the fractional error is larger in the – relation than in the – relation. Even then, the EoSvariation of the relation is less than for . For a system of , the maximum absolute fractional difference is roughly . On the other hand, as one increases the mass for a fixed , the EoS variation becomes larger for the same reason as it does in the – case.
2.3 – Relations
Yet another parametrization of finite size effects in the GW phase of inspiraling NS binaries is through a Taylor expansion of the dimensionless tidal deformability around a reference mass , taking the coefficients in the expansion as tidal parameters:
(17) 
Of course, what enters the GW phase is and . A similar parametrization was suggested e.g. in Refs. [50, 39, 42, 48]. The advantage of choosing this kind of Taylorexpanded parametrization is that the coefficients are constant, and thus, they are the same irrespective of the NS mass and one can increase the measurement accuracy of such tidal parameters with multiple GW detections. The disadvantage is that one might have to include not only but also higher order tidal parameters if the NS mass is not close to ; this introduces correlations between the tidal parameters that could deteriorate the accuracy to which can be measured.
Whether Eq. (17) is a good representation of depends on whether the Taylor expansion converges and the rate at which it does so. The left panel of Fig. 8 shows the relative fractional error between this series representation and the correct , for different truncations , with the MS1 EoS and an expansion mass of . One sees that the more terms that are added to the series, the smaller the error only in the region . Thus, the series converges in this region, but it diverges in the high mass region . This shows that the parameterization is not suitable for NS binaries with masses that are sufficiently different from . The right panel of Fig. 8 shows the maximum absolute fractional error between and using the eleven realistic EoSs in the plane. This contour plot shows the range of for a given with which one can use up to 3 terms in the series to approximate to a given accuracy. For example, the fractional difference is always smaller than for when .
2.3.1 Newtonian Limit
Let us now consider the – relation in the Newtonian limit. First, from Eq. (5), one finds that . Thus, by taking derivatives of this relation with respect to , one finds that
(18) 
where is the Gamma function. for representative values of (corresponding to WFF1, MS1 and mean EoSs in Table 1) and are shown in Table 3, together with the fractional difference from the mean EoS shown in brackets. Observe how the EoS variation increases as one increases .
k  1  2  3 

WFF1  
MS1  
mean 
k  1  2  3 

0.9984  0.9955  0.9923 
2.3.2 Relativistic Relations
We now look at the EoS variation in the – relations in the relativistic regime. We first calculate these relations numerically and then fit them to
(19) 
where and the coefficients for each are given in Table 4. This fit will then allow us to estimate the degree of EoS variability of the numerical results.
The top left panel of Fig. 9 presents the regions spanned by the relations with varying EoSs obtained numerically as a function of . One finds a single curve for a fixed and EoS (see Fig. 1 of [78]) with the single parameter that varies along each curve. The top panel of the figure also shows the fit of Eq. (19), while the bottom panel shows absolute values of the relative fractional difference between the numerical results and the fit. Observe that the difference becomes larger as one increases . The top panel also shows the Newtonian relation of Eq. (18) with (dashed lines). Again, observe that the relations approach the Newtonian ones in the large limit. The right panel of Fig. 9 shows the EoS variation in the – relations for each fixed (analogous to the bottom left panel). For example, when as chosen in [50, 42, 48], the universality holds to for but to and for and 3 respectively.
3 Parameter Estimation
In this section, we estimate the impact of the universal binary Love relations on EoS constraints with GW observations of NS binaries. We obtain this estimate through a Fisher analysis [90, 91], which should provide a rough measure of the accuracy with which one may be able to extract bestfit parameters. For a more robust estimate, in particular for low SNR signals, one needs to carry out a Bayesian analysis [42, 46, 47, 48].
3.1 Preliminary
We begin by reviewing the main ideas behind a Fisher analysis. For stationary and Gaussian noise, and in the large SNR limit, the posterior probability distribution of the template parameter vector with a given signal is approximately given by [92, 93]
(20) 
with the bestfit parameter that maximizes the probability distribution and the prior distribution. The Fisher matrix is defined by
(21) 
where is the waveform template (a model for the response of the instrument due to an impinging GW) and the inner product is defined by
(22) 
with a tilde and * denoting the Fourier transform and complex conjugate operations respectively, and the noise spectral density. We follow [91, 92, 93] and consider a crude prior distribution on certain parameters, given by a Gaussian centered around :
(23) 
The explicit choice of is irrelevant for calculating the statistical error on . Since a product of two Gaussian distributions gives a new Gaussian, one obtains the rootmeansquare of as