Horizon Surface Gravity in Corotating Black Hole Binaries

First, we consider an isolated, stationary black hole of mass $M$, angular momentum $S$, and horizon surface area $A$. Hereafter, we will frequently use the irreducible mass $\mu=\sqrt{A/(16\pi)}$ instead of the horizon surface area $A$ in various formulas.

In a seminal work, Friedman, Uryū and Shibata Fr.al.02 generalized the zeroth and first laws discussed above to spacetimes containing multiple black holes and/or a generic distribution of perfect fluid matter sources with compact support, having in mind the application of their general results to binary systems of black holes and/or neutron stars.³³3See Ref. Ur.al.10 for the generalization of their work to the magnetohydrodynamic case.

Gralla and Le Tiec GrLe.13 have shown how the zeroth and first laws of black hole mechanics can be extended to the case where a Kerr black hole is perturbed by a small orbiting compact object (a “moon”) on the corotating orbit, the unique circular equatorial orbit whose orbital frequency matches the angular velocity of the Kerr horizon. The natural framework for such a problem is that of black hole perturbation theory. Therefore, in this section we consider a Kerr black hole of mass $M$ and spin $S$, perturbed by a corotating point particle of mass $m\ll M$, and we work to linear order in the small mass ratio $q\equiv m/M$.

Next, we consider binary systems of compact objects, modelled as (spinning) point particles with constant masses $m_{a}$. Although the concept of a point particle does not make sense in exact GR GeTr.87—the closest thing to a point particle in GR is a black hole—, in the context of approximation schemes such as PN theory Bl.14, Einstein’s equation can be coupled to a distributional matter source in a meaningful manner, provided that a regularization scheme (e.g. dimensional regularization) is used to handle the divergent self-field of each particle.

Since a point particle has no “spatial” extension, and in particular no horizon, a putative zeroth law of mechanics is meaningless for binary systems of point particles. Still, a first law can be derived by assuming the existence of a global helical Killing field $k^{\alpha}$, and by taking a formal point-particle limit in the generalized first law of Fr.al.02, in the case of two self-gravitating balls of perfect fluid. The resulting first law of binary point-particle mechanics reads Le.al.12

where $z_{a}\equiv-u_{a}^{\alpha}k_{\alpha}$ is the “redshift” of particle $a$, a conserved orbital quantity associated with the helical Killing symmetry, with $u_{a}^{\alpha}$ the four-velocity of particle $a$. Choosing coordinates adapted to that isometry, such that $k^{\alpha}={(\partial_{t})}^{\alpha}+\Omega\,{(\partial_{\phi})}^{\alpha}$ holds everywhere, and not merely in a neighbourhood of infinity, the redshift of particle $a$ is simply $z_{a}=\mathrm{d}\tau_{a}/\mathrm{d}t$, the ratio of the proper times elapsed along the worldline of the particle and along that of a distant static observer. As expected, the first integral associated with the variational formula (19) is

The first law (19) was later recovered, and extended to spinning point particles, using the canonical ADM Hamiltonian framework Bl.al.13. Assuming that the conservative dynamics of a binary system of spinning point particles with masses $m_{a}$ and canonical spins $\mathbf{S}_{a}$ derives from an autonomous Fokker-type two-body Hamiltonian $H_{\text{F}}(\mathbf{x}_{a},\mathbf{p}_{a};m_{a},\mathbf{S}_{a})$, with $\mathbf{x}_{a}$ and $\mathbf{p}_{a}$ the canonical positions and momenta, the first law for circular orbits and spins aligned or anti-aligned with the orbital angular momentum reads⁶⁶6Strictly speaking, the mass appearing in the left-hand side of Eq. (21) is the on-shell value of the local, autonomous Fokker Hamiltonian $H_{\text{F}}(\mathbf{x}_{a},\mathbf{p}_{a};m_{a},\mathbf{S}_{a})$.

Interestingly, while the precession frequencies $\Omega_{a}$ of the spins are given by the partial derivatives of the two-body Hamiltonian with respect to the particle’s spins Da.al.08, it was shown in Ref. Bl.al.13 that (at least to linear order in the spins) the redshifts $z_{a}$ are given by the partial derivatives of the Hamiltonian with respect to the particle’s masses:

Using the canonical ADM framework, the first law (19) for circular orbits was also extended to generic bound (i.e. eccentric) orbits in Refs. Le.15; BlLe.17. These first laws for binary systems of (spinning) point masses, as well as the associated first integrals, have been explictly checked to hold true up to 3PN order included, as well as the logarithmic contributions at the 4PN and 5PN orders.

Motivated by the first law (12) for binary systems of corotating black holes, the authors of Ref. Bl.al.13 suggested modelling each spinning point particle as an equilibrium black hole in a tidal environment, endowed with an “irreducible mass” $\mu_{a}$ and a “proper rotation frequency” $\omega_{a}$. Then, by analogy with the first law (5) for a single, isolated black hole, each spinning point particle with mass $m_{a}$ and spin $S_{a}$ was assumed to obey the variational formula

To obtain explicit expressions for the coefficients $c_{a}=(\partial m_{a}/\partial\mu_{a})_{S_{a}}$ and $\omega_{a}=(\partial m_{a}/\partial S_{a})_{\mu_{a}}$, it was further assumed that each spinning point particle obeys a Christodoulou mass formula of the type $m_{a}^{2}=\mu_{a}^{2}+S_{a}^{2}/(4\mu_{a}^{2})$. Then, for each particle the “response coefficient” $c_{a}$, which is the point-particle analogue of the normalized surface gravity $4\mu\kappa$ of an isolated black hole, can be expressed as a function of e.g. $\mu_{a}$ and $\omega_{a}$, and similarly for the mass $m_{a}$ and spin $S_{a}$. From Eqs. (7) and (8), we immediately find

Substituting for Eq. (23) into the first law (21) for binary systems of spinning point particles yields a formula that involves the variations $\delta\mu_{a}$ and $\delta S_{a}$ of the irreducible masses and the spins, namely $\delta M_{\text{ADM}}=\Omega\,\delta J+\sum_{a}\left[c_{a}z_{a}\,\delta\mu_{a}+(\Omega_{a}-\Omega+z_{a}\omega_{a})\,\delta S_{a}\right]$. By analogy with the first law (12) for binary systems of corotating black holes, the proper rotation state of each spinning particle is constrained to a corotation state if, and only if, the coefficients in front of the variations $\delta S_{a}$ vanish, yielding the corotation conditions (for $a=1,2$)

Physically, this means that the redshifted proper rotation frequency of each spinning point particle, $z_{a}\omega_{a}$, must be equal to the circular orbit frequency $\Omega$, as seen in a frame rotating at the angular rate $\Omega_{a}$ with respect to an inertial frame of reference.

When the corotating condition (25) is imposed, the first law (21) for binary systems of spinning point particles reduces to

Because the first law (12) for corotating black hole binaries can be written in the equivalent form $\delta M_{\text{ADM}}=\Omega\,\delta J+\sum_{a}4\mu_{a}\kappa_{a}\,\delta\mu_{a}$, it appears that the renormalized redshifts $c_{a}z_{a}$ in corotating point particle binaries are analogous to the normalized surface gravities $4\mu_{a}\kappa_{a}$ in corotating black hole binaries. This will be further discussed in the next Sec. II.6. To conclude this section, we point out that the first integral relation associated with the variational first law (26) reads, as expected,

We can now compare the first law (12) for corotating black hole binaries to the first law (26) for corotating point particles binaries, as well as the first integrals (13) and (27). Clearly, if the irreducible masses $\mu_{a}$ are identified for each body, then the surface gravity $\kappa_{a}$ of each corotating black hole is analogous to the redshift $z_{a}$ of each corotating point particle, and more precisely

Recall that the coefficient $c_{a}$ introduced in Eq. (23) is nothing but the normalized horizon surface gravity of an isolated black hole of irreducible mass $\mu_{a}$ and horizon angular velocity $\omega_{a}$, say $4\mu_{a}\bar{\kappa}_{a}$ [see Eqs. (8) and (24c)]. Therefore, the analogy (28) can equivalently be written as

Since the redshift of a particle satisfies $z_{a}\leqslant 1$, the analogy (29) suggests that in corotating black hole binary systems the tidal interaction decreases the surface gravity of each black hole (with respect to its value in isolation). This agrees with the prediction of the perturbative calculation discussed in Sec. II.3.

In the large separation limit, i.e. when $\Omega\to 0$, we know that $z_{a}\to 1$ for each corotating particle, while $\kappa_{a}\to\bar{\kappa}_{a}=1/(4\mu_{a})$ for each corotating black hole. Indeed, in that limit, the spin and proper rotation frequency of each body must vanish to maintain the corotation. To test the soundness of the analogy (28) beyond the large separation limit, we will compute the normalized black hole surface gravities $4\kappa_{a}\mu_{a}$ in Sec. III, as well as the 4PN expansions of the particle’s renormalized redshifts $c_{a}z_{a}$ in Sec. IV, and compare those results, expressed a functions of the circular-orbit frequency $\Omega$, in Sec. VI.

Interestingly, for black hole binaries with large mass ratios, the analogy (29) is supported by a comparison, at the Hamiltonian level, of the expression for the surface gravity of a black hole tidally perturbed by a corotating moon, and that for the redshift of the most massive particle (say body $2$). Indeed, recalling Eqs. (15a) and (22), we find the analogy

where we introduced the Hamiltonian $\tilde{H}\equiv M+zH$, including the mass $M$ of the background hole. The occurence of the moon’s redshift $z$ in $\tilde{H}$ comes from the fact that the Hamiltonian $H$ of a test mass in Kerr is parameterized by the proper time $\tau$, while the Fokker-type two- body Hamiltonian $H_{\text{F}}$ is parameterized by the coordinate time $t$.

Similarly, by comparing the expression (15b) for the perturbed angular velocity of the black hole to the corotation condition (25) with the formulas (22) for the point particle’s proper rotation frequency, we note a formal analogy between $\bar{\omega}$ and $\omega_{2}$, namely

In Secs. IV and V below, we will find that the large mass-ratio limits of the PN expansions of $z_{2}(\Omega)$ and $\omega_{2}(\Omega)$ are in full agreement with the PN expansions of the perturbative expressions for $(\kappa/\bar{\kappa})(\Omega)$ and $\bar{\omega}(\Omega)$, thus supporting the analogies (30) and (31).

On the other hand, still for black hole binaries with large mass ratios, the analogy (29) can be made rigorous for the smaller compact object (say body $1$). Indeed, using a method of matched asymptotic expansions, the author of Ref. Po3.15 recently proved that the horizon surface gravity $\kappa_{1}$ of a small black hole of mass $m$ in the tidal environment of a large black hole of mass $M\gg m$ is closely related to the redshift $z\equiv z_{1}$ of a point mass $m$ moving like a test particle in a certain effective metric satisfying the linearized vacuum Einstein equation. More precisely, to linear order in the small mass ratio $q=m/M$, it can be shown that Po3.15

At this order of approximation, the smaller black hole can be taken to be nonrotating, such that $m=\mu_{1}$, $\bar{\kappa}_{1}=1/(4\mu_{1})$ and $c_{1}=1$. Hence the result (32) proves that, at least to linear order in the mass ratio and for the smaller body $1$, Eq. (28) is an identity, and not a mere analogy. See also the discussion in Ref. Zi.al.16.

The numerical results presented here are similar to those previosly published in Gr.al.02; Gr.10; Ur.al.12. The basic properties of such data are as follows. We employ the 3+1 formalism of GR so that the four-dimensional line element reads

where $N$ is the lapse, $B^{i}$ the shift vector and $\gamma_{ij}$ the three-dimensional spatial metric. Spatial indices are lowered and raised by means of $\gamma_{ij}$ and its inverse $\gamma^{ij}$, so that e.g. $B_{i}\equiv\gamma_{ij}B^{j}$. Circularity of the orbit is enforced by imposing the existence of a helical Killing vector $k^{\alpha}$; recall Eq. (9). From a practical point of view, the equations are written in the corotating frame, such that $k^{\alpha}=\left(\partial_{t}\right)^{\alpha}$ and the time derivatives of the various fields vanish.

Our corotation setup allows us to consider binary black hole systems with an exact helical Killing symmetry, such that the apparent horizon is a Killing horizon, and the surface gravity is well defined over each component of the horizon. This should be contrasted to the recent work of Ref.Zi.al.16, where the authors consider more generic and astrophysically realistic binary configurations, but do not have a well-defined (geometrical and unambiguous) notion of black hole surface gravity.

Far from the system, it is demanded that the spacetime be asymptotically flat. However, as discussed in Sec. II.2, this is not strictly compatible with the existence of a global helical symmetry. Indeed, if a true helical Killing vector existed, the system would have an infinite lifetime. This fact, coupled with the emission of gravitational waves would lead to diverging metric quantities, being the signature of this accumulated emission. There are two possible approaches to cure that problem. Either one can (to some extend) relax the Killing condition, or one needs to find a way to suppress the gravitational radiation. In this work, we opt for the latter solution. There are various ways to proceed. The most widely used is to impose that the solution be spatially conformally flat. The spatial metric $\gamma_{ij}$ is then simply related to the three-dimensional flat metric $f_{ij}$ through a conformal factor $\Psi$, namely

Of course, the price to pay is that not all of Einstein’s equations can be exactly satisfied.

With conformal flatness, the spacetime depends on five fields: the lapse $N$, the conformal factor $\Psi$ and the shift vector $B^{i}$. Those fields obey a subset of Einstein’s equations, namely the trace of the evolution equation (35a), the Hamiltonian constraint (35b) and the momentum constraints (35c). The slicing (i.e. the choice of time) is done via maximal slicing, for which the trace of the extrinsic curvature tensor vanishes. The spacelike coordinates are fixed by the choice of corotating coordinates and the conformal flatness approximation. The resulting system of equations is then

where $D$ denotes the covariant derivative compatible with the three-dimensional Euclidean metric $f_{ij}$, while $A^{ij}\equiv D^{\langle i}B^{j\rangle}/N$ is the conformal extrinsic curvature tensor.

The presence of two corotating black holes is enforced by imposing the following boundary conditions on two spheres corresponding to the apparent horizons:

where $s^{i}$ is the unit normal to each hole (with respect to the metric $f_{ij}$). Equation (36a) is a free choice of the lapse on the boundaries of the physical domain; typically one uses $n_{0}=0.1$. This basically comes from the fact that maximal slicing is a differential gauge choice, defined uniquely up to the boundaries. This was discussed is some length in Refs. Ca.al2.06; GoJa.06. Equation (36b) enforces that the two spheres are apparent horizons, and (36c) that the black holes are corotating. Asymptotic flatness is explicitly recovered by requesting that $N=1$, $\Psi=1$ and $B^{i}=\Omega\left(\partial_{\varphi}\right)^{i}$ at infinity. (Recall that we work in the corotating frame).

For each separation, the constant angular velocity $\Omega$ is determined by the numerical code to ensure the equality of the ADM and Komar masses:

This condition is closely related to relativistic generalizations of the virial theorem GoBo.94; Fr.al.02, and it has been used extensively to compute binary black hole configurations, e.g. in Gr.al.02; Gr.10; Ur.al.12. The five elliptic equations (35) are solved by spectral methods implemented via the Kadath library. Bispherical coordinates are used and the algorithm is the same as the one already described in Refs. Gr.10; UrTs.11. The version of the code used for this work is more general in the sense that it can deal with black holes of different masses. When doing so, one needs to find the location of the rotation axis of the binary. This is done numerically and it is obtained from the fact that the total (ADM) linear momentum of the binary must vanish. The linear momentum can be computed as a surface integral at spatial infinity, and the condition can be rewritten as

where $\vec{e}_{y}$ is the unit vector (with respect to the flat metric) along the $y$-axis and $\mathrm{d}S_{i}$ is the flat surface element at infinity. The various axes are defined as follows : the $x$-direction joins the center of the holes and the $z$-direction is the one of the rotation axis, so that the orbital plane is $z=0$. The $y$-direction is then in the orbital plane, perpendicular to the line joining the holes. Due to the symmetries of the problem, only the $y$ component of the momentum is not trivially zero. It vanishes only when the axis is at the right location on the $x$-axis.

The coordinate radii of the holes are also unknowns. The sequences of corotating binaries are obtained by requesting that the irreducible mass of each hole is held fixed (see Sec. III.3.2 below). The size of the numerical domain is allowed to vary and the numerical code finds the appropriate radii so that the masses have the desired values. More precisely, the irreducible masses are defined by $\mu_{a}\equiv\sqrt{A_{a}/(16\pi)}$, where $A_{a}$ is the surface area of the $a$^th hole, given by the following surface integral over a cross section $S_{a}$ of the corresponding horizon:

Since, by construction, the horizons have a vanishing expansion, their areas $A_{a}$ do not depend on a particular choice of slicing As.al.00. We computed quasi-equilibrium initial data for irreducible mass ratios $q_{\mu}\equiv\mu_{1}/\mu_{2}\in\{1,1/2,1/10\}$.

In order to compute the horizon surface gravity (10) of each hole, we employ a formula given in terms of the 3+1 quantities. A suitable expression can be found, e.g., in Eq. (10.10) of Ref. GoJa.06. There, the normal vector is normalized with respect to the spatial metric $\gamma_{ij}$, whereas the vector $s^{i}$ appearing in our Eqs. (36b) and (36c) is normalized with respect to the flat metric $f_{ij}$. It follows that the two normal vectors differ by a factor of $\Psi^{2}$. Moreover, the term proportional to $\ell^{\alpha}\nabla_{\alpha}\ln N$ in Eq. (10.10) of Ref. GoJa.06 is related to the evolution of the lapse function, and it vanishes for our quasi-stationary configurations. Taking all of that into account, the 3+1 expression for the surface gravity reduces to

Next, we explore the variations of some quantities of physical interest with respect to the numerical resolution. The three-dimensional physical space is split into several (typically a dozen) computational domains. In each domain the number of spectral coefficients is the same with respect to the three space dimensions. We denote by $N$ this number and compute binary configuration for $N=9$, $11$, $13$ and $15$ coefficients.

The first test has to do with the zeroth law of binary black hole mechanics, being the fact that the horizon surface gravity (10) of each hole is constant. To measure this we define the normalized root-mean-square variation (or relative standard deviation) of the surface gravity $\kappa_{a}$ of body $a$ by

At each point one can compute $\kappa_{a}$ by Eq. (40), and the various averages appearing in (41) are taken using all the collocation points located on the horizon of each black hole. Figure 1 shows $\bar{\sigma}(\kappa_{a})$ as a function of $N$ for one configuration with equal masses, one with a mass ratio $q_{\mu}=1/2$ and two with $q_{\mu}=1/10$. All these configurations correspond to separations close to the innermost circular orbit (ICO), namely the configuration of minimum binding energy; see Sec. III.3.2 below. For a mass ratio $q_{\mu}=1/10$ we show an additional configuration with a smaller separation, corresponding to an orbital frequency $\mu\Omega=0.13$, where $\mu\equiv\mu_{1}+\mu_{2}$ is the total irreducible mass. (In that case the ICO is located at $\mu\Omega=0.025$; see Tab. 1). All the curves exhibit the same behavior, being the convergence of $\bar{\sigma}(\kappa_{a})$ to a finite value as the resolution increases. It means that the error on the zeroth law is not coming from a numerical error, since that would depend on the resolution. Rather, it arises from the physical assumptions used when computing the binary configurations, such as our use of the CFC approximation, which implies that only a subset of Einstein’s equations are solved.

Another test is provided by the first integral relation (13) associated with the variational first law (12). For each configuration, we compute the global quantities $M_{\text{ADM}}$, $J$ and $\Omega$, as well as the quasi-local quantities $\mu_{a}$ and $\langle\kappa_{a}\rangle$. Figure 2 shows the violation of this algebraic relation, as measured by the difference

More precisely, Fig. 2 displays the variation of the dimensionless quantity $\varepsilon/\mu$ as a function of resolution, for the same configurations as the ones of Fig. 1. The data for $q_{\mu}=1$ and $1/2$ exhibit a similar behavior, where the error saturates. The convergence with resolution is slightly less clean than for $\bar{\sigma}(\kappa_{a})$, as the curves show a slight increase of the error for $N=15$. This increase is probably due to the difficulty in computing very accurately the global quantities $M_{\text{ADM}}$ and $J$. At this stage it is unclear how the error (42) would behave at even higher resolution, which is unfortunately out of reach of our current computational resources. Nevertheless, from Fig. 2, it is clear that the deviation from the first integral does not converge to zero as the resolution increases, but instead to a small but finite value. From this point of view it is very similar to what is observed for the zeroth law.

The data for $q_{\mu}=1/10$ at a separation close to the ICO (red triangles) exhibits a different behavior where no saturation is visible. The plausible explanation is that, in this case, the deviation from the first integral is so small that the error is still dominated by the effects of a finite resolution. This is consistent with results from Fig. 1, where the values of $\bar{\sigma}(\kappa_{a})$ for $q_{\mu}=1/10$ at the ICO are significantly smaller that in all the other cases depicted. This explanation is also supported by data for a smaller separation (green diamonds). The black holes being closer, the errors induced by the CFC approximation are larger, as confirmed by Fig. 1. In that case, the error on the first law is no longer dominated by resolution effects, such that the curve shows a plateau.

Unless stated otherwise, all the results shown below correspond to the highest resolution, namely $N=15$ for irreducible mass ratios $q_{\mu}=1$ and $1/2$. For $q_{\mu}=1/10$, most of the results were computed while using $N=13$. A few configurations with $N=15$ were also computed to confirm that they were almost indistinguishable from the $N=13$ ones.

For nonspinning (NS) binaries, the relations $z^{\text{ns}}_{a}(\Omega)$ were computed up to 2PN and later to 3PN order in Refs. De.08; Bl.al.10. The leading and next-to-leading order logarithmic contributions at 4PN and 5PN orders, which originate from gravitational-wave tails, were then obtained in Refs. Bl.al2.10; Da.10; Le.al.12. Those results were all derived starting from the definition of the redshift in terms of the spacetime metric $g_{\alpha\beta}(x)$ evaluated at the coordinate location $x^{\alpha}=y_{a}^{\alpha}$ of the particle $a$, namely

Moreover, the nonlogarithmic 4PN terms were also computed, from the known 4PN contribution to the circular-orbit binding energy JaSc.13; BiDa.13 combined with the first law of mechanics (19) for binary systems of nonspinning point masses Le.al.12; Le.15; BlLe.17. All these results are valid for arbitrary mass ratios. Up to 4PN order, we have

where $x\equiv(m\Omega)^{2/3}$ is the small frequency-related PN parameter, with $m\equiv m_{1}+m_{2}$ the total mass, while $\nu\equiv m_{1}m_{2}/m^{2}$ is the symmetric mass ratio and $\Delta\equiv(m_{2}-m_{1})/m=\sqrt{1-4\nu}$ is the reduced mass difference. (We assume that $m_{1}\leqslant m_{2}$.) The expression for $z_{2}^{\text{ns}}$ is easily deduced by setting $\Delta\to-\Delta$ in Eq. (IV.1). Note the occurence of a logarithmic running and of Euler’s constant $\gamma_{\text{E}}$ at 4PN order.

Since our goal is to compute the redshifts for corotating spinning point-particle binaries, we need to account for the spin contributions to the redshifts. At the relative 4PN accuracy, however, it will prove sufficient to include the leading order 1.5PN and the next-to-leading order 2.5PN spin-orbit (SO) contributions. Up to next-to-leading order, the SO contributions to the redshift of body 1 read Bl.al.13