Black Hole Polarimetry I: A Signature of Electromagnetic Energy Extraction

In 1977, Blandford and Znajek showed that the electromagnetic field surrounding a rotating black hole can harvest its spin energy and use it to power a collimated astrophysical jet, such as the one launched from the center of the elliptical galaxy M87. Today, interferometric observations with the Event Horizon Telescope (EHT) are delivering high-resolution, event-horizon-scale, polarimetric images of the supermassive black hole M87* at the jet launching point. These polarimetric images offer an unprecedented window into the electromagnetic field structure around a black hole. In this paper, we show that a simple polarimetric observable -- the phase $\angle\beta_2$ of the second azimuthal Fourier mode of the linear polarization in a near-horizon image -- depends on the sign of the electromagnetic energy flux and therefore provides a direct probe of black hole energy extraction. In Boyer-Lindquist coordinates, the Poynting flux for axisymmetric electromagnetic fields is proportional to the product $B^\phi B^r$. The phase $\angle\beta_2$ likewise depends on the ratio $B^\phi/B^r$, thereby enabling an observer to experimentally determine the direction of electromagnetic energy flow in the near-horizon environment. Data from the 2017 EHT observations of M87* are consistent with electromagnetic energy outflow. Currently envisioned multi-frequency observations of M87* will achieve higher dynamic range and angular resolution, and hence deliver measurements of $\angle\beta_2$ closer to the event horizon as well as better constraints on Faraday rotation. Such observations will enable a definitive test for energy extraction from the black hole M87*.


INTRODUCTION
The Event Horizon Telescope (EHT) has produced resolved images of the hot, magnetized, synchrotronemitting plasma around the supermassive black holes in the Galactic Center (Sgr A*; The Event Horizon Telescope Collaboration et al. 2022) and the elliptical galaxy M87 (M87*; The Event Horizon Telescope Collaboration et al. 2019a) on scales comparable to their respective projected black hole event horizons.The 230 GHz images of both sources feature similar ring-like morphologies.Despite a difference of three orders of magnitude in the black hole mass, both rings have diameters d ≈ 10 GM/Dc 2 ,1 consistent with predictions from analytic and numerical models of low-luminosity accretion flows in the extreme near-horizon regime, where the ef-Corresponding author: Andrew Chael achael@princeton.edufects of strong gravitational lensing and redshift become important (The Event Horizon Telescope Collaboration et al. 2019b, hereafter EHTC V).
Near-horizon polarimetric images carry more information than total-intensity images and more stringently constrain models of the accretion flow and jet-launching region.In particular, EHT images of M87* reveal a helical pattern in the electric vector position angles (EV-PAs) of the linearly polarized intensity around the emission ring (The Event Horizon Telescope Collaboration et al. 2021a, hereafter EHTC VII).While EHT images of M87* in total intensity are consistent with a wide variety of astrophysical models (EHTC V), the polarimetric images strongly prefer models of magnetically arrested (MAD) accretion flows (The Event Horizon Telescope Collaboration et al. 2021b, hereafter EHTC VIII).In MAD flows, the horizon-scale magnetic field is strong, ordered, and dynamically important (Bisnovatyi-Kogan & Ruzmaikin 1974;Narayan et al. 2003).
Strong horizon-scale magnetic fields near black holes naturally produce powerful electromagnetic outflows.Penrose (1969) showed that energy can be extracted from a black hole's spin via particle interactions within the ergosphere.This version of the Penrose process may have astrophysical implications for high-energy emission (e.g., Williams 1995;Schnittman 2015), but the primary astrophysical channel for tapping a black hole's spin energy is widely believed to be the Blandford-Znajek (BZ;Blandford & Znajek 1977) mechanism, whereby electromagnetic energy is extracted from a black hole's spin via magnetic fields that thread its event horizon.
The BZ mechanism is the dominant source of energy outflow in general-relativistic magnetohydrodynamic (GRMHD) simulations of MAD accretion flows around spinning black holes (e.g., Tchekhovskoy et al. 2011;McKinney et al. 2012; EHTC V).The BZ mechanism is thus a natural candidate for the launching power of extragalactic jets (Begelman et al. 1984), including the famous jet from M87, which has been observed in over a century of multi-frequency observations to deliver energy from horizon to galactic scales (e.g., EHT MWL Science Working Group et al. 2021;Lu et al. 2023).GRMHD simulations of MAD flows featuring BZ powered jets provide good matches to both the observed 230 GHz image morphology of M87* near the event horizon (EHTC V; EHTC VIII), as well as the jet power, profile, and core-shift farther downstream (Chael et al. 2019;Mizuno et al. 2021;Cruz-Osorio et al. 2022).
Polarimetric EHT observations were critical for constraining the space of astrophysical models for M87*'s accretion in EHTC VIII, and for providing evidence for strong, ordered, horizon-scale magnetic fields.However, it has not been clear if or how EHT observations can directly probe the electromagnetic energy flow close to the black hole or, more ambitiously, test whether M87's jet is truly powered by the BZ mechanism.
In this paper, we will show that, provided that the sense of rotation of magnetic field lines is known, and provided that Faraday rotation of linear polarization in the emission region or by an external screen is not too severe, the local direction of electromagnetic energy flow is directly mapped onto near-horizon polarimetric images observed by the EHT.In particular, the observable ∠β 2 , a quantity that characterizes the helicity of the spiral of polarization vectors around the black hole image (Palumbo et al. 2020), can be used to infer the direction of the Poynting flux.As we will show, the EHT's polarimetric image of M87* (EHTC VII) is consistent with outflowing electromagnetic energy on scales of ∼ 5 GM/c 2 .While this energy outflow on horizon scales near the jet-launching point suggests that the black hole spin of M87* truly powers the extragalactic jet, current EHT observations cannot yet conclusively rule out that the rotation power of the accretion disk (Blandford & Payne 1982) also plays a significant role.
In the following sections, we develop the connection between ∠β 2 and the direction of electromagnetic energy flux using simple arguments, and we test this connection in both simple analytic models and complex numerical simulations.In section 2, we review the main features of degenerate, stationary, axisymmetric magnetospheres in the Kerr spacetime.We show that electromagnetic energy outflow requires a certain sign of the ratio of magnetic field components B ϕ /B r in Boyer-Lindquist coordinates.We argue that, for synchrotron-emitting plasmas, this ratio is probed by the observed polarization helicity encoded in the observable ∠β 2 .In section 3, we investigate semi-analytic models of the emission close to M87* and show that relativistic effects of aberration, parallel transport, and lensing do not substantially change the relationship between ∠β 2 and the sign of energy outflow.Furthermore, using the Blandford & Znajek (1977) monopole solution to model the magnetic fields close to M87*, we show that there is a clear relationship between black hole spin and ∠β 2 in observed images, which results from larger spins "winding up" the magnetic fields more rapidly and thereby increasing the ratio |B ϕ /B r |.
In section 4, we move from simple analytic models to an analysis of a full library of GRMHD-simulated images.Remarkably, these reveal that the connection between the helicity of the EHT polarization spiral and the direction of electromagnetic energy flow continues to hold even when considering a time-dependent, turbulent, Faraday-rotating plasma.While we will show that published EHT observations are consistent with an outflowing electromagnetic energy flux in M87*, they still do not unambiguously connect the energy flux with the extraction of black hole spin energy.To do this and test the BZ mechanism, we argue in section 5, will require pushing EHT measurements of the linear polarization pattern closer to the apparent boundary of the black hole event horizon, or "inner shadow" (Chael et al. 2021).With ongoing and planned upgrades to the EHT expected to increase its resolution and dynamic range, these tests may become feasible within the next decade.We summarize our findings in section 6.
In the Appendices, we provide a detailed review of properties of degenerate electromagnetic fields around black holes.While these results are not new (and draw heavily on previous reviews such as Gralla & Jacobson 2014), they present a unified treatment of the field structure and energetics around black holes that we use throughout this paper and will reference in later papers in this series on black hole polarimetric images.In particular, in Paper II (Lupsasca et al. in prep), we will use these results to quantify the expected near-horizon polarization pattern for field lines that thread the event horizon, providing a direct future observational test of black hole energy extraction.In Paper III (Wong et al. in prep), we will develop simplified analytic models for the dependence of the near-horizon polarization on black hole spin, providing a physical interpretation of existing GRMHD simulation results and a path toward quantitatively constraining black hole spin using polarimetry.

BLACK HOLE ENERGY EXTRACTION AND
THE SIGN OF ∠β 2

Kerr metric in Boyer-Lindquist coordinates
We use units where G = c = 1.We work in the Kerr spacetime of a black hole with mass M and angular momentum J = aM = a * M 2 in Boyer-Lindquist (BL) coordinates (t, r, θ, ϕ).The metric (Equation C78) is expressed in terms of three poloidal functions ∆ = r 2 − 2M r + a 2 , Σ = r 2 + a 2 cos 2 θ, Π = r 2 + a 2 2 − a 2 ∆ sin 2 θ. (1) The metric determinant is √ −g = Σ sin θ. (2) The event horizons are located at radii where ∆ = 0; in particular, the outer Kerr event horizon lies at a radius Normal observers in BL coordinates have a four-velocity η µ = (−α, 0, 0, 0), where α = 1/ −g tt = ∆Σ/Π is the lapse function.For a ̸ = 0, BL normal observers have a nonzero angular velocity At the event horizon, the angular velocity of the normal observer equals that of the horizon: (5) 2.2.Stationary, axisymmetric, and degenerate electromagnetic fields We consider electromagnetic fields F µν in Kerr that are both degenerate and magnetically dominated: As shown in Appendix A and B, for such fields there exists an infinite family of timelike frames u µ in which the electric field vanishes: e µ = F µν u ν = 0. 2 In particular, electromagnetic fields obeying the equations of ideal GRMHD are constrained to be degenerate and magnetically dominated.In GRMHD, the fluid fourvelocity picks out a unique frame u µ from among the family of frames in which e µ vanishes.Fields obeying 2 This family of timelike frames can be parameterized by the Lorentz factor γ ≥ γ ⊥ relative to a normal observer.The frame with minimal Lorentz factor γ ⊥ = 1 − E 2 /B 2 −1/2 has a fourvelocity that is perpendicular to the normal observer's magnetic field.See Appendix A.
the equations of general-relativistic force-free electrodynamics (GRFFE) are also required to be degenerate, though without additional constraints they can evolve into states in which magnetic domination does not hold (F 2 < 0).The time-space components of the Maxwell tensor ⋆F define the "lab-frame" magnetic field: B i = ⋆F i0 . 3In Appendix C and D, we review (see, e.g., Phinney 1983;Gralla & Jacobson 2014) how degenerate electromagnetic fields F that are also stationary (∂ t F → 0) and axisymmetric (∂ ϕ F → 0) are completely characterized by three quantities: the magnetic flux function (or poloidal potential) ψ(r, θ) = A ϕ (r, θ), the poloidal current I(ψ), and the field-line angular velocity Ω F (ψ).In particular, the lab-frame magnetic field B i in BL coordinates is determined by the coordinate-independent flux ψ and current I as For stationary, axisymmetric, degenerate fields in BL coordinates (Equation 7), the radial electromagnetic energy flux density is (see subsection C.4) In BL coordinates, B ϕ diverges on the horizon as 1/∆, but the product ∆B ϕ , and hence the energy flux density J r E , remain finite.In particular, Appendix D shows how by considering the field in Kerr-Schild coordinates, which are regular on the horizon, one can derive the Znajek condition (Znajek 1977) that fixes the current I on the horizon as a function of Ω F (see also MacDonald & Thorne 1982, Thorne et al. 1986).We can also express the Znajek condition as a condition on the ratio ∆B ϕ /B r in BL coordinates: Hence, the radial Poynting flux through the horizon is From Equation 10, it is clear that the energy extraction from a black hole is maximized when the field-line angular velocity is Ω F = 1 2 Ω H (Blandford & Znajek 1977).This result holds for all values of the black hole spin a.The total outward electromagnetic energy flux through a sphere of radius r is ) is consistent with a ratio B ϕ /B r < 0 in the emission region.Assuming that the black hole spin is aligned with ẑ, that the magnetic field is axisymmetric and degenerate, and that the magnetic field lines co-rotate with the black hole, this field structure produces an outward Poynting flux by Equation 8.
Assuming that Ω F = Ω H /2 and expanding in powers of Ω F , plugging Equation 10 into Equation 11gives recovering the familiar result that the total energy extracted from a spinning black hole by electromagnetic fields is proportional to the square of Ω H and the square of the horizon magnetic flux Φ H = [B r √ −g] r=r+ dθdϕ.
For the BZ monopole solution, the proportionality constant is k = 1/24π. 4With a modified constant k to reflect differences in the field geometry, the scaling in Equation 12 is observed in force-free (e.g., Tchekhovskoy et al. 2010) and GRMHD (e.g., Tchekhovskoy et al. 2012;Narayan et al. 2022) simulations of black hole jets, even up to large values of the black hole spin.By inspection of Equation 8, it is clear that if the field lines co-rotate with the black hole (Ω F > 0), then outward energy flow J r E > 0 requires a ratio of the BL lab-frame fields B ϕ /B r < 0. This preferred orientation of the azimuthal field should manifest in images of synchrotron-emitting electrons close to the event horizon, such as those obtained by the EHT of M87* at 230 GHz (EHTC VII; EHTC VIII). 5 4 Tchekhovskoy et al. (2010) use Φ to represent the flux through one hemisphere of the horizon, in which case k = 1/6π. 5While B ϕ and B r are both coordinate-dependent quantities, both the total energy flux and the observed polarization pattern from the near-horizon accretion flow are coordinate-invariant.Thus, while we motivate and derive the connection between the observed linear polarization and the direction of energy flux using BL coordinates, our main conclusions are frame-independent.

Polarization and energy flux in near-horizon images
The direction of the polarization vector produced locally by synchrotron emission is perpendicular to its wavevector and to the magnetic field.In the absence of the complicating effects of Faraday rotation and parallel transport through curved spacetime, the pattern of linear polarization in a resolved image can be used to directly infer the magnetic field geometry around a black hole.We consider these important complications in more realistic models in section 3 and section 4.
Figure 1 shows a cartoon picture of the synchrotron polarization pattern for several simple axisymmetric magnetic field configurations viewed face-on.Here, we have assumed that the emitting matter is in the equatorial plane and that there is no vertical magnetic field, B θ = 0.In Figure 1, we orient the +ẑ axis into the page, corresponding to an observer inclination θ o ≈ π.We further assume that the spin of the black hole is oriented along +ẑ into the plane of the sky, matching the spin orientation inferred for M87* by the EHT (EHTC V).
As in EHTC VIII, we use the β 2 statistic (Palumbo et al. 2020) to quantify the structure of near-horizon linear polarization in EHT images.β 2 is the second azimuthal Fourier mode of the linear polarization image: In Equation 13, I is the total intensity brightness, P = Q + iU is the complex linearly polarized brightness, and (ρ,φ) are polar coordinates in the image plane.For narrow, axisymmetric rings of polarized emission viewed face-on, ∠β 2 is equal to twice the value of the EVPA at φ = 0.6 Neglecting all relativistic effects and assuming a faceon observer and a narrow axisymmetric emission ring of radius r, the observed value of ∠β 2 can be simply related to the BL field components in the emission region by In the cartoon picture of Figure 1, where the +ẑ axis points into the page, observed ∠β 2 values in the range (−π, 0) correspond to fields with B ϕ /B r < 0.
If we assume that the black hole spin is oriented along +ẑ and that magnetic field lines co-rotate with the black hole (Ω F > 0), then negative values of ∠β 2 correspond to electromagnetic energy outflow for axisymmetric, stationary, degenerate fields.By contrast, ∠β 2 values in the range (0, π) correspond to fields with B ϕ /B r > 0 and would therefore represent a field configuration with energy flowing into the black hole.
In this simple picture, we can thus constrain the sign of the near-horizon electromagnetic energy flux J r E by measuring ∠β 2 , assuming that we know the signs of the inclination cos θ o and field-line angular velocity Ω F .In the next subsection, we turn to more general orientations of the black hole spin, field-line angular velocity, and observer inclination.

General relationship between ∠β 2 and J E
In a black hole magnetosphere, the black hole spin, magnetic field lines, and emitting matter all have angular momenta that can in principle take independent orientations.Throughout this work, we assume all of these are either aligned or anti-aligned with the +ẑ coordinate axis (θ = 0) in BL coordinates.We thus use signed scalars to represent the black hole spin a, field-line angular velocity Ω F , and fluid angular velocity v ϕ = u ϕ /u t .Each scalar is either positive or negative depending on whether it is aligned (positive) or anti-aligned (negative) with the +ẑ-axis.The observer at infinity can view the system at any inclination angle in the range 0 ≤ θ o ≤ π.
In principle, we could follow a convention in which the black hole spin is always positive by aligning the +ẑ axis of our coordinate system with the black hole spin.In GRMHD simulations, however, the community's convention for simulating retrograde accretion flows is to use an initial equilibrium torus far from the black hole with a positive angular velocity v ϕ > 0, and to set the black hole spin along the −ẑ direction, a < 0. For this reason, in this work, we also treat retrograde accretion flows in both analytic models and simulations by setting a < 0.
From Equation 8, for arbitrary signs of black hole spin and field-line angular velocity, it is always true for stationary axisymmetric magnetic fields that Furthermore, a rule of thumb for the sign of ∠β 2 that accounts for arbitrary observer inclination is Combining Equation 15and Equation 16 yields an approximate relation valid for all inclinations θ o and values of Ω F : 2.6.Application to M87* In order to match the M87* total-intensity image asymmetry, EHTC V found that the fluid angular velocity in the emission region must be directed into the plane of the sky.Furthermore, EHTC V and Wong et al. (2021) found that in nearly all GRMHD simulations, the material producing the 230 GHz emission observed by the EHT corotates with both the black hole spin and the magnetic field lines (see also Figure 4).
That is, for prograde accretion flows (a > 0), both v ϕ > 0 and Ω F > 0 in the emitting region.In order for the image to match the observed brightness asymmetry, the observer inclination must then be θ o ≈ π.Electromagnetic energy outflow in this case requires B ϕ /B r < 0 (Equation 8), which by Figure 1 is observable as −π < ∠β 2 < 0.
For retrograde accretion models (a < 0), instead both v ϕ < 0 and Ω F < 0 in the emitting region.To match the brightness asymmetry in this case requires θ o ≈ 0. For there to be energy outflow in this case, the azimuthalto-radial field ratio in BL coordinates must be positive, B ϕ /B r > 0. However, because θ o ≈ 0 in this case, outflowing electromagnetic energy still corresponds to an observable ∠β 2 in the range −π < ∠β 2 < 0.
In summary, we do not know either the sign of cos θ o or Ω F in M87* independently.From EHTC V, however, we do know that sign (Ω F cos θ o ) = −1.Thus, we can use the value of ∠β 2 measured by the EHT to constrain the direction of the near-horizon electromagnetic energy flow in the idealized picture of Figure 1.The observed range of ∠β 2 for M87* from the EHT is −163 deg < ∠β 2 < −129 deg (EHTC VII).We thus find from Equation 17 that the EHT results are consistent with an outward Poynting flux in M87*'s black hole magnetosphere.
This correspondence between the observed ∠β 2 with the direction of Poynting flux in M87*'s magnetosphere is complicated by several effects.First, relativistic parallel transport (from the curved spacetime close to the horizon) and aberration effects (from the rapidly rotating emitting plasma) may result in observed polarization directions that are not perpendicular to the projected magnetic fields in the emission region.We explore these effects in analytic models in section 3. Second, realistic electromagnetic fields around supermassive black holes are not stationary or axisymmetric, and propagation effects through the turbulent plasma (especially Faraday rotation) may significantly alter the emitted polarization direction.We explore these effects in full GRMHD simulation images of M87* in section 4.

∠β 2 AND ENERGY FLUX IN SEMI-ANALYTIC MODELS
In this section, we examine whether the straightforward connection between the observed ∠β 2 in polarized images of synchrotron radiation, the ratio B ϕ /B r , and the direction of electromagnetic energy flux we established in section 2 persists in semi-analytic models of the near-horizon emission from M87*.
We restrict ourselves to axisymmetric, equatorial emission models in this section.We describe the degenerate electromagnetic field in these simple models using the same variables common in GRMHD simulations (e.g., Gammie et al. 2003).That is, we use the fluid-frame four-velocity u µ and the "lab-frame" magnetic field As a toy model of a degenerate electromagnetic field configuration in the Kerr spacetime with outward Poynting flux, we adopt the perturbative Blandford-Znajek monopole solution (Blandford & Znajek 1977) for the potential ψ and current I (see subsection E.5).We then determine B r and B ϕ by Equation 7; B θ vanishes in the equatorial plane.While the BZ monopole solution does not capture the detailed structure of a black hole magnetosphere as seen in GRMHD, GRFFE, or particle-in-cell simulations, it has been shown to match certain time-averaged properties of these simulations well (e.g., McKinney & Gammie 2004;Parfrey et al. 2019).
Real black hole magnetospheres may have a significant vertical magnetic field component in the equatorial plane (Blandford & Globus 2022); the paraboloidal BZ solution may provide a better match to the near-horizon magnetic field structure in these systems (e.g., McKinney & Narayan 2007;Penna et al. 2013).We will explore polarized images of the paraboloidal solution in a future work.For now, we take the BZ monopole as a toy model that captures both an outgoing energy flux and a field ratio B ϕ /B r that depends on black hole spin, similar to what is observed in GRMHD simulations (e.g., Palumbo et al. 2020;Emami et al. 2023).
For the fluid velocity u µ , we explore two different models.First, we use the unique field-perpendicular or driftframe velocity u µ ⊥ determined by the full BZ solution for both the lab-frame magnetic field B i = ⋆F i0 and electric field E i = F 0i (see Equation D137, Appendix B).The drift-frame velocity in BL coordinates produces an outflow in the equatorial plane; while we can include arbitrary boosts along the magnetic field and not change the observed polarization pattern, this family of fourvelocities is not fully representative of black hole accretion flows.
Consequently, we also explore a family of parametrized equatorial inflow models for u µ from Cárdenas-Avendaño et al. (2023).This velocity model is described by three parameters (ξ, β r , β ϕ ) taking values between 0 and 1.The parameter ξ represents the ratio of the fluid angular momentum to the Keplerian value, while β r and β ϕ represent the relative magnitude of the Keplerian and infall velocities, respectively.
In addition to the parameters (ξ, β r , β ϕ ), we can also choose the (sub)Keplerian velocity to be either prograde or retrograde with respect to the black hole spin.We summarize this four-velocity model in Appendix F.
After fixing the magnetic field and velocity profile, we model the 230 GHz emitting region as a narrow ring of constant rest-frame emissivity J(r eq ) in the equatorial plane centered at r ring = 2r + .We generate images from this model using the kgeo code (Chael 2023).We present more details of the model construction in Appendix G.
We assume uniform priors on the velocity model parameters (ξ, β r , β ϕ ) and generate model images covering 250,000 combinations of these parameters for each black hole spin a * ∈ [±0.1, ±0.3, ±0.5, ±0.7, ±0.9],where negative spins correspond to retrograde flows.We fix the observer inclination θ o = 163 deg for prograde flows (θ o = 17 deg for retrograde flows), as measured for M87* (Mertens et al. 2016).We also generate a single image from each model using the outflowing drift-frame velocity u µ ⊥ .We then blur the images to the EHT resolution of 20 µas using a circular Gaussian kernel and compute ∠β 2 using (Equation 13) We plot histograms of the results for ∠β 2 over the velocity model parameter space in Figure 3.The star symbols show the results for the drift-frame BZ velocity u µ ⊥ for each black hole spin.With the exception of a few retrograde, low-spin (a * = −0.1)models, nearly all models in this parameter space have −π < ∠β 2 < 0, corresponding to energy outflow in the simple picture of Figure 1.Unlike our arguments in section 2, these models include effects from parallel transport, aberration and nonzero inclination.They demonstrate that the sign of ∠β 2 is a robust probe of energy extraction from black holes, relatively immune to uncertainty in the kinematics of the emitting plasma.
We also show blurred EHT-resolution images of three extreme cases corresponding to fully prograde and retrograde Keplerian orbiting material at radius r = 2r + for a spin |a * | = 0.5 and the outflowing BZ drift velocity solution.While the images from a uniform emission ring for these two velocity models are quite different (in particular, the overall Stokes I image asymmetry in the retrograde model is substantially reduced from the prograde model), both have ∠β 2 values less than zero.
There is a clear spin-dependence to the value of ∠β 2 in Figure 2. In particular, higher-spin black holes produce values of ∠β 2 closer to zero.In the BZ monopole solution, higher-spin black holes "wind up" magnetic field lines closer to the black hole.As a result, they produce more radial polarization patterns when observed nearly face-on.This same spin-dependence is present in GRMHD simulation images (Palumbo et al. 2020) (see also Figure 3).We will discuss the comparison of analytic models for this spin dependence with GRMHD simulations and EHT observations in more detail in Paper III.

∠β 2 AND ENERGY FLUX IN GRMHD SIMULATIONS
In this section, we investigate snapshots from magnetically arrested (MAD) GRMHD simulations of M87*.Unlike the analytic models in section 3, GRMHD simulations feature complex, non-axisymmetric, timedependent structure in their magnetic fields and fluid velocities.Furthermore, the radiative transfer from GRMHD simulations fully accounts for a non-uniform synchrotron emitting region, self-absorption, and Faraday rotation and conversion of the emitted polarization along the photon trajectory.We find that the picture presented in Figure 1 and Figure 2 remains valid even when subject to these astrophysically important complications.In particular, when observing MAD simulations of M87* at 230 GHz with the black hole spin oriented into the plane of the sky, nearly all simulation snapshots have ∠β 2 < 0, consistent with outward electromagnetic energy flux, as in the simple picture of Figure 1.
We use nine MAD GRMHD simulations from Narayan et al. ( 2022 We generate 230 GHz images from these simulations using the GR radiative transfer code ipole (Noble et al. 2007;Mościbrodzka & Gammie 2018) following the parameters discussed in Appendix H.In Figure 3, we show the distributions of β 2 in the complex plane and ∠β 2 for all snapshot images in the library, colored by the black hole spin.The distributions for a given spin cover a wide range of simulation snapshots and electron-to-ion temperature ratio R high (Mościbrodzka et al. 2016).Nearly all snapshots have −π < ∠β 2 < 0, matching the intuition from section 2 that negative ∠β 2 corresponds to outward electromagnetic energy flow in M87*, assuming that the black hole spin vector is into the plane of the sky.
As in the semi-analytic models of M87* using the BZ monopole (Figure 2), there is also a clear spindependence among the distributions of ∠β 2 in the GRMHD images (Figure 3).Higher-spin simulations produce ∠β 2 values closer to zero (see also Palumbo et al. 2020, EHTC V, Jia et al. 2022 for examples of this same trend in other simulation image libraries).This trend of ∠β 2 with spin in GRMHD simulations is a result of more rapidly spinning black holes "wind- ing up" the magnetic field lines more rapidly, producing more toroidal fields in the emission region (Emami et al. 2023).The black hole spin in GRMHD simulations of MAD accretion flows thus significantly alters the magnetic field structure of the flow in the 230 GHz emission region, with the same qualitative dependence on spin as in the BZ monopole model.In Paper III, we will investigate in more detail the connection between ∠β 2 and spin in both simple models and GRMHD simulations.
It was not immediately obvious that the relationship between the sign of ∠β 2 and the direction of electromagnetic energy flux should be as robust in full GRMHD simulations as it is in the cartoon picture of section 2 or the simple models of section 3.For one, the instantaneous dynamics in GRMHD simulations are nei-ther time-stationary nor axisymmetric, as assumed in section 2 and section 3. Perhaps more surprisingly, the simulations also frequently have significant Faraday rotation in the emission region (Ricarte et al. 2020, EHTC VIII).A significant degree of internal Faraday rotation is in fact necessary to sufficiently de-polarize 230 GHz simulation images to match EHT observations (EHTC VIII), and one might expect Faraday rotation to shift ∠β 2 substantially from its "intrinsic" value without Faraday rotation.Nonetheless, we find that Faraday rotation in the MAD GRMHD simulation images of M87* used here depolarize the 230 GHz images but do not produce an overall EVPA rotation severe enough to change the sign of ∠β 2 .

DISCUSSION
We have established that, assuming we know the sign of the magnetic field line angular velocity Ω F , the helicity of linear polarization in near-horizon synchrotron radiation around a black hole-quantified by the polarimetric Fourier mode ∠β 2 -is a probe of the direction of electromagnetic energy flow in the system.We tested this relationship in both simple axisymmetric models and in images from full 3D, turbulent GRMHD simulations.The latter include non-axisymmetric and nonequatorial emission regions, relativistic parallel transport, aberration, and Faraday rotation.Remarkably, none of the complicating effects present in GRMHD simulations changed the qualitative insight of Figure 1.For M87*, negative values of ∠β 2 such as those observed by the EHT correspond to an outward electromagnetic energy flux on spatial scales of ∼ 5GM/c 2 .
Here, we probe the relationship between ∠β 2 , the magnetic field structure, the electromagnetic energy flow, and spin more deeply in the magnetically arrested GRMHD simulations from section 4. In Figure 4, we show time-and azimuth-averaged quantities (see Appendix I) from the simulations for three values of the black hole spin: a * = −0.5 (left column), a * = 0 (middle column), and a * = +0.5 (right column).In all panels, the cyan contour shows the surface where the magnetization σ = b 2 /ρ = 1 and the white contours show surfaces of constant potential ψ ≡ A ϕ , corresponding to magnetic field lines in the axisymmetrized data.Direct synchrotron emission observed by the EHT can arise in the equatorial plane or at higher latitudes, but it typically originates at characteristic radii in the range r emis ≈ 3 − 5 r g (EHTC V).
The top row of Figure 4 shows the sign of the ratio B ϕ /B r in BL coordinates in the averaged data from the three simulations.The middle row shows the angular frequency of field lines Ω F (in units of t −1 g ), and the bottom row shows the outward Poynting flux √ −gJ r E .In the bottom row, each panel is normalized independently relative to its maximum value, as the overall magnitude of the electromagnetic energy flux in the a * = 0 simulation is far smaller than in the two simulations with nonzero spin.
In the prograde simulation average, the field ratio B ϕ /B r < 0 everywhere.In the retrograde simulation, B ϕ /B r > 0 inside and somewhat exterior to the σ = 1 contour where the field corotates with the black hole.At larger radii, B ϕ /B r < 0 and the field rotates with the disk in the opposite sense to the black hole.Both panels reinforce the physical insight of Equation 8; outflowing Poynting flux requires the field to be wound up such that sign(Ω F B r B ϕ ) < 0 (Equation 15).Both spinning simulations have field-line angular velocity Ω F that are much larger on field lines connected to the black hole horizon than on field lines that are unconnected to the black hole.The value of |Ω F | on horizon-penetrating field lines is close to the value predicted by Blandford & Znajek (1977): for |a * | = 0.5.The field lines that are connected to the black hole are also where the outward electromagnetic energy flux √ −gJ r E is strongest.It is this outward electromagnetic energy flux, extracted from the black hole spin energy, that powers the large-scale jets in these simulations.
The behavior of B ϕ /B r is more complicated in the a * = 0 model.In this simulation, electromagnetic energy flows outward at large radii on field lines disconnected from the black hole via the Blandford & Payne (1982) mechanism; that is, outflow in this simulation is driven by accretion power rather than black hole spin energy.At smaller radii and on field lines connected to the black hole horizon, electromagnetic energy flows into the black hole in the a * = 0 simulation.In order to have electromagnetic energy outflow (inflow), the product Ω F B r B ϕ must be negative (positive).Close inspection reveals that these relationships are always satisfied in the averaged data from the a * = 0 simulation.On field lines disconnected from the black hole, B ϕ /B r < 0 and Ω F > 0, producing a weak energy outflow.On field lines close to and connected to the black hole horizon, the ratio B ϕ /B r and the angular velocity Ω F do not have fixed signs, but the total product Ω F B r B ϕ is always positive.In individual snapshots from the a * = 0 simulation, the field ratio and angular velocity are less ordered and much more turbulent than in the a * = ±0.5 simulations, which largely resemble their time-averaged structure.
In Figure 5, we investigate the behavior of ∠β 2 as a function of radius in these three simulations.Here, we fix R high = 1 and consider the time-averaged 230 GHz images from the simulations, computing ∠β 2 from the time-average image at native resolution in annuli of width 0.2 µas.When considering ∠β 2 as a function of radius in the a * = −0.5 simulation, there is an obvious sign change at an image radius of 7r g ≈ 25 µas, corresponding to where the fields transition from counterrotating to corotating with the black hole in Figure 4.The image-averaged value of β 2 in this case is dominated by emission close to the black hole where the magnetic field lines corotate; as a result, the majority of retrograde snapshots still have image-averaged ∠β 2 < 0 in Figure 3. Changes in the sign of ∠β 2 with radius can be used to diagnose changes in the sign of B ϕ /B r and Ω F as the flow begins to corotate near the horizon (Ricarte et al. 2022).
All three simulations show interesting behavior in ∠β 2 in the range of image radii corresponding to the n = 1 photon ring, where there are contributions from photons that have executed a half-orbit or more around the black hole (Johnson et al. 2020).Shifts in ∠β 2 in the n = 1 photon ring from parallel transport and radiative transfer effects (Himwich et al. 2020;Jiménez-Rosales et al. 2021;Palumbo & Wong 2022) may be detectable in EHT targets on long 345 GHz baselines (Palumbo et al. 2023).Furthermore, in Figure 5 there is a clear radial trend in ∠β 2 in all three simulations as the image radius ρ approaches the value corresponding to the directly lensed image of the event horizon (the "inner shadow"; Chael et al. 2021).Close to the horizon, the ratio B ϕ /B r changes rapidly with radius.The precise behavior of ∠β 2 with radius close to the inner shadow depends strongly on black hole spin.Interestingly, in both simulations with nonzero black hole spin, ∠β 2 trends to zero or positive values as the image radius approaches the inner shadow.In the a * = 0 simulation, ∠β 2 becomes more negative, likely due to the necessarily inwardflowing energy flux very close to the horizon.In Paper II and Paper III, we will explore the radial dependence of ∠β 2 and its value at the inner shadow in analytic models and GRMHD images in detail.
Figure 4 indicates that the electromagnetic energy outflow that powers the jet in MAD simulations originates via the BZ mechanism and is concentrated on field lines that thread the black hole event horizon.In this paper, we have argued that the observed ∠β 2 from EHT images of M87* also indicates electromagnetic energy outflow.However, current observations are not conclusive as to whether or not the observed 230 GHz synchrotron emission originates on field lines that thread the event horizon (thereby extracting energy via the BZ mechanism) or if they only thread the accretion disk at larger radii, extracting energy from the accretion disk's rotation (as in the a * = 0 simulation in Figure 4).
Future EHT observations with more sites and at 345 GHz frequency will have higher resolution and significantly more dynamic range than the pioneering 2017 observations of M87* (Doeleman et al. 2019;Raymond et al. 2021).By imaging the faint linear polarization signal closer to the event horizon, these observations may be able to conclusively determine if the observed emission originates on horizon-threading field lines that extract black hole spin energy.With more work on calibrating simulation and analytic models, and accounting fully for the effects of Faraday rotation, these observations have the potential to conclusively determine whether or not the extragalactic jet from M87 is powered by the spin of the M87* supermassive black hole.
6. CONCLUSION In this paper, we have investigated the link between the structure of resolved polarimetric images of synchrotron radiation near a black hole event horizon and the direction of electromagnetic energy flux in the black hole magnetosphere.We have shown that: • The sign of a polarimetric observable from EHT images of supermassive black holes-∠β 2 , which quantifies the helicity of the observed linear polarization-can be used to infer the sign of the ratio B ϕ /B r in simple axisymmetric models.
• Provided the orientation of the angular velocity of magnetic field lines around a black hole is known, the sign of the radial electromagnetic energy flux J r E is determined by the sign of the ratio B ϕ /B r .• Full GRMHD simulations including Faraday effects and non-axisymmetric fields (Figure 3) show the same overall trends in ∠β 2 as we expect from idealized analytic arguments (Figure 1) and see in simple semi-analytic models (Figure 2).These results provide strong support for connecting ∠β 2 and the direction of energy flux in real polarized images of black holes.
• If the emission comes from a magnetic field line threading the horizon, then the sign of ∠β 2 and its dependence on image radius measures whether energy is flowing into or being extracted from the black hole itself.
• Applied to published EHT observations of M87* (EHTC VII; EHTC VIII), these conclusions are consistent with electromagnetic energy flowing outward in M87*'s magnetosphere.
• As first noted in Palumbo et al. (2020), GRMHD values of ∠β 2 show a strong dependence on black hole spin, with larger spins featuring more azimuthal magnetic fields in the synchrotronemitting region.This trend (although not the exact values of ∠β 2 ) is well-reproduced in images from the Blandford-Znajek monopole model (Figure 2).Larger spins "wrap up" the magnetic field to be more azimuthal in the 230 GHz emission region.This physical intuition and correspondence between models suggests a path forward for robustly measuring black hole spin with β 2 .
• It is unclear from current EHT images whether or not the observed emission lies on field lines that thread the black hole event horizon and thus extract energy directly from the black hole spin.Future observations with an expanded EHT array will have more resolution and dynamic range and should be able to observe polarized emission from close to the projected event horizon, or "inner shadow" (Chael et al. 2021).These observations will be critical for distinguishing whether the outward energy flow observed by the current EHT is a result of the Blandford-Znajek process (Blandford & Znajek 1977) or whether it is powered by the accretion flow somewhat farther out from the event horizon (Blandford & Payne 1982).
Software: kgeo (Chael 2023), KORAL (Sądowski et al. 2013(Sądowski et al. , 2014)) These appendices review key properties of degenerate electromagnetic fields in 4D spacetime, and more particularly in the background of a Kerr black hole.Each appendix begins with a detailed summary of its results.
First, in Appendix A, we define degenerate 2-forms and describe equivalent ways of characterizing them.Our discussion is mathematical, but the results that we derive are then shown in Appendix B to be relevant for electromagnetic fields in ideal magnetohydrodynamics (GRMHD) and force-free electrodynamics (FFE).
In Appendix C, we focus on stationary, axisymmetric electromagnetic fields in the Kerr spacetime.We give explicit forms for such fields in both Boyer-Lindquist and Kerr-Schild coordinates.We also derive their energy and angular momentum fluxes in both coordinate systems.
Then, these two threads join in Appendix D, where we specialize to degenerate, stationary, axisymmetric fields in Kerr.We review how such fields may be completely described by three functions in the poloidal (r, θ) plane: the magnetic flux function ψ(r, θ), the field-line current I(ψ), and the field-line angular velocity Ω(ψ).We derive the Znajek (1977) condition for regularity of these fields on the horizon, and we describe its implications for the horizon fluxes of energy and angular momentum.Our treatment covers similar ground as Carter (e.g., 1979); Gralla & Jacobson (e.g., 2014).
Next, Appendix E reviews the known solutions to FFE in the Kerr background, with a particular emphasis on the Blandford & Znajek (1977) monopole solution that we use in section 3.In Appendix F, we review the fluid velocity model that we used to produce the simulated images shown in section 3; Appendix G explains our procedure for generating these images in detail.Finally, Appendix H describes the GRMHD simulations that we used in section 4, and our strategy for timeand azimuth-averaging the GRMHD simulation data is described in Appendix I.

A. DEGENERATE FORMS
In this appendix, we first review general 2-forms in 4 dimensions before specializing to degenerate 2-forms, for which we present several characterizations and then prove their equivalence (sections A.1 and A.2).
We show that if a degenerate 2-form is magnetically dominated-the case of astrophysical relevance-then it must admit a timelike vector in its kernel (section A.3), and that if it is closed-as in electromagnetism-then it can be described in terms of two scalar Euler potentials (section A.4).We also introduce an "electromagnetic" decomposition for general 2-forms, which simplifies for degenerate fields (section A.5), and an associated local frame (section A.6).
At this stage in our discussion, these statements are purely mathematical, but their physical significance will become clearer in Appendix B, where we will examine degenerate 2-forms that also obey Maxwell's equations.
A.1.General 2-forms in 4 dimensions Consider a 4-dimensional spacetime M with metric g, and let v • w = g µν v µ w µ denote its inner product.The Levi-Civita tensor is defined as , where g = det g µν is the metric determinant and [µναβ] denotes the completely antisymmetric symbol, which is defined to be ±1 according to whether µναβ is an even or odd permutation of the coordinates {0, 1, 2, 3}, and zero otherwise. 7The volume form is the contraction which, by definition, is the 2-form such that9,10 , where brackets denote antisymmetrized indices, and The same manipulations as in (A2) reveal that .
8 The wedge product A ∧ B of a p-form A with a q-form B is the (p + q)-form given by the antisymmetrized tensor product Finally, a direct computation shows that 11 By ( A6)-(A7), the degeneracy of F is also equivalent to which implies that F is degenerate if and only if ⋆F is.
A 2-form F has a nontrivial kernel if there exists some A 2-form F is simple if it is a wedge product F = v∧w of two 1-forms v and w, or explicitly, This decomposition is far from unique, since As a result of this freedom, we can assume without loss of generality that v and w are orthogonal. 12 We now show that these three properties are in fact all pointwise equivalent, that is, that they all imply each other provided that we work locally at a fixed point in spacetime, treating F µν as a matrix rather than a tensor field. 13Our treatment mirrors the discussion of degeneracy in Gralla & Jacobson (2014).
First, suppose that F is simple.Then F = v ∧ w for some 1-forms v and w.As such, F ∧F = v∧w∧v∧w = 0 (by antisymmetry of the wedge product) and hence, F is degenerate.Conversely, suppose that F is degenerate, so that F [µν F αβ] = 0.There must nonetheless exist at least two nonzero vectors x and y such that f 2 ≡ F µν x µ y ν ̸ = 0 (or else, F = 0 identically).From (A5), it follows that Letting f v ν ≡ x µ F µν and f w ν ≡ y µ F µν , this shows that is manifestly orthogonal to w = w ′ and v∧w = F by construction. 13The equivalence almost holds at the level of fields, except that degeneracy only implies simplicity locally: there is no guarantee that the vectors x and y in (A11), which are defined in each tangent space separately, can be glued into smooth vector fields across spacetime; §(3.5.35) of Penrose & Rindler (1984) presents a counter-example that illustrates this phenomenon.
Hence, (A10) holds and we conclude that F is simple.This proves that degeneracy ⇔ (local) simplicity.
Next, suppose that F ̸ = 0 is simple, so that F = v ∧ w for some (orthogonal) 1-forms v and w.Then at every spacetime point p, the vectors v and w define a 2-plane in the tangent space T p M ∼ = R 4 .Every point on this codimension-2 surface is intersected by a perpendicular 2-plane spanned by vectors x and y, such that the set of vectors (v, w, x, y) forms an orthogonal basis of R 4 .14Thus, x and as such, it follows from (A10) that x µ F µν = y µ F µν = 0. Hence, F admits a nontrivial (2-dimensional) kernel spanned by x and y.Conversely, suppose that F admits a nontrivial kernel.Then by assumption, det F µν = 0, and so it follows from (A7) that F µν (⋆F ) µν = 0.By (A6), this implies that F is degenerate and therefore simple.Thus, a nontrivial kernel ⇔ (local) simplicity.This concludes our proof.
We now argue that if F is a magnetically dominated, degenerate 2-form, then there exists a timelike vector u such that u µ F µν = 0.By the preceding discussion, since F is degenerate, it is also (locally) simple, and so at every point in spacetime it can be written as F = v∧w for some orthogonal 1-forms v and w.Then by (A10), Since F 2 > 0, the signs of v 2 and w 2 must be identical.Moreover, this sign must be positive, so that v and w are both spacelike (since they must be orthogonal and there is only one independent timelike direction in spacetime).
As a result, the nontrivial kernel of F , which is spanned by x and y, must admit a timelike vector u (some linear combination ax+by) such that u µ F µν = 0.By a suitable rescaling, we may always assume u to have unit norm.Conversely, if F is a degenerate 2-form with a timelike vector u in its kernel, then it must be magnetically dominated.Indeed, given that there is a unique timelike direction t in spacetime, if t lies in the kernel of F = v∧w (which is spanned by x and y), then t cannot lie within the perpendicular 2-plane spanned by the orthogonal vectors v and w, which must therefore be spacelike.By (A13), it then follows that F 2 > 0. This shows that for a degenerate F , magnetic domination (F 2 > 0) ⇔ the existence of a (unit-norm) timelike vector u in its kernel (u µ F µν = 0). 16Evidently, a degenerate F always has a spacelike vector in its kernel.

A.4. Euler potentials for closed, degenerate 2-forms
Given any two scalars λ 1 and λ 2 , the simple 2-form F = dλ 1 ∧ dλ 2 is manifestly degenerate and also closed: dF = 0 because d 2 = 0. 17 Conversely, if a degenerate 2-form F is also closed, then there (locally) exist scalars λ 1 and λ 2 -known as Euler potentials 18 -such that (A14) Gralla & Jacobson (2014) present multiple proofs of this fact in their §3.2.As discussed therein, the pair of Euler potentials (λ 1 , λ 2 ) is not unique, as there are infinitely many other pairs , but the intersections of hypersurfaces or constant λ 1 and λ 2 are well-defined (i.e., independent of the choice of pair).
If F solves the Maxwell equations (B24), then these 2-dimensional surfaces are called field sheets, and if F is magnetically dominated, then these sheets are timelike.
The inverse relations provide a general decomposition of F and its dual ⋆F in terms of three vectors u, e and b, in agreement with Baumgarte & Shapiro (2003) and (only when u 2 = −1) with ( 4)-(5) of McKinney (2006).Physically, if F is an electromagnetic field strength and u is a unit-norm timelike vector, then the electric and magnetic fields in the local frame of an observer with 4-velocity u are given by e and b, respectively.As such, we will call e and b the electric and magnetic fields, even when F does not obey the Maxwell equations (B24).
In terms of e and b, the scalar invariants F 2 = F µν F µν and (⋆F ) 2 = (⋆F ) µν (⋆F ) µν take the form 20 which is manifestly independent of the choice of u in the decomposition (A16).By (A6)-(A7), the degeneracy of F is equivalent to the frame-invariant property 18 In plasma physics, they are also known as Clebsch coordinates.
If F is both degenerate and magnetically dominated, then (as we have shown) there must exist a unit-norm timelike vector u in its kernel, so that in the frame of an observer with 4-velocity u, the electric field e µ = u ν F µν vanishes.In other words, a 2-form F satisfies both (A18) and (A19) if and only if there exists a vector u such that e = 0, in which case (A16) simplifies to 21 making it manifest that ⋆F = b ∧ u is a simple 2-form.
We reiterate that the choice of u is not unique.We will explicitly construct all such 4-velocities in (B40) below.

A.6. Local frames for degenerate 2-forms
The antisymmetry of F µν and (⋆F ) µν guarantees that the projections (A15) satisfy u • e = u • b = 0.Moreover, if F is degenerate, then (A18) also holds, and so the vectors u, e and b are mutually orthogonal.In that case, the following vector fields define an orthogonal frame (provided that u is not in the kernel of F , so e ̸ = 0): Using (A16), we can compute the projections If u is timelike (u 2 < 0), then e, b and z are spacelike, and we can form an orthonormal frame with frame fields where B. DEGENERATE ELECTROMAGNETISM In this appendix, we now consider degenerate 2-forms that also describe physical electromagnetic fields.Such 2-forms F must obey the Maxwell equations, which in terms of the current J sourcing the field take the form 23 We give expressions for the electric and magnetic fields in the frames of the "normal observer," of the "lab," and of a fluid.We also construct the unique one-parameter family of 4-velocities for which the electric field vanishes, and we show explicitly that velocity components parallel to the electric field do not enter the stress-energy tensor.

B.1. Electromagnetic stress-energy tensor
A 2-form F that obeys the Maxwell equations (B24) describes an electromagnetic field with stress-energy The change in the energy-momentum is equal to the force exerted, so the Lorentz force density is24 where the last step follows from (B24) and (B25).The general projection (A15)-(A16) decomposes T EM µν as 22,25 where u µ z ν +u ν z µ = − (u µ ϵ νρκλ + u ν ϵ µρκλ ) u ρ e κ b λ , and is a projection operator onto hypersurfaces normal to u.

B.2. Application to GRMHD
A prime example of degenerate electromagnetism is general-relativistic magnetohydrodynamics (GRMHD).In GRMHD, u is the 4-velocity of the plasma sourcing the electromagnetic field strength F , and e and b are the electric and magnetic fields in the fluid rest frame.In the ideal GRMHD approximation, the plasma is assumed to be a perfect conductor.As a consequence, the electric field in its rest frame is completely screened (assuming that there is sufficient free charge): this is the "ideal MHD condition" e µ = u ν F µν = 0 (Gammie et al. 2003).By the preceding discussion, e = 0 implies that F is both degenerate (e•b = 0•b = 0) and magnetically dominated (b 2 > 0 = e 2 ), so that it can be decomposed as in (A20).
The fluid 4-velocity u is not the only timelike vector field in the kernel of F .In fact, as we will explicitly show in (B40), there are infinitely many such vectors u ′ (each with its own associated magnetic field b ′ ), allowing for infinitely many decompositions of When F is degenerate and magnetically dominated (as in GRMHD), its kernel contains a unit-norm, timelike u (the fluid 4-velocity).Then in the rest frame of u, e = 0 and the stress-energy tensor (B25) takes the simple form

B.3. Normal observer and lab frame
Given a time coordinate t and spatial coordinates x i , the metric ds 2 = g µν dx µ dx ν decomposes as where the lapse α and shift vector β i are defined by Here, Latin indices run only over the spatial coordinates, and we reserve Greek indices for spacetime coordinates.
The normal observer is defined to have the 4-velocity which is manifestly timelike with unit norm, η • η = −1.
Numerical GRMHD codes often work with the normal observer's electric and magnetic fields (Noble et al. 2006) where we have now replaced the Greek spacetime index ν in (A15) by a Latin spatial index i since by construction, E and B have vanishing time components.Occasionally, instead of working in the frame of the normal observer (B32), GRMHD codes use as primitive variables which may be regarded as the electromagnetic fields in the "lab frame" defined by the (non-normalized) vector In flat spacetime, the lab frame becomes normal (ζ = η), and E i = E i and B i = B i reduce to the usual electric and magnetic fields.In asymptotically flat spacetimes, the lab frame describes an observer "at rest at infinity"

B.4. Relation between fields in different frames
We now assume that F is degenerate and magnetically dominated, and let u denote a unit-norm, timelike vector in its kernel.We can then use (A20) to relate the electric and magnetic fields e and b in the frame of u, given by (A15), to the lab-frame fields E and B given by (B34): In terms of the projection tensor (B28), which in this case is simply P µν = g µν + u µ u ν , this can be recast as where ζ • u = −u t , while γ denotes the Lorentz factor of the flow relative to the normal observer frame, B.5. Explicit timelike frames for magnetically dominated, degenerate fields We showed that if F is degenerate and magnetically dominated, then there must exist a unit-norm, timelike vector u in its kernel.In fact, this vector is not unique.
We now explicitly construct all possible such vectors.The answer takes the general form26 where the subscript (γ) indicates that the Lorentz factor (B39) of this frame relative to the normal observer is γ.
To prove this, we work in the orthonormal frame (A23) associated with the normal observer, with frame fields The general form of a vector u expanded in this basis is This condition is necessary but not sufficient, since27 is generically nonzero.Thus, to ensure that u really lies in the kernel of F , we must also demand that This then leaves us with the general linear combination Lastly, we impose the normalization condition u•u = −1.
Since the frame (B41) is orthonormal, this amounts to Solving this equation for γ results in Since (u • Y ) 2 ≥ 0 and B 2 > E 2 (magnetic domination), this implies that there is a minimum Lorentz factor γ 0 : At last, rewriting (B48) as gives us the most general unit-norm timelike vector u in the kernel of F , parameterized by its Lorentz factor γ relative to the normal observer frame (B41) and a sign: This concludes the derivation of (B40).When γ = γ 0 , this expression agrees with (17) of McKinney (2006).

B.6. Field-perpendicular and field-parallel velocities
The 4-velocity decomposition (B42) can be recast as where γṽ = (u Z is a purely spatial 3-velocity.In fact, ṽ lies fully within the spatial slices normal to the timelike vector η, and it is the spacelike vector obtained by projecting u onto these surfaces: Indeed, ṽt = 0 by definition (as γ = αu t and η t = 1/α) and since η Typically, numerical GRMHD codes (e.g., McKinney & Gammie 2004) do not evolve the 4-velocity u of the flow directly, but rather its spatial projection (B52), which by (B53) obeys the nice relation If F is degenerate and magnetically dominated (as in GRMHD), then its associated timelike flow must lie in the family (B50), and hence ṽ admits a decomposition ṽ = ṽ⊥ + ṽ∥ , (B55) where the (magnetic) field-perpendicular velocity is 28 while the (magnetic) field-parallel velocity is Since ṽ∥ vanishes when γ = γ 0 , (B53) implies that so that the purely field-perpendicular flow (with no fieldparallel component) has the minimal Lorentz factor Conversely, when the Lorentz factor diverges (γ → ∞), the field-parallel velocity is maximized: We can thus parameterize the spatial velocity (B55) as Lastly, we note-as in (18) of McKinney (2006)-that the coordinate 3-velocity is For instance, when i = ϕ is an azimuthal coordinate, the flow (B40) has angular velocity Ω (γ) ≡ u ϕ (γ) /u t (γ) .

B.7. Representations of the stress-tensor
When F is degenerate and magnetically dominated, its energy-momentum tensor (B29) can be written as 28 With εijk the antisymmetric symbol, ṽi where u (γ) is any one of the unit-norm timelike vectors (B50) in its kernel, with associated magnetic field (B38): In particular, the norm of this magnetic field is Among all the flows (B50), one of them is preferred given our choice of coordinate system, because it has minimal Lorentz factor γ 0 = γ ⊥ relative to our normal observer, and it has no component aligned with the magnetic field (so ṽ∥ = 0).Since and the resulting representation of stress-energy tensor, is manifestly independent of the flow velocity parallel to the magnetic field, depending only on its perpendicular velocity (whereas u (γ) = γ(u ⊥ /γ ⊥ + ṽ∥ ) in general.)

B.8. Application to FFE
Force-free electrodynamics, or FFE, is a limit of ideal MHD in which the plasma becomes dilute and the local energy is dominated by the electromagnetic energy, In that case, the Bianchi identity for the Riemann tensor ∇ µ T µν = 0 ensures the conservation of electromagnetic energy: ∇ µ T EM µν ≈ 0. 29 This approximation is known as the force-free condition, since by (B26), it implies that the Lorentz force density vanishes.That is, f ν = 0, or This implies that the current is always perpendicular to the electric field e µ = F µν u ν , for any 4-velocity u: 30 In particular, in the frame (B41), J • X = J • E = 0, so 29 The force-free condition can hold even when (B69) does not. 30Thus, a force-free field obeys both ρ ⃗ To be more explicit, first note that we can write where the first term vanished by (B70) and the last by (B71).Hence, we can decompose the current (B72) into field-parallel and field-perpendicular currents which are explicitly given by and have norms Thus, J ∥ is always spacelike, while J ⊥ is timelike if and only if F is magnetically dominated.We also note that where we used (B39), so we can interpret J ⊥ as a charge density J t /u t flowing allowing u ⊥ .
In GRMHD, the field strength F is always guaranteed to be both degenerate and magnetically dominated by the ideal MHD condition, which requires F to have a timelike vector in its kernel: the flow 4-velocity u.
In FFE, there is no longer such a flow u.Nevertheless, the force-free condition (B70) still requires F to have a nontrivial kernel (containing J), which implies that F is degenerate.However, since J need not be timelike, 31 F is not automatically magnetically dominated.Instead, magnetic domination is a separate assumption, which in fact must be imposed in order to ensure well-posedness of the initial value problem, that is, to ensure that the evolution equations are hyperbolic (rather than elliptic, as in the electrically dominated case where the charges making up the plasma are accelerated away).For further discussion, see Komissarov (2002) or Gralla & Jacobson (2014) and other references therein. 31In fact, J is often spacelike, as in the split monopole solution: physically, the plasma consists of two oppositely charged species (electrons and positrons, say) with timelike velocities u + and u − producing a spacelike net current J ∝ u + − u − .

C. STATIONARY AND AXISYMMETRIC ELECTROMAGNETIC FIELDS ON KERR C.1. Kerr spacetime in Boyer-Lindquist coordinates
The Kerr metric in Boyer-Lindquist coordinates is The outer/inner event horizons r ± are the zeros of ∆, and the angular velocity of the (outer) event horizon is The metric determinant is The lapse α and shift vector β i defined in (B31) are and β r = β θ = 0, where we introduced The Kerr geometry is symmetric under time translations and rotations about the spin axis.These two isometries are respectively generated by the Killing vector fields which leave the metric invariant: letting L ξ denote a Lie derivative along ξ, K and R obey the Killing equation32 An object with 4-momentum p has energy −p t = −p • K and spin angular momentum p ϕ = p • R.

C.2. Kerr spacetime in Kerr-Schild coordinates
The Boyer-Lindquist coordinates x µ = (t, r, θ, ϕ) and the Kerr-Schild coordinates xµ = ( t, r, θ, φ) share the same poloidal coordinates (r, θ), whereas their toroidal coordinates are related by shifts which cancel the divergence of the Boyer-Lindquist line element (C78) as ∆ → 0. This is the main advantage of the Kerr-Schild coordinates: they remove the spurious coordinate singularities across the event horizons (C79).Using (C87), one can convert the line element (C78) to Kerr-Schild coordinates, or find the Jacobian for the coordinate transformation x µ → xµ .For instance, in Kerr-Schild coordinates, the Killing vectors (C84) are while the Boyer-Lindquist normal observer becomes which is no longer normal with respect to Kerr-Schild coordinates; instead, the Kerr-Schild normal observer is where the Kerr-Schild lapse ᾱ and shift vector βi are and βθ = β φ = 0.This normal observer η has zero angular momentum (η • R = η φ = 0) and is thus also a ZAMO.Moreover, this ZAMO also has zero angular velocity, since η φ/ ηt = − β φ = 0.However, this ZAMO also has negative radial velocity ηr /η t = − βr < 0, and is therefore infalling.Lastly, the Kerr-Schild "lab frame" is defined by the (non-normalized) vector ζ = − d t ̸ = ζ.

C.3. Stationary and axisymmetric fields on Kerr
A Kerr electromagnetic field configuration is a 2-form F that obeys the Maxwell equations (B24) in the metric background (C78).In particular, F is necessarily closed (dF = 0) and must therefore be (locally) exact, so that F = dA for some gauge potential 1-form A = A µ dx µ .
Such an electromagnetic field is stationary if L K F = 0 and axisymmetric if L R F = 0.A field that is stationary and axisymmetric always admits a 1-form potential A whose components A µ (r, θ) are independent of (t, ϕ).
. This is the integrability condition for the linear system L K d λ = dλ K , L R d λ = dλ R , which therefore admits a simultaneous solution λ.Thus, the gauge-transformed potential The electric and magnetic field components (B34) in the Boyer-Lindquist lab frame take the explicit form For future reference, we also record the inverse relations The field F = dA is not automatically degenerate, since (A6) does not generically vanish: In the Kerr-Schild coordinates xµ , a stationary and axisymmetric gauge potential In terms of the associated electromagnetic field strength F = dA = Fµν dx µ dx ν with Kerr-Schild components Fµν = ∇µ Āν − ∇ν Āµ , the local electric and magnetic fields (A15) in the Kerr-Schild lab frame of ζ = − d t are Using (C95), these are related to the fields (C92) in the Boyer-Lindquist lab frame by where we introduced a convenient quantity The inverse transformation is 34 C.4.Electromagnetic fluxes and (conserved) currents Given a Killing vector field ξ, its associated Noether current in the field configuration F is defined in terms of the associated electromagnetic stress-tensor (B25) as By construction, such a current obeys 35 where f is the Lorentz force density (B26).Therefore, the Noether current J ξ associated to the isometry along ξ is conserved if and only if f • ξ = 0.This is always the case in FFE, where the Lorentz force f ν = J µ F µν vanishes identically by the force-free condition (B70).It can also happen whenever dissipation vanishes along the symmetry direction (that is, when F µν J µ ξ ν = 0).The Kerr geometry (C78) admits two Killing vector fields (C84), with associated Noether currents 36 C102) corresponding to the flows of electromagnetic energy and angular momentum, respectively.Here and henceforth, we suppress the subscript 'EM' on T and need no longer track the ordering of its indices, as they are symmetric.
The electromagnetic energy and angular momentum of a field F on Kerr are conserved if and only if the currents (C102) are conserved, which by (C101) happens if and only if f t = f ϕ = 0. Otherwise, the field exerts a force (or a force is exerted on it) and it loses (or gains) energy (if f t ̸ = 0) and/or angular momentum (if Regardless of whether the currents J E,L are conserved, they always measure the flow of electromagnetic energy and angular momentum.For fields that are stationary and axisymmetric, these flows only exhibit nontrivial behavior in the poloidal directions r and θ.To see this, let P denote a curve in the poloidal plane (r, θ), and let P ×S 1 denote the 2-surface that is obtained by revolving P around the spin axis (that is, in the direction R = ∂ ϕ ). 34This agrees with (50)-( 52) of McKinney & Gammie (2004) given that their ω (defined in (26) therein) is our Ω in (C98). 35Since T µν EM is symmetric, the term T µν EM ∇µξν = 2T µν EM ∇ (µ ξ ν) missing in the first equation vanishes by the Killing equation. 32 36 The signs here are conventional and match the ADM formulas.
Finally, let S = P × S 1 × ∆t be the extension of this surface by a time interval ∆t along K = ∂ t .The energy and angular momentum fluxes through P × S 1 are then where the integral over the 3-surface S is performed over the 3-currents that are Hodge-dual to (C102), 37 By the same manipulations as in (A2), on S we have 38 The ϕ integral is trivial by axisymmetry, so the outward fluxes of energy and angular momentum through P are An explicit computation-using (C93) to eliminate the gauge potential in favor of the electric and magnetic field components-results in the flux densities We emphasize that we have not yet assumed that F is degenerate (we will do so only in Appendix D below).
37 By definition, 9 the Hodge dual of a p-form in η in n dimensions has components where by (A2), ϵ tϕrθ = −ϵ tϕθr = √ −g.Positive orientation has dϕ to the right.

C.5. Horizon fluxes
We will be particularly interested in the radial fluxes through spheres of constant radius r, which correspond to circular poloidal curves P. In such cases, it follows from (C106) that the (outward) radial fluxes are 39

Ėr
with the radial flux densities J r E , J θ E given in (C107).Naively, J r E and J θ E both appear to vanish identically at the event horizon (C79), where ∆ = 0.However, since the Boyer-Lindquist coordinates are singular across the horizon, the electric and magnetic field components (C92) can diverge at r = r + even if F is in fact regular there (so that the horizon fluxes are actually finite).
This issue is remedied by working with the Kerr-Schild lab-frame electric and magnetic fields Ēi and Bi defined in (C96), which are always regular across the horizon and are related to the singular Boyer-Lindquist fields (C92) by ( C97)-(C99).Using (C99), the horizon radial flux densities of energy and angular momentum take the manifestly regular Kerr-Schild form 40 where Ω H is the angular velocity (C80) of the horizon and Ω + is the quantity (C98) evaluated at r = r + .
D. STATIONARY, AXISYMMETRIC, AND DEGENERATE FIELDS ON KERR D.1.General form of the field strength A Kerr electromagnetic field F is necessarily closed.If F is also degenerate, then as discussed above (A14), there exist (infinitely many) pairs of Euler potentials (λ 1 , λ 2 ) such that F = dλ 1 ∧ dλ 2 (recall section A.4).Moreover, if F is also stationary and axisymmetric, then as shown in (61) of Gralla & Jacobson (2014), one may always choose a pair of Euler potentials of the form 41 Technically, the derivation of these potentials assumes that ∂ ϕ • F ̸ = 0, excluding all-toroidal magnetic fields. 4239 Our Ėr agrees with (32) of McKinney & Gammie (2004). 40This agrees with (34) of McKinney & Gammie (2004). 41The product rule If ∂ ϕ • F = 0, the Euler potentials can be taken to be of the form (D15) in Gralla & Jacobson (2014), namely λ 1 = χ(r, θ) and λ 2 = χ 2 (r, θ) + t.Since these can be obtained as a limit of (D110),43 it follows that the potentials (D110) are in fact still general.
Next, we define a function I(r, θ) as in (64) of Gralla & Jacobson (2014): Explicitly, this function takes the form Eliminating ψ 2 in favor of I, we may then express F as We conclude that an electromagnetic field F on Kerr that is stationary, axisymmetric, and degenerate must be of the general form (D113).In particular, such a field can be fully specified by three poloidal functions ψ(r, θ), I(r, θ), and Ω(ψ).These functions all admit intuitive physical interpretations, which we will describe in the next sections following §7 of Gralla & Jacobson (2014).
To do so, we first note that F = dψ ∧ dλ 2 = d (ψ dλ 2 ), so we may take the associated gauge potential to be Then a direct comparison of (D113) against F = dA with components A = A µ (r, θ) dx µ quickly shows that ψ(r, θ) = A ϕ (r, θ), (D115a) The MHD version of this decomposition, as well as the energy fluxes and boundary conditions discussed in the next sections, were analyzed in detail by Phinney (1983).

D.2. Poloidal field lines and magnetic flux function ψ
First, recall that the 2-surfaces of constant λ 1 and λ 2 define field sheets, and that these sheets are timelike provided that F = dλ 1 ∧ dλ 2 is magnetically dominated (section A.4).In that case, these field sheets may be viewed as the world sheets of magnetic field lines, which are defined with respect to a given time coordinate t as the intersection of a field sheet with a hypersurface of constant t (a spatial slice).
For fields that are also stationary and axisymmetric, we can further define poloidal field lines as the level sets of λ 1 = ψ(r, θ).These are projections in the poloidal plane (r, θ) of the magnetic field lines, whose bending in the azimuthal direction is controlled by ψ 2 (r, θ).Thus, ψ may be viewed as a coordinate on poloidal field lines.
At this stage, we note that the spin axis of the black hole must always be a poloidal field line-that is, ψ(r, 0) and ψ(r, π) must be constant-in order to ensure that the gauge potential (D114) stays regular on the spin axis, where dλ 2 diverges (because dϕ is singular at the poles).It is clear from (D115a) that by a suitable gauge transformation, we can arrange for ψ to vanish on axis, a choice that we will always assume from now on.Another physical interpretation of ψ is as a magnetic flux function.Consider a loop of revolution obtained by rotating a poloidal point (r, θ) in the azimuthal direction ϕ at a fixed time t.The magnetic flux through this loop is the integral of F = dA over any 2-surface S enclosed by this loop ∂S.By (D114) and Stokes' theorem, this is (Note that dλ 2 = dϕ on ∂S.)Hence, 2πψ(r, θ) measures the magnetic flux through the loop of revolution at (r, θ).
We note that mathematically, this argument relies on the regularity of A = ψ dλ 2 everywhere on S, including at the poles, which requires the imposition of (D116).Physically, this condition ensures that the magnetic flux vanishes as the loop ∂S shrinks to zero, whereupon it only encloses the spin axis.

D.3. Field-line angular velocity Ω(ψ)
Since the poloidal function Ω(ψ) is a function of ψ(r, θ) only, it may be viewed as a function on poloidal field lines.By (D110), it may be interpreted as the angular velocity of field lines as they rotate about the spin axis.
The fact that every field line of fixed ψ rotates at the same rate is encapsulated in (D115b), which is known as Ferraro's law of isorotation.The first equality in (D115b) is consistent with our earlier definition (C98), which is now extended by the second equality, which is required to ensure that (C94) vanishes and hence that F is degenerate.
The surfaces where the norm of dϕ − Ω(ψ) dt vanishes are known as light surfaces.They correspond to critical surfaces where the field lines rotate at the speed of light.D.4.Field-line (conserved) current I(ψ) The quantity I(r, θ) is known as the poloidal current for two reasons.First, it completely controls the poloidal components of the Maxwell current J µ = ∇ ν F µν : Second, if the Noether currents (C102) associated with flows of electromagnetic energy and angular momentum are conserved, then f t = f ϕ = 0.A direct computation reveals that this is equivalent to dψ ∧ dI = 0, that is, In other words, I(ψ) is the conserved current flowing along the poloidal field line labeled by ψ.

D.5. Electromagnetic fluxes
For the field (D113), the flux densities (C107) are just as can be obtained either via direct computation or, as a useful consistency check, by plugging in (D135) below.
As a result (recall from (C81) that √ −g = Σ sin θ), the energy and angular momentum fluxes (C106) become44 These fluxes vanish when integrated along a contour of fixed ψ, which makes it clear that the electromagnetic energy and angular momentum both flow along poloidal field lines.Moreover, when the Noether currents (C102) are conserved, so that (D119) holds, Ω and I are both functions of ψ only and so these fluxes become ordinary 1-dimensional integrals, which makes it clear that energy and angular momentum are conserved along field lines.
D.6.Znajek condition for regularity at the horizon Under the coordinate transformation to Kerr-Schild coordinates (C87), the general stationary, axisymmetric, and degenerate electromagnetic field (D113) becomes Since these coordinates are regular, the field F is regular at the (future) horizon (where ∆ → 0) if and only if Z vanishes at r = r + .Thus, the Kerr field (D113) is regular on the (future) horizon if and only if it obeys the Znajek condition (Znajek 1977) on the horizon,45 We emphasize that every quantity in this expression is to be evaluated at r = r + , where √ Π = r 2 + + a 2 = 2M r + .Using the Znajek condition to evaluate the fluxes (D121) on the horizon results in ( 126)-( 127) of Gralla & Jacobson (2014).Explicitly, we obtain the expressions which are consistent with (C108)-(C109).There, we had to use Kerr-Schild coordinates to obtain regular horizon flux densities.Here, we directly derived regular expressions by working with the variables (D115), but to do so we had to apply the Znajek condition (D123), whose derivation also involved Kerr-Schild coordinates.

D.7. Field-line angular velocity from Kerr symmetry
Besides the Killing vectors (C84), the Kerr metric (C78) also possesses an antisymmetric rank-2 tensor Y that obeys the Killing-Yano equation ∇ (λ Y µ)ν = 0: The existence of this Killing 2-form underlies many of the special properties of the Kerr geometry, including: the separability of its wave equation, the integrability of its geodesics, and the existence of the Penrose-Walker constant κ, which reduces to simple algebra the problem of parallel transport of polarization along null rays.
Here, we show how to use this higher-rank symmetry of Kerr to project out the field-line angular velocity Ω(ψ) from the full 2-form (D113).To do so, we first define Raising indices, we obtain the explicit vector fields Then a direct calculation reveals that ) which shows that Ω(ψ) can be extracted from F using only the Killing symmetries of Kerr.This also confirms that Ω(ψ) is an invariant scalar.
It is also possible to extract Ω(ψ) from ⋆F using only the Killing symmetries.Define four projections Then a direct computation shows that and likewise, that To the best of our knowledge, these expressions for Ω(ψ) as scalar symmetry projections are new in the literature.
For completeness, we also note that the current I can be similarly projected out, for instance as (D134) D.8.Boyer-Lindquist electric and magnetic fields Using (C81) and (D115), we can re-express the electric and magnetic fields in the Boyer-Lindquist lab frame (C92) in terms of (ψ, Ω, I) as with ω the angular velocity (C86) of the Boyer-Lindquist ZAMO.We can also write the E i in terms of the B i as If the field (D113) is magnetically dominated, then we showed that its kernel contains a family u (γ) of timelike 4-velocities (B40), which are conveniently parameterized by their Lorentz factor γ = −η • u (γ) relative to the normal observer.The one with minimal Lorentz boost (B59) is perpendicular to the lab-frame magnetic field, ) while the other ones all have a nonvanishing field-parallel component.Here, we chose to express (B40) in terms of the lab-frame fields (B34) rather than the normal fields (B33). 26Next, recall the decomposition (B51): (D138)

E. STATIONARY AND AXISYMMETRIC FORCE-FREE ELECTRODYNAMICS IN KERR
As we showed in section B.8, a force-free field is always degenerate.Hence, stationary, axisymmetric, force-free fields on Kerr must necessarily be of the form (D113).
The force-free condition (B70) requires the Lorentz force density f ν to vanish.This imposes two constraints f t = f ϕ = 0 that are satisfied if and only if the currents (C102) are conserved (or equivalently, if (D119) holds).Thus, a force-free field on Kerr is entirely specified by three poloidal functions: the magnetic flux ψ(r, θ), and the field-line angular velocity Ω(ψ) and current I(ψ).The two other constraints f r = f θ = 0 are equivalent to each other and to the (Grad-Shafranov) stream equation, which we write as in (3.4) of Camilloni et al. (2022):

E.1. Exact analytic force-free fields in Kerr
In the Kerr spacetime, only two exact stationary and axisymmetric force-free fields are known analytically.The first, due to Menon & Dermer (2007), allows for the magnetic flux to be any invertible function f (θ), with inverse θ(ψ) = f (−1) (ψ).As such, the field lines are purely radial, and so it is possible to satisfy (D116).
The angular velocity of field lines is always the same regardless of the specific choice of f , while the field-line current is (up to a sign and constant) In order to avoid an unphysical line current along the spin axis, we must set I 0 = 0 and f ′ (θ) = 0 at the poles.
In that case, the field becomes null, The flux densities (C102) are Besides not being magnetically dominated, this solution also fails to describe energy extraction from the black hole, since J r E < 0 for any choice of f .There exists a non-stationary, nonaxisymmetric generalization of this solution with the same properties (Brennan et al. 2013).
The second exact solution is due to Menon (2015).It is "dual" to the first, in the sense that it has the opposite (perpendicular) field geometry, with field lines that are circular rather than radial.The magnetic flux of this dual solution can be any invertible function g(r), with inverse r(ψ) = g (−1) (ψ).The angular velocity of field lines is always the same regardless of g, ) while the field-line current is (up to a sign and constant) This solution is not physically admissible as it cannot satisfy (D116), so we will not investigate it any further.Since neither of these exact solutions is realistic, we must turn to perturbative force-free solutions to study electromagnetic energy extraction from Kerr black holes.

E.2. Perturbing the force-free equations in small spin
We follow Blandford & Znajek (1977) and solve the force-free equations on Kerr perturbatively in slow spin.That is, we expand the stream equation (E149) in small |a| ≪ M , as reviewed in §4 of Gralla et al. (2016).
At leading order, a = 0 and so we are in Schwarzschild.In that case, we ought to have no rotation and therefore no current.Hence, to order O a 0 , we expect Ω = I = 0, and as such, ψ(r, θ) ought to correspond to a vacuum Maxwell field in a Schwarzschild background.
This motivates the perturbative expansion ) valid to subleading order O(a) in spin.Here, ψ 0 , Ω 1 , and I 1 are all independent of the spin a, while the scaling of the error terms is fixed by the Znajek condition (D123).The leading-order magnetic flux function ψ 0 must solve (E149) with a = Ω = I = 0, that is, the linear equation Its solutions are vacuum Maxwell fields in Schwarzschild.They are reviewed in Appendix B of Gralla et al. (2016).At subleading order O(a), the stream equation ( E149) is identically satisfied by the perturbative ansatz (E159): the only requirement imposed by FFE is that Ω 1 and I 1 be functions of ψ 0 alone.The nontrivial content of the full Kerr stream equation (E149) only kicks in at the higher order O a 2 in perturbation theory, which we do not review here as the details become quite complicated; see, e.g., Armas et al. (2020) for an in-depth treatment.
Thus, to solve the force-free equations to linear order in a, one needs a vacuum Schwarzschild flux ψ 0 together with two additional relations to determine Ω 1 and I 1 . 46 One such relation is universal and follows from the Znajek condition (D123), which must always hold at the horizon r + = 2M + O(a): to leading-order in a, we have As for a second relation, Blandford & Znajek (1977) had the key insight to demand that the ansatz (E159) match onto an exact force-free solution in Minkowski spacetime at radii r ≫ M , where the Kerr geometry becomes flat.
46 If one only wishes to solve the force-free equations to linear order in a, then in principle one can freely choose to impose any second relation between ψ 0 , Ω 1 , and I 1 .However, almost all the resulting fields will be spurious linearized solutions that do not correspond to any family of exact solutions, and hence cannot be extended to solutions at higher order in the spin perturbation.
This reproduces the flat-space parabolic field as r → ∞, is regular on the Schwarzschild horizon r + = 2M , and moreover respects the boundary condition (D116).This matching procedure may seem rather haphazard, but it can be approached more systematically.Such an approach is indispensable for the hyperbolic field, whose extension to Schwarzschild is much more complicated.
The key idea is that all the flat-space force-free fields have magnetic flux functions ψ 0 that are also solutions to the vacuum Maxwell equation in flat spacetime (since one can choose Ω = I = 0), which is (E160) with M = 0.
Since (E160) is separable, one can canonically identify its mode solutions when M = 0 (the flat-space modes) with its mode solutions when M ̸ = 0 (in Schwarzschild).
These sets of modes and their canonical identification are described in Appendix B of Gralla et al. (2016).As an application, the vacuum magnetic flux function for the flat-space hyperbolic field, corresponding to (E168) with Ω = I = 0, is explicitly given in (B.18) therein, and its canonical extension to a vacuum Maxwell field in Schwarzschild is explicitly given in (B.22)-(B.25).

E.5. Perturbative force-free monopole field on Kerr
We now have all the ingredients needed to solve for the perturbative force-free solution (E159) with a radial field-line geometry: the Blandford-Znajek monopole.
The leading-order magnetic flux function ψ 0 (r, θ) is given in (E171), so only the field-line angular velocity Ω 1 (ψ 0 ) and current I(ψ 0 ) remain to be determined.The perturbative Znajek condition (E161) yields one relation between ψ 0 , Ω 1 and I 1 , which is technically imposed at the horizon but must in fact hold everywhere, as both Ω 1 and I 1 are only functions of ψ 0 .Matching onto the flat-space solution (E163) yields a second relation between ψ 0 , Ω 1 and I 1 , which is technically imposed at r → ∞ but must again hold everywhere, as both Ω 1 and I 1 are only functions of ψ 0 .Equating these relations lets us solve for ) in terms of which I 1 (ψ 0 ) is then given by (E175).Therefore, to leading order in perturbation theory, the force-free monopole solution is approximately ψ(r, θ) = 1 − cos θ, (E177a) Although we will not reproduce the details, it is possible to continue this perturbative expansion to higher orders (Armas et al. 2020).For instance, at the next order, where R(x) is a complicated function given in (4.49) of Armas et al. (2020), and such that 72 R(2) = 6π 2 − 49.
E.6.Perturbative force-free parabolic field on Kerr Next, we solve for the perturbative force-free solution (E159) on Kerr with a parabolic field-line geometry.
The leading-order magnetic flux function ψ 0 (r, θ) is given in (E173), so only the field-line angular velocity Ω 1 (ψ 0 ) and current I(ψ 0 ) remain to be determined.The perturbative Znajek condition (E161) yields the first needed relation between ψ 0 , Ω 1 and I 1 (which, again, is only imposed at the horizon, but must hold everywhere as both Ω 1 and I 1 are only functions of ψ 0 ): Here, W k (z) is the k th branch of Lambert's W-function, defined such that w = W k (z) solves z = we w .49For real x and y, the equation x = ye y only admits solutions if x ≥ −e −1 , in which case the two branches W 0 and W −1 suffice to represent them: if x ≥ 0, then there is a unique solution y = W 0 (x), whereas if −e −1 ≤ x < 0, then there is a second solution y = W −1 (x).
To obtain a second relation between ψ 0 , Ω 1 and I 1 , we match onto the flat-space solution (E164) with arbitrary field-line velocity Ω(ψ) and field-line current (E164c): where in the last step, we expanded to leading order in Ω ∼ O(a).The resulting truncation is consistent with the perturbative ansatz (E159), which assumes that Ω and I are both linear in the small spin 0 < |a| ≪ M .Since to leading order in Ω, the magnetic flux (E164b) is ψ 0 (u) ≈ u, we thus obtain a second relation This condition is only imposed at r → ∞, but must hold everywhere as both Ω 1 and I 1 are only functions of ψ 0 .Equating this expression to (E179a) lets us solve for50 Therefore, to leading order in perturbation theory, the force-free parabolic solution is approximately ψ(r, θ) = r (1 − cos θ) − 4M (1 − log 2) (E183a) + 2M (1 + cos θ) [1 − log (1 + cos θ)] , We note that our derivation differs from the treatment by Blandford & Znajek (1977): while their (7.4)agrees with our (E181), their (7.2) is only the restriction to the horizon of our more general (E179a).For this reason, we are able to find the field-line angular velocity (E182) everywhere, while they only obtain its restriction to the horizon, Ω + (θ) = Ω(ψ(r + , θ)), which they give in (7.5):For the perturbative force-free solution (E159) on Kerr with a hyperbolic field-line geometry, the above steps are repeated in detail in §4 of Gralla et al. (2016).Since ũϕ = Ωũ t , this completes the definition of ũ.
For most values of the parameters (ξ, β r , β ϕ ), the flow ũ is not geodesic.However, if β r = β ϕ = 0, then ũ reduces to the radial inflow ū, which is ξ-independent and geodesic, corresponding to an equatorial inflow with zero angular momentum. 52Likewise, if β r = β ϕ = 1, then ũ reduces to the circular flow ũ(ξ), which is sub-Keplerian and hence non-geodesic for ξ < 1, but reduces to geodesic circular Keplerian motion when ξ = 1.This appendix describes in detail the construction of the semi-analytic models used to generate the simulated images of M87* in section 3.As discussed therein, in all models, we take the lab-frame magnetic field B i = ⋆F i0 to be that of the Blandford-Znajek monopole solution (E178), and then we specify some fluid 4-velocity u µ to complete the definition of ⋆F µν .Since by degeneracy, the electric field is assumed to vanish in the fluid frame, ) where in this formula, the "fluid-frame" magnetic field is As mentioned in section 3, we use two different models for u µ .First, we adopt the drift-frame velocity (D137) that is determined by the full Blandford-Znajek solution for both B i = ⋆F i0 and E i = F 0i .Given a fixed ⋆F µν , (D137) provides the unique frame u µ ⊥ in which the local electric field vanishes (e µ = 0) and u • B = 0.As we will show in Paper II, replacing u µ ⊥ by a general member of the family (B40)-which differs from u µ ⊥ by a boost along the magnetic field-would not affect the direction of linear polarization for synchrotron radiation.
In order to survey different models of fluid inflow, we also explore the family of parameterized accretion flow models for u µ introduced by Cárdenas-Avendaño et al. (2023) and reproduced in Appendix F. That is, we fix the lab-frame magnetic field components B i = ⋆F i0 to their values in the Blandford-Znajek monopole solution (E178), and then we freely prescribe some choice of u µ drawn from Appendix F. Generically, such a choice is inconsistent with the allowed family of flows (B40) and therefore modifies the spatial components ⋆F ij (and hence the electric field E i = F 0i ) according to (G194).
Although commonly used in semi-analytic modeling of near-horizon flows, this approach suffers from a major drawback: it typically results in a field F µν that does not satisfy Maxwell's equations (B24).In particular, it is the Bianchi identity dF = 0, or equivalently, the "induction equation" ∇ µ (⋆F ) µν = 0, that generically fails to hold. 52This explains why Ω equals the ZAMO angular velocity (C86).
More precisely, while the induction equation's time component (the "no-monopoles" constraint ∇ i B i = 0) is always satisfied, on the other hand, its spatial part will typically not be.Thus, if one wishes to regard the resulting tensor F µν as a valid electromagnetic field, it can no longer be viewed as a stationary configuration, but rather as the initial state of a time-evolving field.Finally, having specified B i and u µ everywhere in the equatorial plane of the black hole, we need only define the emissivity J(r eq ) for emission received at 230 GHz.We consider a thin ring centered at r ring = 2r + and use the same emissivity as in (12) of Gralla et al. (2020): with a ring width of σ = 0.3.

G.2. Ray tracing and parallel transport
For a given spin a, viewing angle θ o , magnetic field B i , and 4-velocity u µ , we generate polarized synchrotron images of an equatorial emission ring using the public code kgeo (Chael 2023). 53We parameterize the image plane with coordinates (α, β) defined by Bardeen et al. (1972) such that the β axis is aligned with the black hole spin projected onto the plane perpendicular to the line of sight.kgeo analytically computes null geodesics by tracing backwards from each point on the observer screen using the formalism of Gralla & Lupsasca (2020).To determine the EVPA of the received emission, kgeo also solves for the parallel transport along each ray of the polarization vector f µ , which is initially taken to point in a direction locally perpendicular to both the magnetic field and photon momentum at the source.
In particular, for axisymmetric, equatorial emission models, kgeo solves for the radius r eq (α, β; n) where the null geodesic that reaches the observer at coordinate (α, β) makes its (n + 1) th pass through the equatorial plane θ = π/2.This radius, together with the geodesic's conserved specific angular momentum λ = −α sin θ o and Carter constant η = (α 2 − a 2 ) cos 2 θ o + β 2 , determines everywhere along the trajectory the local momentum k µ of a photon with energy E loaded onto the ray: )

Figure 1 .
Figure 1.Schematic diagram illustrating the connection between different simple equatorial axisymmetric magnetic field configurations around a black hole and the polarimetric observable ∠β2.No relativistic or propagation effects (e.g., parallel transport, aberration, Faraday rotation) are included.The +ẑ axis is into the page, and the local polarization direction is given by ẑ × ⃗ B. The range of values obtained by the EHT for M87* (−163 deg < ∠β2 < −129 deg; EHTC VIII) is consistent with a ratio B ϕ /B r < 0 in the emission region.Assuming that the black hole spin is aligned with ẑ, that the magnetic field is axisymmetric and degenerate, and that the magnetic field lines co-rotate with the black hole, this field structure produces an outward Poynting flux by Equation8.
Figure 2. (Left) Histograms of ∠β2 covering different velocity prescriptions in axisymmetric, equatorial models of emission around M87*.In these models, the emission radius is fixed to r = 2r+ and the lab-frame magnetic field B i given by the BZ monopole solution.For each value of black hole spin, the histograms cover uniform priors in the velocity parameters ξ, βr, β ϕ in the range [0, 1] for both prograde (a * > 0) and retrograde (a * < 0) orbits.The observer inclination θo is 163 deg relative to the spin axis for a * > 0 and 17 deg for a * < 0. The stars on the x axis indicate the values of ∠β2 for the same fields B i assuming the self-consistent BZ drift-frame velocity.The range of ∠β2 values consistent with EHT observations of M87* in 2017 (EHTC VII) is indicated with the blue band.(Right) Representative images of the equatorial BZ monopole model for |a * | = 0.5 with a prograde Keplerian velocity model, a retrograde Keplerian velocity model, and the outflowing drift-frame velocity from the complete BZ monopole solution.
) performed using the code KORAL(Sądowski et al. 2013(Sądowski et al. , 2014)).These simulations consider both prograde and retrograde accretion disks at five values of the dimensionless black hole spin a * ∈ [0, ±0.3, ±0.5, ±0.7, ±0.9].The simulations where |a * | > 0 have substantial outward Poynting flux driven by BZ energy extraction (see Figure 4 of Narayan et al. 2022).The a * = 0 simulation also has an outward Poynting flux, which is launched by the rotational energy of the accretion flow rather than the black hole spin (Blandford & Payne 1982).The total energy in the outflow from the a * = 0 simulation is much smaller in than in the simulations with |a * | > 0 (Narayan et al. 2022).

Figure 3 .
Figure 3. (Left column) Histograms of snapshot values of β2 in the complex plane for images from nine MAD GRMHD simulations raytraced with M87* parameters.Prograde models (a * > 0) are shown in the top row, and retrograde models (a * < 0) in the bottom row.The histograms include all snapshot images at different simulation time slices between 50, 000 tg and 100, 000 tg. Simulation images were blurred to the EHT resolution of 20 µas before computing β2.The histograms cover six different models for the low-magnetization ion-to-electron temperature ratio R high .(Right column) Histograms of ∠β2 from the same simulation images.With a few exceptions for predominantly low-spin, retrograde simulation snapshots, nearly all MAD GRMHD snapshots have ∠β2 < 0, corresponding to energy outflow in the simple picture of Figure 1.The range of ∠β2 values consistent with EHT observations of M87* in 2017 is shown in the blue band (EHTC VII).

Figure 4 .
Figure 4. (Top) The sign of B ϕ /B r for three MAD GRMHD simulations averaged in time and azimuth.From left to right, the GRMHD simulations correspond to a retrograde accretion disk around a spinning black hole (a * = −0.5), a Schwarzschild disk (a * = 0), and a prograde disk (a * = +0.5).(Middle) The angular frequency of the magnetic field lines ΩF in the same simulations in natural units of t −1 g .(Bottom) The averaged radial Poynting flux √ −gJ r E in BL coordinates in the same simulations.In the bottom panel, all three panels are independently normalized and are not plotted in the same units; the energy flux in the a * = 0 simulation (driven by the accretion disk) is an order of magnitude smaller than in the a * = ±0.5 simulations (driven by black hole spin).In all plots, the cyan contour indicates the surface where the magnetization σ = 1.White contours show surfaces of constant vector potential ψ ≡ A ϕ .

Figure 5 .
Figure 5. (Left column) ∠β2 as a function of image radius in time-averaged images of M87* from MAD GRMHD simulations.The top row shows the time-averaged image for a * = −0.5, the middle row for a * = 0, and the bottom row for a * = 0.5.All images use R high = 1 to set the electron temperature in radiative transfer.The vertical magenta line shows the outermost radius of the "inner shadow," or the lensed image of the equatorial event horizon; the cyan band shows a range of radii corresponding to the n = 1 photon ring.(Right column) The time-averaged images used in computing ∠β2 in the left panels.Images are displayed in a gamma color scale and are not blurred to EHT resolution.The black hole shadow/critical curve is indicated by the cyan contour, and the inner shadow with the magenta contour.The image-average value of ∠β2 is displayed on each panel.A version of the top row of this plot was first produced in Ricarte et al. (2022), Figure 8.

F
briefly summarize the model for accretion flow introduced in Cárdenas-Avendaño et al. (2023)full derivations are provided in Appendix B therein.We consider an equatorial 4-velocity of the form ũ = ũt ∂ t + ũr ∂ r + ũϕ ∂ ϕ .
G. SEMI-ANALYTIC MODEL DESCRIPTION G.1.Magnetic field, fluid velocity, and emissivity