# Gravitational memory in binary black hole mergers

###### Abstract

In addition to the dominant oscillatory gravitational wave signals
produced during binary inspirals, a non-oscillatory component arises
from the nonlinear “memory” effect, sourced by the emitted
gravitational radiation. The memory grows significantly during the
late inspiral and merger, modifying the signal by an almost
step-function profile, and making it difficult to model by
approximate methods. We use numerical evolutions of binary black
holes to evaluate the nonlinear memory during late-inspiral, merger
and ringdown. We identify two main components of the signal: the
monotonically growing portion corresponding to the memory, and an
oscillatory part which sets in roughly at the time of merger and is
due to the black hole ringdown. Counter-intuitively, the ringdown
is most prominent for models with the lowest total spin.
Thus, the case of maximally spinning black holes *anti-aligned* to the
orbital angular momentum exhibits the highest
signal-to-noise (SNR) for interferometric detectors. The largest memory
offset, however, occurs for highly spinning black holes, with an estimated
value of in the maximally spinning
case. These results are central to determining the detectability of
nonlinear memory through gravitational wave interferometers and
pulsar timing array measurements.

###### Subject headings:

black hole physics — gravitational waves## 1. Introduction

Accelerating massive bodies generate gravitational waves (GWs), a potentially rich source of astrophysical information which a number of large experiments have been designed to mine over the next decade. The principle sources for the interferometric detectors (LIGO, Virgo, LISA) are orbiting bodies, for which the dominant signal is oscillatory. However, it has long been known that GWs will also contain non-oscillatory features, resulting from a net change in time derivatives of the multiple moments of the system (Zel’Dovich & Polnarev, 1974; Braginskii & Thorne, 1987). A nonlinear contribution also results from the interaction of the waves with themselves (Payne, 1983; Christodoulou, 1991; Blanchet & Damour, 1992). Evaluating these modes involves an integral in time, requiring knowledge of the entire past history of the spacetime. As such, they have been given the name “memory” or “hereditary” components of the GW.

The nonlinear memory is sourced by the emitted GWs. For a relativistic binary inspiral, the typical profile is of a slow growth over time which sees a rapid increase during the late inspiral and merger, and reaches a constant value as the merger remnant rings down and ceases to emit. The non-oscillatory nature of the memory suggests that it will not be a dominant feature in interferometric detectors, for which drifts in the background metric are factored out. However, the step induced during merger may be observable in pulsar timing measurements, as has been noted in a recent set of papers (Pshirkov et al., 2009; van Haasteren & Levin, 2009; Seto, 2009).

Post-Newtonian (PN) calculations provide an accurate estimate of the memory to within a few orbits of the merger (Kennefick, 1994), however fail to model the important merger phase where the effect is largest. In a series of papers, Favata (2009a, c, b, 2010) has included an estimate of the merger contribution for binary black holes (BHs) by using an effective-one-body (EOB) model which was tuned using the results of numerical simulations. Recent progress in numerical models of BH spacetimes provides the opportunity to measure the effect directly, though as outlined in Favata (2009c), due to the low amplitude of the memory, and systematic error in standard techniques of wave measurement in numerical relativity, the problem is a challenging one.

In this paper, we apply a newly developed evolution code (Pollney et al., 2009b, a) and wave extraction techniques (Reisswig et al., 2009a, 2010) to carry out fully relativistic binary BH simulations to directly measure the nonlinear GW memory through inspiral, merger and ringdown in a gauge-invariant and mathematically unambiguous notion. We concern ourselves with equal-mass, spinning binaries. For these models, the memory is contained in the spin-weighted spherical harmonic components of the GW strain, . These modes, previously studied in the context of head-on collisions (Anninos et al., 1995), exhibit two dominant effects. The first is a non-oscillatory term which rises as the inspiral progresses, and is associated with the memory. During the merger, an oscillatory signal is superposed onto these modes, induced by the ringdown of the BH remnant. After the ringdown has subsided, the modes are offset from their original value. We find that for the dominant memory mode, this offset is largest in the case of the merger of aligned, maximally spinning BHs, as might be expected given that these are the strongest gravitational emitters (Reisswig et al., 2009b). Interestingly, we observe the strongest oscillatory ringdown signal in the case of lowest total spin (maximally spinning BHs anti-aligned with the orbital angular momentum). As such, and as opposed to the observability of the dominant mode, these models produce a more visible mode to interferometric detectors. We conclude by estimating the strain offset which would be visible to a pulsar timing array.

## 2. Numerical methods

We integrate the BSSNOK formulation of the vacuum Einstein equations (see, e.g., Alcubierre 2008) numerically, using the Llama evolution code (Pollney et al., 2009a, b) within the Cactus framework (Goodale et al., 2003). The evolution equations are discretized via finite differences in space on a grid composed of multiple patches with locally adapted coordinates in the wave zone. Adaptive mesh-refinement is implemented via the Carpet grid driver (Schnetter et al., 2004) on the central patch surrounding the BHs. Importantly for the measurement of small non-oscillatory features such as the memory, the artificial grid outer boundary is causally disconnected from the wave measurements.

We determine initial data for binary BHs by the conformal puncture data method (Brandt & Brügmann, 1997), and evolve them using the “moving puncture” approach (Baker et al., 2006; Campanelli et al., 2006) whereby the choice of gauge (Alcubierre et al., 2003) prevents the spacetime slicing from encountering the curvature singularity (Hannam et al., 2007).

Crucial for determining the memory, we compute the GW signal of the source using the method of Cauchy-characteristic extraction (CCE), which computes unambiguous and coordinate invariant signals at future null infinity, , corresponding to a detector far removed from the source (Bishop et al., 1999; Babiuc et al., 2005; Winicour, 2005; Reisswig et al., 2009a, 2010). By this technique, metric data is collected on a world-tube at finite radius, and a formulation of the full Einstein equations along null hypersurfaces is used to transport the signal to .

We measure two quantities that encode the gravitational signal at . The first is the Newman-Penrose scalar , the Weyl curvature component with slowest falloff in asymptotically flat spacetimes (Newman & Penrose, 1962). Alternatively, we measure the Bondi “news” (Bondi et al., 1962; Sachs, 1962), defined by

(1) |

where is a derivative operator defined on an outgoing null geodesic, , and is the Newman-Penrose shear scalar (Newman & Penrose, 1962). These scalars, and , can be evaluated directly from the local curvature on the sphere at and in the naturally defined gauge of an asymptotic observer (Bondi et al., 1962; Sachs, 1962). To calculate the GW strain which is observed by a detector, these quantities need to be numerically integrated in time:

(2) |

We expand in terms of a basis of spin-2 spherical harmonics . The dominant non-oscillatory hereditary component is contained in the mode,

(3) |

where is the sphere at and is the distance to the source.

Since the simulations are necessarily of finite length, the time integration of , Eq. (2), leaves a constant to be determined, corresponding to the initial value of the signal at the start time of the simulation. For purely oscillatory signals (e.g., ) we can fit this by adjusting the post ringdown value to zero, or by some averaging procedure over multiple inspiral wavelengths. This is not possible in the case of the monotonically growing memory. Instead, we approximate an initial value of the mode by matching to a PN calculation, using the 3PN expansion for derived in Favata (2009c). For spinning binaries we incorporate 2.5PN spin contributions to the frequency evolution (Blanchet et al., 2006). The initial offset of the numerical waveforms is performed by shifting its amplitude so that it fits against the PN models, over an interval from to . For instance, the PN solution for in the non-spinning case is plotted as a dotted line in the lower panel of Figure 1. In Table 1, we report both the overall amplitude offset of the ringdown signal in the mode matched to the PN inspiral, as well as the offset from the point which is computed entirely from numerical data.

## 3. Characteristics of the memory modes

In Figure 1, we show the mode during late inspiral and merger for the case of an equal-mass, non-spinning binary using data from the model presented in Pollney et al. (2009b) and evaluated at in Reisswig et al. (2009a). The upper panel shows the development of the component of the news, , whereas the lower shows its integral, . The time coordinate places at the position of the amplitude peak of the dominant mode (not shown), corresponding roughly to the time of merger. Two main features are present. The first is a slow monotonic growth of the mode during the inspiral and plunge phase. This corresponds to the non-oscillatory memory contribution. At , there is a local maximum in , and shortly afterward a prominent oscillatory signal sets in. This signal represents the binary merger and subsequent BH ringdown. It is not strictly an aspect of the GW “memory”, as it is not dependent on the entire past history of the spacetime, however it is superposed on the growing memory signal.

The lower plot shows the strain, determined via Eq. 2. The superposition of the memory and ringdown signals results in a notable kink in the waveform near , and exponentially decaying oscillation before reaching the steady-state amplitude of

(4) |

We measure approximately 4-th order convergence in when comparing three different grid resolutions for the same model, and we estimate an error on the order of in . The much larger error of in is an upper limit of the estimate on the total error composed of integration error (2) and error in the PN fit for all of the binary models considered in this paper.

Perturbative results provide a prediction for the observed quasi-normal modes (QNM) of (Berti et al., 2009). We perform a least-squares fit of the numerical results to a function of the form

(5) |

over the ring-down portion, , with fitting parameters, , and , , , and arrive at a frequency

(6) |

using the measured mass of the remnant of (Pollney et al., 2009b). Berti et al. (2009) provide a tabulated value of for the mode of a BH with the measured dimensionless spin . Thus, we find a difference of

(7) |

representing an error of less than . We list the final spin and the associated QNM frequencies of the binary systems for all simulations performed in Table 1.

Finally we have also examined modes for and . For this model, the overall amplitudes of the final memory are , and . However, these modes are difficult to distinguish from the numerical error, and we conclude that the higher order modes contribute less than to the total memory.

## 4. Memory from spinning binaries

We have performed simulations using the code infrastructure described in Section 2 in the two-dimensional parameter space of non-precessing equal-mass binaries with spin that has also been studied in Reisswig et al. (2009b); Rezzolla et al. (2008). We have focused on initial data configurations along the two main axes in the space of individual dimensionless BH spins aligned with the orbital angular momentum. The first set of models have equal but opposite spins, ; The second have identical spins, . The model details are listed in Table 1.

Figure 2 plots the evolution of the modes for the two sequences. The upper panel shows models for which the spins of the individual BHs are anti-aligned, , so that the spacetimes have the same total angular momentum. The nearly identical lines are consistent with previous observations that the modes also do not vary appreciably along this direction of parameter space (Vaishnav et al., 2007; Reisswig et al., 2009b), as well as the radiated energies (Reisswig et al., 2009b), at least not to within the error-bars as stated in the latter reference.

The lower plot shows the modes for configurations with equal and aligned spins, . The amplitude of the memory offset increases with higher total spin. This is expected as the memory is sourced by the emitted GWs which are known to be more energetic in the higher spin models (Reisswig et al., 2009b; Lousto et al., 2009). The total offsets, , are listed in Table 1.

We find that the amplitude of the ringdown is largest for models with lowest total spin, . In these cases, spin-orbit coupling increases the final angle of impact so that the merger is more nearly head-on. The ringdown is suppressed in the high spin model. Fits to the QNM frequencies are listed in Table 1 and show good agreement with perturbative results. However, for simulations with the weakness of the ringdown contributes to larger variation in the estimates, particularly if the fitting interval is varied, though still within of the analytic values.

The memory is well estimated by a quartic polynomial in the total spin ,

(8) |

determined by a least-squares fit to the measured data. Figure 3 plots as a function of spin for the aligned models, , as well as the leading order multipole moment estimate resulting from the radiated energy (see Favata (2009a), Eq. (5))

(9) |

with

(10) |

The measured radiated energy, , can also be estimated using the quadratic fit developed in Reisswig et al. (2009b), verified here with the plotted data points corresponding to the new and more accurate simulations (see Table 1). The radiated energy provides a reasonably good estimate of the memory offset, though with a slight under-prediction. The measured memory values also appear to have a somewhat stronger dependence on the total spin. Extrapolating to the extremal case, we estimate a maximum in the case of aligned spins. Fits to for generic spin configurations (Lousto et al., 2009) suggest that this will indeed be the maximum value for arbitrary spins. Additional simulation in the high-spin regime would be required to establish the accuracy of the extrapolation.

() | |||||||||

## 5. Detectability of the memory modes

The modes are small (on the order of ) compared to the dominant components of the GW signal. However, they exhibit both non-oscillatory growth and ringdown features, which present some interesting aspects for detection and identification of these modes. We discuss the prospects for observation in both interferometers, and pulsar timing arrays, in the following sections.

### 5.1. Interferometers

The signal-to-noise ratio (SNR) within a given detector, , in GW searches based on matched filtering is given by Flanagan & Hughes (1998):

(11) |

where is the sensitivity curve of the instrument. Restricting the integral to the mode, we can determine its contribution to the overall SNR, testing against the proposed advanced-LIGO (Shoemaker, 2010) and LISA noise curves (LISA-wiki, 2010).

Table 1 lists the maximum SNR attained over a range of masses for a given model. For advanced-LIGO, the optimal total mass is (independent of the model and invariant under the distance), i.e. a binary, for which we list the SNR at a reference distance of (the SNR scales linearly with the distance). The LISA results refer to an optimal (redshifted) total mass of at a luminosity distance corresponding to a redshift , where models suggest at least one expected merger event per year (Sesana et al., 2005). By assuming a minimum SNR of for detection, it will be possible to observe the mode in LISA out to the given distance. They are unlikely to be visible with significant event rates within the planned advanced ground-based detectors, though future generation experiments may have prospects.

For both detectors, and as opposed to the dominant mode (Reisswig et al., 2009b), we note that the highest SNR occurs for the cases with spins anti-aligned to the orbital angular momentum which is due to their stronger ringdown. In the aligned cases, the weaker ringdown causes a much flatter high-frequency fall-off (Figure 4), and this effect overrides the non-oscillatory step within these detectors.

An important aspect of the matched-filtering algorithm is the construction of templates against which signals are compared. We examine the influence of neglecting the contribution by computing the mismatch between a signal, , containing all of the modes, with one in which the mode is left out. For advanced-LIGO, the mismatch for the model is largest in the mass range , but negligible, with a magnitude of . For the model, the mismatch grows to be on the order of . Thus, although the mode induces a notable offset in the wave signal when viewed in the time domain, it makes little difference to the detection algorithms.

### 5.2. Pulsar timing arrays

An alternative proposal for measuring GWs may be more sensitive to the non-oscillatory step-function nature of the memory. Precise, distant clocks will experience a residual in time-of-arrival measurements which can be associated with GWs (Estabrook & Wahlquist, 1975; Sazhin, 1978; Detweiler, 1979). Several such projects are in development (Hobbs et al., 2009; Manchester, 2008; Jenet et al., 2009; Lazio, 2009). While their principal aim is to examine the stochastic background, recently a number of papers have pointed out the potential for observing non-oscillatory step function burst signals, such as the memory (Pshirkov et al., 2009; van Haasteren & Levin, 2009; Seto, 2009).

Primary candidates for detection are supermassive BH binary mergers. While there is a great deal of ongoing work understanding the nature and evolution of such systems, if there is a tendency for the bodies to align, then the results obtained here are directly applicable. For instance, Dotti et al. (2009) indicate that in gas-rich mergers, the spins are aligned on a short timescale, with spin magnitudes in the range to . Assuming a binary source at a distance of which has aligned spins with a magnitude of , the results of Section 4 lead to an angle-averaged memory offset of:

(12) |

(The corresponding constants for zero and maximally spinning bodies are and , respectively.) The results of Pshirkov et al. (2009) and van Haasteren & Levin (2009) indicate that a burst will be observable over a 10-year observation period at amplitude , suggesting that the merger of a pair of with would be observable at .

## References

- Alcubierre (2008) Alcubierre, M. 2008, Introduction to Numerical Relativity (Oxford, UK: Oxford University Press)
- Alcubierre et al. (2003) Alcubierre, M., Brügmann, B., Diener, P., Koppitz, M., Pollney, D., Seidel, E., & Takahashi, R. 2003, Phys. Rev. D, 67, 084023
- Anninos et al. (1995) Anninos, P., Hobill, D., Seidel, E., Smarr, L., & Suen, W.-M. 1995, Phys. Rev. D, 52, 2044
- Babiuc et al. (2005) Babiuc, M., Szilágyi, B., Hawke, I., & Zlochower, Y. 2005, Class. Quantum Grav., 22, 5089
- Baker et al. (2006) Baker, J. G., Centrella, J., Choi, D.-I., Koppitz, M., & van Meter, J. 2006, Phys. Rev. Lett., 96, 111102
- Berti et al. (2009) Berti, E., Cardoso, V., & Starinets, A. O. 2009, Class. Quant. Grav., 26, 163001
- Bishop et al. (1999) Bishop, N., Isaacson, R., Gómez, R., Lehner, L., Szilágyi, B., & Winicour, J. 1999, in Black Holes, Gravitational Radiation and the Universe, ed. B. Iyer & B. Bhawal (Kluwer, Dordrecht, The Neterlands), 393
- Blanchet et al. (2006) Blanchet, L., Buonanno, A., & Faye, G. 2006, Phys. Rev., D74, 104034
- Blanchet & Damour (1992) Blanchet, L., & Damour, T. 1992, Phys. Rev., D46, 4304
- Bondi et al. (1962) Bondi, H., van der Burg, M. G. J., & Metzner, A. W. K. 1962, Proc. R. Soc. London, A269, 21
- Braginskii & Thorne (1987) Braginskii, V. B., & Thorne, K. S. 1987, Nature, 327, 123
- Brandt & Brügmann (1997) Brandt, S., & Brügmann, B. 1997, Phys. Rev. Lett., 78, 3606
- Campanelli et al. (2006) Campanelli, M., Lousto, C. O., Marronetti, P., & Zlochower, Y. 2006, Phys. Rev. Lett., 96, 111101
- Christodoulou (1991) Christodoulou, D. 1991, Phys. Rev. Lett., 67, 1486
- Detweiler (1979) Detweiler, S. 1979, Astrophys. J., 234, 1100
- Dotti et al. (2009) Dotti, M., Volonteri, M., Perego, A., Colpi, M., Ruszkowski, M., & Haardt, F. 2009, Mon. Not. R. astr. Soc., 1795
- Estabrook & Wahlquist (1975) Estabrook, F. B., & Wahlquist, H. D. 1975, General Relativity and Gravitation, 6, 439
- Favata (2009a) Favata, M. 2009a, J. Phys. Conf. Ser., 154, 012043
- Favata (2009b) —. 2009b, Astrophys. J., 696, L159
- Favata (2009c) —. 2009c, Phys. Rev., D80, 024002
- Favata (2010) —. 2010, arXiv:1003.3486
- Flanagan & Hughes (1998) Flanagan, E. E., & Hughes, S. 1998, Phys. Rev. D, 57, 4535
- Goodale et al. (2003) Goodale, T., Allen, G., Lanfermann, G., Massó, J., Radke, T., Seidel, E., & Shalf, J. 2003, in Vector and Parallel Processing – VECPAR’2002, 5th International Conference, Lecture Notes in Computer Science (Berlin: Springer)
- Hannam et al. (2007) Hannam, M., Husa, S., Pollney, D., Brugmann, B., & O’Murchadha, N. 2007, Phys. Rev. Lett., 99, 241102
- Hobbs et al. (2009) Hobbs, G., et al. 2009, arXiv:0911.5206
- Jenet et al. (2009) Jenet, F., et al. 2009, arXiv:0909.1058
- Kennefick (1994) Kennefick, D. 1994, Phys. Rev., D50, 3587
- Lazio (2009) Lazio, J. 2009, arXiv:0910.0632
- LISA-wiki (2010) LISA-wiki. 2010, LISA parameter estimation wiki, http://www.tapir.caltech.edu/dokuwiki/lisape:home
- Lousto et al. (2009) Lousto, C. O., Campanelli, M., & Zlochower, Y. 2009, arXiv:0904.3541
- Manchester (2008) Manchester, R. N. 2008, AIP Conf. Proc., 983, 584
- Newman & Penrose (1962) Newman, E. T., & Penrose, R. 1962, J. Math. Phys., 3, 566, erratum in J. Math. Phys. 4, 998 (1963)
- Payne (1983) Payne, P. N. 1983, Phys. Rev., D28, 1894
- Pollney et al. (2009a) Pollney, D., Reisswig, C., Dorband, N., Schnetter, E., & Diener, P. 2009a, Phys. Rev., D80, 121502
- Pollney et al. (2009b) Pollney, D., Reisswig, C., Schnetter, E., Dorband, N., & Diener, P. 2009b, arXiv:0910.3803
- Pshirkov et al. (2009) Pshirkov, M. S., Baskaran, D., & Postnov, K. A. 2009, arXiv:0909.0742
- Reisswig et al. (2009a) Reisswig, C., Bishop, N. T., Pollney, D., & Szilagyi, B. 2009a, Phys. Rev. Lett., 103, 221101
- Reisswig et al. (2010) —. 2010, Class. Quant. Grav., 27, 075014
- Reisswig et al. (2009b) Reisswig, C., Husa, S., Rezzolla, L., Dorband, E. N., Pollney, D., & Seiler, J. 2009b, Phys. Rev. D, 80, 124026
- Rezzolla et al. (2008) Rezzolla, L., et al. 2008, Astrophys. J, 679, 1422
- Sachs (1962) Sachs, R. 1962, Proc. Roy. Soc. London, A270, 103
- Sazhin (1978) Sazhin, M. V. 1978, Vestn. Mosk. Univ. Fiz. Astron., 19N1, 118
- Schnetter et al. (2004) Schnetter, E., Hawley, S. H., & Hawke, I. 2004, Class. Quantum Grav., 21, 1465
- Sesana et al. (2005) Sesana, A., Haardt, F., Madau, P., & Volonteri, M. 2005, Astrophys. J., 623, 23
- Seto (2009) Seto, N. 2009, arXiv:0909.1379
- Shoemaker (2010) Shoemaker, D. 2010, Advanced LIGO anticipated sensitivity curves, Tech. Rep. LIGO-T0900288-v3, LIGO Scientific Collaboration
- Vaishnav et al. (2007) Vaishnav, B., Hinder, I., Herrmann, F., & Shoemaker, D. 2007, Phys. Rev., D76, 084020
- van Haasteren & Levin (2009) van Haasteren, R., & Levin, Y. 2009, arXiv:0909.0954
- Winicour (2005) Winicour, J. 2005, Living Rev. Relativ., 8, 10, [Online article]
- Zel’Dovich & Polnarev (1974) Zel’Dovich, Y. B., & Polnarev, A. G. 1974, Soviet Astronomy, 18, 17