Comparison of post-Newtonian templates for compact binary inspiral signals
in gravitational-wave detectors
Alessandra Buonanno,
1,
*
Bala R. Iyer,
2,3,†
Evan Ochsner,
1,‡
Yi Pan,
1,x
and B. S. Sathyaprakash
3,k
1
Maryland Center for Fundamental Physics, Department of Physics, University of Maryland, College Park, Maryland 20742, USA
2
Raman Research Institute, Bangalore, 560 080, India
3
School of Physics and Astronomy, Cardiff University, 5, The Parade, Cardiff, United Kingdom, CF24 3YB
(Received 6 July 2009; published 29 October 2009)
The two-body dynamics in general relativity has been solved perturbatively using the post-Newtonian
(PN) approximation. The evolution of the orbital phase and the emitted gravitational radiation are now
known to a rather high order up to O ðv
8
Þ, v being the characteristic velocity of the binary. The orbital
evolution, however, cannot be specified uniquely due to the inherent freedom in the choice of parameter
used in the PN expansion, as well as the method pursued in solving the relevant differential equations. The
goal of this paper is to determine the (dis)agreement between different PN waveform families in the
context of initial and advanced gravitational-wave detectors. The waveforms employed in our analysis are
those that are currently used by Initial LIGO/Virgo, that is, the time-domain PN models TaylorT1,
TaylorT2, TaylorT3, the Fourier-domain representation TaylorF2 (or stationary phase approximant), and
the effective-one-body model, and two more recent models, TaylorT4 and TaylorEt. For these models we
examine their overlaps with one another for a number of different binaries at 2PN, 3PN, and 3.5PN orders
to quantify their differences. We then study the overlaps of these families with the prototype effective-one-
body family, currently used by Initial LIGO, calibrated to numerical-relativity simulations to help us
decide whether there exist preferred families, in terms of detectability and computational cost, that are the
most appropriate as search templates. We conclude that as long as the total mass remains less than a
certain upper limit M
crit
, all template families at 3.5PN order (except TaylorT3 and TaylorEt) are equally
good for the purpose of detection. The value of M
crit
is found to be 12M
for Initial, Enhanced, and
Advanced LIGO. From a purely computational point of view, we recommend that 3.5PN TaylorF2 be used
below M
crit
and that the effective-one-body model calibrated to numerical-relativity simulations be used
for total binary mass M>M
crit
.
DOI: 10.1103/PhysRevD.80.084043 PACS numbers: 04.30.Db, 04.25.Nx, 04.80.Nn, 95.55.Ym
I. INTRODUCTION
Sensitivity of several interferometric gravitational-wave
detectors has either already reached, or is close to, the
design goals that were set more than a decade ago [1–7].
Upgrades that are currently underway and planned for the
next four to five years will see their sensitivity improve by
factors of a few to an order of magnitude [8]. Coalescing
binaries consisting of neutron stars and/or black holes are
probably the most promising sources for a first direct
detection of gravitational waves. At current sensitivities,
initial interferometers are capable of detecting binary neu-
tron star inspirals at distances up to 30 Mpc, the range
increasing to 60 Mpc for enhanced detectors (circa
2009–2011) and 450 Mpc for advanced detectors (circa
2014+). Binary black holes or a mixed system consisting of
a neutron star and a black hole can be detected to a far
greater distance depending on the total mass and the mass
ratio. For example, a ð10 þ 10ÞM
binary can be detected
out to distances of 160 Mpc by initial detectors and
2200 Mpc by advanced detectors [9].
The range of interferometric detectors for coalescing
binaries is computed by assuming that one can pull the
signal out of noise by matched filtering. This in turn means
that one is able to follow the phasing of gravitational waves
typically to within a fraction of a cycle over the duration of
the signal in band. The reason for this optimism comes
from the fact that one knows the phase evolution of the
signal to a high order in post-Newtonian (PN) formalism
[10]. Several authors have assessed whether the accuracy
with which the formalism provides the waveforms is good
enough for the purpose of detection and parameter estima-
tion [11–26]. The problem, as we shall see below, is
complicated since the PN approximation does not lead to
a unique model of the phase evolution. Moreover, though
PN results are good up to mildly relativistic velocities, the
standard PN approximants become less and less accurate in
the strongly relativistic regime as one approaches the last
stable orbit (LSO). Resummation methods [15] and, in
particular, the effective-one-body (EOB) [27–29] exten-
sions of the PN approximants are needed for analytical
treatments close to and beyond the LSO.
*
buonanno@umd.edu
†
bri@rri.res.in
‡
evano@umd.edu
x
ypan@umd.edu
k
B.Sathyaprakash@astro.cf.ac.uk
PHYSICAL REVIEW D 80, 084043 (2009)
1550-7998=2009=80(8)=084043(24) 084043-1 Ó 2009 The American Physical Society
The success in numerical-relativity simulations of bi-
nary black holes [30–34] now provides results for gravita-
tional waveforms that can be compared to standard PN
results and other resummed extensions. On the one hand,
the analytical PN results for the inspiral phase of the
evolution are needed to calibrate and interpret the
numerical-relativity waveforms of coalescence and
merger. On the other hand, the numerical-relativity results
extend the analytical approximations beyond the inspiral
phase and provide the important coalescence and merger
phases, producing the strongest signals that are crucial for
the detection of binary black holes. However, numerical
simulations are still computationally expensive and time-
consuming, and presently only a small region of the pa-
rameter space can be explored. Even in the foreseeable
future, numerical relativity may not be able to handle the
tens of thousands of cycles that are expected from highly
asymmetric systems (e.g., a neutron star falling into an
intermediate-mass black hole of 100M
) or low-mass
symmetric systems (e.g., a binary neutron star).
Analytical models that smoothly go from the inspiral
through coalescense to quasinormal ringing would be
needed, and this has led to phenomenological templates
[35–37] and EOB waveforms [36,38–46]. In particular, the
recent, improved EOB models [45,46], which also incor-
porate a multiplicative decomposition of the multipolar
waveform into several physically motivated factors supple-
mented by a suitable hybridization (using test particle
results) [47], and an improved treatment of nonquasicircu-
lar corrections, show evidence of remarkable success in
modeling accurately the numerical-relativity waveforms
for different mass ratios.
The emphasis of this work is different. Recently, there
have been investigations [48] on the ability of various
standard families of PN templates to detect a specific signal
model, TaylorEt [49–51], and the often-used TaylorF2 to
detect a complete numerical-relativity signal including
merger and ringdown [36,37]. Reference [48] modeled
the signal by the TaylorEt approximant at 3.5PN order
and looked at the effectualness and systematic biases in
the estimation of mass parameters for TaylorT1, TaylorT4,
and TaylorF2 templates in the LIGO and Virgo detectors. It
also looked into the possibility of improving the effectual-
ness by using unphysical values of beyond the maximum
value of 0.25. It was found that the overlaps of a TaylorEt
signal with the TaylorT1, TaylorT4, and TaylorF2 tem-
plates are smaller than 0.97 and involved for equal-mass
systems a large bias in the total mass. For unequal-mass
systems higher overlaps can be obtained at the cost of a
large bias in mass and symmetric mass ratio and can be
further improved by unphysical values of >0:25. The
templates are more unfaithful with increasing total mass.
To detect optimally the complete numerical-relativity sig-
nal, including merger and ringdown, Ref. [36] suggested
the possibility of using the TaylorF2 template bank with a
frequency cutoff f
c
larger than the usual upper cutoff (i.e.,
the Schwarzschild LSO) and closer to the fundamental
quasinormal mode frequency of the final black hole.
Moreover, they proposed to further improve this family
by allowing either for unphysical values of or for the
inclusion of a pseudo 4PN (p4PN) coefficient in the tem-
plate phase, calibrated to the numerical simulations.
Reference [37] extended the results of Ref. [36] to more
accurate numerical waveforms, found that 3.5PN templates
are nearly always better and rarely significantly worse
than the 2PN templates, and proposed simple analytical
frequency cutoffs for both Initial and Advanced LIGO—
for example, for Initial LIGO they recommended a strategy
using p4PN templates for M 35M
and 3.5PN templates
with unphysical values of for larger masses. However, we
notice that there is no reason for changing the template
bank above 35M
. Reference [37] could have used the
p4PN templates over the entire mass region, if they had
not employed in their analysis the p4PN coefficient used in
Ref. [36], but had calibrated it to the highly accurate
waveforms used in their paper.
1
In this work our primary focus is on binary systems
dominated by early inspiral and on a critical study of the
variety of approximants that describe this. Towards this
end, in this paper we will provide a sufficiently exhaustive
comparison of different PN models of adiabatic inspiral for
an illustrative variety of different systems and quantify
how (dis)similar they are for the purpose of detection.
The choice of the PN models used in this paper is moti-
vated by the fact that they are available in the LIGO
Algorithms Library (LAL), and some of them have been
used in the searches by Initial LIGO. We also compare all
these PN models with one fiducial EOB model calibrated
to numerical-relativity simulations [40] to delineate the
range of mass values where one must definitely go beyond
the inspiral-dominated PN models to a more complete
description including plunge and coalescence. The choice
of this fiducial, preliminary EOB model is only motivated
by the fact that it is the EOB model available in LAL and it
is currently used for searches by Initial LIGO. It will be
improved in the future using the recent results in
Refs. [45,46]. We will conclude that for total masses below
a certain upper limit M
crit
, all template families at 3.5PN
order (except for TaylorT3 and TaylorEt) are equally good
for the purpose of detection. M
crit
is found to be 12M
for Initial, Enhanced, and Advanced LIGO. Based solely
on computational costs, we recommend that 3.5PN
TaylorF2 be used below M
crit
and that EOB calibrated to
numerical-relativity simulations be used for total binary
mass M>M
crit
.
1
We computed that the p4PN coefficient calibrated to the
highly accurate waveforms used in Ref. [37]isY ¼ 3714,
instead of Y ¼ 3923 as found in Ref. [36].
BUONANNO, IYER, OCHSNER, PAN, AND SATHYAPRAKASH PHYSICAL REVIEW D 80, 084043 (2009)
084043-2
The paper is organized as follows. In Sec. II we sum-
marize the present status of the PN approximation. In
Sec. III we recapitulate for completeness the main PN
approximants and try to provide a ready reckoner for the
equations describing them and the relevant initial and
termination conditions. In Sec. IV we discuss the fre-
quency evolution in each of these models. In Sec. V we
discuss overlaps and the maximization used in this work.
Sections VI and VII present the results of our analysis
related to the effectualness, while Sec. VIII summarizes
the results related to the faithfulness. In Sec. IX we sum-
marize our main conclusions. Readers who are interested
in the main results of the paper and want to avoid technical
details could skip Secs. II, III, IV, and V, read the main
results of Secs. VI, VII, and VIII, and mainly focus on
Sec. IX.
II. CURRENT STATUS OF POST-NEWTONIAN
APPROXIMATION
Post-Newtonian approximation computes the evolution
of the orbital phase ðtÞ of a compact binary as a pertur-
bative expansion in a small parameter, typically taken as
v ¼ðMFÞ
1=3
(characteristic velocity in the binary), or
x ¼ v
2
, although other variants exist. Here M is the total
mass of the binary and F the gravitational-wave frequency.
In the adiabatic approximation, and for the restricted
waveform in which case the gravitational-wave phase is
twice the orbital phase, the theory allows the phasing to be
specified by a pair of differential equations
_
ðtÞ¼v
3
=M,
_
v ¼F ðvÞ=E
0
ðvÞ, where M is the total mass of the
system, F its gravitational-wave luminosity, and E
0
ðvÞ is
the derivative of the binding energy with respect to v.
Different PN families arise because one can choose to treat
the ratio F ðvÞ=E
0
ðvÞ differently starting formally from the
same PN order inputs [18]. For instance, one can retain the
PN expansions of the luminosity F ðvÞ and E
0
ðvÞ as they
appear (the so-called TaylorT1 model), or expand the
rational polynomial F ðvÞ=E
0
ðvÞ in v to consistent PN
order (the TaylorT4 model), recast as a pair of parametric
equations ðvÞ and tðvÞ (the TaylorT2 model), or the
phasing could be written as an explicit function of time
ðtÞ (the TaylorT3 model). These different representations
are made possible because one is dealing with a perturba-
tive series. Therefore, one is at liberty to ‘‘resum’’ or
‘‘reexpand’’ the series in any way one wishes (as long as
one keeps terms to the correct order in the perturbation
expansions), or even retain the expression as the quotient
of two polynomials and treat them numerically. There is
also the freedom of writing the series in a different vari-
able, say (suitably adimensional) E (the so-called TaylorEt
model).
In addition to these models, there have been efforts to
extend the evolution of a binary beyond what is naturally
prescribed by the PN formalism. Let us briefly discuss two
reasons why the PN evolution cannot be used all the way
up to the merger of the two bodies. PN evolution is based
on the so-called adiabatic approximation, according to
which the fractional change in the orbital frequency F
orb
over each orbital period is negligibly small, i.e.
_
F
orb
=F
2
orb
1. This assumption is valid during most of
the evolution, but begins to fail as the system approaches
the LSO where f
LSO
¼ð6
3=2
MÞ
1
. In some cases, the
frequency evolution stops being monotonic and
_
f changes
from being positive to negative well before reaching the
LSO—an indication of the breakdown of the
approximation.
From the viewpoint of maximizing the detection poten-
tial, one is also interested in going beyond the inspiral
phase. The merger and ringdown phases of the evolution,
when the luminosity is greatest, cannot be modeled by
standard PN approximation. The use of resummation tech-
niques more than a decade ago was followed by the con-
struction of the EOB model [27–29], which has
analytically provided the plunge, merger, and ringdown
phases of the binary evolution. As mentioned before,
more recently, these models have been calibrated to
numerical-relativity simulations [36,38–46]. We now
have a very reliable EOB model that can be used to model
the merger dynamics.
An astronomical binary is characterized by a large num-
ber of parameters, some of which are intrinsic to the system
(e.g., the masses and spins of the component stars and the
changing eccentricity of the orbit) and others that are
extrinsic (e.g., source location and orientation relative to
the detector). In this paper we will worry about only the
detection problem. Furthermore, we will assume that a
coincident detection strategy will be followed so that we
do not have to worry about the angular parameters such as
the direction to the source, the wave’s polarization, etc. If
binaries start their lives when their separation r is far larger
compared to their gravitational radius (i.e., r GM=c
2
),
by the time they enter the sensitivity band of ground-based
detectors, any initial eccentricity would have been lost due
to gravitational radiation reaction, which tends to circular-
ize
2
a binary [59,60]. Therefore, we shall consider only
systems that are on a quasicircular inspiraling orbit. We
shall also neglect spins, which means that we have to
worry, in reality, about only the two masses of the compo-
nent bodies.
Our goal is to explore how (dis)similar the different
waveform families are. We do this by computing the
(normalized) cross correlation between signals and tem-
plates, maximized either only over the extrinsic parameters
2
Though this assumption is justified for the prototypical bi-
naries we focus on in this work, there exist credible astrophysical
scenarios that lead to inspiral signals from binaries with non-
negligible eccentricity in the sensitive detector bandwidth. A
more involved treatment is then called for and available. See e.g.
[52–58].
COMPARISON OF POST-NEWTONIAN TEMPLATES FOR ... PHYSICAL REVIEW D 80, 084043 (2009)
084043-3
of the templates (faithfulness) or over the intrinsic and
extrinsic parameters of the templates (effectualness), the
noise power spectral density (PSD) of the detector serving
as a weighting factor in the computation of the correlation
(see Sec. V). Our conclusions, therefore, will depend on
the masses of the compact stars as well as the detector that
we hope to observe the signal with.
The overlaps (i.e., the normalized cross correlation
maximized over various parameters and weighted by
the noise power spectral density) we shall compute are
sensitive to the shape of the noise spectral density of a
detector and not on how deep that sensitivity is. Now,
the upgrade from initial to advanced interferometers will
see improvements in sensitivity not only at a given
frequency but over a larger band. Therefore, the agreement
between different PN models will be sensitive to the
noise spectral density that is used in the inner product.
Thus, we will compare the PN families using power spec-
tral densities of initial and advanced interferometric
detectors.
We end this brief overview with the following observa-
tion. As mentioned earlier, following all present
gravitational-wave data analysis pipelines, this paper
works only in the restricted wave approximation. This
approximation assumes the waveform amplitude to be
Newtonian and thus includes only the leading second
harmonic of the orbital phase. Higher PN order amplitude
terms bring in harmonics of the orbital phase other than
the dominant one at twice the orbital frequency. Their
effects can be significant [61,62], especially close to
merger [45], and they need to be carefully included in
future work.
III. THE PN APPROXIMANTS
For the convenience of the reader, in this section, we
recapitulate the basic formulas for the different PN fami-
lies from Refs. [18,19]. While comparing the expressions
below to those in Refs. [18,19], recall ¼1987=3080
[63,64] and ¼11831=9240 [65,66]. In addition to the
evolution equations, we shall also provide initial and final
conditions. From the perspective of a data analyst, the
initial condition is simply a starting frequency F
0
and
phase
0
, which can be translated, with the help of evolu-
tion equations, as conditions on the relevant variables. We
shall also give explicit expressions for the evolution of the
gravitational-wave frequency, namely,
_
F dF=dt,or
more precisely, the dimensionless quantity
_
FF
2
,in
Sec. IV, where they will be used to study the rate at which
the binary coalesces in different PN families, which will
help us understand the qualitative difference between
them. The contents of this section should act as a single
point of resource for anyone who is interested in imple-
menting the waveforms for the purpose of data analysis and
other applications.
The basic inputs for all families are the PN expressions
for the conserved 3PN energy (per unit total mass)
[63,64,67–70] E
3
ðvÞ and 3.5PN energy flux [65,66,71–
73] F
3:5
ðvÞ,
E
3
ðvÞ¼
1
2
v
2
1
3
4
þ
1
12
v
2
27
8
19
8
þ
1
24
2
v
4
675
64
34 445
576
205
96
2
þ
155
96
2
þ
35
5184
3
v
6
; (3.1)
F
3:5
ðvÞ¼
32
5
2
v
10
1
1247
336
þ
35
12
v
2
þ 4v
3
44 711
9072
9271
504
65
18
2
v
4
8191
672
þ
583
24
v
5
þ
6 643 739 519
69 854 400
þ
16
3
2
1712
105
þ
41
48
2
134 543
7776
94 403
3024
2
775
324
3
856
105
logð16v
2
Þ
v
6
16 285
504
214 745
1728
193 385
3024
2
v
7
;
(3.2)
where ¼ 0:577 216 . . . is the Euler constant. In the adia-
batic approximation one assumes that the orbit evolves
slowly so that the fractional change in the orbital velocity
! over an orbital period is negligibly small. That is,
!
!
1, or, equivalently,
_!
!
2
1. In this approximation, one
expects the luminosity in gravitational waves to come
from the change in orbital energy averaged over a period.
For circular orbits this means one can use the energy
balance equation F ¼dE=dt where E ¼ ME.
In the adiabatic approximation one can write an equation
for the evolution of any of the binary parameters. For
instance, the evolution of the orbital separation rðtÞ can
be written as
_
rðtÞ¼
_
E=ðdE=drÞ¼F =ðdE=drÞ. Together
with Kepler’s law, the energy balance equation can be used
to obtain the evolution of the orbital phase
3
:
d
dt
v
3
M
¼ 0; (3.3a)
dv
dt
þ
F ðvÞ
ME
0
ðvÞ
¼ 0; (3.3b)
or, equivalently,
3
Recall that the gravitational-wave phase is twice the orbital
phase for the restricted waveform and leads to differences in
factors of 2 between the equations here for the orbital phase and
those in [18] for the gravitational-wave phase.
BUONANNO, IYER, OCHSNER, PAN, AND SATHYAPRAKASH PHYSICAL REVIEW D 80, 084043 (2009)
084043-4
tðvÞ¼t
ref
þ M
Z
v
ref
v
dv
E
0
ðvÞ
F ðvÞ
; (3.4a)
ðvÞ¼
ref
þ
Z
v
ref
v
dvv
3
E
0
ðvÞ
F ðvÞ
; (3.4b)
where t
ref
and
ref
are integration constants and v
ref
is an
arbitrary reference velocity.
A. TaylorT1
The TaylorT1 approximant refers to the choice corre-
sponding to leaving the PN expansions of the luminosity
F ðvÞ and E
0
ðvÞ as they appear in Eq. (3.3) as a ratio of
polynomials and solving the differential equations numeri-
cally,
d
ðT1Þ
dt
v
3
M
¼ 0; (3.5a)
dv
dt
þ
F ðvÞ
ME
0
ðvÞ
¼ 0: (3.5b)
In the above v v
ðT1Þ
, but for the sake of notational
simplicity we write only v; from the context the meaning
should be clear. In the formulas of this section, and in the
sections that follow, the expressions for F ðvÞ [EðvÞ] are to
be truncated at relative PN orders 2[2], 3[3], and 3.5[3] to
obtain 2PN [18,74–76], 3PN, and 3.5PN [19,65,73] tem-
plate or signal models, respectively.
To see how to set up initial conditions, refer to Eq. (3.4).
Let the initial gravitational-wave frequency be F
0
or,
equivalently, let the initial velocity be v
0
¼ðMF
0
Þ
1=3
.
One normally chooses t ¼ 0 at v ¼ v
0
. This can be
achieved by choosing v
ref
¼ v
0
and t
ref
¼ 0 in Eq. (3.4).
The initial phase
ref
is chosen to be either 0 or =2 in
order to construct two orthogonal templates (see Sec. VA
for details).
B. TaylorT4
TaylorT4 was proposed in Ref. [23] and investigated in
Refs. [33,38,77], thus many years after the other approx-
imants discussed in this paper were proposed (with the
exception of TaylorEt, which is even more recent).
However, it is a straightforward extension of TaylorT1,
and at 3.5PN order, by coincidence, it is in better agree-
ment with numerical simulations of the inspiral phase
[33,36,38,41,43,50,77]. The approximant is obtained by
expanding the ratio of the polynomials F ðvÞ=E
0
ðvÞ to the
consistent PN order. The equation for v
ðT4Þ
ðtÞvðtÞ at
3.5PN order reads
dv
dt
¼
32
5
M
v
9
1
743
336
þ
11
4
v
2
þ 4v
3
þ
34 103
18 144
þ
13 661
2016
þ
59
18
2
v
4
4159
672
þ
189
8
v
5
þ
16 447 322 263
139 708 800
þ
16
3
2
1712
105
þ
451
48
2
56 198 689
217 728
þ
541
896
2
5605
2592
3
856
105
logð16v
2
Þ
v
6
4415
4032
358 675
6048
91 495
1512
2
v
7
: (3.6)
The orbital phase
ðT4Þ
is determined, as in the case of
TaylorT1, by Eq. (3.3a), and the numerical solution of
Eqs. (3.3a) and (3.6) yields the TaylorT4 approximant.
Note that although TaylorT1 and TaylorT4 are perturba-
tively equivalent, the evolution of the phase can be quite
different in these two approximations. The asymptotic
structures of the approximants are also quite different:
While
_
v can have a pole (although not necessarily in the
region of interest) when using Eq. (3.5b), no pole is pos-
sible when Eq. (3.6) is used. Differences of this kind can, in
principle, mean that the various PN families give different
phasing of the orbit. The hope is that when the PN order up
to which the approximation is known is large, then the
difference between the various PN families becomes
negligible.
Setting up the initial conditions for TaylorT4 is the same
as in the case of TaylorT1.
C. TaylorT2
TaylorT2 is based on the second form of the phasing
relations, Eq. (3.4). Expanding the ratio of the polynomials
F ðvÞ=E
0
ðvÞ in these equations to consistent PN order and
integrating them, one obtains a pair of parametric equa-
tions for ðvÞ and t ðv Þ, the TaylorT2 model.
ðT2Þ
n=2
ðvÞ¼
ðT2Þ
ref
þ
v
N
ðvÞ
X
n
k¼0
^
v
k
v
k
; (3.7a)
t
ðT2Þ
n=2
ðvÞ¼t
ðT2Þ
ref
þ t
v
N
ðvÞ
X
n
k¼0
^
t
v
k
v
k
: (3.7b)
Of all the models considered in this study, TaylorT2 is
computationally the most expensive. This is because the
phase evolution involves solving a pair of transcendental
equations, which is very time-consuming.
COMPARISON OF POST-NEWTONIAN TEMPLATES FOR ... PHYSICAL REVIEW D 80, 084043 (2009)
084043-5
