Superradiance of a charged scalar field coupled to the Einstein-Maxwell equations

Superradiance of a charged scalar field coupled to the Einstein-Maxwell equations

Olaf Baake olaf.baake@aei.mpg.de Max Planck Institute for Gravitational Physics (Albert Einstein Institute), Am Mühlenberg 1, 14476 Potsdam, Germany    Oliver Rinne oliver.rinne@aei.mpg.de Hochschule für Technik und Wirtschaft Berlin, Treskowallee 8, 10318 Berlin, Germany Max Planck Institute for Gravitational Physics (Albert Einstein Institute), Am Mühlenberg 1, 14476 Potsdam, Germany
October 26, 2016
Abstract

We consider the Einstein-Maxwell-Klein-Gordon equations for a spherically symmetric scalar field scattering off a Reissner-Nordström black hole in asymptotically flat spacetime. The equations are solved numerically using a hyperboloidal evolution scheme. For suitable frequencies of the initial data, superradiance is observed, leading to a substantial decrease of mass and charge of the black hole. We also derive a Bondi mass loss formula using the Kodama vector field and investigate the late-time decay of the scalar field.

I Introduction

Even though particles can never cross the event horizon from the inside of a black hole (at least classically), it is nevertheless possible to extract energy from a rotating or charged black hole. For particles this can occur via the Penrose process Penrose (1969). The field-theoretic analogue is superradiance, which involves the scattering of scalar, electromagnetic or gravitational waves off a black hole. For a comprehensive recent review article of this field of research see Brito et al. (2015).

In general relativity, superradiance is mostly studied for rotating black holes (the Kerr family of solutions) and is based on linear perturbation theory Press and Teukolsky (1973). It should be stressed that this involves a mode analysis. There has been some controversy as to whether the predicted amplification factors are still meaningful for realistic wave packets. Csizmadia et al. Csizmadia et al. (2013) reported an almost perfect reflection for a wave packet whose initial frequency content lies entirely in the superradiant regime. Until recently little was known about superradiance in the nonlinear case. East et al. East et al. (2014) presented fully nonlinear numerical evolutions of the vacuum Einstein equations for gravitational waves scattering off a rotating black hole and confirmed the existence of superradiance.

Due to the lack of symmetries such simulations are very demanding. A simpler model is obtained by considering a charged scalar field scattering off a charged (Reissner-Nordström) black hole in spherical symmetry. This 1+1-dimensional problem can be tackled with modest computational resources. The test-field case (fixed black hole background spacetime) was studied by Di Menza and Nicholas Di Menza and Nicolas (2015). The authors computed the frequency-dependent energy gain by constructing a conserved flux. Such an approach is not feasible when the matter fields are coupled to the Einstein equations and the mass and charge of the black hole change during the superradiant scattering. Instead it is the changes in those quantities themselves that need to be monitored, as in East et al. (2014) and in the present paper. Numerical evolutions of the coupled Einstein-Maxwell-Klein-Gordon system were first carried out by Hod and Piran Hod and Piran (55) and more recently e.g. by Torres and Alcubierre Torres and Alcubierre (2014), however with a focus on gravitational collapse rather than superradiance. Recently there has been an increasing interest in superradiance in the context of the anti-de Sitter (AdS)–conformal field theory (CFT) correspondence. An example is the numerical study of a superradiant instability of Reissner-Nordström-AdS black holes in Bosch et al. (2016). Due to the timelike infinity one effectively has reflective boundary conditions that lead to an unbounded amplification of the superradiant waves. A similar situation, namely a Reissner-Nordström black hole enclosed in a cavity, was investigated in Sanchis-Gual et al. (2016). To our best knowledge the present paper is the first study of superradiance in the coupled Einstein-Maxwell-Klein-Gordon system in asymptotically flat spacetimes.

Most numerical studies of asymptotically flat spacetimes truncate the computational domain at a finite distance from the black hole, where boundary conditions must be imposed. (This approach was taken e.g. in Di Menza and Nicolas (2015).) Spurious reflections of outgoing waves must be avoided on the relevant time scales. Designing such boundary conditions in the nonlinear case is highly non-trivial. A more elegant solution is to incorporate future null infinity in the computational domain, e.g. by foliating spacetime into hyperboloidal surfaces that may be compactified towards future null infinity. This also enables us to evaluate the Bondi mass and total charge at in a straightforward way. We use a conformal 3+1 decomposition of the Einstein equations on hypersurfaces of constant mean curvature (CMC) developed in Moncrief and Rinne (2009). Such CMC surfaces have the additional advantage that they are able to penetrate black hole horizons, so that an excision boundary may be placed just inside the horizon, where all characteristics leave the computational domain. In Rinne (2010) this approach was first implemented for the vacuum axisymmetric Einstein equations. In Rinne and Moncrief (2013) we included matter sources and studied late-time power-law tails of matter fields in spherical symmetry. The present formulation is based on Rinne and Moncrief (2013).

This paper is organized as follows. In Sec. II we describe our formulation of the Einstein-Maxwell-Klein-Gordon equations and the hyperboloidal evolution scheme. We also investigate mass conservation with the help of the Kodama vector field and derive a Bondi mass loss formula. Details on redundant evolution equations and their regularity at and on alternative electromagnetic gauge conditions are deferred to Appendices A and B. Section III is concerned with the numerical evolution of the system. We describe our numerical methods, construction of initial data, and various notions of charge and mass used to analyze the results. After performing a number of code tests, we present our main results on the occurrence and amount of superradiance depending on the chosen parameters. A separate subsection is devoted to the analysis of the late-time decay of the scalar field. We conclude in Sec. IV.

Ii Formulation and theoretical analysis

In this section we describe our formulation of the Einstein-Maxwell-Klein-Gordon equations and their reduction to spherical symmetry. Our gauge choices for the Einstein and Maxwell equations are explained, in particular our use of constant-mean-curvature slices extending to future null infinity. In Sec. II.2 we investigate the Kodama vector field, construct a conserved flux and relate it to the Hawking mass. This enables us to derive a Bondi mass loss formula.

ii.1 Field equations

We consider a massive charged complex scalar field minimally coupled to the Einstein-Maxwell equations. The action is given by

(1)

where is the spacetime metric, its volume element and the scalar curvature. We use geometric and Gaussian units. In terms of the vector potential , the Maxwell field strength tensor is

(2)

The gauge-covariant derivative is defined as

(3)

where is the covariant derivative of . The mass and charge of the scalar field are and respectively.

We perform a conformal transformation

(4)

where the conformal factor at . In spherical symmetry we may write the conformal metric in isotropic coordinates as

(5)

with . We consider an ADM decomposition Arnowitt et al. (1962) with respect to the time coordinate . Constant mean curvature (CMC) slicing is used; i.e., the mean curvature of the slices is a spacetime constant . The tracefree part of the ADM momentum (cf. Moncrief and Rinne (2009)) has only one degree of freedom in spherical symmetry, which we take to be . The gravitational field is thus described by the four variables and , which are functions of and only. In the following we use an overdot to denote -derivatives and a prime to denote -derivatives.

Preserving the isotropic spatial coordinate condition and the CMC slicing condition under the time evolution yields

(6)
(7)

In spherical symmetry, the geometry is completely determined by the Einstein constraint equations,

(8)
(9)

The source terms , and in (7)–(9) are components of the conformally rescaled energy-momentum tensor and are defined for our matter model below in (19)–(21). There are redundant evolution equations for and , given in Appendix A, that we monitor during the evolution in order to check the accuracy of our code.

The Maxwell equations are conformally invariant and hence we may define the fields in terms of the conformal metric , indicated by tildes in the following. In spherical symmetry the vector potential may be written as

(10)

and the field strength tensor

(11)

We may impose one gauge condition on the vector potential, and the one we choose is temporal gauge

(12)

Alternative gauge conditions are discussed in Appendix B.

The definition (2) of the field strength tensor in terms of the vector potential implies an evolution equation for :

(13)

The Maxwell equations imply an evolution equation for the electric field

(14)

and the Gauss constraint

(15)

The source terms and are given below in (22)–(23). We normally use (15) to solve for and monitor (14) during the evolution but we have checked that the opposite choice gives identical results within numerical accuracy.

We introduce a conformally rescaled scalar field

(16)

and write the Klein-Gordon equation in first-order form (in time) by introducing a new variable :

(17)
(18)

The somewhat non-standard definition of given by (17) is used in order to avoid time derivatives of the lapse and shift in its evolution equation (18), cf. Rinne and Moncrief (2013). In deriving (18), we have re-expressed the four-dimensional scalar curvature in terms of the trace of the energy-momentum tensor using the Einstein equations, which produces terms quadratic in the scalar field. This is required in order to remove a term containing a negative power of the conformal factor, which would be formally singular at . In Rinne and Moncrief (2013) a different approach based on conformal rather than minimal coupling of the scalar field is taken.

Numerically, we have found it to be advantageous to also perform a first-order reduction in space by introducing a new variable , and we evolve the real and imaginary part of the complex scalar field variables separately.

The source terms appearing in the slicing condition (7) and the Einstein constraint equations (8)–(9) are

(19)
(20)
(21)

The source terms appearing in the Maxwell equations (14)–(15) are

(22)
(23)

ii.2 Mass conservation

In the dynamical spacetimes we study, there is no timelike Killing vector field. Nevertheless, in spherical symmetry one can construct a preferred timelike vector field Kodama (1980); Rácz (2006); Csizmadia and Rácz (2010). For a general spherically symmetric metric of the form

(24)

with a two-dimensional Lorentzian metric (indices ranging over ), this Kodama vector field is given by

(25)

where is the volume element of . For our form of the metric (4)–(5), we obtain

(26)
(27)
(28)

The Kodama vector field has the remarkable properties

(29)

and

(30)

Together with Einstein’s equations, the latter implies that the current

(31)

is conserved, . In our formulation we obtain

(32)
(33)
(34)

(The source term is defined in (72).)

The Kodama vector field is closely related to the Misner-Sharp or Hawking mass Csizmadia and Rácz (2010)

(35)

Its derivatives turn out to be

(36)

so that the integral of the Kodama flux is essentially given by the difference of the Hawking masses at the ends:

(37)
(38)

Taking the limit of (38) at , where the Hawking mass coincides with the Bondi mass , and using the regularity conditions stated in Appendix A, a lengthy calculation results in the Bondi mass loss formula (see also Scholtz and Holka (2014))

(39)

which is manifestly non-positive.

Iii Numerical evolution

In this section we present our numerical evolutions of the system derived in the previous section. In Sec. III.1 we briefly describe the numerical methods we use. The initial data for a Reissner-Nordström black hole with scalar field perturbation are constructed in Sec. III.2. In Sec. III.3, various notions of charge and mass are introduced, which are needed to evaluate our numerical evolutions. In Sec. III.4 we perform mass/charge conservation and convergence tests in order to check the accuracy of our code. Sec. III.5 contains our main results: evolutions for various choices of parameters are presented and evaluated with regard to the existence and the amount of superradiance. Finally in Sec. III.6, we investigate the late-time behavior of the scalar field and compare with known perturbative results on quasi-normal modes and power-law tails.

iii.1 Numerical method

We discretize the equations in space using fourth-order finite differences. A mapping of the radial coordinate with an adjustable parameter Rinne and Moncrief (2013) is used in order to provide more resolution where it is needed, especially near the black hole horizon where the fields typically develop steep gradients. The outermost grid point is placed at , which we choose to correspond to . Typical resolutions used for the simulations in this paper range from to radial grid points.

Following the method of lines, the evolution equations (13), (17) and (18) are first discretized in space and then integrated forward in time using a fourth-order Runge-Kutta method with sixth-order Kreiss-Oliger dissipation Kreiss and Oliger (1973). At each time step, the ODEs (6)–(9) and (15) are solved using a Newton-Raphson method, at each iteration solving the resulting linear system using a direct band-diagonal solver. (Alternatively, as mentioned in Sec. II.1, we may replace the Gauss constraint (15) with the evolution equation (14) for the electric field, which yields identical results within numerical accuracy.)

Our treatment of the boundaries follows Rinne and Moncrief (2013). We place an inner excision boundary just inside the apparent horizon of the black hole, whose location is determined at each time step as the outermost zero of the expansion of outgoing null rays,

(40)

One-sided finite differences are used at the excision boundary. Since this boundary lies inside the black hole, all characteristics leave the domain and hence no boundary conditions are required for the evolution equations. We choose to freeze the conformal lapse at the inner boundary, yielding a Dirichlet boundary condition for the slicing condition (7). Inner Dirichlet boundary conditions for the Einstein constraint equations (8)–(9) and the Gauss constraint (15) are obtained by evolving , and there according to their evolution equations (70)–(71) and (14).

The outer boundary is an outflow boundary and hence we use one-sided differences there as well, with no boundary conditions for the evolution equations. The conformal lapse is set to at , which ensures that our time coordinate coincides with Bondi time Rinne and Moncrief (2013). Outer Dirichlet boundary conditions on , and follow from the regularity conditions at , equations (73), (76) and (80).

In all our evolutions the value of the mean curvature is taken to be .

The code has been written in and the figures produced with Python, making use of the NumPy, SciPy and matplotlib extensions.

iii.2 Initial data

We choose initial data that is close to the Reissner-Nordström spacetime. First the geometry variables (, , and ) are set to coincide with this solution, then the initial data for the scalar field are specified and finally the constraints and elliptic gauge conditions are re-solved for the geometry variables.

The Reissner-Nordström spacetime in (uncompactified) CMC coordinates is given by Brill et al. (1980)

(42)

with

(43)
(44)

where (mass), (mean curvature) and are constants. We transform the radial coordinate to a new radial coordinate by demanding that the spatial metric be manifestly conformally flat. For convenience we work with due to its finite range, which yields the ODE

(45)

with

(46)
(47)

and we choose to correspond to (). We obtain the -coordinate of the horizon by numerically integrating

(48)

where , and

(49)

are radii of the outer and inner horizon in Schwarzschild coordinates. We place the excision boundary just inside the outer horizon: , where typical values of are between 0.8 and 0.9. The ODE (45) is then solved numerically on the interval with the initial condition .

The geometry variables can be expressed in terms of the numerically determined function as

(50)
(51)
(52)
(53)

For the scalar field we choose initial data that are supported sufficiently far outside the black hole, where the background is almost flat. We take an ingoing solution of the wave equation on a Minkowski background:

(54)

This is then expressed in terms of our CMC coordinates , which are related to the standard Minkowski coordinates via

(55)
(56)

Now we choose such that the wave packet is localized sufficiently far outside the black hole (we use , which corresponds to ) and take as initial data for the scalar field.

The initial data for the electromagnetic field have to correspond to the Reissner-Nordström spacetime. Together with the gauge condition (12) we set

(57)

Alternative gauge conditions (and our reasons for not using them) are discussed in Appendix B.

iii.3 Charge and mass

As in Torres and Alcubierre (2014) we introduce the charge inside a sphere of radius :

(58)

which using the Gauss constraint (15) can be written as

(59)

Let us first consider a static black hole. Its irreducible mass is given by the area of the event horizon :

(60)

The mass of the black hole can now be calculated as

(61)

These formulas are valid if spacetime is static in a neighbourhood of the horizon. This is the case for the final state of the evolution as well as for the initial data, which consist of a scalar field pulse supported far away from the horizon. During the evolution, we still compute and according to the above expressions, but evaluated at the apparent horizon. The second equality in (60) still holds if h refers to the apparent horizon, but this differs slightly from the event horizon during the dynamical phase of the evolution. Hence when we plot (60) and (61) as functions of time, they cannot strictly be interpreted as the irreducible mass and black hole mass during the dynamical phase but its initial and final values are correct.

Analogously to the rotational energy in East et al. (2014), we define the charge energy of the black hole as

(62)

During superradiant scattering, both the charge energy and the mass of the black hole decrease:

(63)

Part of the charge energy is carried away by the wave, part of it increases the irreducible mass and hence is no longer extractable. Following East et al. (2014) we define the efficiency of this process as

(64)

Alternatively, we may compare the energy of the initial scalar field pulse

(65)

with the total energy radiated away at infinity,

(66)

From these quantities we define an alternative efficiency

(67)

iii.4 Code tests

The first test of our code consists in the evolution of the Reissner-Nordström background spacetime, without any scalar field. On the timescale relevant for the scalar field scattering experiments presented later (), mass (Bondi, black hole and irreducible) is conserved to a relative error of and charge (at future null infinity and at the horizon) is conserved to a relative error of . Both are negligible compared with the observed changes of charge and mass during scalar field scattering.

Next we perform a convergence test for an evolution with scalar field. Figure 1 shows the residual of the evolution equation (14) for , which we do not impose actively as we solve the Gauss constraint (15) for . (To compute the time derivative in the evolution equation numerically, we use fourth-order finite differences, using data from five subsequent time levels.)

Figure 1: Convergence test for three different resolutions: (solid red), (dashed blue) and (noisy green) grid points. Shown is the residual of the evolution equation (14) of the electric field as a function of time. The parameters used here are , , , , , , and .

The observed decrease of the residual for successively doubled resolutions is close to the expected value of for a fourth-order accurate finite-difference method. Similar plots are obtained for the residuals of the evolution equations (70) and (71) for and .

iii.5 Numerical results

For most of the evolutions shown here, we choose the scalar field mass to vanish, . The scalar field charge is taken to be . The parameters of the Reissner-Nordström background solution are and . (The relatively large value of is chosen so that there is enough room for superradiance to occur, as the maximum superradiant frequency is proportional to , see Eq. (68) below.) For the initial data parameters in (54) we first choose , , , and . The corresponding evolution is shown in Fig. 2. It shows the “normal” behavior we expect for a scalar field that mainly falls into the black hole, while only a small amount escapes to infinity: the black hole mass increases by about to almost the initial value of the Bondi mass, while the Bondi mass only decreases by a small amount; similar behavior is seen in the evolution of the charges at the horizon and at . The energy radiated away at infinity is much smaller than the energy of the initial scalar field pulse .

Figure 2: Non-superradiant evolution (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

Next, we change the frequency to , leaving all the other parameters unchanged. The resulting evolution (Fig. 3) is markedly different: now the black hole mass decreases by about (and its final value agrees with the final Bondi mass)—a clear indication of superradiance. The horizon charge decreases by about , the total charge at is nearly halved. This time the energy radiated away at infinity is larger than the energy of the initial scalar field pulse .

Figure 3: Superradiant evolution (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

For this choice of parameters, the superradiant efficiency defined in (64) is shown as a function of in Fig. 4. The evolution in Fig. 3 corresponds to the maximum efficiency of . This is comparable to the maximum efficiency reported by East et al. East et al. (2014) for superradiant scattering of gravitational waves off a rotating black hole, . There is a sharp decrease of towards zero at , beyond which there is no superradiance.

The alternative definition (67) of the superradiant efficiency is also plotted in Fig. 3. Compared with , the maximum of is shifted towards lower frequencies, with a peak value of . The cut-off frequency is virtually the same for both definitions of the efficiency. The quantity shows an almost linear decrease as .

Comparing our findings with perturbative results for monochromatic waves is not straightforward because we work with wave packets of finite extent, the black hole mass and charge change drastically during the scattering process, and the amount of this change depends on the frequency. From monochromatic perturbation theory one would expect superradiance to vanish as , and the upper cutoff frequency should be at Brito et al. (2015); Di Menza and Nicolas (2015)

(68)

For the initial black hole parameters in our simulation this results in , quite close to the observed cutoff in the efficiency in Fig. 4.

Figure 4: Superradiant efficiency (Eq. (64), solid red) and (Eq. (67), dashed blue) as a function of for , , , , , , . Shown are results from several simulations (dots) and a cubic spline interpolant.

In Fig. 5 we report an anomalous case where the black hole charge decreases even though the black hole mass increases by a small amount (about ) and hence there is no superradiance. The scalar field charge is still positive in this evolution; note however that the charge current density (23) may still be negative depending on the form of the scalar field , which is indeed what we observe.

Figure 5: Non-superradiant evolution with decreasing charge (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

Similar behavior (decreasing charge, increasing mass of the black hole) is found generically for negative .

While so far the charge-to-mass ratio of the Reissner-Nordström black hole has been moderate (), Fig. 6 shows the evolution of a near-extreme black hole with . This case is superradiant with efficiency , .

Figure 6: Superradiant evolution for a near-extreme black hole (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

Finally, we present evolutions with non-zero scalar field mass. In Fig. 7 we choose ; the other parameters are the same as in the previous evolutions. This case is superradiant with efficiency , . For (Fig. 8) we obtain a non-superradiant evolution.

Figure 7: Superradiant evolution with scalar field mass (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

Figure 8: Non-superradiant evolution with scalar field mass (, , , , , , , ). Top: black hole mass (solid red), Bondi mass (dashed blue) and irreducible mass (dotted green). Bottom: horizon charge (solid red) and total charge at future null infinity (dashed blue).

As shown in Di Menza and Nicolas (2015) in the test field case, there is no superradiance if , which evaluates to for our choice of parameters. Our results are consistent with this but note that the particular value of at which superradiance ceases to exist depends on the frequency , and can be lower than this bound (in our simulation already was non-superradiant).

For comparison with Di Menza and Nicolas (2015) we also compute the (physical) radius of the effective ergosphere (defined therein as the region where the potential becomes negative). In the case this is , in the case . In both cases the initial data are concentrated around with width so the initial data lie outside the effective ergosphere.

iii.6 Quasi-normal modes and tails

At late times we observe a decay of the scalar field that is well described by an exponential ringdown followed by a power-law tail (Fig. 9; here we return to the massless case ).

Figure 9: Absolute value of the scalar field at the horizon (upper red curve) and at future null infinity (lower blue curve) as functions of time in a log-log plot. The predicted exponential decay during the quasi-normal mode phase (dashed line on the left) and the power-law decay during the tail phase (dashed lines on the right) are shown for comparison. The same parameters as in Fig. 3 are used (, , , , , , , ).

Quasi-normal modes of charged scalar fields on a Reissner-Nordström (and more generally Kerr-Newman) background were computed in Hod (2012); Konoplya and Zhidenko (2013). The relevant regime for the simulation shown in Fig. 9 is that of large , in which the authors find asymptotically with the fundamental () quasi-normal mode frequency

(69)

(In our case there is no rotation.) The imaginary part of (dashed line on the left in Fig. 9) agrees roughly with our simulation. The real part predicted from (69) is and the value measured from the simulation is . (Higher-order corrections in computed in Richartz and Giugno (2014) are negligible in our case.)

A detailed analysis of late-time power-law tails of a charged scalar field on a Reissner-Nordström background was carried out by Hod and Piran in Hod and Piran (1998). In the regime their results imply for spherical perturbations () that at and at the horizon. This agrees reasonably well with our simulation (Fig. 9).

Iv Conclusions

In this paper we investigated superradiance of a spherically symmetric charged scalar field scattering off a Reissner-Nordström black hole in asymptotically flat spacetime. Unlike in previous studies, we solve the fully coupled Einstein-Maxwell-Klein-Gordon system. A hyperboloidal evolution scheme on hypersurfaces of constant mean curvature is used, which is ideally suited to black hole scattering experiments as these slices extend smoothly from inside the horizon to future null infinity.

Our main result is that for sufficiently low frequency of the initial data, superradiance occurs and leads to substantial losses of mass and charge of the black hole. The maximum superradiant efficiency we observed, defined to be the change in black hole mass divided by the change in charge energy of the black hole, was . If the efficiency is defined instead by comparing the total Bondi mass loss with the energy of the initial ingoing scalar field pulse, is obtained. In the non-superradiant regime (increasing black hole mass) we found a somewhat anomalous case in which the black hole charge decreased even though the scalar field charge was positive. Our superradiant evolutions include near-extremal black holes () and nonzero scalar field mass (). It is clear that the massive case is challenging numerically due to the terms in the scalar field evolution equation (18), which are formally singular at (and which vanish analytically only because the massive scalar field falls off faster than any power of at Winicour (1988)). While the evolutions presented here are stable, we do observe a numerical instability for scalar field masses in the range that we have not been able to cure.

We also analyzed the late-time decay of the massless scalar field and found approximate agreement with known perturbative results on quasi-normal modes and power-law tails.

A nice by-product on the theoretical side is the construction of a conserved current using the Kodama vector field, which we related to the Hawking mass as in Csizmadia and Rácz (2010). This enabled us to derive a Bondi mass loss formula.

Acknowledgements.
We are grateful to Lars Andersson, Piotr Bizoń, Maciej Maliborski, Vincent Moncrief, and particularly István Rácz and Claudio Paganini for helpful discussions. This research was supported by grant RI 2246/2 from the German Research Foundation (DFG) and a Heisenberg Fellowship to O.R.

Appendix A Einstein evolution equations and regularity at future null infinity

Here we provide the (redundant) Einstein evolution equations that we do not enforce actively but monitor during the evolution. The evolution equation for the conformal factor is

(70)

and the traceless momentum obeys

(71)

where

(72)