Electromagnetic Emission from a Binary Black Hole Merger Remnant in Plasma: Field Alignment and Plasma Temperature

Comparable-mass black-hole mergers generically result in moderate to highly spinning holes, whose spacetime curvature will significantly affect nearby matter in observable ways. We investigate how the moderate spin of a post-merger Kerr black hole immersed in a plasma with initially uniform density and uniform magnetic field affects potentially observable accretion rates and energy fluxes. Varying the initial specific internal energy of the plasma over two decades, we find very little change in steady-state mass accretion rate or Poynting luminosity, except at the lowest internal energies, where fluxes do not exhibit steady-state behavior during the simulation timescale. Fixing the internal energy and varying the initial fixed magnetic-field amplitude and orientation, we find that the steady-state Poynting luminosity depends strongly on the initial field angle with respect to the black hole spin axis, while the matter accretion rate is more stable until the field angle exceeds $\sim 45\degree$. The proto-jet formed along the black hole spin-axis conforms to a thin, elongated cylinder near the hole, while aligning with the asymptotic magnetic field at large distances.


I. INTRODUCTION
Black holes are the unique end-point of massive stellar evolution in our current understanding of stellar astrophysics, informed by Einstein's general relativity (GR). They are also the inevitable result of the merger of highmass neutron stars, as well as of black holes of all masses, including the supermassive ones believed to reside at the center of most galaxies. Most supermassive black holes are expected to have significant spin through accretion [1]; even during the course of the merger of initially nonspinning holes, enough orbital angular momentum is retained to produce a final hole with a dimensionless spin of ∼ 0.69. This spin angular momentum produces an azimuthal distortion of the surrounding nearby spacetime ("frame-dragging"), acting to concentrate magnetic fields and potentially produce strong steady-state electromagnetic characteristics. Spinning supermassive black holes power active galactic nuclei (AGN) [2,3], with a radio jet likely powered by the black hole's spin, mediated by polodial magnetic field lines pinned to a surrounding ac-cretion disk [4].
We are particularly interested in the scenario of black holes newly formed after merger in a potentially complicated gas-rich environment. Here we expect a transition from whatever was happening through the merger toward a new quasi-steady state centered on the newly formed black hole. Generally the black hole may form in environment characterized by a larger-scale magnetic field, perhaps anchored in the poloidal component of an accretion disk that had surrounded the premerger binary. For some configurations of spinning black hole mergers, the final black hole may be misaligned from the core axis of the broader accretion disk and its poloidal field. It can be particularly valuable to understand the general physics driving jet formation in this kind of environment, which may be robust against detailed variations in the turbulent local environment of the black hole at the point of merger. Toward this we consider here a simplified scenario involving a black hole in an asymptotically uniform magnetic field within a structureless uniform distribution of plasma. Though we are primarily interested in this as a generic model for an accreting binary just after merger, we can also recognize this as a generalization of spherical Bondi accretion with a magnetic field and a spinning black hole.
In the context of the post-merger scenario, we build on [5], hereafter referred to as Paper I, where we used the tools of numerical relativity and ideal general-relativistic magnetohydrodynamics (GRMHD) to investigate how the merger of a comparable-mass black hole binary affect a surrounding plasma. While other studies explore circumbinary accretion disks [6][7][8][9][10], Paper I took a deliberately simplistic approach to initial conditions, in order to elucidate the effects of the merger with minimal assumptions about likely matter configurations. In particular, we explored how an equal-mass black hole binary merger affects an initially uniform density plasma with uniform magnetic field parallel to the binary's orbital angular momentum vector. A parameter survey was performed varying the initial plasma β −1 parameter, to measure the system's dependence on the relative strength of the magnetic field.
We found that the time-development of Poynting luminosity, which may drive jet-like emissions, is relatively insensitive to aspects of the initial configuration. In particular, over a significant range of initial β −1 , the central magnetic field strength is effectively regulated by the gas flow to yield a Poynting luminosity of 10 45 − 10 46 ρ -13 M 8 2 erg s −1 , with the binary black-hole mass M scaled to M 8 ≡ M/(10 8 M ) and ambient density ρ -13 ≡ ρ/(10 −13 g cm −3 ). We also calculated the direct plasma synchrotron emissions processed through geodesic ray-tracing. Despite lensing effects and dynamics, we found the observed synchrotron flux varies little leading up to merger.
Here we extend the results of Paper I, paying special focus to the plasma dynamics near the remnant black hole. In particular, we concentrate on the initial magnetic field angle relative to the remnant black hole's spin axis, and on the initial temperature of the plasma.
Great uncertainty still exists about the environment immediately around supermassive black holes. Even for Sgr A * , the closest, most-studied black hole in the universe, the best estimates for plasma temperature, density, and magnetic field strength vary by orders of magnitude [11][12][13][14]. Similarly, the low-density gas around M87 also appears to be described by a radiatively inefficient accretion flow, but produces powerful radio jets on enormous galactic scales [15][16][17][18]. Therefore we acknowledge that the parameters used in this paper represent only a small region of the potential astrophysical parameter space, but we will show that these idealized conditions still provide valuable insight into some of the fundamental questions about the behavior of magnetized accretion flows around supermassive black holes.
We begin by describing the numerical methods we used in our simulations in Section II. In Section III, we provide the details of the parameter space survey, and in Section IV, we present results from all configurations consid-ered. Subsection IV A considers how varying the initial specific internal energy (a proxy for temperature) affects bulk behavior such as mass accretion rates and Poynting luminosity. Subsection IV B considers how varying the initial magnetic field direction affects bulk behavior such as mass accretion rates and Poynting luminosity. In Section IV C, we more closely investigate the nature of the "proto-jet" region that develops in the vicinity of rotating spacetimes, introducing several measures to help quantify the jet features. Throughout our paper, unless otherwise noted, we use geometrized units where G = c = 1, and Greek (Latin) indices are space-time (space) indices.

II. METHODS
In Paper I, the aforementioned black hole binary simulations in an initially uniform plasma were carried out using the "moving puncture" formalism [19,20], with a simultaneous evolution of the space-time metric and MHD fields, using the McLachlan [21,22] implementation of the BSSNOK equations [23][24][25] for the former, and the IllinoisGRMHD [26,27] implementation of the conservative GRMHD equations (see e.g. [28]).
For the new simulations presented here, since we consider the post-merger end-state of the system, it is more computationally efficient to use a fixed Kerr background with mass and spin appropriate to the spacetime after the merger of an equal-mass, nonspinning binary with initial ADM mass M = 1: M = 0.97, a/M = 0.69 [19,20,29,30]. There are still infinitely many ways to express such a spacetime as a metric; we choose the horizon-penetrating "Kerr-Schild" slicing used by [31,32], as implemented by NRPy+ [33,34]. This slicing has the advantage of placing the horizon at a fixed constant radial coordinate value, identical to that of the better-known Boyer-Lindquist slicing: r hor = r + = M + √ M 2 − a 2 . It has been used for IllinoisGRMHD evolutions of Fishbone-Moncrief initial conditions as part of the community Event Horizon Telescope comparison project [35] and validated within the NRPy+ infrastructure to satisfy the ADM constraints. While the numerical simulations of each configuration use the Einstein Toolkit's Cartesian mesh-refinement driver called Carpet [36,37], for the purpose of post-processing data analysis with GRMHD analysis, a suite of Python-based tools [38], we regularly interpolate MHD fields to a spherical polar grid. Additionally, avoiding the spacetime evolution greatly reduces the computational cost of our studies.
Our primary diagnostic is again the EM (Poynting) luminosity: where S r is the radial component of the relativistic Poynting vector, expressed in terms of the fluid four-velocity u a and magnetic four-vector b a : Rather than rely on the EinsteinToolkit Multipole thorn's output of the spherical harmonic (l, m) components of L Poynt , as we did in [5], we output the 3D MHD field data onto the aforementioned spherical coordinate mesh with uniform sampling in each of the coordinate directions (r, θ, φ) and perform the analysis in post-processing. This output procedure used first-order Lagrange polynomial interpolation, as supplied by the EinsteinToolkit.
We apply the same post-processing suite to estimate the rate of accretion of fluid into the black hole as well (while in [5] we used the Outflow code module in the Einstein Toolkit). In particular, we calculate the flux of fluid across the event horizon S via: where D ≡ ραu 0 is the Lorentz-weighted fluid density, and σ i is the ordinary (flat-space) directed surface element of the horizon. Armed with these, the Poynting efficiency η EM can be computed from

III. NEW ILLINOISGRMHD RUNS
Our investigations began with a canonical case, 'KS', of a single Kerr black hole in Kerr-Schild coordinates surrounded by plasma with uniform density and pressure, initially satisfying a polytropic equation of state P = κρ Γ , with Γ = 4/3, appropriate for a radiationpressure-dominated gas. The plasma is threaded by a uniform-magnitude magnetic field oriented parallel to the hole's spin axis (k). The initial fluid density, pressure, and magnetic field strength are the same as those used in Paper I's canonical configuration, yielding a fluid that is everywhere magnetically sub-dominant, with β −1 = P mag /P gas = 0.025. The canonical configuration's pressure is dominated by the radiation P rad = (a/3)T 4 , implying a temperature of T = 2.906 × 10 5 K. Working from this canonical case, we carried out two suites of simulations at moderate resolutions.
In the first suite, we kept the magnetic field oriented parallel to the spin axis, but varied the initial polytropic coefficient κ in the uniform plasma, and thus the uniform specific internal energy , and hence the gas temperature of the plasma. These configurations are presented in Table I.
In the second suite, we kept the initial canonical temperature fixed, and varied the angle θ B between the initial TABLE I. Initial field configurations for the canonical case and temperature-varied simulations. v Alf is the Alfvén speed (Eq. B2); cs is the fluid sound speed (Eq. B1); the temperature T is deduced assuming a radiation-dominated gas. For analysis, the fields representing fluid density ρ, fluid pressure p, fluid three-velocity v i , and magnetic field B i were interpolated onto an evenly spaced sphericalpolar grid of size r ∈ [0.35M, 150M ], θ ∈ [0, π], and φ ∈ [0, 2π], with N r = 450, N θ = 50, N φ = 50; hence ∆r ≈ M/3, ∆θ ≈ π/50, ∆φ ≈ 2π/50. The interpolation method used was a simple first-order Lagrange polynomial interpolation scheme, supplied by EinsteinToolkit's Carpet mesh-refinement driver [39].

A. Dependence on plasma temperature
In Paper I, we investigated the dependence of the Poynting luminosity on initial density and magnetic field strength while holding fixed the initial specific internal energy 0 . As noted then, the luminosity should satisfy the scaling relation where σ 0 ≡ b 2 0 /(2ρ 0 ) is the initial ratio of magnetic to rest-mass energy density, and F (t/M ; 0 , σ 0 ) is a dimen-sionless function of time. Paper I primarily addressed the σ 0 dependence of F (t), while leaving 0 fixed.
Here we investigate the variation in infall rate and Poynting luminosity with 0 , which serves as a proxy for the initial plasma temperature T 0 . For a Gamma-law gas, where we have assumed Γ = 4/3. The set of configurations are presented in Table I. In analogy to our varying of magnetic field strength in Paper I, here we vary the polytropic constant κ over two orders of magnitude; as a result the initial gas pressure P 0 and specific internal energy 0 also vary over two orders of magnitude, while the temperature varies by roughly a factor of three. In Fig. 1 we show the Poynting luminosity L Poynt during the evolution of each of the initial temperature configurations listed in Table I. We see that the luminosity generally takes longer to settle down with higher 0 , due to the lower Alfvén speed in these cases. After settling, however, the late-stage luminosity shows little variation with 0 , except for the case of the very lowest 0 . This case shows extremely high luminosity, which shows no sign of settling down over the simulation time.  Table I. Dotted lines indicate data from a higher-resolution run.
To complement the Poynting luminosity, in Fig. 2 we show the accretion rate over time of the same configurations. Again, the lowest-temperature case, 0 = 0.06, shows the least stable behavior.
The settling-down time for these configurations also depends on temperature, being later for the highertemperature cases. To investigate the steady state for each configuration, we use a time-average value for each configuration, from a common starting time of t = 2, 000M until the end of the available data. In Fig. 3, we show the resulting Poynting luminosity L Poynt (top Accretion rate |Ṁ | as a function of time for the temperature configurations listed in Table I. Color labels are the same as for Fig. 1. panel), accretion rateṀ (middle panel), and resulting efficiency η EM (Eq. 4) (bottom panel). For each configuration, the "error bars" shown are simply the standard deviation over the time window.
It is noticeable that both the Poynting luminosity L Poynt and mass accretion rateṀ are highest for the lowest values of κ, and hence fluid temperature, though subject to greater variations in time. There also appears to be a shallow local minimum in L Poynt around κ = 0.1, and inṀ around κ = 0.5; the combination of these yields a minimum in efficiency η EM around κ = 0.3, close to our canonical case. However, as this is a shallow minimum, the efficiency is around 20% over most of our temperature range.

B. Dependence on magnetic field orientation
Here we investigate the effect of varying the angle θ B between the Kerr spin vector a and the initial orientation of the uniform magnetic field B. In practice, we fix the former -a = ak -and vary the latter. However we demonstrate in Appendix C that we achieve equivalent results when fixing the field direction and varying a instead.
In Fig. 4 we show the late-time state of the magnetic field integral curves passing near the central black hole, for the KS B45deg configuration. The black hole has not only twisted and concentrated the field, but has tilted it toward the spin axis (z direction), but only out to a radius r < ∼ 30M . In Fig. 5 we show the time-development of the Poynting luminosity L Poynt during the evolution of each of the initial magnetic-field orientations θ B . It is clear that the "post-settling" luminosity has a strong dependence on θ B .  Looking at the late-time (t > ∼ 1, 500M ) behavior of the systems, in Fig. 6 we plot L Poynt as a function of initial inclination angle θ B . We also show a fit (dashed red line) of these L Poynt data points to a functional form quadratic in the cosine of θ B , similar to that seen by [40] in the force-free limit. Our results seem to show a flatter behavior at low and high θ B , captured better by a hyperbolic tangent dependence on θ B (solid blue line), but cannot rule out the cos 2 θ B scaling. It is entirely possible Steady-state (t > 1, 500M ) Poynting luminosity as a function of field alignment angle θB. The luminosity is calculated as a "late-time average" value in each casethe average value for all t > 1, 500M . Error bars show the RMS deviation from the time-average values, beginning at t = 1, 500M . The solid (blue) and dashed (red) curves are best-fit results from assuming a hyperbolic tangent or cosinesquared dependence on θB, respectively. that the inclusion of MHD and matter (as opposed to the force-free scenario) introduces additional physics scaling that lead to the step-function behavior.
As we can see in Fig. 4, even at late times, the magnetic field lines are only oriented toward the BH spin axis relatively close to the hole itself, remaining substantially along its initial direction further out. We can try to quantify the transition region from the BH's "sphere of influence" by examining the Poynting luminosity over a set of extraction spheres. In Fig. 7, we show the integrand in Eq. 1 -essentially the Poynting vector, weighted by the local area measure -as a function of (θ, φ) for R ∈ {20M, 30M, 40M, 50M } for the KS B45deg configuration. We see that the angular location (i.e. "point in the sky") of peak contribution moves with extraction radius; we also see that the tube seems to contract in angle. We will attempt to quantify these observations in Sec. IV C.
In Fig. 8, we show the rate of mass loss into the Kerr horizon,Ṁ (Eq. 3) during the evolution of each of the initial magnetic-field orientations θ B . Again, the accretion rates for different θ B show little variation until t ≈ 300M . Even at late times, the different configurations'Ṁ deviate by only around 50%, with the highest rates associated with the greatest deviation of the initial magnetic field angle. As with the Poynting luminosity, we can produce a time-averaged accretion for the steady state (t > 1, 500M ) of each configuration. This is presented in Fig. 9. Viewed in this way, we see that the steady-state accretion rate is relatively constant for 0°≤ θ B < ∼ 40°, dropping off steeply for larger θ B .
In Fig. 10, we plot the resulting efficiency (Eq. 4). Dominated by the field orientation, it shows levels of ∼ 25% for small θ B , dropping an order of magnitude for θ B > ∼ 40°.

C. Features of Proto-Jet
In studies of black-hole neutron star mergers, [41] identify an "incipient, magnetized jet" as an "unbound, collimated, mildly relativistic outflow (Lorentz factor of ∼ 1.2), which is at least partially magnetically dominated". Informally, we identify a "proto-jet" as a magnetically dominated region showing concentrated twisting of magnetic field lines, and strong localized Poynting flux. [42][43][44] We use the term "proto-jet" here, because while it shows intense winding of magnetic fields in a tra- ditional jet-like funnel region, the net fluid flow in this region is inward, with low Lorentz factor. In this subsection, we attempt to clarify this definition by studying more carefully the nature of the magnetic fields and Poynting vector at late times.
Beam size. To evaluate the Poynting luminosity at a radius R, we integrate the Poynting vector over a coordinate sphere at that radius (Eq. 1). Figure 7 shows the distribution of the integrand (that is, the Poynting vector weighted by the local angular Jacobian) over the sphere, with contours showing regions containing 10%, 50%, and 90% of the Poynting flux. We estimate the size of the proto-jet by calculating the solid angle subtended by the 50% contours. At late times, we can plot this solid angle as a function of extraction radius. In Fig. 11, we show the width of the beam in the northern hemisphere as measured from the 50% contours for each of the configurations. We see that the solid angular width is smaller for the larger field initial inclination angles θ B . Moreover, the width generally decreases with radius, especially for configurations with θ B > ∼ 15°(upper panel). For extraction radii R ∈ {20M, 40M }, the falloff in angular width is approximately 1/R 2 -fast enough to keep the jet's absolute cross-sectional area roughly constant, or "pencil-like" (lower panel) Beam shape. As can be seen from Fig. 7, the crosssectional shape of the beam deviates strongly from circular when the magnetic field is misaligned with the black hole spin. We present in Fig. 12 the beam shape as represented by the 50% contour for a range of field alignments, measured at R = 30M .
For an aligned field, the beam cross-section is annular at all extraction radii, as the magnetic field drops to zero on the axis due to symmetry. Here we see that the beam shape becomes steadily less annular with increasing θ B . Simultaneously, the overall luminosity decreases, and the beam weakens, becoming harder to distinguish Beam position. We present in Fig. 13 the positions of the center of the proto-jet for each configuration, showing how it varies with extraction radius. To avoid highfrequency variations, at each extraction radius R, we decomposed the Poynting vector over the sphere into (real) spherical harmonics up to = 2: The center positions are then the maxima of this smoothed functional form. While the jet positions are properly given as a pair of angles (θ, φ), we find it easier to display as a pair of Cartesian-like projected coordi- nates X ≡ sin θ cos φ, Y ≡ sin θ sin φ, so that the hole's spin direction lies unambiguously at the origin in each panel.
From the figure, we can see that all configurations have jet directions that approach the asymptotic initial magnetic field direction at large R (denoted by × in the figure). As we move inward along each configuration's curve, we see twisting of the jet direction around the origin (that is, the BH spin axis). For initial inclination angles θ B between 0°(i.e. parallel to the spin axis) and ∼ 60°, the jet direction approaches the spin axis for small R. For larger θ B , the jet's direction stops short of the pole. Jet Strength. We noted at the start of this subsection that our "proto-jet" has not yet been demonstrated to produce ultra-relativistic particle speeds. In particular, as in Paper I, fluid inflow in the jet region is both subrelativistic and inward-pointing. While analyzing the aftermath of a BHNS merger, [41] encounter a similar situation; they point out, however, that strong magnetic dominance in the asymptotic jet region is expected to lead to much higher Lorentz factors: Γ ∼ b 2 /2ρ [45].
In our case, the peak energy ratio drops to below ∼ 5 outside a few horizon radii, implying that actual relativistic jet conditions may not be reached for the fluid particles present. This can be misleading, as the MHD fluid is ion-dominated, and unlikely to be the source of significant high-energy EM emission. If a mechanism is present to seed the magnetically pressure-dominated region with electrons or electron-positron pairs, these can be expected to experience much greater accelerations, leading perhaps to jet-like electromagnetic emission.

V. DISCUSSION
In this paper, we have extended the work of [5] (Paper I), focusing on the steady-state behavior of plasma around a post-merger Kerr black hole. While Paper I featured merging equal-mass nonspinning black-hole binary systems, with a spacetime dynamically simulated via the moving puncture formalism, here we concentrated on the end-state of such a merger, a single spinning Kerr black hole with mass M final = 0.97M , and dimensionless spin a/M final = 0.69. Since the spacetime is quiescent after merger, we used a fixed Kerr-Schild metric representation of the spacetime to reduce the computational load, and simplify post-simulation analysis. Nevertheless, the initial MHD fields are fully dynamical, and our canonical MHD configuration is that of Table I of Paper I. We have concentrated on the Poynting luminosity L Poynt , mass accretion rateṀ , and resulting "Poynting efficiency" η EM of different initial plasma configurations.
First we investigated the dependence of L Poynt andṀ on the specific internal energy of the plasma for fixed fluid density. We find that higher-temperature configurations take longer to settle down to a steady state, due to the lower Alfvén speed in these cases. After reaching a steady state, it is noticeable that both the Poynting luminosity L Poynt and mass accretion rateṀ are highest for the lowest temperatures, though subject to greater variations in time. There also appear to be shallow local minima in both steady-state L Poynt andṀ around moderate temperatures; the combination of these yields a Poynting efficiency η EM of around 20% over most of our temperature range, with a shallow minimum close to our canonical configuration.
Returning to the canonical configuration, we also studied the result of varying the angle θ B between the asymptotic magnetic field and the Kerr hole's spin direction. We found that L Poynt falls swiftly with θ B , consistent with expectations from force-free MHD. We find the mass accretion rate is less sensitive to θ B until around 40°, dropping steeply thereafter.
Finally, we investigated the form of the "proto-jet" formed by the spacetime dragging of plasma and magnetic field lines. We showed how the jet width varies with radius, and how the jet direction moves continuously from aligning with the black-hole spin axis for small R( < ∼ 5M ) to adopting the direction of the asymptotic magnetic field for R > ∼ 30M . Having summarized the main results of our investigations, we may ask what implications these have for the astrophysical question of electromagnetic counterparts of black-hole mergers. As noted above, we see little variation in either mass-accretion rate or Poynting luminosity over a broad range of specific internal energies around our canonical value. We do see a strong dependence on the angle θ B between the black hole's spin and the global magnetic field, with luminosity dropping quickly as the misalignment angle increases. Additionally, while the proto-jet is aligned with the black hole spin in the strong-gravity region (within a few Schwarzschild radii of the horizon), it soon relaxes to lie parallel to the initial asymptotic magnetic field direction. Expectations of jet alignment are ambiguous in the absence of surrounding matter: should they align with the hole's spin or with the asymptotic magnetic field [40,46]? Our results indicate a transition between the two states, with the hole's spin's influence declining rapidly with distance. Assuming that the asymptotic magnetic field is seeded by the plasma, our result here agrees qualitatively with GMRHD explorations by [47] of tilted-disk simulations of highly spinning black holes, where the jet aligns with the black hole spin at r = 4M , but with the disk's angular momentum at r = 40M . This latter observation raises the question of what we should assume for the shape, strength, and orientation of the external magnetic field. This is a complicated question, beyond the scope of this paper, which we have neglected in favor of a survey over orientations, regardless of cause.
In astrophysical units, the steady-state luminosity for our canonical case is consistent with the "peak" luminosity of Paper I, but drops steeply as the angle θ B between spin and asymptotic magnetic field increases, to about 10% of its maximum. The rate of drop-off in θ B is consistent with, the cos(θ B ) 2 expectations from forcefree models, or with a slightly steeper step-function. We can combine the results of Paper I with the θ B dependence seen in Fig. 6 to obtain a more general expression for the steady-state Poynting luminosity at arbitrary spin inclination angle θ B to the asymptotic magnetic field direction: where H(θ B ) is a function that captures the smooth steplike behavior observed in Fig. 6.
In performing the studies presented here, we have fulfilled some of the additional investigations outlined in the discussion of Paper I, focusing on the bulk behavior of MHD fields around the post-merger Kerr black hole. Since the work here was performed on a post-merger stationary Kerr background spacetime, a full radiation transport treatment of the resulting MHD fields using, e.g. the Pandurata code [48] could not be expected to produce novel results including an EM signature of the merger process, and we did not attempt it here.
Meanwhile, other extensions of Paper I are being carried out in parallel to this work, looking at the effect of significant spins on the pre-merger black holes [49]. With these in place, we anticipate turning our attention to rotationally supported matter distributions and magnetic field configurations, and to more realistic inclusion of radiation transport with the fully dynamical merging binary metric, using new developments in Pandurata. The form of the background Kerr metric used for the evolutions is one frequently called "Kerr-Schild" by accretion-disk theorists [32]:
As the radial and polar coordinates here are unchanged from that of the original Boyer-Lindquist form, the horizon is still a coordinate sphere defined by ∆(r) = 0: r + = M + √ M 2 − a 2 . total matter content of the domain (and hence the MHD fluid density ρ) is coupled to the black hole mass M . With M = 1 in code units, we chose ρ 0 = 10 −6 . With this initial code-units density, the horizon mass of the black hole increased by ∼ 20% over the course of 1000M of evolution. Dropping to ρ 0 = 10 −8 and scaling b 0 accordingly decreased the accretion rate by a factor of ∼ 100, indicating that the accretion roughly scales linearly with the fluid density.
A code-units density of ρ 0 really means ρ 0 M/M 3 . The length unit appearing in the denominator is ∼ 1.5(M/M )km; thus the density in physical units will be ρ 0 M/(3.375 · (M/M ) 3 )km −3 . For our canonical total mass M = 10 8 M , this becomes ρ = ρ 0 10 8 M 3.375 · 10 24 km −3 ≈ 50ρ 0 g · cm −3 Thus our choice of ρ 0 = 10 −6 in code units is equivalent to a physical density of ∼ 5 × 10 −5 g · cm −3 . As this is more than seven orders of magnitude greater than our assumed canonical plasma density (10 −13 g cm −3 ), we conclude that for all configurations considered in the main text, the black hole mass could increase by no more than one part in 10 6 , even with feedback switched on.

Varying B or a direction?
As mentioned in Section IV B, we investigate the θ Bdependence of our results by keeping a ≡ ak, and setting B = B cos θ B . for one representative case, however, we instead chose the spin vector to be oriented as a = a sin π/4î + a cos π/4k, with B = Bk.
In Fig. 14, we show the Poynting luminosity L Poynt for two configurations with θ B = 75°, with either a or B fixed along the z axis. We can see that the two curves track exactly until after the peak, where small differences begin to set in.
In this figure we also demonstrate how the luminosity changes with resolution, and with field orientation in the equatorial plane.