The Gravitational-Wave Discovery Space of Pulsar Timing Arrays
Curt Cutler, Sarah Burke-Spolaor, Michele Vallisneri, Joseph Lazio, Walid Majid
TThe Gravitational-Wave Discovery Space of Pulsar Timing Arrays
Curt Cutler,
1, 2
Sarah Burke-Spolaor, Michele Vallisneri,
1, 2
Joseph Lazio, and Walid Majid Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Dr., Pasadena, CA 91109 Theoretical Astrophysics, California Institute of Technology, Pasadena, California 91125 (Dated: October 30, 2018)Recent years have seen a burgeoning interest in using pulsar timing arrays (PTAs) as gravitational-wave (GW) detectors. To date, that interest has focused mainly on three particularly promisingsource types: supermassive–black-hole binaries, cosmic strings, and the stochastic background fromearly-Universe phase transitions. In this paper, by contrast, our aim is to investigate the PTApotential for discovering unanticipated sources. We derive significant constraints on the availablediscovery space based solely on energetic and statistical considerations: we show that a PTA de-tection of GWs at frequencies above ∼ × − Hz would either be an extraordinary coincidenceor violate “cherished beliefs;” we show that for PTAs GW memory can be more detectable thandirect GWs, and that, as we consider events at ever higher redshift, the memory effect increasinglydominates an event’s total signal-to-noise ratio. The paper includes also a simple analysis of theeffects of pulsar red noise in PTA searches, and a demonstration that the effects of periodic GWsin the ∼ − . –10 − . Hz band would not be degenerate with small errors in standard pulsarparameters (except in a few narrow bands).
PACS numbers: 04.80.Nn, 04.25.dg, 04.30.-w, 95.85.Sz, 97.80.-d, 97.60.Gb
I. INTRODUCTION
The idea of detecting gravitational waves (GWs) bymonitoring the arrival times of radio pulses from neu-tron stars (i.e., by pulsar timing ) was first proposed bySazhin [1] and Detweiler [2]; its modern formulation byHellings and Downs [3] emphasizes the importance ofsearching for correlations in the pulse-timing time devi-ations among an array of intrinsically stable millisecondpulsars. The last few years have seen a strong renewedinterest in these searches, with the formation of three ma-jor pulsar timing programs: the European Pulsar TimingArray (EPTA, [4]), the North American Nanohertz Ob-servatory for Gravitational Waves (NANOGrav, [5]), andthe Australian Parkes Pulsar Timing Array (PPTA, [6]),which have now joined into a global collaboration, theInternational Pulsar Timing Array (IPTA, [7]).The most promising known sources of GWs for PTAsare in-spiraling supermassive black hole binaries (SMB-HBs). Some estimates suggest that these will be detectedby PTAs as soon as ∼ GW ≡ ρ GW /ρ (the ratio of the energy density in GWsto the closure density of the Universe), and is beginningto impact standard theories of hierarchical structure for-mation via constraints on the SMBH merger rate.In this article we explore the discovery potential ofPTAs. Our main motivation is to minimize the risk thatcurrent observing strategies and planned data-analysispipelines artificially preclude the discovery of varioustypes of sources. For instance, most pulsars in PTAsare currently observed with irregular cadences of ∼ ∼ − Hz), where PTAs are par-ticularly sensitive. However a search for GW bursts last-ing (say) 10 s would clearly benefit greatly from coordi-nated timing observations (using a few radio telescopes)that get repeated several times a day. Thus we addressthe following questions: • Is there a strong motivation for increasing the ob-serving cadence to improve our sensitivity to GWswith frequencies ∼ − –10 − Hz? • What constraints can we impose on the PTA dis-covery space based on simple energetic, statistical,and causality arguments?In addressing the first question, an important issuethat arises is whether, even if strong sources exist in thisband, our sensitivity might be degraded by degeneraciesbetween GW effects and small errors in the timing-modelparameters of the monitored pulsars. In addressing thesecond question, we are necessarily retracing some ofthe trails blazed by Zimmermann and Thorne [16] (here-inafter ZT82) in their classic paper, “The gravitationalwaves that bathe the Earth: upper limits based on the-orists’ cherished beliefs.” However there are importantdifferences between our paper and theirs: • ZT82 restricted attention to sources at z ∼ <
3, whilewe consider the case of very high- z sources as well. • Unlike ZT82, we include the “memory effect”among potential observables; its detection turnsout to be especially promising in the high- z case. a r X i v : . [ g r- q c ] S e p • ZT82 restricted attention to GWs in the frequencyrange 10 − < f < Hz (the band of interestfor ground-based and space-based interferometers),while we focus on GWs with f ∼ < − Hz. (How-ever, there are several instances for which the ZT82estimates extend trivially to lower frequency; wewill note these instances in our paper as they arise.)This paper is organized as follows. In Sec. II we de-scribe a simple general framework for thinking about pul-sar timing observations, and we characterize how the de-tection signal-to-noise ratio scales with quantities such asthe number of pulsars surveyed, the timing accuracy pro-vided by each pulsar observation, the observing cadence,and the total observation time. We also briefly reviewpulsar timing noise, with some emphasis on its red noisecomponent. In Sec. III we summarize salient results re-garding PTA searches for SMBHBs and cosmic strings,largely to provide points of comparison with possible un-known GW sources. In Sec. IV we demonstrate that thetiming residual signatures of GWs in the 10 − . –10 − . Hz band are not degenerate with small errors in the pul-sar parameters, except for very narrow frequency bands;had this been otherwise, there would have been littlepoint in considering more fundamental constraints onpossible sources in this band. In Secs. V and VI we in-vestigate what constraints on source strengths arise fromfundamental considerations of energetics, statistics andcausality. In Sec. VII we discuss how are estimates getmodified for highly beamed sources, and for sources inour Galaxy. In Sec. VIII we summarize our conclusions,listing some caveats.Regarding notation, we adopt units in which G = c =1. Also, the signal frequency f , observation time T obs and signal duration T sig all refer to time as measured inthe observer’s frame, at the Solar system barycenter. II. THE PTA SIGNAL-TO-NOISE RATIO FORGW SIGNALS OF KNOWN SHAPEA. Signal-to-noise ratio for white noise signals
In the rest of this paper, we are going to assume anidealized, general scaling law for the detection signal-to-noise ratio (SNR) of an individual GW source, as ob-served by a pulsar timing array: to wit,SNR = M N (cid:28) δt δt (cid:29) , (2.1)where • δt GW is the timing residual due to GWs; • δt noise is the noise in the residuals, which in-cludes contributions from the observatory, frompulse propagation, and from intrinsic pulsar pro-cesses; • (cid:104)· · · (cid:105) denotes the average over all pulsars in thePTA and over all observed pulses; • M is the number of pulsars in the PTA; and • N is the total number of observations for each pul-sar.In what follows, purely for simplicity we will assume that δt rms is roughly the same across PTA pulsars and obser-vations, so we define (cid:28) δt δt (cid:29) = (cid:104) δt (cid:105) δt . (2.2)with δt rms a representative rms value for the noise.The term “timing residual” requires definition: it is thedifference between the time of arrival (TOA) of a trainof pulses observed at the radio telescope and the TOA predicted by the best-fitting timing model for a pulsar.This deterministic model includes parameters (such asthe sky position of and distance to the pulsar) that affectthe propagation of signals to the observatory, as well asparameters (such as the pulsar period and its derivativesand, if needed, orbital elements for pulsars in binaries)that describe the intrinsic time evolution of the pulsar’semission.The pulses from millisecond pulsars are usually tooweak to be observed individually, so the TOAs refer to integrated pulse profiles obtained by “folding” the outputof radiometers with the putative pulsar period over ob-servations with durations of tens of minutes to an hour.Typically, such pulsar timing observations are repeatedat intervals of two to four weeks, yielding sparse data sets;however, the individual observations are often run quasi-simultaneously at multiple receiving frequencies (typi-cally one hour to two days apart, since the feeds need tobe switched), yielding a set of TOAs at the same epoch .See [17, 18] and references therein for more detail.In analogy with other applications in GW data anal-ysis [19], our scaling for the SNR can be motivated byconsidering a ratio of likelihoods : namely, the likelihoodof the residual data r i (with i indexing both epochs andpulsars) under the hypothesis that r i = g i + n i , with g i describing a GW signal of known shape, and n i denot-ing noise; and the likelihood of the residuals under thenoise-only hypothesis r i = n i . For Gaussian noise, whenthe GW signal is really present, the likelihood ratio isexp { g i ( C − ) ij g j / n i ( C − ) ij g j } (2.3)(summations implied), where C ij = (cid:104) n i n j (cid:105) is thevariance-covariance matrix for the noise. The first termin the exponent, which depends only on the GWs, isidentified as SNR /
2, while the second term is a ran-dom variable with mean zero and variance (over noiserealizations) equal to SNR . This can be proved, e.g., byconsidering that Gaussian noise with covariance C canbe written as √ C ¯ n , with √ C √ C T = C the Choleskydecomposition of C , and with ¯ n a vector of uncorrelated,zero-mean/unit-variance Gaussian variables. Then (cid:104) ( n i ( C − ) ij g j )( n l ( C − ) lm g m ) (cid:105) = ( C − ) ij g j √ C ki (cid:104) ¯ n k ¯ n p (cid:105)√ C pl ( C − ) lm g m = ( C − ) ij C il ( C − ) lm g j g m = ( C − ) jm g j g m . Equation (2.1) follows immediately under the (strong)assumption that noise is uncorrelated and homogenousamong pulsars and epochs, so that it can be repre-sented by ( C − ) ij = δ ij /δt . We are assuming alsothat the sampling of pulsars and epochs in the datasetis sufficiently broad and non-pathological that (cid:80) i g i (cid:39) M N (cid:104) δt (cid:105) ; that is, that the sampling can effectivelyperform an average over time and pulsar sky position. Ifthe noise is uncorrelated (i.e., white ), but not homoge-neous, Eq. (2.1) still stands, provided that δt can betaken to represent a suitable averaged noise.Under these assumptions, Eq. (2.1) is remarkable inthat the actual form of the signal to be detected appearsonly through its variance (cid:104) δt (cid:105) , and that the structureof observations appears only through their overall num-ber M × N and rms noise δt . By contrast, one mayhave imagined that detecting (say) quasi-sinusoidal sig-nals of high frequency f GW would require rapid-cadenceobservations spaced by ∆ t ∼ < /f GW , according to theNyquist theorem. However, that theorem is a statementabout the reconstruction of the whole of a function on thebasis of a set of regularly spaced samples, but it does notapply to our case—computing the likelihood that a signalof known shape is present in the data [20]. In effect, weare checking that the measured data are consistent withour postulated signal: for uncorrelated noise, it does notmatter when we check, but only how many times we doit. B. Relaxing the assumption of white noise
There are two important considerations that challengeour assumption of white, uncorrelated noise.The first is that the residuals include a stochastic con-tribution due to the over-fitting of noise (and possibleGWs) at the time of deriving the timing model. We dis-cuss this further in Sec. IV, where we show empiricallythat the detection of quasi-sinusoidal signals at most fre-quencies would not be affected. From a formal stand-point, van Haasteren and colleagues [21] show that itpossible to marginalize the likelihood over timing-modelparameter errors δξ by replacing the inverse covariancein Eq. (2.3) with C − − C − M ( M T C − M ) − M T C − ,where M is the design matrix for the timing model fit, sothat the extra contribution to the residuals has the form M δξ (A similar strategy of “projecting out” parametererrors was employed earlier by Cutler and Harms [22],in the context of removing residual noise from slightly incorrect GW foreground subtraction). For uncorrelatednoise, Eq. (2.1) is modified only by restricting the com-putation of (cid:104) δt (cid:105) to the GW components that are notabsorbed away by the timing model (and this is indeedwhat we investigate in Sec. IV).The second important consideration (and for whichthe GW frequency does matter) is the impact of cor-related noise. The physically interesting case here is thatof long-term correlations, which generate red noise thatis stronger at low frequencies. To understand the impactof red noise, we study a toy model in which the N obser-vations are organized in P “clumps” of Q TOAs taken atnearby times (with N = P × Q ), and where noise con-sists of two components: uncorrelated noise with vari-ance σ and noise with variance κ that is completelycorrelated within clumps, and completely uncorrelatedbetween clumps. (We use κ since κ ´ oκκινoζ is Greek for“red.”) We consider a single pulsar, although the gener-alization to more is trivial.The resulting C has the structure C = σ I + κ P (cid:88) i =1 O i , (2.4)where each O i is a matrix that has ones for every com-ponent corresponding to a combination of samples in thesame burst, and zeros everywhere else. Each O i can alsobe written as u i u Ti , where u i is a vector that has ones forthe components in clump i , and zeros everywhere else.From the block structure of C and the Woodbury lemma[23], it follows that C − = σ − I − σ − P + κ − /σ − P (cid:88) i =1 O i . (2.5)If the characteristic frequency of the GW signal is“slower” than the timescale of a clump (i.e., the time overwhich the Q samples in a clump are collected), then thesum (cid:80) i g Ti O i g i (cid:39) P Q (cid:104) δt (cid:105) , because the same valueof g is being summed over and over in each burst. Itfollows thatSNR = (cid:104) δt (cid:105) P Qσ + Qκ = (cid:104) δt (cid:105) Pσ /Q + κ ; (2.6)that is, the repeated observations in each clump averageout the uncorrelated component of noise (as ∝ / √ Q ),but not its correlated part. Increasing the number ofobservations in a clump provides diminishing returns as σ /Q → κ .Let us follow the other branch of our derivation: if thecharacteristic frequency of the GW signal is “faster” thanthe timescale of the clumps, then, barring special coin-cidences, (cid:80) i g Ti O i g i (cid:39) P Q (cid:104) δt (cid:105) , and SNR reduces(modulo an O [1 /Q ] correction) to the general expression(2.1), with N = P Q . C. Noise characteristics inferred fromobservational data
In this section, we consider the characteristics of noisefor real pulsars. Namely, to what extent is our analysisapplicable to timing residuals from actual PTAs?For the radiometer noise due to thermal effects in thereceiving system, the assumption of no correlations (i.e.,“white”) is well justified: for observations over a ra-dio frequency bandwidth ∆ ν , the correlation timescaleis (∆ ν ) − , so this noise contribution is effectively uncor-related in time. Further, from thermodynamic considera-tions, the assumption of Gaussianity is also well justified.Pulsars can show correlated, red spectrum fluctuationsin their TOAs, and Cordes and Shannon [24] present asummary of various effects, ranging from intrinsic spinfluctuations to magnetospheric and propagation effects;see also [25]. These effects have spectral densities ∝ f − x ,with x typically > >
4. On timescales ∼ f ∼ − . Hz), the residuals appear to bedominated by white components ([5, 26]; see also Figs.10 and 11 of [6] for a visual representation of noise ef-fects in PPTA pulsars). Even if σ ≈ κ at frequencies ∼ − . Hz, at higher frequencies ( ∼ > − Hz), the vari-ance from white processes will exceed that of any redprocesses with relatively shallow spectra ( x ≈
1) by afactor of approximately 15; for red processes with steeperspectra ( x ≈ Q [Eq. (2.6)].For more general observation schemes and red-noise pro-cesses, we may think of the number of clumps P as T obs /T red , where T obs is the total duration of observation,and T red is the correlation timescale of the most signifi-cant red-noise process; then Q (cid:39) N ( T red /T obs ). For GWsignals with frequency ∼ < /T red , our toy model wouldthen suggest thatSNR = (cid:104) δt (cid:105) σ /N + κ ( T red /T obs ) ; (2.7)that is, the 1 / √ N averaging of noise becomes limited byred noise once N ∼ ( σ /κ ) × ( T obs /T red )—an interestingscaling in its own right. For GW signals with frequencies ∼ > /T red , the simpler scaling (2.1) applies.In the remainder of this paper, we neglect the effects ofred noise in the scaling of SNR and assume the expressionof Eq. (2.1). Our assumption is correct because one ormore of the following circumstances will be true (or true enough ) in practice: • The characteristic GW frequency of interest willbe greater than 1 /T red for the most significant red-noise component. • For a majority of the pulsars in the PTA, the white-noise variance will exceed that of the most domi-nant red-noise process for the time scales of inter-est. • The number of observations will not saturate theaveraging of white noise with respect to sub-dominant red noise (i.e., in the “clump” picture, σ /Q > κ ). III. BRIEF REVIEW OF PROSPECTS FOR PTASEARCHES OF SUPERMASSIVE–BLACK-HOLEBINARIES AND COSMIC STRINGS
Here we collect a few salient points concerning PTAsearches for SMBHBs and cosmic strings, mostly to pro-vide points of comparison with the hypothetical sourceswe consider in the next sections. We refer the reader tothe literature cited below for more details.
A. The detectability of GWs from supermassiveblack hole binaries
When two galaxies merge, the SMBHs at their centersare brought together by tidal friction from the surround-ing stars and gas. It seems likely that their separationeventually shrinks to the point at which gravitationalradiation emission dominates the inspiral, and the twoSMBHs eventually coalesce [27]. The GWs from all in-spiraling SMBHBs in the observable Universe contributeto a stochastic background of GWs with characteristicamplitude h c ∼ h rms √ f given by h c ≈ A ( f /f ) − β (3.1)in the PTA band, where β ≈ / A is predicted tobe in the range 5 × − –5 × − for f = 10 − Hz[8, 14, 28, 29]. Depending on the actual A , the first PTAdetection of GWs is expected between 2016 and 2020 [30].The background is expected to be dominated by bina-ries with chirp masses M c ≡ ( m m ) / ( m + m ) − / ∼ M (cid:12) at z ∼ <
2. At frequencies above f ≈ − Hz,sources are sparse enough that the central limit theoremdoes not apply, so the distribution is significantly non-Gaussian and a few brightest sources would appear abovethe background. Thus, the first PTA discovery could ei-ther be an individual strong (and possibly nearby) source,or the full background.
B. The detectability of GWs from cosmic strings
There are several mechanisms by which an observablenetwork of cosmic (super)strings could have formed inthe early Universe [31]. Simulations have shown thatstring networks rapidly approach an attractor: the dis-tribution of straight strings and loops in a Hubble volumebecomes independent of initial conditions. The networkproperties do depend on two fundamental parameters:the string tension µ and the string reconnection proba-bility p . The size of string loops at their birth shouldin principle be derivable from µ and p , but the studiesare difficult and different simulations have produced verydifferent answers. Therefore most astrophysical analysestoday assume that the size of loops at their birth can beparametrized as α H − ( z ), where H − ( z ) is the Hubblescale when the loop is “born,” and where α is treatedas a third unknown parameter. We refer the reader to[31, 32] for nice reviews. To make matters more com-plicated, Polchinski has argued that the distribution ofloop size at birth is actually bimodal, with both rela-tively large and small loops being produced at the sameepoch [33]. Regarding the string tension µ , physicallymotivated values range over at least six orders of magni-tude: 10 − ∼ < µ ∼ < − .Once formed, string loops oscillate and therefore loseenergy and shrink due to GW emission. These wavesform a stochastic GW background. In addition to thisapproximately Gaussian background, the cusps and kinksthat form on the string loops emit highly beamed GWbursts [34, 35] Depending on the string parameters, PTAscould discover the stochastic background, the individualbursts, both or neither. The current limit on Ω GW ( f )from pulsar timing is Ω GW ( f ∼ − ) ∼ < × − [4],corresponding to a limit on the string tension of Gµ ≤ . × − . C. Current constraints on Ω GW ( f ) As mentioned, the current limit on Ω GW ( f ) from pul-sar timing is Ω GW ( f ∼ − ) ∼ < × − [4]. By com-parison, the limit from first-generation ground-based in-terferometers is Ω GW ( f ∼
100 Hz) < . × − [36].From Big Bang nucleosynthesis, we know also that anyGW stochastic background that existed already whenthe Universe was three minutes old satisfies Ω GW < . × − today [37]. Combined measurements of CMBangular power spectra (which are sensitive to lensing bya stochastic GW background) with matter power spec-tra also yield Ω GW ∼ < − today, but this method issensitive to any GWs produced before recombination at z ≈ H measurements gives the limit Ω GW ∼ < × − [39]. IV. SPECTRAL ABSORPTION EFFECTSFROM PULSAR TIMING-MODEL FITTING
The best knowledge of pulsar parameters comes fromthe iterative observation and refinement of a timingmodel, which predicts the times of arrival of all the pulsesas a function of all relevant parameters, such as the pe-riod and period derivatives of the pulsar’s intrinsic spin;the position, proper motion, and parallax of the pulsar;and possibly parameters that describe the motion of thepulsar in a binary system. Depending on the cadence andtotal time of observation, and on the shape and durationof the GWs, the effects of the GWs on pulse arrival times may correlate with the effects of changing the pulsar pa-rameters, so the GW power may be partly on entirelyabsorbed by the parameter-fitting process (see, e.,g., thestudy of the effect of a GW background on pulsar timingparameter estimation [40]).As a specific study of this effect, here we investigate theabsorption of sinusoidal GWs to demonstrate frequency-dependent signal loss to pulsar parameter fitting. To dothis, we use the
Tempo2 software suite [41] to simulatea set of timing residuals for pulsar J0613-0200 [6], asobserved with the Parkes observatory. We generate oneTOA every other day for T obs = 1,000 days, at a randomtime compatible with the pulsar being visible from theobservatory, and we add a white-noise component withrms amplitude of 100 ns. Into these simulated residualswe inject sinusoidal GWs from a circular SMBHB binarylocated at RA = Dec = 0, varying the GW frequency f between 10 − . and 10 − . Hz (corresponding to GWperiods of ∼ h + = h × = 10 − f , so that the SNR is fixed.For each GW frequency, we measure the power spectraldensity of the relevant frequency component before thetiming-model fit and after seven different types of fit: a fitagainst the full set of parameters and individual fits forpulsar frequency, frequency derivative, position, propermotion, parallax, and binary period. Figure 1 shows theratio of the power spectral densities in each case, as afunction of the source GW frequency. In effect, we areshowing the absorption spectrum of sinusoidal GWs, asfiltered by the timing-model fit. Above 10 − Hz, ∼ > f = 1 / year (corresponding to pulsar posi-tion/proper motion) and 1 / / (2 days) = 10 − . Hz),and to the sidereal day (at 1 / (23 .
934 hr) (cid:39) − . Hz).The former can be avoided with higher-cadence or irreg-ular observations, but the latter reflects the limitationsof using a single observatory, which can only observe asource while it is above the horizon. In our simulation wehave chosen random observation times within the windowof coverage, but more structured observing cycles can en-gender even deeper features. By contrast, this feature canbe avoided for a polar target that never sets. log(GW frequency, Hz)
Full FitFrequencyFrequency Derivative PS D r a t i o PositionProper MotionParallaxBinary Parameters
FIG. 1: GW power absorbed by fitting for various pulsar parameters as a function of GW frequency, for pulsar J0613 − To summarize, our example study suggests that fora majority of PTAs a high-frequency GW signal will bewell preserved through the standard timing-model fittingprocess, save for narrow features at roughly the observ-ing cadence and the sidereal day. GWs at frequenciesclose to either (1 year) − ) or (6 months) − will be signif-icantly impacted, as will GWs with periods approachingthe longest-duration pulsar observations. V. DISCOVERY SPACE FOR SOURCES IN THELOW-REDSHIFT UNIVERSE, z ∼ < O (1) In this section we begin to characterize the PTA dis-covery space for the case of sources in the low-redshiftUniverse, by which we mean z ∼ < O (1). We imaginethat there is some heretofore undiscovered GW source,and we ask what it would take for it to be detectablevia pulsar timing. We consider separately the case of modeled signals (for sources already conceived by the-orists, so that a parameterized waveform model can beused in a matched-filtering search), the case of unmodeledbursts , and the case of the gravitational memory effectfrom modeled sources.We will assume that the GW sources are distributedisotropically and that we do not occupy a preferred lo-cation in space and time with respect to them—that is, we assume that the Earth is not improbably close (spa-tially) to one of the sources, and that the sources havebeen emitting GWs for a significant fraction of the last10 years.We parametrize our projections in terms of the en-ergy density Ω GW . Because we consider sources in thelow-redshift Universe, in what follows we ignore redshifteffects. Nevertheless our results at z ∼ A. Discovery space for modeled GW signals in thelow-redshift Universe
As we established in Sec. II, the SNR modeled GWsignals as observed by a PTA isSNR = (cid:104) δt (cid:105) δt M N = (cid:104) δt (cid:105) δt M p min { T sig , T obs } , (5.1)where (cid:104) δt (cid:105) and δt are the mean-square-averagedtiming residuals due to GWs and measurement/pulsarnoise; M is the number of pulsars in the array; N is thenumber of times each pulsar is observed, which we rewritein terms of the cadence of observation p (e.g., 1/day), thetotal duration of observation T obs (e.g., 3 years), and thetypical duration of the GW signal T sig .For a sinusoidal GW signal of frequency f and rms am-plitude at Earth h = (cid:113) h + h × , the root-mean-squaretiming residual averages to¯ δt GW ≡ (cid:113) (cid:104) δt (cid:105) = 14 √ π hf (cid:39) hf . (5.2)Furthermore, the average rate at which the sources radi-ates energy in GWs is ˙ E = ( π / h f d (cid:39) h f d [19,Eq. (1.160)], where d is the distance to the source, and G = c = 1 (as we will set throughout). The GW energydensity from source of this kind isΩ GW ρ (cid:39) ( ˙ ET sig )( R τ ) , (5.3)where R is the spacetime rate–density of sources, and τ ∼ yr is the current age of the Universe. Approxi-mating the closure density ρ = 3 H / π as τ − /
10 (since τ (cid:39) H − ) and rewriting R ≡ ( V R T R ) − in terms of afiducial volume V R and the total event rate T R in that vol-ume, we can re-express the expected GW-induced timingresidual as¯ δt GW (cid:39) f − d − (cid:18) Ω GW V R T R τ T sig (cid:19) / . (5.4)We estimate the distance to the closest source thatwould be observed over time T obs by setting43 π d max { T obs , T sig } R = 1 (5.5)(where the maximum accounts for the persistence of mul-tiple emitting sources if T sig > T obs ), whence d near (cid:39) (cid:20) π V R T R T sig min { , T sig /T obs } (cid:21) / . (5.6)Folding together all the results of this section, we obtainthe corresponding, largest SNR that would be observedasSNR (cid:39) − Ω GW f τ M p T obs δt (cid:20) V R T R T sig min { , T sig T obs } (cid:21) / (cid:39) × − Ω GW f τ M p T obs δt d near . (5.7)We would now like to determine how large an SNR near we could expect for a given Ω GW , and for given observa-tional parameters M , p , T obs , and δt . This amounts To derive Eq. (5.2) we compute the Estabrook–Wahquist [42]fractional Doppler response (for the pulsar “Earth term” alone)to a sinusoidal GWs given by h + ( t )+ ih × ( t ) = ( h/ √
2) exp 2 πift ,take the antiderivative to obtain the corresponding pulse-time de-lay, square and average over time, sky position, and polarizationangle. to maximizing SNR near with respect to the GW-sourceparameters V R , T R , and T sig ; since these appear togetherin d near , we obtain the largest possible SNR near by setting d near = τ , the Hubble distance. We dare not place theGW source farther, since we are considering the “local”Universe and neglecting redshift effects.Note that the scaling SNR ∝ d near of Eq. (5.7)seems counterintuitive, since we would naively think ofthe strongest sources as the closest. However, while thesquared GW strain h at the Earth scales as 1 /d , it alsoscales with the total energy ∆ E that is emitted by eachsource, and that is “available” to each source given a fixedΩ GW ; this ∆ E increases with decreasing source density,and is proportional to d . This surprising intermediateresult was already shown in ZT82 [16].We can now plug in fiducial values for the observationalparameters (as well as τ = 3 × s), arriving atmax { SNR } ∼ < (cid:18) f − Hz (cid:19) − (cid:20) Ω GW − (cid:21) / × obs . ∼ < . (cid:18) f − Hz (cid:19) − (cid:20) Ω GW − (cid:21) / × obs ., (5.8)where obs . = (cid:20) δt rms − s (cid:21) − (cid:20) M p T obs (cid:21) / . (5.9)While we derived these constraints for the case of small z , we shall see below that they become even stronger forhigh- z sources.The fiducial values for f and Ω GW in the second row ofEq. (5.8) are motivated by our original question, whetherPTA searches should be extended to frequencies as highas ∼ − Hz. The current upper limit (from structureformation) on the energy density of hot dark matter isΩ
HDM ∼ < . × − (at 95% confidence) ; this limit ap-plies also to Ω GW . Our conclusion is that a PTA detec-tion of GWs at frequencies above ∼ × − Hz shouldbe considered very unlikely on fundamental grounds.
B. Discovery space for unmodeled GW bursts inthe low-redshift Universe
Quite simply, a burst is a signal with T sig ∼ /f .Since it contains only ∼ T sig , wecan still adjust R so that d near , as defined in Eq. (5.6),is equal to τ .For instance, if T sig ∼ < T obs and T obs = 10 s, this re-quires one burst every 10 s within a Hubble volume. Sofor this rate, the instantaneous SNR is the same as givenin Eq. (5.8) for modeled signals. This seems promising,because since bursts require no model for their detec-tion, they could potentially reveal phenomena that no-body ever thought of. At the same time, their detectionwould require the utmost care in excluding instrumentaland astrophysical artifacts. C. Discovery space for GW memory in thelow-redshift Universe
GWs with memory (for a recent review see [43]) causea permanent deformation – a “memory” of the passage ofthe waves – in the configuration of an idealized GW de-tector. They are emitted by systems with unbound com-ponents (linear memory), and by generic
GW sourcesbecause of the contribution of the energy–momentumof their “standard” GWs to the changing radiative mo-ments of the source (nonlinear memory). Several authorshave discussed the detectability of the GW memory ef-fect by PTAs for known source types, especially mergingmassive–black-hole binaries [12, 13, 44, 45]. Here we con-sider the effect from the point of view of the PTA discov-ery space, and again we ask in which region of parameterspace PTAs could discover previously unimagined sourcesby way of their GW memory.For a source at distance d from Earth, which emits atotal energy of ∆ E in GWs, the amplitude of the memoryeffect is [12] h mem ∼ α √ Ed , (5.10)where α < T sig (cid:28) T obs .Then we can approximate the “turn on” of the memoryeffect as a step function, and the effect on any pulsar isto create a timing residual that grows linearly in time: δt GW ∼ θ ( t − t ) h mem t, (5.11)where the memory passes over the Earth at time t .In any single pulsar, a linear-in-time residual can beinterpreted simply as a glitch causing an instantaneouschange in the pulsar frequency. However, all pulsars inthe PTA would show such apparent glitches at the sametime, with relative amplitudes following a simple pat-tern on the sky [12] determined by four parameters (thesky-location angles and two amplitudes that specify thetransverse–trace-free part of the metric), so in principlethe detection problem is well posed. The correspondingPTA SNR is [12]SNR mem ∼ h mem T obs δt rms ( M p T obs ) / , (5.12)where the factor 1 /
20 accounts for the facts that δt GW will typically be zero for a significant fraction of T obs , and that a large part of the effect will be absorbed in thepulsars’ timing models (and especially by the fitting oftheir periods and period derivatives) [12]. Note that GWmemory effect is essentially a low-frequency effect: SNRcan build up precisely because memory remains constant,but non-zero, for a sizable fraction of T obs . Thus thereis no particular advantage to high-cadence timing mea-surements.We can now derive how large an SNR we may expectfor detecting GW memory for a given Ω GW and for givenobservational parameters. As above, we relate the energydensity in GWs to the energy emitted in GW bursts,Ω GW ∼
10 ∆
E R τ ; (5.13)we then combine Eqs. (5.10), (5.12), and set d = d near =(4 πR T obs / − / , to obtainSNR mem , near ∼ α
300 Ω GW τ R − / T / ( M p T obs ) / δt rms . (5.14)Again, for fixed Ω GW we maximize SNR mem , near by tak-ing R to be as small as possible, subject to the constraintthat d near < τ , leading tomax { SNR mem } (cid:39) α
500 Ω GW τ T ( M p T obs ) / δt rms (cid:39) α (cid:20) Ω GW − (cid:21)(cid:20) T obs s (cid:21) × obs . (5.15)Comparing Eqs. (5.8) and (5.15), we see that–depending on the values of Ω GW and f –the memory effectfrom a burst could be much more detectable than its di-rect waves. More generally, comparing SNR mem with the direct SNR for the same source, as given by Eqs. (5.1)and (5.2), we find:SNR mem
SNR dir = 1 / / h mem T obs h/f (cid:18) M p T obs
M p T sig (cid:19) / = 1 / / π α √ h f d T sig T obs (cid:18) M p T obs
M p T sig (cid:19) / = 11 / π α √ dir δt rms T − T d ( M p T obs ) / (cid:39) α SNR dir (cid:20) T sig s (cid:21) − (cid:20) T obs s (cid:21) (cid:20) dτ (cid:21)(cid:2) obs . (cid:3) − (5.16)where in the second row we have used the fact that h mem (cid:39) ( α/ √ E/d ) and ∆ E = ( π / h f d × T sig ; in the third row we have substituted SNR dir =(1 / h/f ) δt − ( M pT sig ) / and replaced f with 1 /T sig ,as appropriate for a burst signal. Since SNR dir scales as h dir while SNR mem scales as h , the memory effect dom-inates for a sufficiently strong signal. VI. DISCOVERY SPACE AT HIGH REDSHIFT
In the previous section we have considered sources atsmall z , neglecting cosmological effects. We now turnto sources in the early Universe, at z (cid:29)
1. Again, wewill assume that the sources are isotropically distributedand that the Earth does not have a preferred location inspacetime with respect to them. The especially interest-ing cases are GW memory, which we discuss first, andunmodeled bursts.We begin by collecting a few useful formulas. Let t ≡ (cid:82) a − ( τ ) d τ be the conformal time coordinate, interms of which the (spatially flat) Robertson–Walkermetric becomes ds = a ( t ) (cid:2) − dt + dx + dy + dz (cid:3) . (6.1)We find it useful to define the high- z epoch into theradiation-dominated era for z (cid:28) z eq and the matter-dominated era for z (cid:29) z eq , where z eq ≈ ,
200 (theredshift at which the energies of matter and radiationwere equal). Then we can approximate a ( τ ) ∝ τ / for τ < τ eq and a ( τ ) ∝ τ / for τ > τ eq (of course, we nowknow that the Universe is dark-energy, rather than mat-ter dominated for z ∼ < .
7, but we neglect this correctionin keeping with the back-of-the-envelope spirit of this pa-per).We use the subscript “0” to refer to present Universe(e.g., τ ∼ years is the present age of the Uni-verse), and we choose our spatial coordinates so that a ≡ a ( τ ) = 1. Then t ( z ) = (cid:40) (1 + z eq ) (cid:0) τ / τ / ( z ) − τ eq (cid:1) z < z eq , (1 + z eq ) (cid:0) τ / τ / ( z ) (cid:1) z > z eq , . (6.2)and in particular, t (cid:39) (1 + z eq ) (cid:0) τ / τ / (cid:1) , (6.3)and therefore t t ( z ) (cid:39) (cid:40) (1 + z ) / z < z eq , (1 + z )(1 + z eq ) − / z > z eq , . (6.4)Now consider GW bursts produced at z (cid:29)
1. The sizeof the particle horizon at redshift z is ∼ t ( z ) in co-movingcoordinates, and so the number of such particle-horizonvolumes within our horizon volume today is ∼ [ t /t ( z )] .Let B be the average number of GW bursts coming fromeach horizon volume [ t ( z )] . Let the energy (as measuredat z ) of a typical burst ∆ E ( z ); by today that energy hasbeen redshifted to ∆ E = ∆ E z / (1 + z ). The total energytoday, within a Hubble volume, from all such bursts atredshift z is ∆ E B [ t /t ( z )] , and it satisfies∆ E B [ t /t ( z )] ∼ <
110 Ω GW τ . (6.5)We write “ ∼ < ” instead of “ (cid:39) ” because there could be othersignificant sources for Ω GW , besides this early-Universecontribution. A. Discovery space for GW memory from sourcesat high z The generalization of Eq. (5.10) to sources at arbitrary z is h mem ∼ α √ E ( z )(1 + z ) D L , (6.6)where ∆ E ( z ) is the locally measured energy loss and D L is the luminosity distance to the source (this followsfrom the propagation of GW-like perturbations in theRobertson–Walker spacetime [19] and from the defini-tion of D L ). The energy carried by those emitted wavestoday is ∆ E = ∆ E ( z ) / (1 + z ), while for high z we have D L ≈ τ (1 + z ). Thus we have h mem (cid:39) α E (1 + z ) τ . (6.7)It is instructive to determine the high- z version of Eq.(5.16) for the ratio SNR mem / SNR dir . The only change inthe derivation is the replacement d → τ (1 + z ), leadingto:SNR mem SNR dir (cid:39) × α SNR dir × (cid:20) z (cid:21)(cid:20) T sig s (cid:21) − (cid:20) T obs s (cid:21) (cid:20) dτ (cid:21)(cid:2) obs . (cid:3) − (6.8)By combining Eqs. (6.4), (6.5), and (6.7), we can con-strain SNR mem given B and Ω GW : h mem ∼ < α
80 Ω GW B × (cid:40) z ) / (cid:28) z (cid:28) z eq , (1+ z eq ) / z ) z (cid:29) z eq ; (6.9)the corresponding SNR follows from Eq. (5.12). Wewant to have a high probability of seeing one such signalwithin the observation time T obs . Since the local ratecan be shown to be R ∼ π ( B/τ )[ t /t ( z )] . Imposing R T obs ∼ > { SNR mem } (cid:39) α
125 (1 + z ) Ω GW τ T ( M p T obs ) / δt rms (cid:39) α (cid:20) z (cid:21)(cid:20) Ω GW − (cid:21)(cid:20) T obs s (cid:21) × obs ., (6.10)a factor of order (1 + z ) larger than the limit we derivedin Eq. (5.15) for sources at z ∼ <
1. We regard this as apromising result, since current constraints on Ω GW stillleave a great deal of room for possible discovery. Briefly, this can be shown by using Eq. (10) of [22], ap-proximating the term 4 π ( a r ) ≡ π (cid:0) a ( t − t ( z )) (cid:1) by4 π ( a t ) ≡ π ( τ ) and using ˙ n ( z )( dτ /dz )∆ z = ˙ n ( z )∆ τ =( B/τ )( t /t ( z )) . B. Discovery space for unmodeled GW bursts athigh z We now examine the prospects for detecting a GWburst from high z . The total energy emitted by such asource is∆ E ( z ) = ∆ E (1 + z ) (cid:39) π h f T sig D L ; (6.11)using Eq. (6.5) and D L ∼ τ (1 + z ), we then have h ∼ < × − Ω GW B ( f τ ) − ( f T sig ) − × (cid:40) (1 + z ) − / (cid:28) z (cid:28) z eq , (1 / z eq ) / (1 + z ) − z (cid:29) z eq . (6.12)Again, a high probability of observing a signal constrainsthe rate R according to R max { T sig , T obs } ∼ >
1, leading tomax { SNR dir } ∼ < (cid:20) Ω GW z (cid:21) / (cid:0) M p T obs (cid:1) / ( f δt rms )( f τ ) ≈ (cid:20) f − Hz (cid:21) − (cid:20) Ω GW − (1 + z ) (cid:21) / × obs . (6.13)This is basically the same limit we found for the largest-SNR burst at z <
1, but multiplied by the factor (1 + z ) − / . VII. CORRECTIONS FOR BEAMING ANDFOR GALACTIC SOURCES
So far our estimates of signal strengths have implic-itly assumed that the radiation is not strongly beamed.We have also implicitly assumed that detectable PTAsources will be extra-Galactic. In this section we brieflyshow how our estimates get modified if one drops theseassumptions. Both these issues were addressed by ZT82[16], but here we extend their considerations to large z . A. Modifications for highly beamed radiation
Assume that the GW energy is beamed into solid angle4 πF . To see how max { SNR } for “direct” radiation scaleswith F , we will take Ω GW and the total radiated energyto be fixed, which together imply a fixed rate density. Forthe case z ∼ <
1, we can approximate space as Euclidean,so the distance d to the closest source beaming in ourdirection scales as d ∝ F − / ; the observed h scales as h ∝ F − / /d ; and altogether h ∝ F − / . We see thatthe effect of beaming on max { SNR } is extremely weak;for instance, a beaming factor F = 10 − yields only afactor ∼ z ∼ < z (cid:29)
1, to account for beaming, on the right handside of Eq. (6.11) we would replace ∆ E with ∆ E /F .However the condition R T obs ∼ > R F T obs ∼ >
1, which leads to ∆ E ∝ Ω GW B − F . Thusthe F factors cancel, and beaming has basically no effecton max { SNR } for high- z sources. Note that our low- z and high- z upper limits, Eqs. (5.8) and (6.13) respec-tively, have slightly different character: for the formerwe maximize the SNR from the nearest detected source,for the latter we fix z and therefore luminosity distanceunder the constraint of detecting at least one source dur-ing the experiment.What about memory? The effect of beaming is negli-gible, since the memory component of GW strain is notbeamed, even when direct waves are. The dominant ef-fect is that the parameter α changes by a factor of orderone compared to the case of quadrupole emission. B. Modifications for Galactic sources
Throughout Secs. V and VI we have assumed that theEarth does not occupy a preferred location in the Uni-verse. However the Earth lies in the Galaxy; how mightthat modify our results? For sources at low z , universe,we showed in Sec. V A that, for fixed Ω GW , detectionSNR is maximized for sources whose event rate is onceper T obs in a Hubble volume. For a Galactic source to beobservable, this rate must increase to once per T obs perMilky-Way-like galaxy, or ∼ times greater. To main-tain the same Ω GW , the energy ∆ E radiated per eventmust decrease by a factor 10 . (We must also assumethat the Galaxy can sustain such a rate of events.) Onthe other hand, the distance to the extra-Galactic sourceis ∼ ∼
10 kpc for a randomly locatedGalactic source. For the direct radiation, h ∝ ∆ E / /d ,so the ratiomax { SNR
Galdir } max { SNR z ∼ } ∼ − / ∼ , (7.1)as was first shown by ZT82 [16]. Thus, besides beingintrinsically less plausible, putative Galactic sources in-crease max { SNR dir } by only an order of magnitude, com-pared to the z ∼ GW ∼ − (say), these putative Galacticexplosions would have to release ∼ M (cid:12) in GW energyroughly every ∼ h ∝ ∆ E/d , so we may estimatea ratiomax { SNR
Galmem } max { SNR z ∼ } ∼ − / ∼ − . . (7.2)1Finally, we note that if we had focused on sources inthe Local Group instead of just the Milky Way, the eventrate for sources outside the Milky Way would be domi-nated by Andromeda. Since Andromeda has roughly thesame mass as the Milky Way but is ∼
100 times furtheraway than our Galactic Center, the strongest such eventswould be ∼
100 times weaker than Galactic events.
VIII. CONCLUSIONS AND CAVEATS
In the paper we have constrained and characterized theGW discovery space of PTAs on the basis of energeticand statistical considerations alone. In Secs. V and VIwe showed that a PTA detection of GWs at frequenciesabove ∼ × − Hz would either be an extraordinary co-incidence, or have extraordinary implications; this effectresults from an analysis of fundamental constraints onpossible sources across the PTA sensitivity range, ratherthan deficiencies in PTA detection itself. We showed alsothat GW memory can be more detectable than directGWs, and that memory increasingly dominates the to-tal SNR of an event for sources at higher and higherredshifts; indeed, GW memory from high- z sources rep-resents a large discovery space for PTAs.Although we assumed modest beaming in our esti-mates, in Sec. VII we argued that even extreme beamingwould have a minor impact on detection SNRs. Similarly,although we assumed that the strongest GW sources dur-ing PTA observation would be extragalactic, our con-straint on max { SNR } rises only by a factor ∼
10 forGalactic sources. Throughout the paper we adopted anSNR scaling law valid for white pulsar noise; in Sec. IIwe explained, on the basis of toy model and of the ob-servational characterization of pulsar noise, why this wasappropriate.In Sec. IV we demonstrated how to properly incorpo-rate the effects of red noise in PTA searches, and wedemonstrated that the effects of periodic GWs between ∼ − . and 10 − . Hz band would not be degeneratewith small errors in the standard pulsar parameters, ex-cept in a few very narrow bands.Theoretical upper limits are akin to no-go theorems,and the authors are well aware that the history of thelatter in physics is replete with examples of results that,while strictly correct, turned out to be misleading be-cause their assumptions were overly restrictive. For thisreason, our chief motivation in doing this research was not to rule out possibilities, but to uncover promising butneglected areas of search space. With this in mind, wenow recall some of the assumptions that we have made,and point out some of the ways that Nature could beside-stepping them. • In this paper we assumed that the Earth is not ina preferred location in the Universe. In Sec. VII weconsidered the case in which relevant GW sourcesare clustered in galaxies, but still assumed that theEarth is not in some preferred location within theMilky Way. • Even if the Earth does not occupy a preferred loca-tion with respect to relevant GW sources, some mil-lisecond pulsars might do so. For instance, if two ormore pulsars are located in a globular cluster that also contains a BH binary with masses ∼ > M (cid:12) the correlated timing residuals due to the binary’sGWs impinging on the pulsars could well be de-tectable (see, e.g., [46]). • In this paper we assumed that at any redshift z there are no structures (such as phase-transitionbubbles) that are significantly larger than the con-temporaneous horizon size t ( z ). This is a reason-able way to incorporate causality constraints forprocesses that are not correlated on super-horizonscales to begin with, but it certainly does not holdfor all cases: for instance, inflation would imprintcorrelations on much larger scales. So a priori therecould arise strong GW sources that violate this as-sumption.It might be worthwhile to try to come up with reason-able physical scenarios that violate one or more of ourassumptions. Acknowledgments
CC gratefully acknowledges support from NSF GrantPHY-1068881. MV is grateful for support from the JPLRTD program. This work was carried out at the JetPropulsion Laboratory, California Institute of Technol-ogy, under contract to the National Aeronautics andSpace Administration. Copyright 2013 California Insti-tute of Technology. [1] M. V. Sazhin, Soviet Astronomy , 36 (1978), ADS.[2] S. Detweiler, Astrophys. J. , 1100 (1979), ADS.[3] R. W. Hellings and G. S. Downs, Astrophys. J. Lett. ,L39 (1983), ADS.[4] R. van Haasteren, Y. Levin, G. H. Janssen, K. Lazaridis,M. Kramer, B. W. Stappers, G. Desvignes, M. B. Purver,A. G. Lyne, R. D. Ferdman, et al., MNRAS , 3117 (2011), arXiv:1103.0576.[5] P. B. Demorest, R. D. Ferdman, M. E. Gonzalez, D. Nice,S. Ransom, I. H. Stairs, Z. Arzoumanian, A. Brazier,S. Burke-Spolaor, S. J. Chamberlin, et al., Astrophys. J. , 94 (2013), arXiv:1201.6641.[6] R. N. Manchester, G. Hobbs, M. Bailes, W. A. Coles,W. van Straten, M. J. Keith, R. M. Shannon, N. D. R. Bhat, A. Brown, S. G. Burke-Spolaor, et al., Pub. Astron.Soc. of Australia , e017 (2013), arXiv:1210.6130.[7] G. Hobbs, A. Archibald, Z. Arzoumanian, D. Backer,M. Bailes, N. D. R. Bhat, M. Burgay, S. Burke-Spolaor,D. Champion, I. Cognard, et al., Classical and QuantumGravity , 084013 (2010), arXiv:0911.5206.[8] A. Sesana, MNRAS , L1 (2013), arXiv:1211.5375.[9] A. Sesana, A. Vecchio, and M. Volonteri, MNRAS ,2255 (2009), arXiv:0809.3412.[10] D. R. B. Yardley, G. B. Hobbs, F. A. Jenet, J. P. W.Verbiest, Z. L. Wen, R. N. Manchester, W. A. Coles,W. van Straten, M. Bailes, N. D. R. Bhat, et al., MNRAS , 669 (2010), arXiv:1005.1667.[11] A. Sesana, ArXiv e-prints (2013), arXiv:1307.4086.[12] R. van Haasteren and Y. Levin, MNRAS , 2372(2010), arXiv:0909.0954.[13] J. M. Cordes and F. A. Jenet, Astrophys. J. , 54(2012), ADS.[14] A. H. Jaffe and D. C. Backer, Astrophys. J. , 616(2003), arXiv:astro-ph/0210148.[15] F. A. Jenet, G. B. Hobbs, K. J. Lee, and R. N. Manch-ester, Astrophys. J. Lett. , L123 (2005), arXiv:astro-ph/0504458.[16] M. Zimmermann and K. Thorne, in Essays in GeneralRelativity, a Festschrift for Abraham Taub , edited byF. J. Tipler (Academic Press, New York, 1980), pp. 139–155.[17] D. R. Lorimer, Living Reviews in Relativity , 8 (2008),arXiv:0811.0762.[18] I. H. Stairs, Living Reviews in Relativity , 5 (2003),arXiv:astro-ph/0307536.[19] M. Maggiore, Gravitational Waves: Volume 1: Theoryand Experiments , Gravitational Waves (Oxford Univer-sity Press, 2008).[20] G. L. Bretthorst, in
American Institute of Physics Con-ference Series (2001), vol. 567 of
American Institute ofPhysics Conference Series , pp. 1–28.[21] R. van Haasteren, Y. Levin, P. McDonald, and T. Lu,MNRAS , 1005 (2009), arXiv:0809.0791.[22] C. Cutler and J. Harms, Phys. Rev. D , 042001 (2006),arXiv:gr-qc/0511092.[23] W. W. Hager, SIAM Review , pp. 221 (1989),ISSN 00361445, URL .[24] J. M. Cordes and R. M. Shannon, ArXiv e-prints (2010),arXiv:1010.3785.[25] R. M. Shannon and J. M. Cordes, Astrophys. J. ,1607 (2010), arXiv:1010.4794. [26] J. P. W. Verbiest, M. Bailes, W. A. Coles, G. B. Hobbs,W. van Straten, D. J. Champion, F. A. Jenet, R. N.Manchester, N. D. R. Bhat, J. M. Sarkissian, et al., MN-RAS , 951 (2009), arXiv:0908.0244.[27] M. C. Begelman, R. D. Blandford, and M. J. Rees, Na-ture (London) , 307 (1980), ADS.[28] J. S. B. Wyithe and A. Loeb, Astrophys. J. , 691(2003), arXiv:astro-ph/0211556.[29] S. T. McWilliams, J. P. Ostriker, and F. Pretorius, ArXive-prints (2012), arXiv:1211.4590.[30] X. Siemens, J. Ellis, F. Jenet, and J. D. Romano, ArXive-prints (2013), arXiv:1305.3196.[31] J. Polchinski, ArXiv e-prints (2007), arXiv:0707.0888.[32] B. Allen, in Some Topics on General Relativity andGravitational Radiation , edited by J. A. Miralles, J. A.Morales, and D. Saez (1997), p. 3.[33] F. Dubath, J. Polchinski, and J. V. Rocha, Phys. Rev. D , 123528 (2008), arXiv:0711.0994.[34] T. Damour and A. Vilenkin, Phys. Rev. D , 064008(2001), arXiv:gr-qc/0104026.[35] T. Damour and A. Vilenkin, Phys. Rev. D , 063510(2005), arXiv:hep-th/0410222.[36] B. P. Abbott, R. Abbott, F. Acernese, R. Adhikari,P. Ajith, B. Allen, G. Allen, M. Alshourbagy, R. S. Amin,S. B. Anderson, et al., Nature (London) , 990 (2009),arXiv:0910.5772.[37] L. A. Boyle and A. Buonanno, Phys. Rev. D , 043531(2008), arXiv:0708.2279.[38] T. L. Smith, E. Pierpaoli, and M. Kamionkowski, Phys.Rev. Lett. , 021301 (2006), arXiv:astro-ph/0603144.[39] M. Archidiacono, S. Hannestad, A. Mirizzi, G. Raf-felt, and Y. Y. Y. Wong, ArXiv e-prints (2013),arXiv:1307.0615.[40] J. A. Ellis, M. A. McLaughlin, and J. P. W. Verbiest,MNRAS , 2318 (2011), arXiv:1107.4014.[41] G. B. Hobbs, R. T. Edwards, and R. N. Manchester,MNRAS , 655 (2006), arXiv:astro-ph/0603381.[42] F. B. Estabrook and H. D. Wahlquist, General Relativityand Gravitation , 439 (1975), ADS.[43] M. Favata, Classical and Quantum Gravity , 084036(2010), arXiv:1003.3486.[44] N. Seto, MNRAS , L38 (2009), arXiv:0909.1379.[45] M. S. Pshirkov, D. Baskaran, and K. A. Postnov, MN-RAS , 417 (2010), arXiv:0909.0742.[46] F. A. Jenet, T. Creighton, and A. Lommen, Astrophys.J. Lett.627