Binary Black Hole Waveform Extraction at Null Infinity
In this work, we present a work in progress towards an efficient and economical computational module which interfaces between Cauchy and characteristic evolution codes. Our goal is to provide a standardized waveform extraction tool for the numerical relativity community which will allow CCE to be readily applied to a generic Cauchy code. The tool provides a means of unambiguous comparison between the waveforms generated by evolution codes based upon different formulations of the Einstein equations and different numerical approximation.
pacs:04.20Ex, 04.25Dm, 04.25Nx, 04.70Bw
The study of gravitational waves will brighten unexplored features of our universe that are otherwise invisible to conventional astronomy and will increase our knowledge about the very nature of time and space (1). Gravitational wave detectors are already operating, and results from the first LIGO and Virgo collaboration were recently published in Nature (2). The signal predicted by numerical relativity will provide a template bank used for filtering the noise, indispensable to the success of gravitational wave detectors such as LIGO, Virgo, and LISA. The current sensitivity levels of the detectors will be improved substantially in next-generation detection estimated by 2015. Although existing simulations are sufficiently accurate for populating the parameter space in current searches of ground-based detectors, the new generation of advanced detectors will be 10-15 times more sensitive by 2015.
Ideally, the emitted gravitational wave signature should be extracted at spatial or null infinity. However, most present codes impose artificial, finite outer boundaries and are performing the waveform extraction at finite radius. This method introduces systematic errors, especially for higher modes, which is the main obstacle in reaching the desired accuracy. With the exception of Pretorius (3), who uses coordinates that compactify spatial infinity, all the other codes use a computational domain with a finite outer boundary and Sommerfeld-like approximate outer boundary conditions must be imposed, which introduce errors in the computation of gravitational waves. The choice of proper boundary conditions is complicated by gauge freedom and constraint preservation (4). The emitted gravitational wave signature is calculated at finite distance, using either the Newman-Penrose Weyl scalar (5), or the odd and even parity functions , in the Zerilli-Moncrief formalism (6). The strain of the wave used in detection is obtained performing one time integration from the Zerilli-Moncrief multipoles, and two time integrations from the Newman-Penrose curvature. The waveform is affected by gauge ambiguities which are magnified by the integration (7).
Cauchy-characteristic extraction (CCE) (8), which is one of the pieces of the CCM strategy (9), offers a means to avoid the error introduced by extraction at a finite world-tube. In CCE, the inner world-tube data supplied by the Cauchy evolution is used as boundary data for a characteristic evolution to future null infinity , where the waveform can be unambiguously computed by geometric methods. This characteristic initial-boundary value problem based upon a timelike world-tube (10) has been implemented as a mature evolution code, the PITT null code (11); (12), which incorporates a Penrose compactification of the space-time. By itself, CCE does not use the characteristic evolution to inject outer boundary data for the Cauchy evolution, which can be a source of instability in full CCM.
The PITT code has been tested to be second order convergent in a wide range of testbeds extending from the perturbative regime (13) to highly nonlinear single black hole spacetimes (12). However, in cases which require high resolution, such as the inspiral of matter into a black hole, the error in CCE has been a troublesome factor in the postprocessing phase (14). This has motivated a recent project (15) to increase the accuracy of the PITT code. Other results achieved with previous versions of the PITT have been recently reported (16); (17). Recently, the code underwent major improvements and corrections to previous versions, to improve accuracy and convergence (18).
Here we test this improved version of CCE on a realistic application involving a Cauchy evolution of the inspiral and merger of two equal mass non-spinning black holes. We use the same code specifications described in (15) except that the accuracy of angular derivatives has been increased to a 4th order finite difference approximation. The results presented here are a work in progress towards our goal to develop CCE as a reliable and accurate waveform extraction tool for the numerical relativity community. This paper addresses the first two objectives:
To create a robust and flexible interface between a binary black-hole Cauchy evolution code and a characteristic code for wave extraction at infinity.
To prove the robustness of the interface by performing precise computations of gravitational waveforms at infinity from binary black-hole, using this Cauchy-characteristic extraction approach.
We construct an interface that takes the Cartesian data from a Cauchy evolution and converts it into boundary data on a spherical grid for the characteristic evolution. The data are evolved to future null infinity, where it is used to compute the gravitational waveform. The flexibility of the interface is due to implementation of a spectral decomposition of data. This implementation has been tested with a realistic application involving a binary black-hole inspiral. In Sec. II we review the formalism underlying CCE, including enough details of the patching, evolution and extraction, to make clear the difficulties underlying the calculation of an accurate waveform at . In Sec. III we briefly describe the initial data for Cauchy and characteristic evolution. In Sec. IV, we present the details of the CCE interface which allows the data from a Cauchy evolution to be used as boundary data on an inner worltube for a characteristic evolution to , where the waveform is extracted. In Sec. V, we test the CCE interface by extracting the waveform from a Cauchy evolution of a binary black-hole inspiral and merger, and by comparing it to the waveform obtained by other standard method in current practice.
Ii The CCE Formalism
ii.1 Cauchy-Characteristic Patching
Characteristic data are provided by the Cauchy evolution on a world-tube , free initial data being given on the initial null hypersurface , which sets the metric on the entire initial cone (fig. 1).
The metric data from a Cauchy evolution are interpolated onto a timelike inner world-tube to extract the boundary data for the characteristic evolution. The characteristic evolution is embedded into the Cauchy evolution and is extending to future null infinity , where the waveform can be unambiguously computed using the geometric methods developed by Bondi et al (19), Sachs (20) and Penrose (21). The extraction process involves carrying out the complicated Jacobian transformation between the Cartesian coordinates used in the Cauchy evolution and the spherical null coordinates used in the characteristic evolution (the full details are given in (22).)
ii.2 Characteristic Evolution
The characteristic formalism is based upon a family of outgoing null hypersurfaces, emanating from some inner world-tube, which extend to infinity where they foliate into spherical slices. We let label these hypersurfaces, be angular coordinates which label the null rays and be a surface area coordinate. (fig. 2).
where is the Bondi-Sachs conformal 2-metric with . The code introduces an auxiliary unit sphere metric , with associated complex dyad satisfying . For a general Bondi-Sachs metric, the full nonlinear is uniquely determined by the dyad component , since the other dyad component is constrained by the determinant condition . The spherically symmetric case characterized by . We introduce the spin-weighted fields and , where
are also introduced to eliminate all second angular derivatives. In certain applications this has been found to give rise to increased accuracy by suppressing short wavelength error (24).
In this formalism, the Einstein equations decompose into hypersurface equations, evolution equations and conservation conditions on the inner world-tube. As described in more detail in (25); (26), the hypersurface equations take the form
is the curvature scalar of the 2-metric . Those equations have a hierarchical structure in such that the right hand sides, e..g. only depend upon previous variables and their derivatives intrinsic to the hypersurface.
The evolution equation takes the form
where, , , , and are nonlinear terms which vanish for spherical symmetry. Expressions for these terms as complex spin-weighted fields and a discussion of the conservation conditions are given in (8).
The characteristic Einstein equations are evolved in a domain between an inner radial boundary at the interior world-tube, and an outer boundary at future null infinity. The characteristic evolution code implements this formalism as an explicit finite difference scheme, based upon the compactified radial coordinate
so that at . Here is a parameter based upon the extraction world-tube, which in the CCE module is chosen as the radius of the extraction world-tube, as determined by in terms of the Cartesian coordiates used in the Cauchy evolution code. The boundary data for , , , , and on the world-tube supply the integration constants for a radial numerical integration of the hypersurface Einstein equations. The finite difference scheme for integrating the hypersurface and evolution equations is based on the marching equation for a spherically symmetric scalar field :
where the point N is the ”new” point in the evolution scheme, and is defined by the spherically symmetric version of the Bondi-Sachs metric given above. The evolution scheme in the full gravitational case used to determine the metric at the next point on the null hypersurfaces is modeled after this example (see (12); (24) for details).
ii.3 Gravitational Radiation Calculation
The theoretical derivation of the waveform at infinity is carried out in terms of an inverse surface-area coordinate , where at . In the resulting coordinates, the physical space-time metric (1) has the conformal compactification , where is smooth at and takes the form (10)
As described in (15), the Bondi news function and the Newman-Penrose Weyl tensor component which describe the waveform are both determined by the asymptotic limit at of the tensor field
constructed from the leading coefficients in an expansion of the metric in powers of
where and are the 2-dimensional curvature scalar and covariant derivative associated with .
The expansion coefficients , , and (all functions of and ) completely determine the radiation field. Before the gravitational radiation is calculated from the metric in the neighborhood of , it is necessary to determine the conformal factor relating to a unit sphere metric , i.e. to an inertial conformal Bondi frame (10) satisfying
The news function is directly computed by the code in terms of the computational coordinates , as opposed to the inertial coordinates on corresponding to an idealized distant observatory. The transformation to inertial coordinates proceeds first by introducing the conformally rescaled metric in which the cross-sections of have unit sphere geometry, in accord with (18). Then the rescaled null vector is the generator of time translations on , i.e. . The inertial coordinates thus satisfy the propagation equations
where in terms of the computational coordinates. The inertial coordinates are obtained by integrating (19), thus establishing a second pair of stereographic grid patches corresponding to . Then the news function is transformed into .
The Bondi news function is given by (20),
where is the Lie derivative with respect to . The Newman-Penrose Weyl tensor component is given by (21)
In the inertial Bondi coordinates, the expression for the news function (20) reduces to the simple form
and (21) reduces to the single term
This is related to the expression for the news function in inertial Bondi coordinates by
Equation (24) holds true in the linearized approximation of the Einstein equations. In the nonlinear case, the full expression for news and must be used in the code. This introduces additional challenges to numerical accuracy due to high order angular derivatives of and large number of terms.
Iii Initial Data
iii.1 Initial Cauchy Data
For the Cauchy evolution we used the LazEv code (27); (28) along with the Cactus framework cac () and Carpet (30) mesh refinement driver. LazEV is an eighth-order-accurate finite-difference code based upon the Baumgarte-Shapiro-Shibata-Nakamura (BSSN) formulation (31); (32) of Einstein’s equations, which deals with the internal singularities by the moving puncture approach (27); (33). Our simulation used 9 levels of refinement with finest resolution of , and outer Cauchy boundary at . The initial data consisted of a close quasicircular black-hole binary with orbital frequency , leading to more than a complete orbit before merger (see (34)). We output the metric data on the extraction world-tube every .
iii.2 Initial Characteristic data
The initial data for the characteristic evolution consist of the values of on the initial hypersurface . One way of supressing incoming radiation in the data would be to set the Newman-Penrose Weyl tensor component on the initial null hypersurface. For a perturbation of the Schwarzschild metric, this condition implies no incoming radiation in the linearized approximation. However, in order to avoid shocks arising from incompatibility with the Cauchy data on the extraction world tube (with given by 10), we also need to require that and are continuous. In the linearized approximation, the condition that implies that . The combination of those requirements leads to , which would imply that at . For technical simplicity we avoid this complication by initializing according to
which matches the Cauchy data and the derivatives at and is consistent with asymptotic flatness. Since this choice of vanishes at infinity, the initial slice of has a unit sphere metric so that the conformal factor has the simple initialization .
Iv Computational interface
We have designed an interface that takes Cartesian grid data from a Cauchy evolution and converts it into boundary data for characteristic evolution on a spherical grid extending to . We treat each component of the Cauchy metric as a scalar function in the Cartesian coordinates which are used in the evolution.
In order to make the interface as flexible as possible for future development as a community tool for waveform extraction, we have based it upon a spectral decomposition of the Cauchy data in the region between two world tubes or radii and , where is the Cartesian coordinate radius. Then at a given time , we decompose in terms of Tchebychev polynomials of the second kind and spherical harmonics , where are related to in the standard way. The Tchebychev polynomials are conventionally defined as functions on the interval . Here we map them to the interval by the transformation
where the extraction shell thickness is determined by the number of Tchebychev polynomials used. (In tests of binary black holes with mass M we use a relatively small range , a larger value of would be needed if the range were expanded). Thus, for , we expand
For the applications to waveform extraction given in this paper, it is sufficient to consider , where , and , where . The coefficients then allow us to reconstruct a spherical harmonic decomposition of each component of the Cauchy metric on the extraction world-tube , i.e.
This decomposition is carried out at a sequence of Cauchy time steps , where is chosen to be much smaller than the characteristic time scale of the problem but, for purposes of economy, larger than the time step used for the Cauchy evolution. A fifth-order polynomial interpolation is carried out locally over the to provide characteristic boundary data at any time in analytic form.
The extraction module also requires the derivatives and at the extraction world-tube. The -derivative is constructed by a fourth-order-accurate finite-difference stencil using the surrounding Cauchy times . The -derivative is obtained analytically, at each time level , by differentiation of the Tchebychev polynomials.
The spherical harmonic interpolator from the Cartesian to the spherical coordinates is part of the extraction module, but its resolution is controlled by the Cauchy evolution.
The stereographic coordinates used to label the outgoing null rays in the Bondi metric are matched to the spherical coordinates induced by the Cartesian Cauchy coordinates on the extraction worldtube by a standard transformation, using the conventions in (23). The value of the surface-area coordinate in the Bondi-Sachs metric is obtained on the extraction world-tube from the 2-determinant of the Cartesian metric on the surfaces . As a result the radius of the Bondi coordinate on the extraction world-tube. The metric has to be calculated at a common value of the surface coordinate , because the original Cauchy extraction was at constant R. In order to make this calculation possible, the transformation from Cartesian coordinates to Bondi-Sachs coordinates is carried out via an intermediate Sachs coordinate system (20) where is an affine parameter along the outgoing null rays. The affine freedom allows us to set on the extraction world-tube . After carrying out the Jacobian transformation from to , the Cartesian metric and its first derivatives at the extraction world-tube provide a first-order Taylor expansion in (about ) of the null metric in Sachs coordinates. The corresponding Taylor expansion of the metric in Bondi-Sachs coordinates then follows from the computed value of and at , which are obtained from the 2-determinant of the Cartesian metric. In order to obtain a first-order Taylor expansion for the Bondi metric variable , the hypersurface equation (4) must be used to evaluate at the extraction world-tube. All other metric variables are then initialized consistent with second order accuracy. Taylor expansions are also needed to start up the radial integration equations for the auxiliary variables (3) used to convert angular derivatives to first-order form. These expansions are obtained from applying the -operator to the Taylor expansion of the underlying metric. This is a complicated process because the operator intrinsic to the extraction world-tube is not the same as the operator intrinsic to the Bondi spheres (see (18) for a discussion). The low order intermediate Taylor expansions limit the accuracy of the result. A new approach that avoids entirely the use of the Taylor expansion and gives better accuracy is presented in (18). In the original approach, used for this paper, the resulting Taylor expansion of the evolution variables is used to fill the points of the Bond-Sachs grid to start the integration of the characteristic hypersurface and evolution equations (4) - (9). The integration proceeds from the extraction world tube to on a radial grid based upon the compactified -coordinate (10).
Domain of dependence considerations place a constraint between the characteristic time step and the size of the characteristic grid analogous to the CFL condition for the Cauchy evolution. For an estimate, consider the Minkowski space case with the conformally rescaled metric
where the unit sphere metric takes the form
The past light cone is determined by
For typical characteristic grid parameters, , the resulting restriction is
For a Cauchy simulation of a binary black-hole system of total mass with timestep (sufficient to describe the typical frequencies of a binary system), (31) leads to
for the choice of characteristic timestep . The corresponding number of radial gridpoints must roughly satisfy . This places no limit of practical concern on the resolution of the characteristic evolution even for the small extraction radius . Thus, for purposes of CCE, there are no demanding CFL restrictions.
The interface was debugged and calibrated using the analytic Schwarzschild metric in Kerr-Schild coordinates ,
We present results for the characteristic extracted waveform either in terms of , related to the Bondi news by in the linearized regime, or, when comparing to the perturbative waveform, in terms of the Newman-Penrose component . The relationship between the Cauchy and the characteristic waveforms is: . We decompose the signal in spherical harmonic modes but, for illustrative purposes, we concentrate on the dominant and sub-dominant modes. The Cauchy data were given at the extraction radii . The relationship between the Cauchy radius and the characteristic world-tube radius is: . The characteristic extraction module was run with the following specifications: angular gridpoints = radial gridpoints = , and timestep , where . The test was run until , using -order accurate angular derivatives, on stereographic patches with circular boundaries and angular dissipation (see (15) for details on how the angular dissipation is added to the evolution equation 9). The results are shown for the highest resolution. Table 1 gives the convergence rates for the world-tube variables obtained with a small extraction radius at a time corresponding to the peak of the signal (). The rates are given for the real and the imaginary part. All quantities are very close to second order convergent, including , which is the term which determines the waveform.
We are reporting only first-order convergence rates at future null infinity for the Bondi News and the Weyl complonent , but the error is relatively small ( during the late inspiral). From an extraction point of view, these errors are smaller than the error in the Cauchy code data and are of little concern. The data are not convergent at early time when high frequencies dominate the error. For a thorough analysis of the causes for the first-order accurate results and major improvements to the code see (18).
Figure 4 compares the imaginary and real parts of the mode of the Cauchy with the complex conjugate of the mode of the characteristic . We obtain very good amplitude match. Also, we observe improved phase agreement as the extraction radius is increased (from to ), because the phase error in is reduced with increased extraction radius.
Figure 5 compares the imaginary and real parts of the mode of the Cauchy with the complex conjugate of the mode of the characteristic . Here we see two effects. First the improved phase agreement as , but also an attenuation of the amplitude due to dissipation of higher-order modes. Also, the noise is apparent for the mode.
Figure 6 compares the amplitudes and the phases between the absolute value of the mode of the Cauchy extracted at , and the absolute value of the mode of the characteristic for the same extraction radius, at the highest resolution (). The difference in amplitude is relatively small, maximum of the Cauchy amplitude in the wave zone.
Figure 7 compares the real part of the modes of the characteristic extracted at three different extraction radii: , and . The waveform extracted at has the biggest amplitude, and a very small attenuation of the signal with the radius is observed.
We have presented here a method for interfacing outer boundary data from a Cauchy evolution with inner boundary data for a characteristic evolution so that the waveform can be accurately extracted at infinity. We have demonstrated how the PITT null code can be interfaced with the LazEv code, which is a finite-difference BSSN code, to produced calibrated waveforms from a binary black-hole inspiral. The extraction interface has been implemented as a thorn in the Einstein Computational Toolkit ein () In this paper we are reporting only preliminary results (see (18) improvements). Although we are aware of the deficiencies in the characteristic waveform extraction tool presented here, there is pressing interest from several numerical relativity groups to apply the tool to extract waveforms from binary black-hole inspirals.
- B.S. Sathyaprakash, and B. F. Shultz, Living Rev. Rel., 12 2 (2009), arXiv:0903.0338
- The LIGO Scientific Collaboration and the Virgo Collaboration, Nature, Vol. 460 990 (2009)
- F. Pretorius, Class. Quant. Grav, 22 425 (2005), gr-qc/0407110
- L. Lehner, and O. M. Moreschi, Phys. Rev. D, 76 124040 (2007), arXiv:0706.1319
- E. T. Newman, and R. Penrose, J. Math. Phys., 3 566 (1962)
- V. Moncrief, Ann. Phys, 88 323-342 (1974)
- L. Lindblom, B. J. Owen, D. A. Brown,
- N.T. Bishop, R. Gómez, L. Lehner, and J. Winicour, Phys. Rev. D, 54 6153 (1996).
- Winicour, J. (2005), Living Rev. Relativity 8, 10.
- L.A. Tamburino and J. Winicour, Phys. Rev. 150 1039, 1966.
- R.A. Isaacson, J.S. Welling and J. Winicour, J. Math. Phys. 24 1824 (1983).
- N. T. Bishop, R. Gómez, L. Lehner, M. Maharaj and J. Winicour, Phys. Rev. D,56 6298 (1997).
- M. Babiuc, B. Szilágyi, I. Hawke, and Y. Zlochower, Class. Quantum Grav. 22 5089 (2005)
- N. T. Bishop, R. Gómez, S. Husa, L. Lehner, J. Winicour Phys.Rev. D 68, 084015 (2003)
- M.C. Babiuc, N.T. Bishop, B. Szilágyi and J. Winicour, Phys. Rev. D, 79 084011 (2009).
- C. Reisswig, N. T. Bishop, D. Pollney abd B. Szilágyi, Phys. Rev. Lett., 95, 221101 (2009).
- C. Reisswig, N. T. Bishop, D. Pollney abd B. Szilágyi, Class. Quantum Grav. 27, 075014 (2010).
- M.C. Babiuc, B. Szilágyi, J. Winicour and Y. Zlochower [arXiv:1011.4223]
- H. Bondi, M.J.G. van der Burg and A.W.K. Metzner, Proc. R. Soc. A 269 21 (1962).
- R.K. Sachs, Proc. R. Soc. A 270 103 (1962).
- R. Penrose, Phys. Rev. Letters, 10 66 (1963).
- N.T. Bishop, R. Gómez, L. Lehner, B. Szilágyi, J. Winicour and R. A. Isaacson, Black Holes, Gravitational Radiation and the Universe, Kluwer Academic Publishers, Dordrecht, 1998
- R. Gómez, L. Lehner, P. Papadopoulos and J. Winicour, Class. Quantum Grav. 14 977, 1997.
- R. Gómez, Phys. Rev. D 64, 1–8 (2001).
- J. Winicour, J. Math. Phys. 24 1193 (1983).
- J. Winicour, J. Math. Phys. 25 2506 (1984).
- M. Campanelli, C. O. Lousto, P. Marronetti and Y. Zlochower, Phys. Rev. Lett. 96, 111101 (2006) [arXiv:gr-qc/0511048].
- Y. Zlochower, J. G. Baker, M. Campanelli and C. O. Lousto, Phys. Rev. D 72, 024021 (2005) [arXiv:gr-qc/0505055].
- Cactus Computational Toolkit home page: http://www.cactuscode.org.
- E. Schnetter, S. H. Hawley and I. Hawke, Class. Quant. Grav. 21, 1465 (2004) [arXiv:gr-qc/0310042].
- M. Shibata and T. Nakamura, Phys. Rev. D 52, 5428 (1995).
- T. Baumgarte and S. L. Shapiro, Phys. Rev. D 59, 024007 (1999).
- J. G. Baker, J. Centrella, D.-I. Choi, M. Koppitz, and J. van Meter, Phys.Rev.Lett. 96, 111102 (2006).
- M. Campanelli, C. O. Lousto and Y. Zlochower, Phys. Rev. D 74, 041501 (2006) [arXiv:gr-qc/0604012].
- Einstein Computational Toolkit home page: http://einsteintoolkit.org.