当前位置:首页 >> >>

Dynamics of perturbations of rotating black holes

Dynamics of perturbations of rotating black holes
William Krivan(1,2) , Pablo Laguna(1) , Philippos Papadopoulos(1,3) , and Nils Andersson(4)
Department of Astronomy & Astrophysics and Center for Gravitational Physics & Geometry Penn State University, University Park, PA 16802, USA (2) Institut f¨ ur Astronomie und Astrophysik Universit¨ at T¨ ubingen, D-72076 T¨ ubingen, Germany (3) Max-Planck Institut f¨ ur Gravitationsphysik D-14473 Potsdam, Germany (3) Department of Physics Washington University, St Louis, MO 63130, USA (February 7, 2008)

arXiv:gr-qc/9702048v1 24 Feb 1997

We present a numerical study of the time evolution of perturbations of rotating black holes. The solutions are obtained by integrating the Teukolsky equation written as a ?rst-order in time, coupled system of equations, in a form that explicitly captures its hyperbolic structure. We address the numerical di?culties of solving the equation in its original form. We follow the propagation of generic initial data through the burst, quasinormal ringing and power-law tail phases. In particular, we calculate the e?ects due to the rotation of the black hole on the scattering of incident gravitational wave pulses. These results may help explain how the angular momentum of the black hole a?ects the gravitational waves that are generated during the ?nal stages of black hole coalescence. 04.30.Nk

Typeset using REVTEX 1


Numerical relativity, i.e. the numerical construction of solutions to the Einstein equations, is a rapidly advancing ?eld. Reliable multi-dimensional simulations of astrophysically relevant scenarios should become possible in the near future. The systematic exploration of the non-linear content of Einstein’s theory is of major importance for this accelerated progress. In the pursuit of numerically solving non-linear systems, intrinsic accuracy tests, such as convergence, help establish the reliability of the solution and its proximity to the continuum limit. However, recurrent experiences suggest that tests, extrinsic to the discretization process, are also valuable tools for assessing the reliability of complicated numerical algorithms. One versatile tool, providing a number of such extrinsic tests, is the investigation of the linear regime of non-linear systems. In the appropriate limit, numerical solutions of non-linear systems should agree with the considerably simpler perturbative solutions. This testing approach is bound to have fundamental importance to numerical relativity as the ?eld progresses. Perturbative methods enjoy renewed attention in their own right. The emphasis now is on the application of the standard perturbative methodology to particular realistic systems that are thought to be of interest to the emerging ?eld of gravitational wave astronomy. Perhaps the most important one among those systems is the collision of inspiraling black hole binaries. The key premise of the perturbative approach to black hole collisions is that, during the last stages of the coalescence (the close limit), the spacetime can be accurately approximated as that of a single perturbed black hole. A particular subclass of such problems is the study of black hole collisions as perturbations of Schwarzschild spacetimes; these are situations (e.g. head-on collisions) for which the end product of the collision is a slowly or non-rotating black hole. Mathematically they are described by the gauge invariant Zerilli function [1]. The numerical solution of the Zerilli equation does not pose di?culties and is routinely employed in estimating aspects of gravitational wave physics in the context of black holes [2]. The generalization of this framework to rotating black holes is physically relevant; it is likely that, during the inspiral collision of black holes, the system will not be able to radiate all of its angular momentum and will leave behind a single, rotating black hole. The work in this paper constitutes an important step towards the goal of generalizing the close-limit treatment of black hole collisions to that of perturbations about a single, rotating black hole. The general framework to achieve this goal requires: (A) the construction of appropriate initial data describing two orbiting black holes, (B) the subtraction of the Kerr background from the initial data, to read o? the initial perturbations, and (C) the evolution of these perturbations. In this paper, we concentrate on the last point, namely the evolution of gravitational perturbations of the Kerr spacetimes. Previously, we investigated the dynamics of scalar ?elds in the background of rotating black holes [3] (from here on referred to as Paper I). We considered both slowly and rapidly rotating black holes. For slowly rotating black holes the background geometry can be treated as a perturbed Schwarzschild spacetime, with the angular momentum per unit mass playing the role of a perturbative parameter. The study was centered on the late-time dynamics of the scalar ?eld and showed that, for rotating black holes of arbitrary angular momentum, the late-time dynamics is dominated by the lowest allowed multipole (l = 0). Here we undertake the development of a numerical scheme for the study of gravitational perturbations of Kerr 2

black holes in the time domain. The main objective of this work is to follow the propagation of generic initial data through the burst, quasinormal ringing and power-law tail phases of the evolution. In particular, we are interested in investigating the e?ects of black hole rotation on the scattering of incident gravitational wave pulses. A characterization of such e?ects provides an indication of the role that rotation plays on signals produced during the ?nal stages of black hole coalescence.

At ?rst instance, a direct derivation of the equations governing the perturbations of Kerr spacetimes is to consider perturbations of the metric. This path, however, leads to gauge-dependent formulations. A theoretically attractive alternative is to examine curvature perturbations. Using the Newman-Penrose formalism, Teukolsky [4,5] derived a master equation governing not only gravitational perturbations (spin weight s = ±2) but scalar, two-component neutrino and electromagnetic ?elds as well. In Boyer-Lindquist coordinates and with the use of the Kinnersley null tetrad [6], this master evolution equation reads ? (r 2 + a2 )2 4Mar M (r 2 ? a2 ) ? a2 sin2 θ ?tt Ψ ? ?tφ Ψ ? 2s r ? + ia cos θ ?t Ψ ? ? ? 1 a2 1 ?φφ Ψ (1) ?θ (sin θ?θ Ψ) + ? + ??s ?r ?s+1 ?r Ψ + sin θ ? sin2 θ a(r ? M ) i cos θ + 2s ?φ Ψ ? s2 cot2 θ ? s Ψ = 0, + ? sin2 θ

where M is the mass of the black hole, a its angular momentum per unit mass, and ? ≡ r 2 ? 2Mr + a2 . For gravitational perturbations, the function Ψ is given in terms of the Weyl tensor tetrad components ψ0 and ψ4 . That is, Ψ = ψ0 for s = +2 and Ψ = ρ?4 ψ4 for s = ?2, with ρ = ?1/(r ? ia cos θ). The Teukolsky equation (1) reduces to the Bardeen-Press equation [7] in the limiting case of non-rotating black holes. For the case s = 0, it yields the equation for a scalar wave propagating in a Kerr background, a system which we studied numerically in Paper I. A result of great importance for perturbation theory was the discovery by Teukolsky that, when viewed in the frequency domain, Eq. (1) is separable in the r and θ coordinates (separation of the azimuthal angular dependence is always possible). To our knowledge, most of the work on the dynamics of perturbations of Kerr spacetimes has been performed in the frequency domain (or under the assumption of a harmonic time dependence). This has certainly been the case for studying quasinormal mode frequencies [8–11], wave scattering [12], and the motion of test particles in the Kerr background [13–15]. Here we are interested in the time integration of physical initial data, possibly from the inspiral collision of binary black holes. The data may be analytic approximate solutions to the linearized constraints, or numerical solutions of the non-linear constraints. In principle, one can Fourier transform the initial data and perform the evolution of such data for each frequency. Afterwards, the data are transformed back to the time domain if needed. From the computational point of view, however, this approach is far more expensive than it is to solve the problem in the time domain. For several reasons the number of 3

frequencies that one has to consider is orders of magnitude higher than the number of angular components required to resolve the θ-direction: While the angular directions are bounded, the time direction is only bounded from below by the initial data surface. Focussing on the investigation of quasinormal modes (QNM) and power-law tails, one would ?rst of all expect that one needs quite a ?ne resolution near the ω = 0 point in order to be able to correctly resolve the tails. Similarly, the resolution of the QNM would be very sensitive to the spacing in frequencies. With this in mind, we have chosen to solve the Teukolsky equation in t, r, and θ coordinates. The resulting 2+1 evolution equation is a hyperbolic, linear equation which is quite amenable to numerical treatment, provided suitable coordinates, variables and boundary conditions are chosen. The major numerical di?culty in solving the Teukolsky equation in the time domain arises from the term linear in s, involving the ?rst-order time derivative. Depending on the relative sign between the coe?cients of ?t Ψ and ?tt Ψ, one may view ?t Ψ either as a damping term (when the signs of both coe?cients agree) or an anti-damping term (otherwise). Without the factor 2s, the real part of the coe?cient of ?t Ψ reads M (r 2 ? a2 ) ? r. (2) ? √ In the physically allowed range [r+ , ∞), where r+ ≡ M + M 2 ? a2 represents the event horizon, the function C is monotonic in r , with limr→r+ C = ∞ and limr→∞ C = ?∞. Therefore, there exists a point rc such that C (r < rc ) > 0 and C (r > rc ) < 0. On the other hand, the coe?cient of ?tt Ψ is always negative. Consequently, for s > 0, the term ?t Ψ acts as a (anti-)damping term for (r < rc ) r > rc . The situation reverses for s < 0. Of course, there is no real damping or anti-damping, in terms of energy balance. The precise polynomial dependence of the coe?cient(s) is a?ecting, among other things, the required asymptotic fall-o? behavior of ingoing and outgoing radiation. However, the above consideration suggests that, in numerical work, the presence of such derivative terms could lead to exponentially growing modes. Those unphysical modes are not part of the family of solutions to the Teukolsky equation. Indeed, in our ?rst attempts of solving the Teukolsky equation, such growing modes were present in our simulations. These modes were clearly of numerical origin. The suppression of these instabilities proved to be an interesting and challenging exercise in the construction of numerical algorithms. The two key factors in successfully solving the Teukolsky equations were: ?rst, to rewrite equation (1) in a form that explicitly exhibits the radial characteristic directions, and second, to carefully select the evolution ?eld and its asymptotic behavior. A successful numerical evolution of the Teukolsky equation was achieved by considering the s = ?2 case. This restriction does not diminish in any way the range of initial data that can be evolved. In retrospect, the choice of ψ4 as the evolution ?eld is natural since the radiation content at in?nity is most succinctly described in terms of that ?eld. Yet in our investigation this choice was dictated by the need to obtain numerically acceptable asymptotic limits near the horizon. In order to reconstruct the perturbation of the metric coe?cients from the solution of the Teukolsky equation, however, one would have to evolve the equation for s = 2 as well. In our approach, i.e. evolving ψ4 , it is possible to recover the metric perturbations only at in?nity. In some sense this restricts the information that we get from our calculations. Local information about the perturbed spacetime (about geodesics etc) is not available. C (r ) = 4

For a given spin weight s, outgoing waves correspond to solutions to Eq. (1) with the limiting behaviour [5]
r ? → +∞

lim |Ψs | ? 1/r 2s+1

(3) (4)

r ? →?∞

lim |Ψs | ? 1

at spatial in?nity and the event horizon, respectively. Meanwhile, ingoing waves behave as
r ? → +∞

r ? →?∞

lim |Ψs | ? ??s .

lim |Ψs | ? 1/r

(5) (6)

Here we have used the Kerr tortoise coordinate r ? , which is de?ned by dr ? = r 2 + a2 dr . ? (7)

The convenient properties of the s = ?2 choice can be veri?ed by looking at the asymptotic form of propagating waves near the horizon (r ? → ?∞). For s = ?2, the solutions are bounded for any direction of propagation; in contrast, for s = 2 the ingoing solution diverges as the horizon is approached (? → 0). The asymptotic behavior for s = ?2 at r ? → +∞ can be subsequently ?xed by requiring that outgoing waves have asymptotically bound amplitude also in this limit. For the Teukolsky function Ψ?2 , this is achieved through a rescaling by an appropriate function of r . A convenient and simple choice is r 3 , a factor that is regular at the horizon. Regarding the choice of spatial coordinates, we use the Kerr tortoise coordinate r ? de?ned ? coordinate instead of the Boyerby Eq. (7). As azimuthal coordinate, we use the Kerr φ ? is de?ned by Lindquist coordinate φ. The coordinate transformation to φ ? ≡ dφ + dφ a dr . ? (8)

? was previously used for the scalar ?eld evolution in order to The azimuthal coordinate φ ameliorate coordinate pathologies near the horizon. A discussion of those pathologies and their precise manifestation in the slow rotation limit is given in Paper I.

We can now introduce the following Ansatz for the solution to the Teukolsky equation:
? 3 ?) ≡ eimφ Ψ(t, r ? , θ, φ r Φ(t, r ? , θ) .


It follows that the Teukolsky equation for Φ(t, r ? , θ) has a structure similar to that of Eq. (1), and therefore the same analysis of the damping or anti-damping nature of ?rst-order time derivative terms applies. After a series of unsuccessful numerical experiments with this second-order in time and space equation for Φ, we found that numerical instabilities were suppressed by introducing an auxiliary ?eld Π that converts the Teukolsky equation to a coupled set of ?rst-order equations in space and time. This can be accomplished by de?ning 5

Π ≡ ?t Φ + b ?r? Φ , r 2 + a2 b≡ , Σ Σ2 ≡ (r 2 + a2 )2 ? a2 ? sin2 θ . The resulting form of the Teukolsky equation has the following structure: ?t u + M ?r? u + Lu + Au = 0,

(10a) (10b) (10c)


where u ≡ {ΦR , ΦI , ΠR , ΠI } is the solution vector with (I) and (R) labeling the imaginary and real parts, respectively. The coe?cient matrices M and A are given by: b ? 0 ? M ≡? ? m31 ?m32

0 b m32 m31 0 0 a32 a31


0 0 0 0 ? ? ? , ?b 0 ? 0 ?b


respectively. Finally, L is a matrix operator that contains all the angular derivatives and has the following non-vanishing elements: 0 ? 0 ? L≡? ? l31 0

0 ? 0 ? A≡? ? a31 ?a32

?1 0 a33 ?a34

0 ?1 ? ? ? , a34 ? a33


0 0 0 l31

0 0 0 0

The coe?cients in the above matrices are given explicitly in Appendix A. The ?rst-order system given by Eq. (11) is hyperbolic in the radial direction since the eigenvalues of M are real and there is a complete set of linearly independent eigenvectors. Diagonalization of the matrix M and construction of evolution schemes based on the eigen?elds generated stable evolutions. However, it turned out that this last step was not necessary. Stable evolutions were also achieved using a modi?ed Lax-Wendro? method (see, for example, [16]) applied to Eq. (11) when this equation was rewritten as: ?t u + D ?r? u = S , (15)

0 0? ? ?. 0? 0


with D = diag(b, b, ?b, ?b) and S = ?(M ? D)?r? u ? Lu ? Au. Equation (15) is discretized on a uniform two-dimensional polar grid. Typically we used ? a computational domain of size ?100M ≤ ri ≤ 500M and 0 ≤ θj ≤ π with 0 ≤ i ≤ 8000 and 0 ≤ j ≤ 32. The Lax-Wendro? updating of ?eld values un i given on a timestep tn n+1/2 consisted of two steps. In the ?rst step, intermediate ?eld values ui+1/2 are obtained from ui+1/2 =

δt 1 n 1 n n ui+1 + un un ? un D i ? i ? S i+1/2 , 2 2 δr ? i+1/2 i+1 6


where we have omitted the angular indices j for clarity. The ?elds are j -centered in the angular direction with angular derivatives approximated by the appropriate second-order centered expressions. Radial derivatives in S n i+1/2 are approximated by centered di?erence n expressions of ?eld values at i and i + 1. Similarly, algebraic terms in D n i+1/2 and S i+1/2 are obtained by averaging over the values at i and i + 1. The second and ?nal step is a Staggered leapfrog step:
+1 un = un i i ? δt

1 n+1/2 n+1/2 n+1/2 n+1/2 . ui+1/2 ? ui?1/2 ? S i Di ? δr
n+1/2 n+1/2


Here the centered di?erence expressions and averages for S i and Di are taken from n+1/2 n+1/2 values at ui+1/2 and ui?1/2 . Examination of radial and angular characteristic directions at arbitrary locations shows that the appropriate CFL condition is δt ≤ max{ δr ?, 5 δθ }, where δt is the evolution timestep. Boundary conditions are required at the horizon, at the rotation axis of the black hole and at radial in?nity. The boundary conditions near the horizon are Φ = Π = 0, as follows immediately from the asymptotic behavior of ingoing waves (see Eq.(6)). The distinction of ingoing and outgoing waves near the Kerr event horizon requires a careful de?nition because of the rotational dragging of inertial frames [5]. However, in practice it clearly does not matter if we set both propagation modes to zero. At the axis, only a condition on Φ is required. Depending on the behavior of the angular eigenfunction belonging to the azimuthal number m near the axis, one imposes either Φ = 0 or ?θ Φ = 0 [12]. Finally, the boundary conditions at in?nity are much harder to impose when using a ?nite radial grid. Approximate schemes allow the transmission of waves across the boundary, but the implied truncation of the equations coe?cients interferes strongly with physical features of the evolution such as quasinormal ringing and tails. A clean solution can be achieved via “Cauchy-characteristic matching” (see [17,18] and references therein). However, this approach will not be pursued here. For the problems of immediate interest the computational domain can be made su?ciently large that errors generated at the outer boundary will not a?ect the results. Consequently, we set Φ and Π arbitrarily to zero at the outer boundary, and our results are computed only from information inside the maximal Cauchy development of the initial data surface. The stability of the code was veri?ed with long-time evolutions, of the order of 1000M . Its accuracy was tested using standard convergence tests. The code was found to be second order convergent for evolution times t < 50M . The convergence rate degrades as the total evolution time is increased, but is consistently above 1.3. A note on the nature of the unstable modes discussed previously is of some relevance here. Linear systems with constant coe?cients may exhibit numerical instabilities if the eigenvalues of the update matrix have modulus larger than unity. This phenomenon is easily analyzed with the so-called Von Neumann analysis of the di?erence equation. A similar analysis can be used in equations with variable coe?cients, via local stability analysis, but the method becomes far less conclusive in that case. For the Teukolsky equation, a number of simple, locally stable, discretizations turned out to be unstable for late times. The instabilities depended on the numerical discretization lengths, with increased resolution generally leading to slower growth rates for the unphysical modes. Yet only impractically high resolutions would suppress the instabilities for the long evolution times required for the applications discussed here. 7


The initial data for gravitational perturbations of rotating, asymptotically ?at spacetimes can be broadly divided into two classes: Waves coming in from in?nity and then scattering o? the black hole background, and waves emanating from the black hole as a result of an external excitation of the black hole. The latter case is perhaps the most interesting since it includes the description of the ?nal stages of the binary black hole coalescence within the close-limit approximation. In general, the evolution of both types of initial data consists of three stages, as seen by an observer located away from the hole. During the ?rst stage, the observed signal depends on the structure of the initial pulse and its re?ection by the curvature potential (burst phase). This phase is followed by the exponentially decaying quasinormal ringing of the black hole (quasinormal phase). In the last stage, the wave slowly dies o? as a power-law tail (tail phase). The precise manifestation of the last two phases is dictated by the subtle interference of the respective amplitudes. The complex frequencies of the QNMs for various values of a are well known [8–11] and we can use them as benchmarks of the present code’s accuracy. Similarly, the late-time tail phenomenon is mathematically well understood in spherically symmetric backgrounds, where tails have been predicted and veri?ed [19,20]. For rotating black holes, the background is no longer spherically symmetric, but numerical results for scalar waves in Kerr backgrounds suggest that power-law tails will also be present for spin-2 ?elds [3]. Furthermore, it has recently been argued that these tails should take the same form as for non-rotating black holes [21]. We have used the numerical algorithm described in the previous section to obtain the time evolution of generic perturbations impinging on the rotating black hole. For all the results presented here, we have used ingoing initial data with compact support for ΦR and have set ΦI = 0. All the evolutions clearly showed the exponentially damped oscillations of the QNMs and the subsequent late-time power-law tail. We will now discuss each of these phases in somewhat more detail.
A. Power-law tails

In Paper I we showed that the late time evolution of scalar ?elds in the background of rotating black holes is qualitatively similar to the non-rotating case. Here we extend this result to the physically more interesting evolution of a spin-2 ?eld. Starting with initial data ∝ ?2 Y2m (where s Ylm denotes a spin-weighted spherical harmonic), we studied the late time regime for a/M = 0, 0.25, 0.5, 0.75, 0.9 and m = 0, ±1, ±2. In all cases, we found a power-law tail behavior |Φ| ∝ t?? . For a nonrotating black hole it is known that the power-law exponent should be ? = 2l + 3 [20,21]. Indeed, for a = 0 we could reproduce the theoretical value ? = 7 with an accuracy of 10%. Moreover, our calculations show that the exponents governing the behavior for a = 0 remain similar to the Schwarzschild value. These results are summarized in Table I. We conclude that the power-law tail behavior is basically determined by the dominant asymptotic form of the potential, namely its “Schwarzschild” part. That this would be the case has recently been suggested by analytic arguments [21]. Analytic results (for non-rotating black holes) can also explain the main discrepancy between our numerically determined power-law exponents 8

and the theoretical values. The ?eld will be governed by a pure power law t?(2l+3) only at very late times. In a somewhat earlier regime “higher order” terms will also be signi?cant [21]. The presence of such higher order terms tend to increase the value of a numerically extracted power-law exponent. The results discussed above correspond to initial data that are dominated by the lowest allowed multipole l for the particular value of m that is being considered (l ≥ |m|). When this is not the case, we ?nd that di?erent multipoles become mixed in the evolution. To investigate this, we started from an initial angular distribution given by ?2 Y3m , and then computed the power-law exponent for di?erent values of m. For m = 3 the initial pulse used for the two-dimensional evolution corresponds to the lowest possible value of l and the late time power-law behavior is consequently dictated by ? = 10.05. The situation is di?erent for the case m = 0, where we ?nd that the power-law exponent is given by ? = 7.72. That is, the late-time evolution is dominated by the quadrupole. In contrast, the corresponding Schwarzschild evolution exhibits no such mixing of multipoles, and evolutions for l = 3, m = 0, 3 agree well with the theoretical power-law tail exponent ? = 9. A result similar to this was discussed in Paper I. This result is, however, not too surprising. The mixing of multipoles is due to the rotational dragging of inertial frames. Basically, there are two reasons why di?erent multipoles will be present in the evolution. The ?rst one is “imperfect” initial data: Here we have expressed the initial data in terms of the spin-weighted spherical harmonics s Ylm . The appropriate angular functions that are used when the r and θ components of Eq. (1) are separated are the spin-weighted spheroidal harmonics s Slm (θ; aω ) [5]). Hence, because of the frequency (time) dependence in the s Slm (θ; aω ), it is di?cult to generate initial data for the Kerr problem that represent a pure multipole l for a speci?ed value of m. Furthermore, the rotation of the black hole will lead to an active coupling between di?erent multipoles (for a discussion of the analogous situation for rotating stars, see [22]). So even if we could initiate the evolution with a pure multipole, other multipoles would be generated and may play a role at later times. The mixing of multipoles can have an interesting result. Suppose that the initial data is dominated by a certain multipole l but that a relatively small part of l ? 1 is also present. Each of these multipoles will give rise to power-law tails, but the amplitude of the latter will be much smaller than the ?rst. In this situation the late time ?eld should contain both a tail of form t?(2l+3) and another term that behaves as t?(2l+1) . Since the ?rst term dies o? faster than the second it seems unavoidable, no matter how small the amplitude of the second tail term is, that the very late time evolution of the ?eld will be dominated by the l ? 1 multipole. The corresponding time evolution will simply exhibit a late-time change in the angular behaviour. An example of this was discussed in Paper I.
B. Quasinormal modes

To further test the performance of the Teukolsky code, we performed a series of simulations for di?erent values of the angular momentum parameter a and the azimuthal number m. Then we focussed on the later part of the ringing sequence, which is dominated by the slowest damped mode of oscillation, and read o? the oscillation frequency ω via a standard Fourier transform of the time signal. Table II summarizes results obtained in this way. 9

A comparison of our values for the real part of each QNM frequency with those given by Leaver [9] and Kokkotas [11] for m = 0, a/M = 0, 0.5, 0.9 shows an agreement to better than 0.3%. For m = ?1, our results agree with the values of Detweiler [8] and Leaver (di?erent sign convention) to within 0.5%. Our results for the imaginary parts are accurate to 1%. Although higher multipole evolutions are readily obtained, accurate estimates of the corresponding QNM frequencies require increasingly higher angular resolution. Hence, we have only estimated QNM frequencies for the quadrupole modes. But it should be pointed out that, as was the case for the power-law tails, QNMs corresponding to other multipoles are present in all signals. Interestingly, we ?nd that the numerically extracted QNM frequencies for non-zero m do not depend on the sign of m, i.e. we get the same values for the QNM frequencies from evolutions for e.g. m = ±1. At ?rst sight this is surprising, since it seems to contradict well established results. According to for example Detweiler [8], the imaginary parts of the m = +1 and m = ?1 modes should become quite di?erent as a increases. Fortunately, the answer is simple: The frequencies of both the m and the ?m QNMs are present in a typical evolution. This also explains the existence of two distinct regimes of the quasinormal ringing that can be seen in Fig. 1. In this ?gure we show the ?eld as seen by an observer at r ? = 20M , θ = π/4 for a = 0.99M , M = 0.5, and m = 2. In the early phase of QNM ringing (120M ≤ t ≤ 160M ), the decay of the ?eld is approximately given by |Φ| ∝ e?0.175t , whereas the behavior for late times is governed by a slower decaying mode, |Φ| ∝ e?0.0605t . These values are in good agreement with the theoretical values for l = 2 and m = ±2, as can be seen by comparison with Fig. 1(a) in [8]. One question remains. Why should we expect both these modes to be present in the signal? The answer can be found in the symmetries of the problem. In the Schwarzschild case we know that there will be modes in the ?rst and second quadrant of the complex ω -plane (assuming that the modes have a time-dependence eiωt ). All these modes will contribute to the signal, but since they occur as complex conjugate pairs ωl and ?ωl? (where the asterisk represents complex conjugation) the contribution from one set of modes can be expressed in terms of the other. Hence, the QNM part of the signal can be evaluated using only modes in the ?rst quadrant (for a study of the Schwarzschild problem, see [21]). In the case of Kerr black holes the situation is more complicated. The QNMs no longer occur as complex conjugate pairs. Let us suppose that we work in the frequency domain, and that we have a solution Ψlm (ω, r, θ)eimφ to the Teukolsky equation for given l and m. Then it is easy to show (using the separated equations, see [5]) that a solution to the equation for ?m will follow as Ψl?m (ω, r, θ) = [Ψlm (?ω ? , r, π ? θ)]? . (18) This means that if ω is a QNM corresponding to l and m then the complex conjugate ?ω ? will be a QNM for l and ?m. This symmetry is nicely illustrated in Fig. 3 in [9]. We can use this symmetry to express the contribution from the Kerr QNMs in the second quadrant in terms of modes in the ?rst quadrant (we will describe this procedure in detail elsewhere). If we represent modes in the ?rst quadrant by ωlmn , we can schematically write the QNM contribution to the signal as
M = ΨQN m ln

F (ωlmn ) + G(?ωl??mn ) 10



This explains the presence of both branches of QNMs in our evolutions.
C. Wave propagation

Finally, an interesting direct application of our code is the identi?cation of e?ects of the background rotation on the propagation of signals. The idea would be to use the same initial data for black holes with di?erent rotation rates and see to what extent the rotation of the black hole can be inferred from the emerging signals. A stumbling block for work in this direction, though, is the di?culty to prescribe initial data that represent the “same” initial physical perturbation for any value of the black hole angular momentum. This prescription is possible for data coming in from in?nity, since the e?ects of rotation are negligible for large enough radii. However, it is not clear how to set up initial data (e.g. position of the pulse) in the vicinity of the black hole in such a way that comparison of evolutions for di?erent angular momenta makes sense. For this reason, we addressed only the e?ect of rotation on pulses that were initially far way from the hole. In Fig. 2 we show a series of snapshots in r ? of the evolution of a pulse for m = 2. At t = 0, the initial data consists of a bell-shaped pulse in r ? centered at 75M and propagating inwards. Then the evolution of ΦR is shown at constant θ = π/2 and 0 ≤ t ≤ 200M for four di?erent values of the angular momentum of the black hole: a/M = 0, 0.5, 0.9, 0.99. One immediately sees that an increase in the angular momentum of the black hole has two e?ects on the scattered pulse. Both the amplitude and wavelength of the emerging signal change. By t = 200M the wavelength in the trailing part of the signal for a/M = 0.99 is approximately half as long as for the non-rotating case. This agrees well with the anticipated increase in the dominating QNM frequency [8]. Furthermore, at t = 200M the amplitudes for a/M = 0 and a/M = 0.99 also di?er by ? 25%. The amplitude of the emerging signal generally decreases as we increase a. The angular behaviour for the scattering scenario also shows interesting features. An example can be seen in Fig. 3 where we show the evolution for a/M = 0.99. The initial data is given by ΦR ∝ sin2 θ, centered in r ? at 75M and propagating inwards, but after some time there is a clear transition of the angular distribution to ?2 Y22 ∝ cos4 θ/2. This is the dominant angular behavior during the later evolution. That is, the quadrupole dominates the scattered wave at later times. This situation is similar to the one we discussed for power-law tails.

We have studied the time evolution of perturbations of rotating black holes by numerically integrating the Teukolsky equation written as a system of equations of ?rst-order in time and space. In particular, we investigated the role of the black holes angular momentum on the dynamics during the burst, quasinormal ringing and power-law tail phases. We have reproduced values for the fundamental quasinormal mode frequencies with an accuracy of better than 1%. As for the late-time regime, we found power-law tail behavior analogous to the scalar ?eld case that we had studied previously. These results are in good agreement with what one would expect and establish the reliability and accuracy of our code. We 11

thus have access to a useful numerical laboratory that can be used to further illuminate the detailed physics of rotating black holes. In the future we plan to use this laboratory to investigate several interesting issues. We would, for example, like to establish what e?ect the so-called super-radiance will have on the evolution of initial data. From studies of scattering of monochromatic waves from rotating black holes (see, for example [12]) it is know that frequencies such that ωM < am/2r+ will correspond to a re?ection coe?cient whose magnitude is larger than unity. That is, these frequencies are “super-radiant”. This is the wave-analogue to the well-known Penrose process. In the results presented in this paper the e?ect of super-radiance was not obvious, but it is plausible that it can be unveiled in a more detailed study. For example, it seems that super-radiance should lead to di?erent results for co- and counter-rotating waves in an evolution. That is, if super-radiance is relevant one should see di?erence between evolutions for |m| and ?|m|. Another interesting issue concerns the excitation of quasinormal modes. It is known that the damping of the ?|m| modes will vanish as a increases [9]. This would potentially make the modes of a rapidly rotating black hole ideal for gravitational wave detection. However, in an approximate study Ferrari and Mashhoon [23] have shown that the amplitude of these extremely slowly damped modes will be negligible. In the extreme Kerr limit no energy is expected to be radiated through these modes. This result has proved to be very hard to verify analytically, but it is testable with the present code. We should thus be able to establish whether the extremely slowly damped quasinormal modes are of astrophysical relevance or not. Another outstanding problem in black-hole physics concerns the dynamical stability of the Kerr black hole to perturbations. So far, only mode-stability has been established [24]. To prove complete stability one must ensure that all relevant quantities remain point-wise bounded during an evolution, cf. [25]. It seems plausible that our Teukolsky code can prove useful also in this context. Perhaps the most interesting future application of our numerical Teukolsky code will be to evolve initial data qualitative similar to the late stages of a binary black hole coalescence. The main problem in approximating black hole collisions with perturbations about Kerr spacetimes is the construction of initial data. Work so far on the close limit approximation to construct initial data has heavily depended on having a conformally ?at background geometry. On the other hand, constant time hypersurfaces in Boyer-Lindquist or Kerr coordinates are not conformally ?at. Hence disentangling these e?ects is not as straightforward as for nonrotating black holes. Possible ways of circumventing these obstacles are presently being investigated [26].

We thank K. Kokkotas, H.-P. Nollert, R. Price and J. Pullin for helpful discussions. This work was supported by the Binary Black Hole Grand Challenge Alliance, NSF PHY/ASC 9318152 (ARPA supplemented) and by NSF grants PHY 96-01413, 93-57219 (NYI) to PL. WK was supported by the Deutscher Akademischer Austauschdienst (DAAD) and the Deutsche Forschungsgemeinschaft (DFG), SFB 382. NA is supported by NSF grant PHY 92-2290. 12


In the following, we give the explicit form of the coe?cients in the matrices M , A, and L in Eq. (11). With ?3Mr 2 + Ma2 + r 3 + ra2 , Σ2 r ? (1 + s) ? M (a2 ? r 2 ) s , c2 ≡ ?2 Σ2 2 Mrm + ? s cos θ c3 ≡ 2 a , Σ2 r 2 + a2 , c4 ≡ ?2 a m Σ2 c1 ≡ 2 s the coe?cients of M can be written as m31 ≡ ?b c1 + b m32 ≡ b c3 ? c4 . De?ning ?(?m2 ? 2 cos θ s m ? cos2 θ s2 + sin2 θ s) , c5 ≡ Σ2 sin2 θ (r ? M ) s m a c6 ≡ ?4 , Σ2 the coe?cients of A are given by as well as for investigating asymptotic behavior a31 a32 a33 a34 Finally, with c7 ≡ ? ? , Σ2 ? , Σ2 (A11) (A12) ≡ c5 , ≡ ?c6 , ≡ c1 , ≡ ?c2 . (A10a) (A10b) (A10c) (A10d) (A10e) (A7) (A8) (A9) ?b + c2 ?r ? (A6a) (A6b) (A1) (A2) (A3) (A4) (A5)

c8 ≡ ? cot θ

the only non-vanishing coe?cient of the operator matrix L reads l31 ≡ c7 ? ?2 + c8 . 2 ?θ ?θ 13 (A13)

[1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] F. J. Zerilli, Phys. Rev. D, 2, 2141 (1970) R.H. Price and J. Pullin Phys. Rev. Lett. , 72, 3297 (1994) W. Krivan, P. Laguna, and P. Papadopoulos, Phys. Rev. D, 54, 4728 (1996) S. A. Teukolsky, Phys. Rev. Lett. , 29, 1114 (1972) S. A. Teukolsky, Astrophys. J., 185, 635 (1973) W. Kinnersley, J. Math. Phys., 10, 1195 (1969) J.M. Bardeen and W.H. Press J. Math. Phys., 14, 7 (1973) S. Detweiler, Astrophys. J., 239, 292 (1980) E. W. Leaver, Proc. R. Soc. London, Ser. A 402, 285 (1985) E. Seidel and S. Iyer Phys. Rev. D, 41, 374 (1990). K. D. Kokkotas, Class. Quantum Grav. 8, 2217 (1991) J. A. H. Futterman, F. A. Handler and R. A. Matzner, Scattering from Black Holes (Cambridge Univ. Press 1988) M. Sasaki and T. Nakamura Progr. Theor. Phys. 67, 1788 (1982). Y. Kojima and T. Nakamura Phys. Lett. A96, 335 (1983); A99, 37 (1983) Y. Kojima and T. Nakamura Progr. Theor. Phys. 71, 79 (1984); 72, 494 (1984). W. H. Press, S. A. Teukolsky, W. T. Vetterling, B. P. Flannery: Numerical Recipes in Fortran, Second Edition (Cambridge University Press, Cambridge, England, 1992) R. G? omez, P. Laguna, P. Papadopoulos, and J. Winicour, Phys. Rev. D, 54, 4719 (1996). P. Papadopoulos and P. Laguna, ‘Cauchy-characteristic Evolution of Einstein-KleinGordon Systems: The Black Hole Regime, to appear in Phys. Rev. D (1997) R. H. Price, Phys. Rev. D, 5, 2419 (1972) C. Gundlach, R. H. Price, and J. Pullin, Phys. Rev. D, 49, 883 (1994) N. Andersson, Phys. Rev. D, 55, 468 (1997) Y. Kojima Phys. Rev. D, 46, 4289 (1992) V. Ferrari and B. Mashhoon Phys. Rev. D, 30, 295 (1984) B.F. Whiting, J. Math. Phys., 30, 1301 (1989) B.S. Kay and R.M. Wald, Class. Quantum Grav., 4, 893 (1987) J. Baker, W. Krivan, P. Laguna, H.-P. Nollert, P. Papadopoulos and J. Pullin, in preparation (1997).


TABLE I. The numerically determined power-law exponents ? governing the late-time behavior of |Φ| ∝ t?? measured at r ? = 20M , θ = π/2. The theoretical value for a = 0 is ? = 7. a/M 0 0.25 0.50 0.75 0.90 ? (m = 0) 7.70 7.73 7.75 7.81 7.68 ? (m = 1) 7.69 7.77 7.68 7.79 7.71

TABLE II. The numerically extracted QNM frequencies that dominate the late-time part of the ringing phase for m = 0, ±1. a/M 0 0.25 0.50 0.75 0.90 ωM (m = 0) 0.373 + i 0.0885 0.376 + i 0.0885 0.383 + i 0.0865 0.398 + i 0.0835 0.412 + i 0.0780 ωM (m = |1|) 0.373 + i 0.0885 0.392 + i 0.0885 0.420 + i 0.0860 0.467 + i 0.0795 0.517 + i 0.0695




-0.175 t


|Φ| 1


-0.0605 t

10 0




300 t/M




FIG. 1. Quasinormal ringing for a = 0.99M , M = 0.5, and m = 2. During 120M ≤ t ≤ 160M , the ?eld decays as |Φ| ∝ e?0.175t ; the subsequent evolution is dictated by a slower decaying mode |Φ| ∝ e?0.0605t .


















FIG. 2. Snapshots of the evolution of ΦR for m = 2, at θ = π/2, in the time interval 0 ≤ t ≤ 200M . For clarity, ΦR /α is displayed with α a scale factor: α = (10, 103 , 103 , 103 ) for t = (50, 100, 150, 200) M , respectively. The values used for the angular momentum are a/M = 0, 0.5, 0.9, 0.99. The initial data consist of a bell-shaped pulse centered at 75M and propagating inwards. Di?erences in the pulses do not become appreciable until the pulse has been re?ected at the potential barrier outside the black hole (see t = 100M snapshot). By t = 200M a decrease of ? 25% in the wavelength can be observed. The maximum amplitude of the pulse also decreases with a.


(a) 60 40 ΦR 20 0 -20 0 50 100 r */ M 150 200 0 0.5 1 1.5 2 θ 2.5 3 ΦR 0 -104 -2 x104 0 50 100 r */ M 150 200 0 0.5 1 1.5 2 θ 2 x104 104




(c) 3 x104 2 x104 104 ΦR 0 -104 -2 x104 -3 x104 0 50 100 r */ M 150 200 0 0.5 1 1.5 2 θ 2.5 3 ΦR 3 x104 2 x104 104 0 -104 -2 x104 -3 x104 0 50 100 r */ M 150 200 0 0.5 1 1.5 2 θ




FIG. 3. Evolution for a/M = 0.99 at t = 50M (a), t = 100M (b), t = 150M (c), and t = 200M (d). The initial data consist of a bell-shaped pulse centered at r ? = 75M , θ = π/2 and propagating inwards. During the evolution there is a clear transition of the angular distribution to 4 2 ?2 Y2 ∝ cos θ/2.




All rights reserved Powered by 甜梦文库 9512.net

copyright ©right 2010-2021。