Probing Quadratic Gravity with Binary Inspirals

In this paper, we study gravitational waves generated by binary systems within an extension of General Relativity which is described by the addition of quadratic in curvature tensor terms to the Einstein-Hilbert action. Treating quadratic gravity as an effective theory valid in the low energy/curvature regime, we argue that reliable calculations can be performed in the early inspiral phase, and furthermore, no flux of additional massive waves can be detected. We then compute the massive dipole (-1PN) leading corrections to the post-Newtonian (PN) expansion of the standard waveform. By confronting these theoretical calculations with available experimental data, we constrain both unknown parameters of quadratic gravity to be $0 \leq \gamma \, \lesssim 5.7\cdot 10^{76}$, and $-\frac{\gamma}{4} \leq \beta \, \lesssim - 4.2\cdot 10^{75}$.


Introduction
The observations of gravitational waves from compact binary systems provide a testing ground for various speculative hypotheses regarding extensions of General Relativity (GR) [1] (see also [2; 3], and the references therein). For example, it has been demonstrated in [4] that remarkably strong constraints on noncommutative spacetime follow from the observed signals. An extension of GR obtained by the addition of higher curvature terms to the simplest diffeomorphism-invariant Einstein-Hilbert action is arguably the least speculative modification. In fact, the quantum nature of matter makes such a modification of the Einstein theory of gravitation unavoidable. In particular, one-loop renormalisation in the semi-classical gravity regime leads to the gravitational action with additional quadratic in curvature tensor terms [5]: where κ = 8πG = 1/M 2 P [M P ≈ 2.4 · 10 18 GeV is the reduced Planck mass in natural units], and β, and γ are the dimensionless parameters 1 . Due to the nonrenormalisability of gravity, these parameters cannot be computed theoretically and must be inferred from experiments. Furthermore, higher loops generate an infinite series of higher curvature terms, with an infinite number of a priori free parameters, which makes the theory intractable in full generality 2 . In this situation, a sensible thing to do is to view the theory described by the action (1) as the effective theory applicable in the low energy/curvature regime only (see below for a more quantitative description). Adopting this approach to quadratic gravity (1) in this work we aim to infer constraints on β, and γ parameters from the existing gravitational wave observations [1].
The study of quadratic gravity has a long history [for a recent review see [8], and references therein]. In fact, it has been known for a long time [9] that quadratic gravity by itself while being a renormalisable theory contains, however, a fatal deficiency -propagating negative energy states, which compromises it's stability 3 . More recently, gravitational waves generated by binary inspirals have also been studied. In [19] dynamical quadratic gravity with scalar field dependent parameters has been studied with the main focus on the additional massless scalar waves. In [20] gravitational waveforms were studied in pure scalar quadratic gravity (with γ = 0) with an additional assumption that the black hole binary is surrounded by the shell of a scalar field. In Ref. [21], which is the closest to our present study, a constraint on the spin-2 wave mass (equivalently on γ) were obtained under the assumption that the massive spin-2 waves can reach the detector. In [22; 23] the limits on the conformal versions of quadratic gravity has been discussed. Finally, in [24] post-merger formation of the horizonless 2-2 hole has been considered with the characteristic prediction of gravitational echoes during the ringdown phase. None of these scenarios are feasible within the consistent effective field theory treatment of quadratic gravity adopted in the current work. 1 A possible term proportional to the square of the Riemann tensor, R µνρσ R µνρσ , in (1) can be eliminated using the Gauss-Bonnet identity. In a spacetime with boundaries there is an additional term proportional to the covariant d'Alambertian of the Ricci scalar [6], R, which is not relevant for our discussion and we omit it in (1). 2 In principle, the string effective action must also contain an infinite and well organized series of higher curvature corrections to the leading Einstein gravity [7]. However, any truncation of such action for practical calculations will bring the problems of the effective field theory back. 3 See, e.g., Refs [10][11][12][13][14][15][16][17][18] for some alternative interpretations, however.
The rest of the paper is organised as follows. In Sec. 2 we rewrite the action (1) in an equivalent and more convenient form separating out the massive spin-0, and the massive spin-2 fields, and will present their solutions in the linearised approximation. Given these solutions we compute the leading order corrections to the waveforms generated by a binary inspiral in Sec. 3. Constraints on quadratic gravity based on these calculations are discussed in Sec. 4. The final Sec. 5 is reserved for our conclusions.

Summary of Methodology
It may be instructive to summarise our method and philosophy in slightly less technical language than as follows in the rest of this paper. First, we reinterpret the quadratic gravity degrees of freedom as massive scalar and tensor fields, within the realm of validity of an effective field theory for gravity. Since we expect these extra fields to be small compared to usual GR, we can find field equations and solutions for these fields to lowest order.
The post-Newtonian (PN) expansion scheme provides a method for perturbative calculations of the effects of gravity on bodies moving slow compared to the speed of light. This is certainly true for black hole binaries in their inspiral phase, and provides a useful framework to study modifications of gravity. As a perturbative scheme, the largest corrections to usual GR would occur at lowest PN order. Nominally, if one was to modify gravity, to place constraints on the modifications from the actual GW observations one would need to run a Bayesian analysis, varying all parameters to determine the maximum possible deviations from GR. However, this is computationally expensive, and technically difficult. In place of this sort of analysis, one can alternatively utilise the results of [1], in which they run a Bayesian analysis, varying all GR parameters to determine the maximum possible deviation from GR of any particular order of the waveform. Since this conservative bound involves varying all parameters in the scenario, one would not expect that a modification of gravity could be "absorbed", for example, into a re-scaling of the black hole masses, or any other parameter. This provides us with an order-of-magnitude estimate of any possible modification of GR, which is convenient for our analysis, which is already to lowest-order, and can only aspire to find an order-of-magnitude constraint. One also sees here that the largest observational constraints on possible modifications exist at the lowest PN order, further justifying the decision to calculate corrections to lowest PN order, since this will provide the most robust constraints on any modification.
With the above in mind, we must then calculate the GW signal to lowest order. This is done by first deriving the motion of the two black holes, and in particular their phase as a function of time. In usual GR, the binary loses energy due to GW emissions. In our modification, additional scalar and tensor radiation affect the binary motion as well, which in turn straightforwardly modifies the GW waveform detected at Earth. The calculation of this waveform is conducted in the well-known stationary phase approximation. As explained above, conservative constraints on the possible deviation of this waveform from GR already exist, and these constraints can then be used to directly constrain our modification of gravity.

Massive Scalar and Spin-Gravitons
In quadratic gravity (1) besides the massless spin-2 graviton, the theory contains additional propagating degrees of freedom: a massive scalar, φ, and a massive spin-2 graviton, π µν . It is convenient to explicitly separate these degrees of freedom by rewriting the quadratic action (1) in the form of the Einstein-Hilbert action supplemented by terms describing the extra massive degrees of freedom. This can be achieved by explicitly introducing φ, and π µν fields in the action through the set of Lagrange multipliers which enforce the relations: φ = 1 2 R, and π µν = R µν − 1 4 g µν R. After field redefinitions and keeping terms quadratic in φ, and π µν only, we arrive at the following action [25]: whereR is the Ricci scalar constructed out of the metric to leading order in φ and π µν , and the mass terms are defined as: We require γ ≥ 0, and 4β ≥ −γ in order to avoid tachyonic instabilities. Notice, however, that the massive spin-2 field in (2) has the wrong sign for its kinetic term. To avoid the associated ghost instabilities, and to maintain theoretical consistency in our calculations we need to clarify the domain of validity of the effective quadratic gravity. This can be done by assuming that the quadratic action is a truncated asymptotic series expansion in powers of curvature tensor, where the error introduced by the truncation at a given power is smaller than the next terms in the expansion. Suppose we are interested in a system with a typical size r (e.g., the size of a binary inspiral). Then |R| ∼ 1/r 2 and the requirement that the Einstein-Hilbert term in (1) dominates over the quadratic terms implies M 2 P r 2 /β, M 2 P r 2 /γ > 1, which can be translated into the constraints: m φ,π r > 1. There are few important ramifications of the above restrictions coming from the validity of quadratic truncation of the effective gravity action. First, it is clear from the action (2) that vacuum solutions of the Einstein gravity (e.g., the Schwarzschild black hole solution) is also a solution in quadratic gravity with trivial background massive fields φ = π µν = 0. The stability of these solutions within the domain of applicability of quadratic gravity are then guaranteed, and all the exotic non-perturbative solutions, such as the 2-2 hole solution, should be dismissed [26]. Second, theoretically, reliable calculations can be performed only in the regime with not too big curvature tensors, e.g. in the inspiral phase (this study), or in the ringdown phase. The merger phase is in principle intractable within the effective theory approach. Finally, the constraints m φ,π r > 1, imply that the flux of massive waves cannot be produced during the inspiral phase. This comes from the requirement that the frequency of the waves arriving at the detector, given approximately by an inspiral angular velocity ω ≃ Ω ≈ v/r, must exceed the mass, ω > m φ,π and this requirement is in contradiction with the constraints m φ,π r > 1, since v < 1 during the inspiral phase. Therefore we do not expect massive waves to be produced during the inspiral phase (see also [27]). We would like to stress that these conclusions are conservative and are purely based on the validity of quadratic gravity as an effective field theory. We do not exclude other treatments of quadratic gravity (perhaps supplemented by extra fields) where our constraints are relaxed.
The binary system during the inspiral phase is simply modelled as two point particles with masses m a , and 4-velocities v µ a , a = 1, 2. The action (in redefined fields) reads: Next we consider linear perturbations,h µν =g µν − η µν , φ, and π µν and combining the actions (2), and (5), we obtain the following Euler-Lagrange equations: whereh µν =h µν − 1 2 η µνh α α , and ∂ µh µν = ∂ µ π µν = π µ µ = 0. Eq. (6) describes the usual massless gravity waves produced by the binary inspiral, while (7), and (8) describe massive scalar and spin-2 waves, respectively. Since the massive waves contain scalar polarizations, they result in monopole and dipole radiation in addition to the standard quadrupole radiation already at the lowest order of the post-Newtonian (PN) expansion. Therefore, to constrain two parameters of quadratic gravity it is sufficient to compute corrections to the waveform at -2PN; however, we should note that observational constraints from gravitational wave experiments are currently only available at -1PN order. Higher-PN orders can also be systematically computed following the known formalism (see the review Ref. [28], and references therein) in this framework. The solution forh µν is well known. The solutions of the massive field (7), and (8) can also be straightforwardly obtained using massive, retarded Green's functions: where t r = t − | x − y a (t r )| is the retarded time. As it was expected, these solutions demonstrate that the massive fields mediate additional short-range Yukawa interactions in the binary system. Note that, while the interactions mediated by the massive scalar are attractive, the massive spin-2 generates a repulsive force. Of course, this is due to the fact that the π µν field describes negative energy states. Using these solutions, we are ready now to compute leading order corrections to the inspiral waveforms in quadratic gravity.

Computing Leading Order Corrections to the Inspiral Waveforms
Our aim is to compute the leading order corrections to the phase of gravitational waves generated in the inspiral regime of a compact binary within quadratic gravity. To this end, we start with the covariant conservation equation, where∇ µ is the covariant derivative constructed out of the metricg µν (3), and the energy momentum tensor of the binary is defined as: Picking up ν = i spatial components in (11), and integrating over the spatial volume surrounding only one of the two compact objects, we obtain the modified Newton's force law describing the dynamics of the object: where P i a is the linear momentum of the a-th object, and F i a is the force acting on it due to the massless graviton exchange [29]: From (13), and using our solutions (9, 10) we calculate the accelerations, a 1,2 = d v 1,2 dt , as: where r 12 = | y 1 − y 2 |, andn 12 = ( y 1 − y 2 )/r 12 ( a 2 = a 1 [1 ↔ 2]). As we are interested in relative motion, we introduce convenient parameters: the total mass M = m 1 + m 2 , the reduced mass µ = m 1 m 2 m 1 +m 2 , and the symmetric ratio ν = µ M . The relative acceleration, a = a 1 − a 2 then reads: wheren =n 12 = −n 21 , dropping the subscript.
In the quasi-circular orbit approximation this defines the angular frequency of the binary system: The total energy of the virialised binary is given by: We see that the massive scalar field decreases the energy of the binary, while the massive spin-2 carries out negative energy, and increases the energy of the binary system. The rate of the energy change is convenient to calculate by introducing a dimensionless frequencyrelated parameter: x ≡ (GMΩ) The relative distance in terms of x reads: Using (20) we replace r by x in (18) keeping leading order terms, and taking the derivative, we obtain: The change of the energy of binary must be balanced by the flux of emitted gravitational waves, hence: We have argued previously, that the massive scalar, and the massive spin-2 waves cannot contribute to the radiated flux, therefore F in (22) is given by the standard GR expression: We now have most of the tools we need to calculate the binary phase. First we introduce the dimensionless time variable Θ: where t c is the instant of coalescence. The orbital phase ϕ of the binary system is defined as dϕ/dt = Ω. This in terms of Θ reads: Using (24), and (25) Next, using (21), and (23) the we obtain the differential equation for the phase, which can be easily integrated: This is the desired result which describes the leading corrections to the binary inspiral phase. The exponentially suppressed terms in (28) are due to additional massive scalar, and massive spin-2 degrees of freedom in quadratic gravity. Among those, corrections ∝ x −7/2 are due to the Yukawa suppressed dipole-like contributions, which is totally absent in GR, and can be classified as -1PN corrections. The corrections ∝ x −5/2 come from Yukawa quadrupole radiation and are at 0PN Newtonian order in the post-Newtonian expansion.

Constraints on Quadratic Gravity
Previous constraints on quadratic gravity parameters are based on the modification of the Newtonian inverse square law introduced by the exchange of the massive scalar, and the massive spin-2 mediators. This can be probed in various experiments, including tabletop torsion balance [30], and laser ranging satellite experiments [31; 32]. While significant constraints on the strength of Yukawa modification of the Newtonian potential are obtained in those experiments, which invalidates quadratic gravity, one must keep in mind that that those constraints are obtained only for specific ranges of interactions. For example, the torsion balance experiments excludes quadratic gravity only for 10 −5 eV m π,φ 10 −3 eV [30] (in terms of the parameters, using (4), 2.9 · 10 60 ≤ γ 2.9 · 10 64 , and 4.8 · 10 59 − γ 4 ≤ β 4.8 · 10 63 − γ 4 ), while satellite experiments are sensitive only to the range comparable to the orbital size they monitor: m π,φ ≃ 3 · 10 −14 eV (earth ranging [31]), and m π,φ ≃ 5 · 10 −16 eV (lunar ranging [32]).
On the contrary, bounds on quadratic gravity that can be inferred from the modification of the gravitational waveform is a priori valid for an arbitrary range and thus we expect genuine upper limits, which we will derive now. The gravitational wave strain signal with amplitude A(t), and phase Φ(t) takes the following form: Since the gravitational wave frequency is twice the orbital frequency of the binary system we have Φ = 2ϕ. Next, using (22) we find x as a function of the dimensionless time parameter Θ, and rewrite (28) as an explicit function of time: The gravitational wave forms are analysed in the frequency space, however. Fourier transform of (29) in the stationary phase approximation [33] reads: where the parameter t f is given by the time when the gravitational wave frequency dΦ(t)/dt is equal to the Fourier frequency f . The phase can be explicitly evaluated and written as expansions in PN orders: where t c , and φ c are the time, and phase at coalescence respectively. In particular, we obtain that the highest order corrections to the standard GR phase to be of negative second, and negative first order PN corrections: where δ = 4GM (8πGMf ) −2/3 . The advantage of using gravitational waves to constrain quadratic gravity steams from the fact that the leading correction arises at lower than Newtonian order (see, (28)). The absolute deviation of the -1PN phase has been constrained from observations of gravitational waves [1] to be |δφ −1P N | < 10 −2 , while (to our knowledge) there does not exist at the moment -2PN constraints [1]. From our calculations of -1PN correction (35), we then obtain: We could also write explicitly higher order corrections to the GR phase. However, the lower PN bounds provide tighter constraints on quadratic gravity than higher PN bounds and we do not display higher PN corrections here. Assuming no accidental cancellation between the two terms on the left hand side of (36) and taking the typical values f = 75Hz, and M = 30M ⊙ , the inequality is reduced to: which leads to the bound: m φ,π 7.0 · 10 −12 eV.
At last, using (4), we rewrite (38) as limits on the original dimensionless parameters of quadratic gravity: The bound on β (properly translated) is 4 orders of magnitude stronger than the one anticipated in [20].
As a final remark we note that ideally we would use -2PN correction to the GW waveform to make our constraints, since throughout our analysis we have been taking lowest-order terms only. This means the exact values of the -1PN coefficients given above should actually be modified; for example, the multiplication of a -2PN term with a 1PN term would result in a -1PN term which needs to have been kept for full accuracy. Our calculation is only tractable, however, by remaining at lowest order. One would not expect the order of magnitude of the coefficients to vary greatly in the full solution, so our order-of-magnitude estimation above is safe. Once the -2PN constraints are known, we could instead use those to more accurately place bounds on quadratic gravity.

Conclusions
In this paper, we have studied gravitational wave signals from the inspiral phase of a compact binary system within the framework of quadratic gravity. Quadratic gravity has been considered as a truncated approximation of a theory represented by an action written in terms of the powers of curvature tensors. Such an action necessarily emerges as a result of renormalisation of the Einstein-Hilbert action within the semi-classical gravity at any fixed order of the quantum loop-expansion. Since semi-classical gravity is not renormalisable there are in principle an infinite number of free parameters corresponding to an infinite series of higher curvature terms. Therefore, in the truncated quadratic theory which contains only two parameters, we must make sure that we are working in the regime, where higher than quadratic curvature terms are sub-dominant. This restricts our consideration to inspiral phase of the evolution of the binary system. We have also argued that non-perturbative solutions significantly deviating from the GR solutions (e.g., the 2-2 hole solution) are not valid within our approximation. Furthermore, constraints on the theoretical validity of quadratic gravity, imply that no flux of additional massive radiation is produced during the inspiral phase.
Under these constraints, we have computed the leading corrections to the phase of the GR waveform due to the quadratic terms (28). These corrections comprise of formally -1PN corrections due to the dipole radiation of massive scalar, and scalar polarization of the massive spin-2 fields. Confronting these corrections with the existing data [1], we have extracted bounds (39, 40) on both β, and γ parameters of quadratic gravity. To the best of our knowledge these are the stringent limits available in the literature.
Higher PN order calculations can in principle be performed. However, no improvement of the obtained bounds (38, 39, 40) within quadratic gravity is expected 4 . On the other hand, one can envisage the program which probes higher than quadratic curvature terms by more accurate PN calculations of the binary inspiral. Ultimately, the best probe of higher curvature gravity must come from the signal produced during the merger phase. However, besides complications related with numerical calculations, theoretical control of the validity of finite truncation in merger phase becomes a significant challenge.