Issue 
A&A
Volume 663, July 2022



Article Number  A169  
Number of page(s)  26  
Section  Astrophysical processes  
DOI  https://doi.org/10.1051/00046361/202142847  
Published online  27 July 2022 
Spinning black holes magnetically connected to a Keplerian disk
Magnetosphere, reconnection sheet, particle acceleration, and coronal heating
^{1}
Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France
email: ileyk.elmellah@univgrenoblealpes.fr
^{2}
School of Mathematics, Trinity College Dublin, Dublin 2, Ireland
Received:
6
December
2021
Accepted:
30
March
2022
Context. Accreting black holes (BHs) may be surrounded by a highly magnetized plasma threaded by an organized poloidal magnetic field. Nonthermal flares and powerlaw spectral components at high energy could originate from a hot, collisionless, and nearly forcefree corona. The jets we often observe from these systems are believed to be rotationpowered and magnetically driven.
Aims. We study axisymmetric BH magnetospheres, where a fraction of the magnetic field lines anchored in a surrounding disk are connected to the event horizon of a rotating BH. For different BH spins, we identify the conditions and sites of magnetic reconnection within 30 gravitational radii.
Methods. With the fully general relativistic particleincell code GRZeltron, we solve the timedependent dynamics of the electron–positron pair plasma and of the electromagnetic fields around the BH. The aligned disk is represented by a steady and perfectly conducting plasma in Keplerian rotation, threaded by a dipolar magnetic field.
Results. For prograde disks around Kerr BHs, the topology of the magnetosphere is hybrid. Twisted open magnetic field lines crossing the horizon power a BlandfordZnajek jet, while open field lines with their footpoint beyond a critical distance on the disk could launch a magnetocentrifugal wind. In the innermost regions, coupling magnetic field lines ensure the transfer of significant amounts of angular momentum and energy between the BH and the disk. From the Y point at the intersection of these three regions, a current sheet forms where vivid particle acceleration via magnetic reconnection takes place. We compute the synchrotron images of the current sheet emission.
Conclusions. Our estimates for jet power and BH–disk exchanges match those derived from purely forcefree models. Particles are accelerated at the Y point, which acts as a heat source for the socalled corona. It provides a physically motivated ringshaped source of hard Xrays above the disk for reflection models. Episodic plasmoid ejection might explain millisecond flares observed in Cygnus X1 in the highsoft state, but are too fast to account for daily nonthermal flares from Sgr A^{*}. Particles flowing from the Y point down to the disk could produce a hot spot at the footpoint of the outermost closed magnetic field line.
Key words: acceleration of particles / magnetic reconnection / black hole physics / radiation mechanisms: nonthermal / methods: numerical / relativistic processes
© I. El Mellah et al. 2022
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1. Introduction
Although their masses and luminosities span many orders of magnitude, stellarmass black holes (BHs) and supermassive black holes (SMBHs) share many similarities, which emphasizes how simple the intrinsic structure of a BH is. Due to their high compactness, they can accrete ambient material and become intense sources of light visible up to cosmological distances. The intensity, the spectrum, and the polarimetry of the radiation emitted by the plasma in the immediate surroundings of accreting BHs allow us to determine two main properties of these universal beacons: their mass and their spin.
Multiple lines of evidence indicate that BHs spin. Nonzero spins were measured in Xray binaries that host a BH surrounded by an accretion disk, using Xray reflection spectroscopy and thermal continuum fitting (Reynolds 2021; Bambi et al. 2021). Although these two diagnostics do not always yield consistent values and rely on many assumptions and degreesoffreedom, they suggest that stellarmass BHs can be fast rotators (e.g., GRS 1915+105; Shreeram & Ingram 2020, or Cyg X1, Zhao et al. 2021). Xray reflection spectroscopy applied to active galactic nuclei shows that SMBHs also have intrinsic angular momentum. It is reasonable to hope for a direct measure of the spin of the SMBH at the center of the Milky Way, Sgr A^{*}, via LenseThirring precession of the orbit of a nearby star thanks to the efforts of the GRAVITY Collaboration (GRAVITY Collaboration 2018a). The ubiquity of jets from accreting stellarmass BHs and SMBHs has also been interpreted as indirect evidence of spinning BHs.
Accretion onto Xray active stellarmass BHs is thought to be partly mediated by an optically thick and geometrically thin accretion disk similar to the one first described by Shakura & Sunyaev (1973). In semidetached Xray binaries, where the donor star fills its Roche lobe, an accretion disk forms owing to the high angular momentum of the stellar material provided at the inner Lagrangian point (see, e.g., Frank et al. 1986). In detached highmass Xray binaries (HMXBs), where the donor star is massive, the formation of a disk is still possible provided the UVdriven stellar wind fueling the BH is slow and/or the orbital separation is small (El Mellah et al. 2019;Sen et al. 2021). The properties of the accretion flow around SMBHs are less constrained. The farinfrared and radio images of the immediate vicinity of the SMBH M 87^{*} reported by Event Horizon Telescope Collaboration (2019) are consistent with an underlying optically thin accretion disk where the magnetic field plays a major dynamic role: the magnetically arrested disk (MAD; Event Horizon Telescope Collaboration 2021a). This type of disk has been shown to be a robust bedrock for jet launching à la BlandfordZnajek (Tchekhovskoy et al. 2011). The presence of an accretion disk around M 87^{*} is thus corroborated by the prominent jet emitted by the SMBH (Algaba et al. 2021). For Sgr A^{*}, numerical simulations by Ressler et al. (2020) suggest that a MAD can form around the SMBH via the capture of the stellar wind emitted by WolfRayet stars orbiting within a few parsecs. Finally, the hot spot observed in orbit around Sgr A^{*} by the GRAVITY Collaboration was interpreted as a sign of an underlying disk (GRAVITY Collaboration 2018b).
Accreting BHs are thought to also be surrounded by a diffuse environment often called the corona. In the HMXB Cyg X1, where a late Otype supergiant supplies material to an approximately 21 M_{⊙} BH (MillerJones et al. 2021), a highenergy polarized component between 0.4 and 2MeV was found (Laurent et al. 2011; Rodriguez et al. 2015). It is present in the lowhard state but also in the highsoft state (see Grinberg et al. 2013, for the spectral state classification in Cyg X1). In the highsoft state, radio emission from the jet is much dimmer and was long thought to be absent (Remillard & McClintock 2006; Rushton et al. 2012; Zdziarski et al. 2020), which suggests that jet launching is quenched and is unlikely to be the source of this highenergy emission. Cangemi et al. (2021) showed that this highenergy component was instead compatible with a hybrid thermal–nonthermal Comptonization spectrum associated with a population of relativistic electrons. They also found that the magnetic field these electrons were embedded in could be strong enough to accelerate particles up to relativistic Lorentz factors via magnetic reconnection in the corona (Malzac & Belmont 2009; Poutanen & Vurm 2009). These results were in agreement with those of Del Santo et al. (2013), who derived upper limits for the magnetic field magnitude in the corona between 10^{5} and 10^{7} G. Similarly, the daily nonthermal flares from Sgr A^{*} in the nearinfrared (NIR) have been interpreted as synchrotron emission from nonthermal relativistic electrons within a few dozen gravitational radii (Witzel et al. 2018), consistent with the presence of a hot corona around the SMBH. Witzel et al. (2021) carried out an extensive multiwavelength analysis of variability in Sgr A^{*}, from submillimeter to Xrays, and ascribed the emission to highenergy electrons in a magnetic field of ∼8−13 G, in agreement with other diagnostics (Mościbrodzka et al. 2009; Dexter et al. 2010; Ponti et al. 2017). Noticeably, Ponti et al. (2017) showed that during a flare, the time evolution of spectral properties are coherent, with a drop in the magnetic field, which suggests that magnetic reconnection could be at play. The recent polarization measurements of M 87^{*} performed by Event Horizon Telescope Collaboration (2021b) point to a largescale poloidal magnetic field, with an intensity similar to that of Sgr A^{*}, that is threading a dilute plasma (Event Horizon Telescope Collaboration 2021a).
The aforementioned Xray reflection spectroscopy diagnostics assume a predefined geometry for the corona that plays the role of a source of hard Xrays reprocessed by the underlying disk. In the common lamppost model, the corona is a point source located on the BH spin axis a few gravitational radii above the BH (see, e.g., Kara et al. 2019). However, both around accreting stellarmass BHs (Zdziarski et al. 2021) and SMBHs (Zoghbi et al. 2019), evidence has emerged in favor of an extended corona (Chainakun et al. 2019). The actual geometry of the corona and the energy sources responsible for its formation, however, remain elusive, and so do the sites where the coronal heating takes place.
Configurations where magnetic field lines are all open and thread either the disk or the event horizon have been the object of much attention. In the equatorial plane within the innermost stable circular orbit (ISCO), a highly dynamic current sheet has been shown to form and be the stage of vivid magnetic reconnection (Parfrey et al. 2019; Crinquand et al. 2021; Bransgrove et al. 2021; Ripperda et al. 2022). Alternatively, the magnetosphere could contain closed magnetic field lines that thread the event horizon but with a footpoint on the surrounding disk. This case has been extensively studied by Uzdensky (2004, 2005) and more recently by Yuan et al. (2019a,b) and Mahlmann et al. (2018), who worked in the purely forcefree regime. They found steady solutions where disk–BH coupling field lines in the innermost region coexist with open field lines anchored farther in the disk. The presence of magnetic loops formed in the disk and extending into the corona was first speculated by Galeev et al. (1979), and a stochastic model of their dynamics was designed by Uzdensky & Goodman (2008), but their formation remains elusive. Recent general relativisticmagnetohydrodynamics (GRMHD) simulations by Chashkina et al. (2021) highlighted how these loops could be advected inward in the disk. On the other hand, Parfrey et al. (2015), Yuan et al. (2019b), and Mahlmann et al. (2020) used timevariable forcefree simulations to show that efficient jet launching and dissipation in the corona could be sustained by the advection of successive magnetic loops of varying polarity.
In this paper we study the global axisymmetric structure of a BH magnetosphere where the main physical ingredients are accounted for: a spinning BH and a magnetic field threading both a highly magnetized hot pair–plasma corona and, in the BH equatorial plane, a disk. In the hybrid magnetospheres we focus on, the topology of the magnetic field is not limited to field lines that thread the BH event horizon and extend to infinity. Initially, magnetic field lines are all anchored in the disk and cross the horizon, forming a single magnetic field loop that allows for exchanges of energy and angular momentum between the two components. We do not address the question of the origin of this configuration nor of the accretion of the loop itself; instead, we study the relaxation and time variability of this setup. We solve the dynamics of the plasma and of the electromagnetic fields in the corona specifically using an ab initio particleincell (PIC) approach, suitable to capture dissipative effects in highly magnetized plasmas. The disk is modeled as a steady and perfectly conducting plasma in Keplerian rotation. In light of the major results obtained by previous authors, we believe this hybrid configuration can be used as a fruitful framework before studying more realistic magnetic field distributions on the disk. Our aim is to characterize the topology of the magnetic field and the possible reconnection sites in the corona, depending on the BH spin, on the disk orientation (prograde or retrograde), and on the disk thickness.
In Sect. 2 we present the BH magnetosphere model we rely on and detail the numerical setup we designed to address it with the GRZeltron code. We present the results we obtained in Sect. 3 and discuss the implications for accreting BHs, with an emphasis on Sgr A^{*}, M 87^{*}, and Cyg X1. Finally, we summarize our main conclusions in Sect. 5 and suggest future tracks to be explored.
2. Model
2.1. Kerr black hole magnetosphere
2.1.1. Kerr spacetime metric
We consider a stationary and axisymmetric Kerr spacetime induced by a rotating BH of mass M and dimensionless spin a (Kerr 1963). We used the 3 + 1 formalism in order to have a universal time coordinate t and spatial hypersurfaces of constant t, mapped by the coordinates x^{i} (MacDonald & Thorne 1982). This introduces a class of observers whose world lines are orthogonal to these spatial hypersurfaces, the FIDucial Observers (FIDOs). If we write h_{ij}, the spatial 3metric on the hypersurfaces of constant t, the Kerr space time line element, ds^{2}, can be written in the ADM form (Arnowitt et al. 1962):
$$\begin{array}{c}\hfill \mathrm{d}{s}^{2}={\alpha}^{2}{c}^{2}\mathrm{d}{t}^{2}+{h}_{\mathit{ij}}(\mathrm{d}{x}^{i}+{\mathit{\beta}}^{i}c\mathrm{d}t)(\mathrm{d}{x}^{j}+{\mathit{\beta}}^{j}c\mathrm{d}t),\end{array}$$(1)
where α is the lapse function, β is the shift vector, and c is the speed of light. α encapsulates the information on time dilation from coordinate time to the proper time measured locally by the FIDO, FIDO β is the 3velocity of the FIDO with respect to the coordinate grid. The properties of interest of this metric for our study are the presence of an event horizon of radius ${r}_{\mathrm{H}}={r}_{\mathrm{g}}(1+\sqrt{1{a}^{2}})$ where r_{g} = GM/c^{2} is the gravitational radius and G is the gravitational constant, an ergosphere and an ISCO (Shapiro & Teukolsky 1983). At the event horizon, stationary observers are forced into rotation at an angular speed ω_{H} = ac/2r_{H} that we identify to the angular speed of a rotating BH. Hereafter, we work in KerrSchild spherical coordinates (r, θ, ϕ) in order to avoid the coordinate singularity at the event horizon, which would arise in BoyerLindquist coordinates for instance.
2.1.2. Keplerian disk
In our model, we rely on a simplified representation of a rotating and perfectly conducting disk entirely parametrized by the magnetic flux function, Ψ, in the disk and the disk aspect ratio, ϵ. The latter is assumed to be uniform and is defined as the ratio of the scale height to the distance R = rsinθ to the BH projected on the disk plane. We do not model accretion in our model since there is no radial speed in the disk. The disk rotation axis is assumed to be aligned with the BH spin axis, and the disk is thus assumed to lie in the BH equatorial plane. We accounted for the disk aspect ratio in the angular velocity profile by relying on the following approximate correction (see, e.g., Fromang et al. 2011):
$$\begin{array}{c}\hfill \mathrm{\Omega}(R)={\mathrm{\Omega}}_{\mathrm{K}}\sqrt{1{\u03f5}^{2}},\end{array}$$(2)
where
$$\begin{array}{c}\hfill {\mathrm{\Omega}}_{\mathrm{K}}(R)=\sqrt{\frac{\mathit{GM}}{{r}_{\mathrm{g}}^{3}}}\frac{1}{{(R/{r}_{\mathrm{g}})}^{3/2}+a}\end{array}$$(3)
is the relativistic counterpart of the Keplerian angular speed profile in Newtonian mechanics. In this paper, we assume a dipolar magnetic field in the disk, centered on the BH. The influence of the disk magnetic flux distribution on the magnetic structure of forcefree BH magnetospheres was found to be modest by Yuan et al. (2019a) who considered, in addition to this dipolar case, asymptotically paraboloid and uniform vertical fields. In a followup paper, when exploring a more realistic setup where two successive poloidal loops of opposite polarity are accreted, (Yuan et al. 2019b) found similar results when the inner loop was sufficiently weak with respect to the outer one. The magnetic field lines are frozen into the disk and rotate at the angular speed of their footpoint in the disk midplane. This magnetic flux function enables us to explore a configuration where initially, magnetic field lines threading the disk are all connected to the BH by extending the Ψ function to the whole simulation space (see Sect. 2.2.4).
2.1.3. Corona
Spectral fits of accreting stellarmass BHs and SMBHs motivate our choice to include in our model a hot collisionless electron–positron pair plasma to represent the corona. We studied the time evolution of this corona with the disk kept stationary. Although several mechanisms were invoked to explain its replenishment (e.g., pair cascade, Crinquand et al. 2020, or magnetocentrifugal loading from the disk, Blandford & Payne 1982), we adopt an empirical approach. Independently of the injection mechanism, we assume that the particle density in the corona is both low enough to ensure that the plasma is collisionless, and high enough to be in the forcefree regime. In this framework, electromagnetic forces dominate the dynamics. The plasma flows along the poloidal magnetic field lines and screens any electric component along these lines. The forcefree approximation requires a cold magnetization σ ≫ 1, where σ is defined as the ratio of the magnetic energy to the inertial mass energy of the particles:
$$\begin{array}{c}\hfill \sigma =\frac{{B}^{2}/4\pi}{({\gamma}_{+}{n}_{+}+{\gamma}_{}{n}_{}){m}_{\mathrm{e}}{c}^{2}},\end{array}$$(4)
with B = B the magnitude of the magnetic field measured by the FIDO, γ_{±} (respectively n_{±}) the Lorentz factor (respectively the number density) of the positrons/electrons and m_{e} the mass of the electron. In these conditions, the plasma satisfies the forcefree equation (see, e.g., Uzdensky 2004), and in the steady state each magnetic field line rotates at a fixed angular speed. In order for this to be achieved, this regime requires a net plasma charge density, ρ, above the threshold set by the GoldreichJulian charge density, ρ_{GJ} (Goldreich & Julian 1969):
$$\begin{array}{c}\hfill {\rho}_{\mathrm{GJ}}=\frac{\mathbf{\Omega}\xb7\mathit{B}}{2\pi c},\end{array}$$(5)
where Ω is the angular velocity vector of the magnetic field line. This charge density translates into a net number density n_{GJ} = ρ_{GJ}/e, with e the charge of a positron (e > 0). Hereafter, we assume that σ ≫ 1 and that the plasma multiplicity κ = (n_{+} + n_{−})/n_{GJ} ≳ 1. We detail in Sect. 2.2.6 how we ensure that the number density n = n_{+} + n_{−} of electrons and positrons is high enough to be in a regime where the net charge density ρ can be higher than ρ_{GJ}.
2.2. Numerics
2.2.1. The GRZeltron code
We solve the dynamics of the plasma and of the electromagnetic fields in the whole corona, including in regions where the forcefree approximation breaks down and where significant amounts of electromagnetic energy can be converted into particle kinetic energy. To do so, we use the general relativistic PIC code GRZeltron first introduced in Parfrey et al. (2019) and based on the Zeltron code (Cerutti et al. 2013). Its principle is described in detail in Crinquand et al. (2021), and here we only recall the main aspects.
We work in the 3 + 1 formalism presented in Sect. 2.1.1 and we note B and D the magnetic and electric fields locally measured by the FIDO. These fields are advanced by solving MaxwellFaraday and MaxwellAmpère equations. In these equations intervenes the gridbased current (hereafter, the current) J = qv, where q is the charge of the particles and v their 3velocity on the spatial hypersurface measured in universal coordinate time. MaxwellFaraday and MaxwellAmpère equations also make use of the auxiliary fields H and E obtained from (Komissarov 2004):
$$\begin{array}{cc}& \mathit{H}=\alpha \mathit{B}\mathit{\beta}\times \mathit{D}\hfill \end{array}$$(6)
$$\begin{array}{cc}& \mathit{E}=\alpha \mathit{D}+\mathit{\beta}\times \mathit{B}.\hfill \end{array}$$(7)
We then solved the full equation of motion for the electrons and positrons, including the Lorentz force (see, e.g., Bacchini et al. 2018, 2019; Parfrey et al. 2019). Although available in GRZeltron, we did not include any radiative drag force or cooling induced by synchrotron or Compton emission (see Mehlhaff et al. 2020; Sridhar et al. 2021 for the impact of radiative feedback on particle acceleration). We now proceed to a description of the numerical setup we designed to address the problem described in Sect. 2.1.
2.2.2. Grid
We performed global PIC simulations on a twodimensional (r, θ) grid in KerrSchild coordinates. Axisymmetry around the polar axis is assumed such that all quantities are independent of the ϕ coordinate. Although the grid is twodimensional, all three components of the electromagnetic fields and of the particle velocities are accounted for. Furthermore, our model being symmetric with respect to the BH equatorial plane, we only cover the region from θ_{min} = 0 up to θ_{max} = π/2. The inner edge of the simulation space is set within the horizon, with r_{min} = 0.9 r_{H}, while the outer edge lies at r_{max} = 30 r_{g}, far enough to limit the impact of the outer boundary conditions detailed in Sect. 2.2.3 on the region of interest. To cover the whole simulation space with a uniform cell aspect ratio, we use a logarithmically stretched spacing in r. Unless stated otherwise, the results we present are for a grid with 2048 cells in the r direction and 1120 cells in the θ direction, which ensures a cell aspect ratio close to 1. We argue in Sect. 2.2.6 that given the parameters we explore, this resolution is enough to resolve the (nonrelativistic) plasma skin depth ${\delta}_{\mathrm{e}}=\sqrt{{m}_{\mathrm{e}}{c}^{2}/4\pi n{e}^{2}}$ everywhere, a statement that we verify a posteriori in our simulations (see Appendix A).
2.2.3. Boundary conditions
The inner edge of the simulation space lies within the event horizon (r_{min} < r_{H}) so the dynamics above the event horizon is not impacted by the boundary conditions at r_{min}. At the pole (θ = 0), we apply boundary conditions based on the symmetry properties of the polar and axial electric and magnetic vectors. At the outer edge of the simulation space, we used an absorbing layer from r = 27 r_{g} to r = r_{max}. In this region, resistive terms were added to Maxwell’s equations in order to damp the electromagnetic fields and avoid spurious wave reflection (Cerutti et al. 2015). Without an ambient disk to hold the magnetic field, the initial magnetic field lines would be expelled from the BH, in agreement with the nohair theorem. We introduced a perfectly conducting disk from θ_{d} = π/2 − arctan(ϵ) to θ_{max} = π/2, where the magnetic field is frozen to its initial value (see Sect. 2.2.4). For r > r_{ISCO}, we enforced the magnetic field lines to corotate with the disk at the Keplerian angular speed of their footpoint in the equatorial plane, given by Eq. (2). Conversely, although the equatorial plane remains perfectly conducting up to the event horizon, we enforced a smooth transition from a Keplerian angular speed to zero angular speed at r < r_{ISCO}. In the disk, the electric field, E, is obtained from the assumption of a perfectly conducting environment. Therefore, the electric field in the frame rotating with each disk annulus must vanish, and we have
$$\begin{array}{c}\hfill \mathit{E}+\frac{{\mathit{V}}_{\mathrm{disk}}}{c}\times \mathit{B}=0\end{array}$$(8)
in the disk, with V_{disk} = Ω∂_{ϕ}. H and D in the disk are then deduced from the constitutive relations (6) and (7).
2.2.4. Initial conditions
We started from a dipolar magnetic field in vacuum defined by the following ϕ component of the 4potential:
$$\begin{array}{c}\hfill {A}_{\varphi}={B}_{0}\phantom{\rule{0.166667em}{0ex}}{r}_{\mathrm{g}}^{2}\frac{{sin}^{2}\theta}{r/{r}_{\mathrm{g}}}\xb7\end{array}$$(9)
Constructing the magnetic field from the 4potential enables us to enforce an initially divergencefree field, a feature that is maintained throughout the simulation. In this initial configuration, all the magnetic field lines connect the disk to the event horizon. The further the footpoint of the magnetic field line on the disk, the higher the latitude of the point where this field line intersects the event horizon.
This is the simplest divergencefree setup that connects the disk to the BH, but it introduces a set of possibly nonphysical magnetic field lines in the equatorial plane within the ISCO. Indeed, this region is thought to be orders of magnitude less dense than the disk (see Potter 2021, though they focus on a case where magnetic fields are not dynamically important). In the highly magnetized regime we are interested in, such a density drop means that this region could be forcefree; however, according to the noingrownhair theorem, closed magnetic field lines anchored in the BH can only be maintained if they intersect or loop around a nonforcefree region such as the disk (MacDonald & Thorne 1982; Gralla & Jacobson 2014). The behavior of the plasma in this region, which however vanishes as a → 1, might thus be altered by this initial magnetic field configuration, especially because it is enforced as a boundary condition in the equatorial plane.
2.2.5. Particle injection
The corona is initially empty and needs to be populated in order to match the forcefree conditions described in Sect. 2.1.3. Then, since particles propagate at relativistic speeds and eventually escape from the simulation space, the corona has to be replenished in order to prevent it from emptying within a few tens of r_{g}/c. In this paper we remain agnostic on the physical mechanism that ensures that the corona is not empty and embrace an empirical approach. At each time step, we inject electron/positron pairs in cells where n < 3Ω_{ISCO}B/2πce with a numerical weight given by
$$\begin{array}{c}\hfill w=\frac{{\mathrm{\Omega}}_{\mathrm{ISCO}}B}{2\pi ce},\end{array}$$(10)
where Ω_{ISCO} = Ω_{K}(R=r_{ISCO}) is the Keplerian angular speed in the disk midplane at the ISCO, and B = B is the magnitude of the local magnetic field. The choice of Ω_{ISCO} as a reference angular speed was found to be a good tradeoff to load both magnetic field lines anchored in the disk and in the event horizon with enough plasma, for any BH spin. Given that we start with a dipolar magnetic field, this weight guarantees that the number of particles per cell remains reasonably uniform through the simulation space since the r^{−3} dependence of the magnetic field compensates the cell volume increase from r_{min} to r_{max}. Particles are deleted when they enter the disk, the absorbing layer or leave the simulation space. They do not experience any radiative feedback in these numerical simulations.
2.2.6. Parameters
The main degreesoffreedom of our model are the BH spin (a), the disk aspect ratio (ϵ), and the reference magnitude of the magnetic field (B_{0}). In the simulations, the magnetic field manifests through the ratio of the gravitational radius to the Larmor radius R_{L, 0} = m_{e}c^{2}/eB_{0}. For the BH spin, we explore the case of a nonrotating BH (a = 0), rotating BHs with increasing spins from a = 0.6 to a = 0.99 surrounded by a prograde disk. At lower spin than a = 0.6, critical features expected for a > 0 like the separatrix (see Sect. 3.1.1) extend up to the outer edge of our grid at r_{max} = 30 r_{g}. In order to avoid numerical artifacts due to the outer boundary conditions, we limit our study of prograde disk to BH spins a ≥ 0.6. We also investigate the case of a rotating BH with a spin a = 0.8 surrounded by a counterrotating disk (i.e., a = −0.8). We consider a geometrically thin disk with ϵ = 5% and obtain qualitatively similar results for a thicker disk with ϵ = 30% around a BH with spin a = 0.8.
Realistic values of the dimensionless magnetic field ${\stackrel{\sim}{B}}_{0}={r}_{\mathrm{g}}/{R}_{\mathrm{L},0}$ are out of reach of the current computational capacities. Indeed, in M 87^{*} and Sgr A^{*}, the magnetic field near the event horizon ranges between approximately 1 and 100 G, which corresponds to ${\stackrel{\sim}{B}}_{0}\sim 4\times {10}^{810}$ for Sgr A^{*} and ${\stackrel{\sim}{B}}_{0}\sim 5\times {10}^{1113}$ for M 87^{*}. In Cyg X1, ${\stackrel{\sim}{B}}_{0}$ is comparable since the gravitational radius, smaller by 5 to 7 orders of magnitude, is counterbalanced by a similar increase in the magnetic field magnitude (Del Santo et al. 2013; Cangemi et al. 2021). To remain in the forcefree regime with such high magnetic fields, we would need a much higher particle number density according to Eq. (5): a high magnetic field then sets a stringent upper limit on the skin depth, which becomes so small that it cannot be resolved anymore. Reasoning in terms of number density, we must match the three following conditions in the whole simulation space: (i) the skin depth is resolved (i.e., $n<{n}_{\delta}=\frac{{m}_{\mathrm{e}}{c}^{2}}{4\pi {e}^{2}}/\mathrm{\Delta}{r}^{2}$, where Δr is the radial extent of the cell and where we assumed a cell aspect ratio close to unity); (ii) the magnetization σ ≫ 1 (i.e., n ≪ n_{σ} = B^{2}/4πm_{e}c^{2}); and (iii) the multiplicity κ ≳ 1 (i.e., n ≳ n_{GJ} = Ω_{ISCO}B/2πce).
The main numerical challenge of this setup lies in the combination of the dipolar magnetic field and the high ratio between the radius of the outer and inner edges of the simulation space since the aforementioned number densities vary by several orders of magnitude from r_{min} to r_{max}. In Fig. 1 we show representative radial profiles for n_{δ}, n_{σ}, and n_{GJ} assuming B ∝ r^{−3}. The arrows indicate where the plasma number density must lie with respect to each profile according to the aforementioned conditions. The scaling of each profile is chosen such that the three conditions can be matched for any r (green hatched region), which requires
$$\begin{array}{c}\hfill 2{\left(\frac{{r}_{\mathrm{max}}}{{r}_{\mathrm{g}}}\right)}^{3/2}\ll {\stackrel{\sim}{B}}_{0}<\frac{1}{2}\frac{{r}_{\mathrm{min}}^{3}}{{r}_{\mathrm{g}}\mathrm{\Delta}{r}^{2}(r={r}_{\mathrm{min}})}\xb7\end{array}$$(11)
Fig. 1. Fiducial radial profiles for the critical densities n_{δ}, n_{σ}, and n_{GJ} defined in the text. The numbers are the powerlaw exponents, and the arrows indicate the allowed regions. The green hatched region represents the region where the forcefree regime can be achieved while still resolving the skin depth. 
With the spin, resolution, and radial boundaries we use, this means that we must work with a dimensionless magnitude of the magnetic field that approximately verifies
$$\begin{array}{c}\hfill 320\ll {\stackrel{\sim}{B}}_{0}<190\phantom{\rule{0.166667em}{0ex}}000.\end{array}$$(12)
In what follows, we set ${\stackrel{\sim}{B}}_{0}={10}^{5}$, a value suitable to guarantee that from r_{min} to r_{max}, the forcefree regime is achievable for plasma densities at which the skin depth is resolved by our grid.
3. Results
3.1. Magnetic structure
In Sect. 3.1.1 we first analyze the global structure of the corona in the light of properties derived in the forcefree approximation and verify a posteriori the validity of our injection method. We then detail the properties of each of the three main regions identified in Sects. 3.1.2 and 3.1.3. In Sect. 3.1.5 we describe the current sheet that separates the regions where magnetic field lines are open and anchored either in the disk or in the BH.
3.1.1. Forcefree regions
In Fig. 2 we represented color maps of four different quantities overlaid on top of the magnetic field lines, after the initial setup has relaxed. The BH spin is 0.8, the disk is thin (ϵ = 5%), and the results are qualitatively similar to what is found for other spins ≥0.6. The two thick magnetic field lines correspond to the innermost and outermost coupling field lines anchored in the BH and in the disk, respectively at the ISCO and of magnetic potential A_{ISCO} and at r = R_{S} and of magnetic potential A_{S}. Following Uzdensky (2005), we call the latter field line the separatrix. Space around the BH is subdivided into four distinct regions delimited by those two critical field lines.
Fig. 2. Poloidal magnetic field lines around the BH (solid black lines), with the event horizon (black disk) and the ergosphere (black dashed line), for a fiducial snapshot from a simulation with a = 0.8 and a thin disk (in chartreuse yellow in the equatorial plane). The thicker field lines delimit the region between the ISCO and the separatrix coupling the BH to the disk, visible in the zoomedin view in the bottom panel. The color maps represent the total plasma density, n, in units of n_{GJ} (top left), the mean Lorentz factor of the particles (top right), the toroidal component of the H field (bottom left), and the radial component of the current, J (bottom right). To compensate for spatial dilution, plasma density and current were multiplied by r^{2}. A smoothing Gaussian kernel a few cells wide was applied to the current map. In the plasma density map, the red circle locates the Y point. Distances to the BH on the x and y axes are given in units of r_{g}. The white line in the topright corner stands for the outer light surface. 
In the innermost region, where A > A_{ISCO}, magnetic field lines are anchored in the equatorial plane and do not rotate, in agreement with the boundary condition at their footpoint. As explained in Sect. 2.2.4, it is an unphysical configuration due to the lack of nonforcefree medium in the region where r < r_{ISCO} and θ ∈ [θ_{d}; θ_{max}] (in white at the equator in Fig. 2). In reality, these magnetic field lines would be stretched in the equatorial plane due to the plasma quickly falling toward the BH. It would form a current sheet where the forcefree approximation would break up. In such a case, the magnetic field lines would either close within the horizon (Komissarov 2004) or be fully open (Crinquand et al. 2021; Bransgrove et al. 2021). In the latter case, though, no coupling between the BH and the disk would exist since no magnetic field line would connect the two. In the perspective of a magnetically coupled configuration, it is likely that the lack of a significant vertical magnetic field in the region where r < r_{ISCO} and θ ∈ [θ_{d}; θ_{max}] would lead to a field line anchored at the ISCO that would penetrate the horizon at much lower latitude than in our simulations, resembling more the results obtained by Uzdensky (2005) and Yuan et al. (2019a) who directly solved the relativistic GradShafranov equation. The bloated region at the equator near the event horizon, visible in the zoom in bottom panel in Fig. 2, is thus probably the least accurate feature in these simulations. On the other hand, the three other regions, further described in the next sections, are both physically realistic and numerically robust.
Except in the current sheet, the electric field colinear to the magnetic field, D ⋅ B/B^{2}, remains marginally small, with values typically below a few times 10^{−4}, which is indicative of a corona close to forcefree. All over the grid, the initial conditions of high magnetization and plasma pair multiplicity above 1 described in Sect. 2.2.6 are respected. The skin depth is also resolved everywhere. We plotted in Fig. 3 the histogram of (A_{ϕ}, Ω) values in each cell of the simulation space. The scale of the color map is logarithmic, with a cutoff at its lower end fixed at 1% of the maximum value of the histogram. The local angular speed, Ω, of the magnetic field lines is estimated from Eq. (8):
$$\begin{array}{c}\hfill \mathrm{\Omega}=\frac{{E}_{\theta}}{\sqrt{h}{B}^{r}},\end{array}$$(13)
Fig. 3. Angular speed of the magnetic field lines. Top panel: color map of the timeaveraged angular speed, Ω, with poloidal magnetic field lines overlaid, for a simulation with a = 0.6 and ϵ = 5%. The dashed white line is the outer light surface. Bottom panel: angular speed of magnetic field lines as a function of their magnetic flux function, A_{ϕ}, with the corresponding distance of the line footpoint on the disk at the top. The color map is the histogram of the (A_{ϕ}, Ω) pairs measured in each cell of the simulation space at a fiducial time. The left (respectively right) vertical line stands for the separatrix footpoint (respectively the ISCO). In orange are the Keplerian profiles (Eq. (2)) without (solid) and with (dashed) the smooth cutoff, and half of the BH angular speed (dashdotted). 
with h the determinant of the spatial 3metric. When overlaid on top of poloidal magnetic field lines, the Ω isocontours match (see top panel in Fig. 3), which indicates that each magnetic field line rotates at a fixed angular speed. They depart from each other only in the immediate vicinity of the current sheet, where B^{r} is close to zero. The two vertical lines in the bottom panel in Fig. 3 stand for the outermost magnetic field line (left line) and for the magnetic field line anchored at the ISCO (right line). The top x axis indicates the footpoints on the disk of the magnetic field lines corresponding to a given magnetic potential, A_{ϕ}, according to Eq. (9), and thus, the distance to the BH increases toward the left. At these footpoints, the local Keplerian angular speed profile (2) is enforced (orange solid line), with a cutoff at the ISCO (orange dashed line). The field lines that remain anchored in the disk follow this profile, whether they are open (for A_{ϕ} < A_{S}, lower branch) or coupled to the BH (for A_{ϕ} > A_{S}). The apparent fading at low A_{ϕ} and low Ω is due to the lower number of cells of the mesh at large distance from the BH. The opening of the magnetic field lines for A < A_{ϕ} leads to the set of open field lines anchored in the BH. They span the same range of A_{ϕ} values as their counterparts anchored in the disk but they rotate at an angular speed close to ω_{H}/2, the optimal value for rotational energy extraction from the BH according to the forcefree BlandfordZnajek process (upper branch for A < A_{S} in Fig. 3).
Although the numerical grid we work with is twodimensional, we use the axisymmetric assumption to provide a threedimensional representation of the magnetic field lines in Fig. 4 for a BH spin a = 0.8 and a disk thickness ϵ = 5%. The upper and lower surfaces of the disk are semitransparent to make the field lines below the disk visible. They extend down to the ISCO and the central black sphere is the event horizon of the BH. The orientation of each field line is specified and their color stands for their magnetic potential A_{ϕ}, from dark green for low values (i.e., footpoint of the magnetic field line far on the disk at initial state) to magenta for high values (i.e., footpoint near the ISCO). Near the BH spin axis, we find a set of open magnetic field lines anchored in the event horizon, reminiscent of the configuration considered by Blandford & Znajek (1977). They correspond to the four twisted field lines threading the event horizon near the pole, visible in both panels in Fig. 4 (dark green, light green, blue and yellow field lines). We show in Sect. 3.1.2 that given their properties, these field lines do form the magnetic backbone of a relativistic collimated jet. In the outermost regions, above the disk, there are open magnetic field lines anchored in the disk. They formed after the closed magnetic field lines of the initial dipole split beyond the separatrix footpoint. A subset is visible in the upper and lower panels in Fig. 4 and they have the same colors as the field lines crossing the event horizon since they have the same magnetic potential A_{ϕ} < A_{S}. Finally, the field lines with a magnetic potential between A_{S} and A_{ISCO} remain anchored both in the BH and in the disk, forming a closed region that couples these two components. In the bottom panel in Fig. 4, they are the two innermost magnetic field line, represented in coral and magenta. In between the two regions where the magnetic field lines are open, a current sheet forms. Magnetic reconnection kicks in and plasmoids form at the point where these two regions meet with the closed part of the magnetosphere, called the Y point (red circle in the upper left panel in Fig. 2).
Fig. 4. Threedimensional representation of the magnetic field lines around the BH (central black sphere), with the upper and lower disk surfaces in semitransparent black. Field lines are colorcoded from low magnetic potential (dark green) to high magnetic potential (magenta). Two additional field lines are represented in the zoomedin view in the bottom panel. 
For a geometrically thick disk (ϵ = 30%) such as the one that is believed to be present around M 87^{*}, the fundamental mechanism responsible remains effective. We carried out a simulation of a thick disk in prograde rotation around a BH with spin a = 0.8 and found qualitatively similar results, with the same regions of closed and open magnetic field lines, a Y point at the end of the separatrix and a tearing unstable current sheet. The location of the footpoint of the separatrix at the surface of the disk does not change significantly with respect to simulations with a lower aspect ratio. The plasma density in the current sheet is more dilute and even if plasmoids still form, magnetic reconnection is less vivid with marginally lower rates of electromagnetic energy dissipated. Although these results suggest that our conclusions might hold in a thick disk around a SMBH, it must be acknowledged that the biases introduced by the steady disk are enhanced in this configuration, where the disk has a larger spatial extent. In this case, a magnetic diffusivity could be set in the disk as boundary conditions (Parfrey et al. 2017). Alternatively, a fluid approach where kinetic effects are encapsulated in parametrized magnetic diffusivity coefficients would be insightful to derive a physically accurate magnetic structure in the disk (Ripperda et al. 2020). In what follows, we detail the results we obtain with thin disks only.
3.1.2. Opened magnetic field lines anchored in the black hole
In the polar regions of the BH, magnetic field lines thread the event horizon at their basis and form a flux tube around the spin axis. In some previous studies devoted to magnetic configurations where the BH is connected to a surrounding Keplerian disk, there were no open magnetic field lines threading the horizon (hereafter the jetlaunching region) because of numerical constrains (Uzdensky 2005; Yuan et al. 2019a). On the contrary, we find that the jetlaunching region coexists with magnetic field lines that couple the disk to the BH. Indeed, Uzdensky (2005) rely on numerical relaxation techniques to find a steady solution of the relativistic GradShafranov equation, with both the BH spin a and the footpoint of the separatrix R_{S} as degreesoffreedom (see also Uzdensky 2002). These methods require an initial guess for the magnetic potential A_{ϕ} but due to the lack of nonideal effects in the forcefree framework, the topology of the magnetic field is frozen as A_{ϕ} successively converges toward the solution. The regularity of the solution at the event horizon (Nathanail & Contopoulos 2014) and the fact that the GradShafranov equation is singular on the light surfaces drove Uzdensky (2005), Mahlmann et al. (2018) and Yuan et al. (2019a) to consider solutions where a separatrix subdivides magnetic field lines in two sets: the open ones anchored in the disk, beyond the separatrix, and the ones coupling the BH to the disk, with a footpoint on the disk closer than R_{S}. Our approach does not suffer from this limitation since we solve a more fundamental set of equations where magnetic reconnection is allowed and susceptible to modify the topology of the initially prescribed magnetic field. In the intermediate approach of timedependent forcefree simulations, where the topology of the magnetic field can change due to numerical or physical dissipation, a hybrid closedopen magnetosphere is also found, although with different initial and boundary conditions (Parfrey et al. 2015; Yuan et al. 2019b; Mahlmann et al. 2020).
The open magnetic field lines threading the horizon are dragged by the rotation of the BH, which manifests as a twisting of the magnetic field lines in Fig. 4. A significant toroidal component of the magnetic field on the grid, H_{ϕ}, develops. Since the magnetic moment of the initial dipole is aligned along the same direction as the angular momentum of the BH, H_{ϕ} is positive (respectively negative) in the lower (respectively upper) hemisphere, as visible in the bottomleft panel in Fig. 2. The upper branch for A < A_{S} in Fig. 3 shows that most of these field lines rotate at an angular speed close to ω_{H}/2 within 20%. The field lines near the spin axis (i.e., at low A_{ϕ}) tend to rotate slower than ω_{H}/2. It is probably an effect of the outer edge of the simulation space since these field lines are precisely those that do not cross the outer light surface inside the simulation space.
We now comment on the magnetic flux and the electromagnetic power carried in this region. We represented these quantities measured for different BH spin values in the middle panel in Fig. 5. We computed the magnetic flux Φ_{jet} through a solid angle from the pole up to the separatrix. Since the magnetic field is divergencefree, we can compute this flux through a fraction of a sphere at any distance from the BH. We average the values obtained over radii from r_{H} up to the sphere that intersects the Y point, in order to avoid the turbulent current sheet beyond. In Fig. 5 we compare it to the magnetic flux threading the full upper hemisphere of the BH at the event horizon, Φ_{H}, given by its initial value for a magnetic dipole:
$$\begin{array}{c}\hfill {\mathrm{\Phi}}_{\mathrm{H}}={\displaystyle \u222f}\mathit{B}\xb7\mathrm{d}\mathbf{\Sigma}=2\pi {A}_{\varphi}({r}_{\mathrm{H}},\pi /2)=2\pi {B}_{0}\phantom{\rule{0.166667em}{0ex}}{r}_{\mathrm{g}}^{2}\frac{{r}_{\mathrm{g}}}{{r}_{\mathrm{H}}}\xb7\end{array}$$(14)
Fig. 5. Reduced quantities as a function of the BH spin. Top panel: ISCO (solid green line), event horizon (solid black line), and ergosphere at the equator (dashed black line) as a function of the BH spin. The blue line indicates the approximate location of the footpoint of the separatrix found by Uzdensky (2005), and the blue squares are the values measured in our simulations. Middle panel: magnetic flux Φ_{jet} of the open magnetic field lines threading the event horizon (green, left axis). In blue (right axis) are the jet electromagnetic powers measured (squares) and predicted by Tchekhovskoy et al. (2011) (crosses). Bottom panel: angular momentum (left axis) and energy (right axis) transferred from the BH to the disk per unit time, with and without the forcefree approximation (crosses and squares). 
The measured fraction Φ_{jet}/Φ_{H} carried by the open magnetic field lines threading the BH event horizon increases from 30% at a = 0.6 up to 45% for a = 0.9 and it plateaus beyond. In parallel, we estimated the power carried by the jet, P_{jet}, by computing the flux of the Poynting vector as seen from an observer at infinity:
$$\begin{array}{c}\hfill \mathit{S}=\frac{c}{4\pi}\mathit{E}\times \mathit{H}.\end{array}$$(15)
This flux hardly varies with the distance to the BH, which indicates that electromagnetic energy is seldom dissipated in this region, in agreement with the forcefree regime. We also checked that the flux of kinetic energy of the particles is negligible compared to the Poynting flux in this region: the basis of the jet is essentially devoid of mass and is controlled by the electromagnetic fields. The jet power quickly increases as the spin increases from 0.6 to 0.99. In order to provide an empirical characterization of this evolution, we compared P_{jet} to a formula inspired by the results of Tchekhovskoy et al. (2011):
$$\begin{array}{c}\hfill {P}_{\mathrm{th}}=\frac{k}{4\pi}\frac{{\omega}_{\mathrm{H}}^{2}}{c}{\mathrm{\Phi}}_{\mathrm{jet}}^{2}f\left({\omega}_{\mathrm{H}}\right),\end{array}$$(16)
where f is a correcting factor at large spin found by Tchekhovskoy et al. (2011) and given by
$$\begin{array}{c}\hfill f({\omega}_{\mathrm{H}})=1+1.38{\left(\frac{{\omega}_{\mathrm{H}}{r}_{\mathrm{g}}}{c}\right)}^{2}9.2{\left(\frac{{\omega}_{\mathrm{H}}{r}_{\mathrm{g}}}{c}\right)}^{4}.\end{array}$$(17)
The k factor encompasses the information on the geometry of the magnetic field lines and does not depend on the BH spin. For a split monopole (respectively a parabolic) geometry, k ∼ 0.053 (respectively k ∼ 0.044). Since we consider a different magnetic configuration, we reinjected the magnetic flux Φ_{jet} we measure into Eq. (16) and fitted for the k factor. We obtained the blue crosses in the middle panel in Fig. 5, which correspond to
$$\begin{array}{c}\hfill k=0.072\pm 0.002.\end{array}$$(18)
We retrieved that the correcting factor f is necessary to capture the dependence of P_{jet} on Φ_{H} for spin values a ≥ 0.95.
3.1.3. Opened magnetic field lines anchored in the disk
Among the magnetic field lines anchored in the Keplerian disk, those beyond the separatrix do not remain connected to the BH and open up. The footpoint of the separatrix on the disk is located at a distance R_{S} from the BH, which essentially depends on the BH spin (and weakly on the disk thickness). The existence of this outermost closed magnetic field line at a finite distance around spinning BHs was first highlighted by Uzdensky (2005). The key argument is that in an initial configuration where magnetic field lines all thread the disk and the event horizon, the field lines will slip along the horizon at different angular speed owing to the different locations of their footpoints on the Keplerian disk. This longitudinal shearing of the field lines induces a toroidal component. For a > 0, near enough from the BH pole, this component will induce a magnetic pressure too high to be confined by the tension of the poloidal component at large distances. It leads to the opening up of the magnetic field lines we observe in the relaxation phase of our simulations.
Uzdensky (2005) estimated the maximal extent of the coupling part of the magnetosphere on the disk by looking for a steady solution of the relativistic GradShafranov equation at a given BH spin a and location of the separatrix R_{S}. For each spin, Uzdensky (2005) showed that there is a maximal value of R_{S} beyond which the solver could not converge. We represented the upper limit Uzdensky (2005) found in the upper panel in Fig. 5 (blue line). The gray shaded region corresponds to the region within the event horizon (solid black line). The green line is the inner edge of the Keplerian disk, set at the ISCO, and lies within the ergosphere (black dashed line) for a greater than approximately 0.95. We located the separatrix by using the change in sign of the H_{ϕ} variable, which happens in the current sheet where the magnetic field reconnects (see the bottomleft quarter in Fig. 2). Once the simulations have relaxed, the position of the footpoint remains stable within 5%. The method we use yields robust estimates of R_{S}, which represent lower limits since transient closed magnetic field lines are observed up to a footpoint ∼20% higher than R_{S} (e.g., in Fig. 2 where a closed magnetic field line anchored at ∼4.7 r_{g} on the disk extends up to ∼23 r_{g}).
Finally, we notice that the open magnetic field lines anchored in the disk are very inclined, with a maximum inclination i with respect to the disk plane approximately given by the inclination of the current sheet. It is visible both in Fig. 2 and in the threedimensional representation in Fig. 4. As a varies between 0.6 and 0.99, we find that i remains approximately between 30° and 35°. In the magnetocentrifugal model introduced by Blandford & Payne (1982), the maximum inclination angle to launch an outflow is 60°, which is > i for any a ≥ 0.6. It means that matter from the disk is susceptible to slide out along the magnetic field lines and provide a possible source of matter to ensure forcefree conditions in the corona. Yuan et al. (2019a) studied this effect in great detail and showed that when accounting for relativistic effects, the range of inclinations suitable to launch a disk wind is even wider, especially at higher spins.
3.1.4. Coupling magnetic field lines
As the BH spin increases, the separatrix moves closer to the BH, and so does the ISCO. The net effect is a smaller disk surface connected to the BH at larger spins, but the connected region lies deeper in the gravitational potential and probes stronger magnetic fields. We now quantify ${\dot{L}}_{\mathrm{BH}\mathrm{disk}}$ and Ė_{BH−disk}, the amount of angular momentum and energy exchanged per unit time, respectively. To do so, we used the flux of energy density at infinity given by Eq. (15), S, and the flux of angular momentum in an axisymmetric configuration, S_{L}, given by (Blandford & Znajek 1977; MacDonald & Thorne 1982; Komissarov 2004)
$$\begin{array}{c}\hfill {\mathit{S}}_{\mathrm{L}}=\frac{1}{4\pi}({E}_{\varphi}\mathit{D}{H}_{\varphi}\mathit{B}),\end{array}$$(19)
where positive values denote transfers from the BH to the disk. We integrate these fluxes over both sides of the disk, from the ISCO up to the separatrix, to obtain ${\dot{L}}_{\mathrm{BH}\mathrm{disk}}$ and Ė_{BH−disk}, displayed as green and blue squares, respectively, in the bottom panel in Fig. 5. We compared these values with what would be expected in a fully forcefree configuration, where the infinitesimal angular momentum and energy exchanged per unit time along a magnetic flux tube delimited by field lines of magnetic potential A_{ϕ} and A_{ϕ} + dA_{ϕ} would be given by (Blandford & Znajek 1977)
$$\begin{array}{cc}& \mathrm{d}{\dot{L}}_{\mathrm{ff}}={H}_{\varphi}\mathrm{d}{A}_{\varphi}\hfill \end{array}$$(20)
$$\begin{array}{cc}& \mathrm{d}{\dot{E}}_{\mathrm{ff}}={H}_{\varphi}\mathrm{\Omega}\mathrm{d}{A}_{\varphi}.\hfill \end{array}$$(21)
For each spin, we integrate the quantities H_{ϕ} and H_{ϕ}Ω over the magnetic potential A_{ϕ} on a wedge just above the disk, from the ISCO to the separatrix, and obtain the green and blue crosses in the bottom panel in Fig. 5 for the exchange rates of angular momentum and energy, respectively.
The good match between the quantities obtained from the angular momentum and energy fluxes (15) and (19) and the ones derived in the idealized forcefree case with Eqs. (20) and (21) indicate that only a small fraction of the energy funneled along the coupling magnetic field lines is dissipated between the BH and the disk. Except in the immediate vicinity of the Y point, the closed part of the magnetosphere is essentially forcefree and little electromagnetic energy is dissipated along the magnetic field lines coupling the disk to the event horizon.
We notice that $\dot{L}$ is positive for any a ≥ 0.6, which can be understood in terms of corotation radii. At a ∼ 0.36, the corotation radius is approximately located at the ISCO so for higher spins, the angular speed from the ISCO to the separatrix is lower than the BH angular speed ω_{H}. Thus, the coupling magnetic field lines can only spin down the BH and transfer angular momentum to the disk, which is consistent with the positive values we measure. The rate at which rotational energy is extracted from the BH and deposited onto the disk is comparable to the jet power for a ≥ 0.6.
3.1.5. Current sheet
For any BH spin a ≥ 0.6 and for both disk thickness ϵ = 5% and ϵ = 30%, a current sheet forms between the two regions of opened magnetic field lines. As mentioned in Sect. 3.1.3, the H_{ϕ} component displays a clear change of sign at the current sheet, visible in the bottomleft quarters in Fig. 2. It is in agreement with the toroidal motion of the magnetic field lines, which rotate with the BH and with the prograde disk (see Fig. 4). Due to the opening of the magnetic field lines, the polarity reverses at the current sheet and the rotational drag of the field lines by the ergoregion or by the disk yields a toroidal component of H of opposite sign on both sides of the current sheet. In the global steady state we reach, MaxwellAmpère’s equation shows that a jump in H_{ϕ} along the longitudinal direction corresponds to a local peak in J^{r}. The current sheet hosts a strong current flowing away from the BH visible in the bottomright quarters in Fig. 2. This positive current is essentially carried by positrons and the associated electric circuit closes in the disk. Magnetic reconnection takes place in this current sheet, mostly at the Y point but also further along the sheet, at local X points. With the magnetic field amplitude we can reach (see Sect. 2.2.6), the sheet is resolved with 3 to 10 grid points and has a typical thickness δ ∼ 5 × 10^{−2} r_{g}, thin enough to be tearing unstable (Zenitani & Hoshino 2007). The current sheet breaks up, which leads to the formation of a chain of magnetic islands where particles gather into plasmoids (Loureiro et al. 2007), visible in the upper left quarters in Fig. 2. The distance between plasmoids vary with the BH spin.
We now estimate the reconnection rate by looking at the physical properties of the upstream flow near the current sheet. In Fig. 6 we report the transverse profiles of the components of the positrons and electrons velocities normal to the current sheet. In a simulation with an intermediate BH spin (a = 0.8) and a thin disk (ϵ = 5%), we worked at the innermost X point on the current sheet, formed just after a plasmoid detached from the Y point. The box around the current sheet we consider is 3 r_{g} long (x axis in Fig. 6), centered on the X point and 0.5 r_{g} wide so as to work with averaged quantities and reduce the numerical noise. Figure 6 displays the poloidal component of the particle velocity normal to the current sheet, v_{⊥}, with v_{⊥} > 0 (respectively v_{⊥} < 0) if particles are above the current sheet and move toward it or if they are below the current sheet and move away from it (respectively if they are above the current sheet and move away from it or if they are below the current sheet and move toward it). The current sheet is approximately along a radius (θ ∼ constant) but since the plasma tends to flow along the magnetic field lines, v^{r} ≫ v^{θ} in general and both velocity components contribute to v_{⊥}. The exact shape of the v_{⊥} profiles depends on the inclination angle of the current sheet we use to perform the projection but the glitches at the X point (i.e., at the vertical gray line) subsist for any reasonable inclination. The profiles have a net nonzero offset, which corresponds to a bulk instantaneous transverse motion of the current sheet, possibly due to the drift kink instability (Zenitani & Hoshino 2007; Barkov & Komissarov 2016). Once corrected for this offset, a change of sign at the X point is measured with a step Δv_{⊥}/c ∼ 5% for positrons and Δv_{⊥}/c ∼ 15% for electrons. In the highly magnetized regime where we work, σ ≫ 1 upstream, the Alfven speed is comparable to the speed of light and the reconnection rate β_{rec} is given by
$$\begin{array}{c}\hfill {\beta}_{\mathrm{rec}}\sim \frac{\mathrm{\Delta}{v}_{\perp}/2}{c}\sim 5\%.\end{array}$$(22)
Fig. 6. Profiles along a line transverse to the current sheet and passing through an X point. The components of the drift velocity (solid black line) and of the positron and electron velocities (respectively blue circles and green squares) normal to the current sheet are represented. 
An alternative way to compute the reconnection rate relies on the electromagnetic fields and the associated drift velocity, v_{d}, defined as
$$\begin{array}{c}\hfill \frac{{\mathit{v}}_{\mathrm{d}}}{c}=\frac{\mathit{E}\times \mathit{B}}{{B}^{2}},\end{array}$$(23)
and Γ = v_{d}/c is the Lorentz factor of the plasma bulk motion. In Fig. 6 we show the transverse profile of the component of v_{d} normal to the current sheet component of the drift speed (solid black line). It is coherent with the particle velocity profiles and the step in velocity is intermediate between the values obtained from the positron and electron profiles. We compute β_{rec} ∼ 5%, which is similar to what Crinquand et al. (2021) measured in PIC simulations of magnetic reconnection in the equatorial plane near the event horizon.
These values are of the same order of magnitude as what was computed in PIC simulations of magnetic reconnection in collisionless pair plasma (Kagan et al. 2015; Sironi & Spitkovsky 2014; Werner et al. 2018). They are approximately an order of magnitude higher than the ones derived from GRMHD simulations of BH magnetosphere where the finite volume approach induces a numerical diffusivity at small scales (see Bransgrove et al. 2021, for a comparison of magnetic reconnection in PIC and GRMHD simulations). Furthermore, it must be noticed that this reconnection configuration presents some specific features that differ from idealized situations. First, magnetic reconnection is asymmetric here since the plasma properties on the disk and jet sides (e.g., density, velocity, and multiplicity) differ (see Mbarek et al. 2021, for recent simulations of asymmetric reconnection in the relativistic regime). A persistent signature of this asymmetry manifests in the leftright asymmetry of the transverse profiles in Fig. 6, where the profiles significantly depart from a simple step function. In the poloidal plane, the bulk motion of the plasma is faster above than below the current sheet (see Fig. 6). This shear of the plasmoids can be seen in the topleft quarters in Fig. 2 and can trigger KelvinHelmholtz instability and enhance the fraction of electromagnetic energy dissipated in the process (Sironi et al. 2021). Second, the angular speed of the magnetic field lines threading the plasma varies, with magnetic field lines above the current sheet rotating faster since they cross the event horizon (see Fig. 3). It creates an additional shear in the toroidal direction, in addition to the strong twisting of the magnetic field lines visible in Fig. 4.
3.2. Particle acceleration
3.2.1. Y point
At the furthest point on the separatrix, the outermost closed magnetic field line, the current sheet connects to the closed magnetosphere, forming a Y point. Below, most of the plasma flows to the disk at relativistic speed (see upper right insert of bottom panel in Fig. 2). It hits the disk at the point where r ∼ R_{S}, where the separatrix is anchored. The Y point is located near the outermost light surface, as visible in the upper right quarter of the upper panel in Fig. 2. It moves back and forth along the current sheet as plasmoids form and flow away. In the left panel in Fig. 7, we show the distance d_{Y} between the Y point and the footpoint of the separatrix it belongs to as a function of time, for a = 0.6 (blue squares) and a = 0.99 (green circles). It shows how the separatrix progressively stretches (i.e., when d_{Y} increases) until magnetic reconnection occurs at a point lower on the separatrix and a plasmoid detaches. This lower point becomes the new Y point, which corresponds to the sharp decreases in the left panel in Fig. 7, and the separatrix stretches again, closing the loop. We can see that for higher BH spin values, the extent of this stretching is more limited and the Y point spans a smaller region. A higher spin also means more frequent plasmoid formation and shorter time intervals τ_{Y} between successive tearing of the separatrix near the Y point.
Fig. 7. Plasmoid detachment from the Y point. Left panel: distance between the Y point and the footpoint of the separatrix, d_{Y}, as a function of time for different BH spin values. The origin of time is arbitrary, but over these time lapses the simulations have reached a numerically relaxed state. Right panel: estimates of the time periodicity of the motion of the Y point, τ_{Y}, as a function of the BH spin. 
In the right panel in Fig. 7, we represented estimates of τ_{Y} as a function of the BH spin, with error bars. It illustrates the aforementioned trend, with τ_{Y} decreasing from ∼9 r_{g}/c to less than 5 r_{g}/c for a increasing from 0.6 to 0.99. These values are related to the reconnection rate β_{rec} measured in Sect. 3.1.5 using the upstream particle speed and the electromagnetic fields near an X point in the current sheet. Indeed, the drops in d_{Y} happen each time a new magnetic island is formed, which depends on how fast the magnetic field lines reconnect. The most unstable mode of the tearing instability corresponds to a spacing L between X points of $2\pi \sqrt{3}\delta $ (Zenitani & Hoshino 2007). Consequently, if the tearing instability controls the rate of formation of plasmoids at the Y point, we can estimate a reconnection rate β_{rec, Y} as
$$\begin{array}{c}\hfill {\beta}_{\mathrm{rec},\mathrm{Y}}=\frac{L}{c{\tau}_{\mathrm{Y}}},\end{array}$$(24)
which ranges between 6 and 11% for τ_{Y} from 9 to 4.8 r_{g}/c. The good match with the order of magnitudes derived in Sect. 3.1.5 suggests that the tearing instability is the dominant mechanism at play near the Y point. Also, the decreasing τ_{Y} we measure as the BH spin increases can be interpreted in the light of the spatial extent of the closed magnetosphere: as the spin increases, the Y point moves closer to the BH and the distance L covered by the Y point as the line stretches shrinks, as visible in Fig. 8. In the highly magnetized regime we work in, the reconnection rate hardly varies (see, e.g., Werner et al. 2018) and the characteristic timescale τ_{Y} between the formation of two plasmoids at the Y point thus shrinks.
Fig. 8. Density maps of the locations of the Y point for different BH spin values. The inner solid lines correspond to the event horizons, and the dashed lines are the separatrix. Distances to the BH are given in r_{g}. 
In order to better appreciate the typical coordinates of the Y point above the disk, we represented in Fig. 8 the spatial occurrence rate of the Y point for different spins from 0.6 (purple) to 0.99 (red). The dashed lines stand for the timeaveraged separatrix and the solid lines for the event horizon. We retrieve that the Y point is closer from the separatrix footpoint for higher BH spins and that it moves within a smaller region. Figure 9 shows how the height above the disk y_{Y} and the projected distance to the BH x_{Y} of the Y point decrease when the spin increases (respectively green and blue points), along with the distance r_{Y} to the BH (black points). The error bars represent the 1σ variation in these quantities around their median value, and the dashed lines are the best fit obtained based on power laws (see Appendix B).
Fig. 9. Coordinates (x_{Y}, y_{Y}) and distance, r_{Y}, to the BH center of the Y point as a function of the BH spin. The dashed lines are bestfit power laws (see Appendix B). 
In our axisymmetric framework, this Y point unfolds into a threedimensional ring above and below the disk. Once the quasiperiodic motion of the Y point along the current sheet is accounted for, it manifests in threedimensional maps integrated over durations > τ_{Y} as a section of a cone between two planes normal to the BH spin axis. In Sect. 3.3.2 we compute the images of the synchrotron emission of the particles. We discuss in Sect. 4 the possible links of this geometry with the lamppost model (Ross & Fabian 2005) and compare the physical values of τ_{Y} to the typical duration between flares in Cyg X1, Sgr A^{*} and M 87^{*}.
3.2.2. Energy budget
Magnetic reconnection at the Y point and beyond emphasizes the importance of particle acceleration in these simulations. For the most dynamic simulation, where a = 0.99, we monitor the evolution of the fluxes of electromagnetic and particle kinetic energy, respectively given by
$$\begin{array}{c}\hfill {L}_{\mathrm{EM}}^{\infty}={\displaystyle \int {\int}_{\mathrm{\Sigma}}}\mathit{S}\xb7\mathrm{d}\mathbf{\Sigma}\end{array}$$(25)
and
$$\begin{array}{c}\hfill {L}_{K}^{\infty}={\displaystyle \int {\int}_{\mathrm{\Sigma}}}n{e}_{\infty}\mathit{v}\xb7\mathrm{d}\mathbf{\Sigma},\end{array}$$(26)
with e_{∞} = −u_{t} the energy at infinity of the electrons and positrons. Since there are nonzero longitudinal components of S and v on the wedge surface just above the disk, we need to consider an integration surface Σ that covers two sides in order to account for (i) the fluxes through a fraction of a sphere at constant r, from the pole up to θ_{d}, and (ii) the fluxes through the disk surface at constant θ = θ_{d} from the event horizon up to a given r. The fluxes through the pole tend to zero by symmetry and due to the vanishing surface element. We compute the sum of the fluxes through these 2 surfaces for each r from r_{h} up to r = 26 r_{g}, and for each time snapshot between t = 110 r_{g}/c and t = 150 r_{g}/c, once the initial conditions have relaxed. The results are presented in Fig. 10, which illustrates how the electromagnetic energy is progressively transferred to the particles.
Fig. 10. Dissipation of electromagnetic energy and particle acceleration. Left column: electromagnetic (top panel) and particle kinetic (bottom panel) energy flux through a surface extending up to a radius r (left axis) as a function of time and for a thin disk with a BH spin a = 0.99. Right column: slices at constant time t = 150 r_{g}/c (top panel) and constant radius r = 26 r_{g} (bottom panel) of the fluxes. The dashed lines in the upper panel are the timeaveraged profiles. The solid, dashed, and dotted yellow lines show the event horizon, the distance to the separatrix footpoint, and the distance to the Y point, respectively. 
The spacetime diagrams in the left column show these fluxes for ${L}_{\mathrm{EM}}^{\infty}$ (top panel) and ${L}_{K}^{\infty}$ (bottom panel). For the sake of visualization, we provide two slices in the right column: the radial profiles of the fluxes at t = 150 r_{g}/c (upper right panel) and the fluxes through the integration surface Σ at r = 26 r_{g} as a function of time (lower right panel). In the upper right panel, we also display with dashed lines the timeaveraged radial profiles in order to smooth out the impact of the plasmoids. The solid, dashed and dotted yellow lines are the event horizon, the distance to the separatrix footpoint and the distance to the Y point, respectively. Near the event horizon, the electromagnetic flux is $\sim 9\times {10}^{3}{B}_{0}^{2}{r}_{\text{g}}^{2}c$, which corresponds to approximately half^{1} of the sum between the jet power and the amount of electromagnetic energy deposited by the coupling magnetic field lines per unit time (see Fig. 5). Since this fraction is accounted for due to the integration surface we consider, there is no significant loss of electromagnetic energy before we reach the separatrix footpoint (dashed line). Then, the particles falling onto the disk from the Y point along the separatrix experience a net acceleration due to the electric field parallel to their motion. It explains the slight decrease (respectively increase) in the electromagnetic flux (respectively the particle kinetic flux) visible in the radial profiles between the separatrix footpoint (yellow dashed line) and the Y point (yellow dotted line). Yet, the main dissipation occurs in the current sheet beyond the Y point. The net conversion of electromagnetic energy into particle kinetic energy manifests as a steady decrease in the timeaveraged radial profiles. Interestingly enough, we find that in the plasmoids, the electromagnetic and kinetic fluxes are locally enhanced due to the higher plasma densities, currents and electromagnetic energy they carry. The motion of the plasmoids propagating outward can be followed in the spacetime diagrams. They produce stripes with an initially low bulk Lorentz factor, which progressively increases up to Γ ∼ 2 (i.e., β ∼ 0.75), corresponding to the asymptotic slope of the stripes. In the lower left panel, we see the impact on the plasmoids crossing the r = 26 r_{g} surface and leading to transient increases in the fluxes over a few r_{g}/c.
It must be noticed that in the immediate vicinity of the event horizon, we measure negative net fluxes of particle energy at infinity. This is due to the process, first described by Penrose et al. (1969) and later by Parfrey et al. (2019), that enables energy from a Kerr BH to be extracted when particles of negative energy at infinity cross its event horizon. The effect is largely subdominant though, with ${L}_{K}^{\infty}$ approaching only a tenth to a fifth of its value at r = 26 r_{g}. Most of the energy provided to the particles comes from the electromagnetic field.
3.2.3. Coronal heating
Magnetic reconnection has long be thought to be the main culprit for the heating of the corona above the disk (Galeev et al. 1979). We thus examine the amount of electromagnetic energy dissipated by magnetic reconnection in the current sheet, from the separatrix footpoint up to r = 26 r_{g}. We isolate a narrow stripe of 1 r_{g} thick around the current sheet and compute the volume integral of J ⋅ E within this stripe. While in forcefree regions this quantity is negligible, it becomes significant in the current sheet and leads to the dissipated powers shown in Fig. 11. With the intent of providing a semianalytic interpretation of these values, we confront them to the distance r_{Y} between the BH and the Y point previously given as a function of BH spin in Fig. 9. All the error bars show the 1σ variability of the computed values from a snapshot to another. For information, we plot the corresponding BH spins on the top of the figure.
Fig. 11. Amount of electromagnetic energy dissipated per unit time in the current sheet as a function of the distance of the Y point to the BH (bottom x axis). The top x axis indicates the corresponding BH spin (analog to the black points in Fig. 9). The dashed line stands for the best semianalytic fit (see text). 
All dissipated powers are positive which means that overall, electromagnetic energy is tapped to accelerate particles (as visible also from the trend in the upper right panel in Fig. 10). In order to understand the steep increase in the dissipated power with the BH spin, we can first estimate the amount P_{Y} of electromagnetic energy dissipated per unit time by continuous magnetic reconnection at the Y point. If we assume that the surface through which magnetic reconnection occurs can be written $K{r}_{\text{Y}}^{2}$ with K a multiplicative constant, we have
$$\begin{array}{c}\hfill {P}_{\mathrm{Y}}\sim K{\beta}_{\mathrm{rec}}\frac{{B}_{\mathrm{Y}}^{2}{r}_{\mathrm{Y}}^{2}c}{8\pi},\end{array}$$(27)
with B_{Y} the magnetic field near the Y point. The magnetic reconnection rate β_{rec} does not vary significantly from a = 0.6 to a = 0.99 in the highly magnetized regime we work in. We suppose that ${B}_{\mathrm{Y}}\propto {r}_{\mathrm{Y}}^{p}$ with p > 0 an unknown exponent. We fit P_{Y} as a function of r_{Y} for the fitting parameters p and K and we obtain the dashed line in Fig. 11 which corresponds to sensible values: p ∼ 2.7 ± 0.1 and K ∼ 45 ± 13. Indeed, a purely dipolar decay of the magnetic field would have given p = 3. The lower value hints at the influence of the toroidal component, whose decay with the distance is slower than r^{−3}, in the reconnection process, consistent with the significant twisting of the magnetic field lines in Fig. 4. The fact that K is an order of magnitude higher than unity is probably due to magnetic reconnection in the current sheet beyond the Y point.
We can thus affirm that as the BH spin increases, the Y point gets closer to the BH and the regions probed by the magnetic reconnection sites host a much higher magnetic field, which yields a much stronger dissipation of electromagnetic energy. Since the dissipated power we obtain is similar to the Poynting and particle kinetic energy fluxes computed in Sect. 3.2.2, we conclude that most of the dissipation takes place in the current sheet and that the rest of the simulation space is essentially forcefree.
Although qualitatively valid, the aforementioned interpretation is to be taken with caution since in reality, the dissipated power will strongly depends on the local conditions. It is also why we discarded a detailed discussion on the uncertainties and on the degeneracy between p and K. For instance, P_{Y} will change with the magnetic flux distribution on the disk, which might not be dipolar: the results obtained by Chatterjee et al. (2021) with GRMHD simulations actually suggest a slower decrease of the magnetic field with the radius (∝r^{−1}), which would lower the dependence between P_{Y} and the spin.
3.3. Light emission
3.3.1. Ray tracing and emission process
We now proceed to characterize the appearance and the time variability of the BH magnetospheres we computed so far using a raytracer adapter from the geokerr code (Dexter & Agol 2009) and first used in Crinquand et al. (2021). Once the initial conditions have relaxed, each particle has a certain probability to emit a photon at each time step. The photons are not advanced by GRZeltron but instead, we store the necessary information for raytracing in postprocess among which the location, direction and time at emission, along with the energy of the photon (see Eq. (28) below). When the photons reach r = 1000 r_{g}, we neglect the curvature of spacetime and the raytracing becomes straightforward. We refer the reader to Crinquand et al. (2022) for more information about the generation of synthetic synchrotron images from GRZeltron simulations.
We focus on synchrotron emission so the photons are emitted with a weight, w_{γ}, defined by
$$\begin{array}{c}\hfill {w}_{\gamma}=\frac{2{e}^{4}}{3{m}_{\mathrm{e}}^{2}{c}^{3}}w{\gamma}^{2}{B}_{\perp ,\mathrm{eff}}^{2},\end{array}$$(28)
where w is the weight of the emitting particle, γ is the particle Lorentz factor, and B_{⊥, eff} is the magnitude of the effective magnetic field transverse to the particle motion after correcting for the drift velocity:
$$\begin{array}{c}\hfill {\mathit{B}}_{\perp ,\mathrm{eff}}=\mathit{D}+\frac{{\mathit{v}}_{\mathrm{p}}\times \mathit{B}}{c}\frac{({\mathit{v}}_{\mathrm{p}}\xb7\mathit{D}){\mathit{v}}_{\mathrm{p}}}{{c}^{2}},\end{array}$$(29)
where v_{p} is the particle velocity. Both the corona and the disk are supposed to be optically thin to synchrotron radiation, and as such photons are neither absorbed nor scattered. We account for the photons from the region θ ∈ [π/2; π] by associating to each photon emitted in the simulation space its symmetric with respect to the equatorial plane, with symmetric emission direction. In Fig. 12 we represented the relativistic average synchrotron emission power P_{sync} assuming an isotropic pitchangle distribution. The current sheet displays an enhanced emissivity, especially where plasmoids are located: as is apparent in the upper right panel in Fig. 10, the plasmoids are regions of high plasma density and magnetic energy. On the other hand, the physical accuracy of the contribution from the ambient medium is uncertain since it might be influenced by our injection method. Its relative influence is artificially enhanced in Fig. 12 due to the isotropic pitchangle assumption. In addition, the transverse motion of the particles responsible for synchrotron emission would be significantly damped in environments such as M 87^{*} where synchrotron cooling is significant, which would lower the pitch angle along magnetic field lines and thus the emissivity in the ambient medium (Crinquand et al. 2022). For these reasons, we discard the contribution of the ambient medium to the synchrotron emission and consider exclusively light emission from the current sheet delimited with a white dotted line in Fig. 12. A comparison of the relative flux from both components is beyond the scope of our setup since due to the unrealistically low magnetization we use, the maximum Lorentz factor of the particles is artificially low and the emissivity of the current sheet is underestimated.
Fig. 12. Map of synchrotron emission power for a BH spin a = 0.8, multiplied by r^{2} to highlight the contribution from the current sheet. Photons emitted from the region delimited with with a dotted line are those shown as originating from the current sheet in the synthetic images in Figs. 13 and C.1. 
3.3.2. Synthetic images
We computed synthetic synchrotron images for a variety of BH spin and for different inclinations of the lineofsight with respect to the BH spin axis (see Crinquand et al. 2022, for a detailed analysis of synthetic images). The goal is to determine the apparent morphology of the current sheet and of the plasmoids. To do so, we compute the bolometric flux on a screen at infinity normal to the lineofsight, with 400 × 400 square pixels of equal size and a fieldofview extending from −18 r_{g} to +18 r_{g} in both directions. We limit the shot noise by applying a local median over the closest neighboring pixels. Then, to account for the point spread function of the instruments, we smooth the image with a Gaussian kernel of width one pixel.
Figure 13 is a synthetic image for a BH spin a = 0.8 and an intermediate viewing angle, where the lineofsight is almost aligned with the current sheet visible in Fig. 12. In the innermost regions, we plotted the image of the photon ring (white dashed line) centered on the BH position (white cross). First of all, we identify a bright crescentshaped region near the photon ring. We interpret it as a manifestation of the significant relativistic beaming for this inclination and BH spin. Out of the lensed photon ring, bright arcs are visible with an intensity varying between 1 and 10% of the brighter inner crescent. When observed under different inclinations (see panels in Fig. C.1), it appears that they correspond to plasmoids propagating along the conic current sheet. Due to the moderately relativistic bulk motion of the plasmoids (Γ ∼ 1.2 to 2), the emission is not strongly beamed. The arcshape of the plasmoids is an artifact induced by our axisymmetric assumption. In full threedimensional, the plasmoids would instead appear as tilted flux ropes of finite azimuthal extent moving on the envelope of the coneshaped current sheet (see Sect. 4.1). In addition to these features, for this grazing viewing angle, we get the cumulated photons emitted in the direction of outward motion along the whole coneshaped current sheet. They produce a faint stripe across the whole image, with an intensity that reaches a few percent of the one of the inner crescent. Contrary to the other sharper arcs, this arc is not due to our axisymmetric assumption and it is characteristic of the presence of a coneshaped current sheet around the BH. At high viewing angle inclinations, the emission is dominated by a hourglassshaped region corresponding to the projection of the conic current sheet.
Fig. 13. Synthetic images of synchrotron emission from the current sheet for a simulation with a thin disk in prograde rotation around a Kerr BH with a spin a = 0.8. 
3.4. Retrograde disk
When the directions of the BH and the disk’s angular momentum are opposite, the dynamics of the corona is very different from the case when the disk is in prograde rotation. First of all, the ISCO of a retrograde disk lies at a larger distance from the BH, especially for the a = −0.8 case we study where r_{ISCO} ∼ 8.4 r_{g} (while r_{ISCO} ∼ 2.9 r_{g} for a = 0.8). It limits the extent of the region on the disk where magnetic field loops anchored in the disk could close inside the event horizon. However, the very existence of these loops if the disk is in retrograde rotation is unlikely, even for low absolute spin values. Indeed, no steady coupling magnetic field line can sustain between the disk and the BH since in the ergoregion, the field lines have to rotate in the same direction as the BH. Consequently, the field line is necessarily twisted and a toroidal component quickly grows but contrary to the prograde case, there is no region where the differential rotation is low enough that the field line can remain closed.
Our results are in agreement with this analysis. After the initial relaxation phase, all magnetic field lines are open, except some in the region where r < r_{ISCO} and θ ∈ [θ_{d}; θ_{max}]. We must keep in mind though that in the plunging region, the field lines do not rotate and their physical accuracy is questionable, as explained in Sect. 2.2.4. This numerical limitation would be alleviated in simulations where inward advection of the magnetic field lines’ footpoints is accounted for, beyond the scope of this paper. Elsewhere, the field lines are either anchored in the disk or thread the event horizon. Between these two regions, a current sheet forms but no magnetic reconnection is observed. In Fig. 14 we represent the toroidal component of the magnetic field. The same color scale as in Fig. 2 is used in order to emphasize that in this configuration, H_{ϕ} does not change sign at the current sheet (i.e., where the radial component of the magnetic field reverses). The currents are essentially radial and close in the BH equatorial plane.
Fig. 14. Map of the toroidal component, H_{ϕ}, for a simulation with a = −0.8 (retrograde) and ϵ = 5%. The magnetic field lines are those surrounding the nonreconnecting current sheet and the pair corresponding to the potential at the ISCO. The Keplerian disk is in light green and extends down to r_{ISCO} ∼ 8.4 r_{g}. 
The lack of magnetic reconnection in the current sheet confirms the preponderant role played by the toroidal component of the magnetic field, as suggested in Sect. 3.1.5. In contrast with the prograde configuration, the jump in magnetic field at the current sheet is exclusively produced by the poloidal component B^{∥}, which is much smaller than the toroidal component beyond a few gravitational radii. The direction of the magnetic field lines on both sides of the current sheet are thus too similar to trigger magnetic reconnection. In these conditions, plasmoids, if present, could not have been formed by a purely kinetic mechanism like the one at work in the case of a disk in prograde rotation around the BH. Particles could not be steadily accelerated in the corona by magnetic reconnection and the Y point would not be present. We cannot draw firm conclusions for the retrograde case since in this configuration, the corona could still be activated by the accretion of magnetic field loops of varying polarity (Parfrey et al. 2015). Furthermore, magnetic reconnection could still take place for a retrograde disk but in the equatorial region between the ISCO and the event horizon (Ripperda et al. 2020; Crinquand et al. 2021).
4. Discussion
4.1. Hot spot
GRAVITY Collaboration (2018b) reported on a centroid shift near Sgr A^{*} during two flares that they interpreted as synchrotron emission from a hot spot orbiting around the SMBH. GRAVITY Collaboration (2020) reproduced the astrometric data with a compact hot spot of less than 5 r_{g} located at 9 r_{g} from the SMBH. The origin of the hot spots is still unknown and the uncertainties concerning the orientation of their trajectory with respect to the lineofsight and the very center of the putative orbit prevent us from determining whether it was a structure inside a disk surrounding Sgr A^{*}. In the magnetically coupled configuration considered in this paper, if the BH is surrounded by a disk in prograde rotation, a hot spot could either be produced on the disk or correspond to a plasmoid ejected along the current sheet. In the former case, we identified two mechanisms susceptible to form a hot spot: either a deposit of electromagnetic energy onto the disk between the ISCO and the separatrix footpoint, or a beam of relativistic particles accelerated by magnetic reconnection at the Y point and hitting the disk at the separatrix footpoint.
In Sect. 3.1.4 we show that the electromagnetic energy deposited per unit time via the coupling magnetic field loops is comparable to the jet power. Yuan et al. (2019a) evaluated the impact of this contribution on the radial profile of the disk emission. They showed that for a maximally spinning BH and a separatrix footpoint similar to the one we find selfconsistently, the thermal emissivity is significantly enhanced with respect to the standard model by Novikov & Thorne (1973) up to ∼10 r_{g}, a distance approximately three to four times larger than the location of the separatrix footpoint. Since the electromagnetic energy deposit rates we compute match very well those estimated in the purely forcefree regime, we do not repeat the derivation of thermal emissivity profiles performed in Yuan et al. (2019a).
We can however explore another source of energy deposit onto the disk, the one associated with the particle kinetic energy. In our simulations, we observe an important flow of relativistic particles accelerated by magnetic reconnection at the Y point and flowing down to the separatrix footpoint. In the left panel in Fig. 15, we show in blue the radial profile of the particle kinetic energy deposited per unit time and per radial bin by the particles onto the disk. The solid lines is the median profile over 50 r_{g}/c and the blue shaded region is the 1σ variability. The negative values for the particle kinetic energy at infinity measured near the event horizon (vertical solid line) are a manifestation of the Penrose process, as mentioned in Sect. 3.2.2. Beyond, the profile is close to 0 everywhere except in a region centered on the separatrix footpoint (vertical dashed line). This positive peak corresponds to the particles cascading onto the disk at relativistic speed and we show in the right panel in Fig. 15 the energy deposited per unit time onto the disk obtained from integrating the radial profiles near the separatrix footpoint for each BH spin. After normalization with the amount of electromagnetic energy deposited onto the coupled disk region per unit time, Ė_{BH−disk}, we conclude that this source of energy is approximately two orders of magnitude lower than Ė_{BH−disk} for any spin above 0.6. Although this component seems minor, it must be noticed that it is more localized since energy is essentially deposited at the footpoint of the separatrix.
Fig. 15. Deposit of particle kinetic energy on the disk. Left panel: flux on the disk of particle kinetic energy per unit time per ring of infinitesimal width, dr, as a function of the distance, r, to the BH, for a = 0.8. The solid blue curve is the median profile, and the light blue shaded region is the 1σ variability. Vertical lines are the event horizon (solid) and the separatrix footpoint (dashed). Right panel: flux of particle kinetic energy per unit time near the separatrix footpoint compared to the rate of energy deposited onto the disk, Ė_{BH−disk}, per unit time via the coupling magnetic field lines. 
Last, it is still unclear whether the hot spots observed were moving in a plane containing the SMBH itself. Matsumoto et al. (2020) found that their motion was superKeplerian, which could indicate that they are associated with structures flowing away from an underlying disk. GRAVITY Collaboration (2020) emphasized that their best fit was obtained for near faceon orbits and an outofplane velocity component of 0.15c. Ball et al. (2021) designed a model where synchrotron emission from plasmoids moving on the hull of an outflowing jet reproduces the features associated with the hot spots and in particular the offset between the center of the trajectory of the centroid and the position of the BH. The presence of these plasmoids in GRMHD simulations, produced by magnetic reconnection in current sheets induced by physical (Ripperda et al. 2020) or numerical resistivity (Nathanail et al. 2020), brings strong support to these types of models where the apparent orbital motion of hot spots is the result of their slowly expanding helicoidal trajectory and of our viewing angle. In our PIC framework, these plasmoids are those moving along the current sheet with mildly relativistic speed (v_{d} ∼ 0.2 to 0.5c) and with a rotation speed close from the one of the separatrix footpoint given by Eq. (3), hence the superKeplerian motion. If the BH spin axis is aligned (or antialigned) with the lineofsight, a hot spot formed at the Y point or farther would be at a projected distance x_{Y} ∼ 6 − 10 r_{g} for a BH spin of at least 0.7 (see the dashed blue curve in Fig. 9). For lower spins, the Y point is located too far from the BH and a hot spot could not form within 10 r_{g} via the mechanism we have exhibited. When we compute the outward projected drift Δx over a full orbit, we obtain
$$\begin{array}{c}\hfill \frac{\mathrm{\Delta}x}{{x}_{\mathrm{Y}}}\sim 2\pi \frac{{v}_{\mathrm{d}}}{{v}_{\mathrm{orb}}},\end{array}$$(30)
with i_{cs} the inclination of the current sheet with respect to the equatorial plane and v_{orb} = R_{S}Ω_{K}(r = R_{S}) the azimuthal speed of the hot spot. With the inclinations i_{cs} ∼ 30° we measure, we find a relative shift of two to seven times the initial distance projected on the BH equatorial plane, with lower values for lower BH spins. We notice that these values are marginally consistent with the motion of the main hot spot observed around Sgr A^{*}^{2}, which is first located at a projected distance of ∼8 r_{g} and is found at a projected distance of ∼16 r_{g} after twothird of an apparent orbit. Among the spins we explored, a = 0.7 is thus the favored value provided the spin axis is aligned with the lineofsight. We notice that for this spin, the amount of electromagnetic energy deposited onto the disk via the coupling magnetic field lines per unit time is ${\dot{E}}_{\mathrm{BH}\mathrm{disk}}\sim 7\times {10}^{4}{B}_{0}^{2}{r}_{\mathrm{g}}^{2}c$ (see bottom panel in Fig. 5), comparable to the dissipated electromagnetic power in the current sheet, ${P}_{\text{Y}}\sim 3\times {10}^{4}{B}_{0}^{2}{r}_{\text{g}}^{2}c$ (see Fig. 11).
4.2. Flares and illuminating spectrum
Daily nonthermal flares are observed from Sgr A^{*} in NIR (Eckart et al. 2006) and in Xrays (Baganoff et al. 2003). They are associated with a drop in the magnetic field amplitude by almost an order of magnitude (Ponti et al. 2017), which indicates magnetic reconnection as a possible culprit for particle acceleration and subsequent nonthermal emission. In all models, NIR flares comes from optically thin synchrotron emission (Witzel et al. 2021). The origin of Xray flares though is still a matter of debate. Ponti et al. (2017) ascribe it to optically thin synchrotron emission with a cooling break induced by a highenergy cutoff in the distribution of the underlying nonthermal population of electrons. Alternatively, it has been argued in favor of an inverse Compton (DoddsEden et al. 2009) or a synchrotron selfCompton scattering (Mossoux et al. 2016) origin.
4.2.1. Particle energy distribution and scaling correction
The results we presented in Sect. 3.2.3 show that a significant fraction of the electromagnetic energy is tapped to accelerate particles, with dissipated powers up to 10% of the jet power: acceleration at magnetic reconnection sites produces a nonthermal population of relativistic electrons and positrons.
In Fig. 16 we show the energy spectrum of the electrons and positrons computed for a fiducial snapshot in the relaxed state of a simulation with a = 0.8 and ϵ = 5%. We only account for particles inside the current sheet (as defined by the region in Fig. 12). We see that in addition to a thermal peak at γ of a few, corresponding to the bulk motion of the plasma, a nonthermal component emerges and dominates beyond γ ∼ 40. The cutoff at γ_{c} ∼ 200 indicates that particle acceleration is limited. In our setup, this limitation arises from the magnetization of the plasma around the current sheet, upstream the reconnection region, since in the relativistic regime, we expect γ_{c} to be on the order of σ (Kagan et al. 2018; Werner et al. 2016).
Fig. 16. Distribution of the electrons and positrons in the current sheet as a function of their Lorentz factor. 
Between γ = 40 and γ = 200, we fit a power law N(γ)∝γ^{−p} and find that the best fit is achieved for an exponent p ∼ 0.9 (black dashed line in Fig. 16). In the collisionless, relativistic, and strongly magnetized regime, Werner et al. (2016) showed that particles accelerated by magnetic reconnection without a guide field in a symmetric current sheet have an energy distribution that is close to a power law with an exponent of p ∼ 1.2. This proximity with the slope index we measure, in a more realistic environment, suggests that this component we find in the particle energy distribution corresponds to particles accelerated up to relativistic Lorentz factors by magnetic reconnection in the coneshaped current sheet.
An additional argument in favor of the accuracy of this exponent is that although our dimensionless magnetic fields ${\stackrel{\sim}{B}}_{0}$ are strongly underestimated compared to the typical values around accreting BHs, our magnetization σ is of the right order of magnitude. Indeed, in M 87^{*}, the analysis of the polarized synchrotron emission leads to a (nonrelativistic) magnetization of the plasma in the innermost regions between 1 and up to 10^{4} for a magnetic field of 30 G and a plasma density of 10^{4} cm^{−3} (Event Horizon Telescope Collaboration 2021a). In Sgr A^{*}, the multiwavelength campaign carried by GRAVITY Collaboration (2021) indicates that the density of nonthermal electrons is on the order of 10^{6} cm^{−3} for an ambient magnetic field comparable to M 87^{*}. This would correspond to a magnetization on the order of 10−100. Finally, in Cyg X1, the estimated magnetization ranges from 10^{3} to 10^{5} for a plasma density of 10^{12} cm^{−3} and a magnetic field between 10^{5} and 10^{6} G (Cangemi et al. 2021). These values are coherent with the maximum particle Lorentz factors. In Cyg X1, M 87^{*} and Sgr A^{*}, the Lorentz factors of nonthermal particles responsible for NIR synchrotron emission are generally considered to be ≲10^{3} (Cangemi et al. 2021; Event Horizon Telescope Collaboration 2021a). In Sgr A^{*}, GRAVITY Collaboration (2021) found that in order to account for the spectral energy distribution and time evolution of the Xray and NIR flares, sustained particles acceleration was necessary with particle Lorentz factors as high as a few times 10^{4}. In our simulations, near the Y point, we reach a plasma magnetization on the order of a few times 200 (accounting for a bulk Lorentz factor of 1.5), which corresponds to the cutoff γ_{c} ∼ 200 we find. It is approximately an order of magnitude smaller than observations. Direct comparison to observed values is difficult though since radiative models are based on a onezone approach, but we conclude that our magnetization lies within the range of possible values for M 87^{*} and Sgr A^{*}, and is somewhat lower than the values in Cyg X1. Consequently, we can assume that the slope p of the particle energy distribution beyond γ ∼ 40 does not depend on the values of ${\stackrel{\sim}{B}}_{0}$ and σ in the regime we explored.
Alternatively, Yuan et al. (2019b) showed that significant amounts of the BH rotational energy could be dissipated per unit time in magnetospheres containing smallscale poloidal magnetic field loops. If magnetic confinement is strong enough, the magnetic field lines funneling the jet are disrupted by the kink instability. They suggest that most of the electromagnetic energy is dissipated through this magnetic reconnection process close to the BH and that it could be a significant heat source for the corona. Another candidate mechanism was found by GRMHD simulations in the disk midplane, where magnetic reconnection happens in current sheets located on the edge of magnetic bubbles episodically formed in MADs (Porth et al. 2019; Dexter et al. 2020; Ressler et al. 2021; Ripperda et al. 2022). This mechanism corresponds to the current sheet found in the plunging region by PIC simulations (Crinquand et al. 2021). Recently Scepi et al. (2022) computed the flaring and quiescent Xray emission produced by these regions and retrieved fluxes and spectral shapes analogous to the observed ones.
4.2.2. Illuminating spectrum
The consequences of magnetic reconnection in the BH corona also manifest in terms of irradiation of the underlying disk. Part of the hard Xray emission produced by the particles accelerated in the current sheet illuminates the disk and produces a typical Xray reflection spectrum characterized by, for instance, fluorescence lines (in particular the Iron Kα line at 6.4 keV Reynolds 2021). Xray reflection spectrum models are used to determine the spin of BHs accreting material from a geometrically thin and optically thick disk. However, in these semiempirical models, the geometry and the location of the emission site above the disk where the hard Xray illuminating photons come from is a degreeoffreedom of the problem, and so is the spectrum of this emission.
We computed the spectral energy distribution of the synchrotron photons emitted by the leptons in the current sheet (Fig. 17). The photon frequency unit is given by the reference synchrotron frequency, ν_{0}:
$$\begin{array}{c}\hfill {\nu}_{0}=\frac{3e{B}_{0}}{4\pi {m}_{\mathrm{e}}c},\end{array}$$(31)
Fig. 17. Spectral energy distribution of synchrotron emission from the current sheet. 
which we can express as a function of the dimensionless magnetic field, ${\stackrel{\sim}{B}}_{0}$:
$$\begin{array}{c}\hfill {\nu}_{0}=\frac{3}{4\pi}{\stackrel{\sim}{B}}_{0}{\left(\frac{{r}_{\mathrm{g}}}{c}\right)}^{1}.\end{array}$$(32)
We retrieve the thermal peak at low energy and a nonthermal component at higher energy. We fitted this component for photon frequencies between ν_{0} and 10ν_{0} with a power law ${F}_{\nu}\propto {\nu}^{{p}^{\prime}}$ with p′ the spectral index. The best fit is obtained for a spectral index p′∼ − 0.04 ± 0.03. This component is approximately what we would expect for the synchrotron emission in a uniform magnetic field from a powerlaw distribution of positrons and electrons energies with an exponent 1 − 2p′ = 1.1 ∼ p, so we can be confident that it is produced by the nonthermal population of particles accelerated by magnetic reconnection that we identified in Sect. 4.2.1. It must be acknowledged that these values differ from the one derived in NIR by Ponti et al. (2017) for Sgr A^{*}. They found a slope of the spectrum that corresponds to a spectral index p′∼ − 0.7 (i.e., a photon index of ∼1.7), coherent with a nearly flat spectrum. One of their model ascribes this emission to a nonthermal population of particles with a powerlaw energy distribution of exponent p ∼ 2.2. The less steep particle energy distribution we find might be due to the limited photon frequency extent of the powerlaw component, over only one decade.
We notice that for magnetic field magnitudes B ∼ 100 G and particle Lorentz factors γ = 100 of Sgr A^{*}, the synchrotron frequency γ^{2}ν_{0} would correspond to infrared according to Eq. (31). We would retrieve this frequency range with Eq. (32) for higher and more realistic values of ${\stackrel{\sim}{B}}_{0}$ than the one we use, which is limited by the necessity to spatially resolve the Larmor radius (see Sect. 2.2.6).
In this section we studied the properties of a simulation with a BH spin a = 0.8 surrounded by a thin disk (ϵ = 5%). However, the more general exploration of the properties of the coronal emission illuminating the disk as a function of the BH spin is now possible. Most importantly, this result provides a physically motivated coneshaped geometry and emission spectrum for Xray reflection spectroscopy models. It reduces the number of degreesoffreedom compared to the more empiric lamppost model where a point source is located above the BH on its spin axis.
4.2.3. Timing
We next wanted to determine whether the episodic detachment of plasmoids at the Y point is susceptible to yield flares with the appropriate timing properties. In Sgr A^{*}, the flares occur on an approximately daily basis and last an hour or so. Using a mass of the SMBH of 4 × 10^{6} M_{⊙}, we find that the timescales reported in the left panel in Fig. 7 is on the order of a few minutes. To our knowledge, there is no variability reported with a recurrence time on the order of a few minutes in Sgr A^{*}. Intrinsic variability of this mechanism is thus too fast to account for the flares from Sgr A^{*}, although it could be the origin of the nonthermal population responsible for the flare emission. The triggering of a flaring episode itself is more likely due to changes in the accretion flow properties. Timedependent forcefree simulations of magnetic loop accretion by Parfrey et al. (2015) show that for loops of a few r_{g} and realistic accretion speeds, a periodicity of a few hundred r_{g}/c (i.e., a few days for Sgr A^{*}) appears. Alternatively, in the GRMHD simulations of the MAD regime performed by Chashkina et al. (2021), the magnetic flux and accretion onto the BH spike approximately every few hundred r_{g}/c. This cycle is due to the formation of a transient current sheet in the BH equatorial plane, similar to the one found by Ripperda et al. (2020) with GRMHD simulations and studied by Crinquand et al. (2021) in the PIC framework. More generally, models (Gutiérrez et al. 2020) and numerical simulations (Dexter et al. 2020; Porth et al. 2021; Chatterjee et al. 2021) suggest that the MAD regime is prone to transient magnetic reconnection events. Consequently, if the mechanism we describe in this paper is active during the flares, it could lead to quasiperiodic variability during Sgr A^{*}’s flares on a timescale of a few minutes. During a flare of M 87^{*}, this timescale would translate into a quasiperiodicity of a few days.
In Cyg X1, the plasma dynamics at the Y point might well be responsible for the millisecond flares reported by Gierliński & Zdziarski (2003). They observed short Xray flares during which the flux increases by an order of magnitude, with stronger flares in the highsoft state when the standard accretion disk extends very close to the BH, probably all the way down the ISCO. Beloborodov (2017) suggested that these flares might come from chains of plasmoids formed near the BH. In the highsoft state, the topology of the magnetosphere could be similar to the one we considered, with a coneshaped current sheet forming above and below the disk. Unfortunately, the very short duration of the flares precluded any conclusive results on the spectrum during the flares, although Gierliński & Zdziarski (2003) suggested that a hybrid emission model with a black body disk component and a nonthermal tail from a population of relativistic electrons is possible. Alternatively, Mehlhaff et al. (2021) proposed that these flares might be due to kinetic beaming based on their PIC simulations of inverse Compton scattering by very highenergy particles in the corona.
4.3. Jet power
In Sect. 3.1.2 we computed the electromagnetic power carried by the Poynting vector through the open field lines threading the BH at its poles. We identified this value to the power of the jet, which, downstream of the region we probe, will gain mass and where electromagnetic energy will eventually be converted into particle kinetic energy and radiation. These values can thus be interpreted as upper limits on the unbeamed jet luminosity, independently of the emission mechanism in the jet. In Fig. 5 we show the maximum jet power we obtained for a BH with spin a = 0.99, approximately
$$\begin{array}{c}\hfill {P}_{\mathrm{jet},99}=0.006{B}_{0}^{2}\phantom{\rule{0.166667em}{0ex}}{r}_{\mathrm{g}}^{2}c.\end{array}$$(33)
It would be misleading to use the values derived from polarization measurements or spectral modeling for B_{0}, which is, in our model, the amplitude of the magnetic field at 1 r_{g}. Indeed, the values derived from onezone models are averaged over extended emission regions where the magnetic field amplitude varies quickly. Instead, we consider the maximal value B_{0} can reach by evaluating the dynamical importance of the magnetic field with the dimensionless ratio ${\mathrm{\Phi}}_{\mathrm{H}}/\sqrt{\dot{M}c{r}_{\mathrm{g}}^{2}}$, with Ṁ the mass accretion rate in the innermost region. When this ratio reaches approximately 50, the magnetic field saturates and the accretion proceeds via a MAD (Tchekhovskoy et al. 2011; Porth et al. 2019). In these conditions, the magnetic field in the innermost regions must reach
$$\begin{array}{c}\hfill {B}_{\mathrm{MAD}}=\frac{25}{\pi}\sqrt{\frac{\dot{M}c}{{r}_{\mathrm{g}}^{2}}},\end{array}$$(34)
which is to be understood as a maximum value. Event Horizon Telescope Collaboration (2021a) performed a broad exploration of the parameter space with GRMHD simulations in an attempt to model the polarized synchrotron emission reported from M 87^{*} in Event Horizon Telescope Collaboration (2021b). They showed that MAD were systematically better at matching the observations and found a mass accretion rate of Ṁ ∼ 3 − 20 × 10^{−4} M_{⊙} yr^{−1}, in agreement with the upper limit derived by Kuo et al. (2014) from Faraday rotation and by Russell et al. (2015) from Xray analysis. With the parameters of M 87^{*}, we find a magnetic field B_{MAD} = 200 − 500 G and, using Eq. (33), a near maximal jet power of P_{jet, 99} ∼ 6 × 10^{42} erg s^{−1} to P_{jet, 99} ∼ 3 × 10^{43} erg s^{−1}, in agreement with the jet kinetic power deduced by Lucchini et al. (2019) from observations.
In Cyg X1, the Xray accretion luminosity is on the order of 1% of the Eddington luminosity in the lowhard state, when jet emission is detected. With a conversion efficiency of 10%, it corresponds to a mass accretion rate Ṁ ∼ 5 × 10^{−9} M_{⊙} yr^{−1}. The corresponding MAD magnetic field amplitude is B_{MAD} ∼ 10^{8} G and the maximal jet power is P_{jet, 99} ∼ 10^{37} erg s^{−1}, approximately one order of magnitude higher than what the timeaveraged estimate derived by Russell et al. (2007) from interaction with the interstellar medium. This discrepancy might be due to the fact that the jet is absent during the highsoft state, although Cyg X1 is in such a state for only 10% of its lifetime (Gallo et al. 2005). Alternatively, it could mean that the magnetic field never reaches the MAD saturation value, contrary to M 87^{*}. A magnetic field amplitude in the jetlaunching region of ∼10^{7} G (Cangemi et al. 2021) would give a jet power consistent with observations.
For Sgr A^{*}, the existence of a jet is more elusive. Li et al. (2013) identified a parsecscale jet based on what they interpreted as a radio shock front due to the jet from Sgr A^{*} crossing the local gas. Diagnostics such as Faraday rotation measurements allow us to estimate the mass accretion rate onto the SMBH, which ranges from 10^{−9} to 10^{−7} M_{⊙} yr^{−1} (Bower et al. 2003; Marrone et al. 2007; Wang et al. 2013). The simulations of stellar wind capture performed by (Ressler et al. 2018) reproduce these values, with an extrapolated mass accretion rate at the event horizon on the order of 10^{−8} M_{⊙} yr^{−1}. It would yield a maximum magnetic field B_{MAD} = 600 − 6000 G, several orders of magnitude above the measured values, and a near maximal jet power P_{jet, 99} ∼ 2 × 10^{37} erg s^{−1} to P_{jet, 99} ∼ 2 × 10^{39} erg s^{−1}. This range is compatible with the values deduced from the ram pressure of the jet on the ambient medium (Li et al. 2013), from emission models (Markoff et al. 2007) and from numerical simulations (Yuan et al. 2012). Even in the most optimistic case though, where Sgr A^{*} would host a near maximally spinning SMBH accreting in a MAD regime, the putative jet for a mass accretion rate of 10^{−8} M_{⊙} yr^{−1} is orders of magnitude less powerful than in M 87^{*}. We also notice that a BlandfordZnajek jet powered by a MAD magnetic field at this mass accretion rate would not be enough to fuel the Fermi bubbles observed in γrays (Su et al. 2010) or the Xray bubbles recently discovered by eROSITA (Predehl et al. 2020).
4.4. Spinning up torques from the BH
We want to evaluate the impact of the deposit of angular momentum from the BH onto a disk in prograde rotation. We do so by extracting a characteristic timescale, τ_{L}, over which the angular momentum of the coupled section of the disk, between the ISCO and the footpoint of the separatrix, changes significantly due to the angular momentum supplied by the BH via the coupling magnetic field lines. An estimate for the angular momentum L_{disk} of the coupled disk section is given by
$$\begin{array}{c}\hfill {L}_{\mathrm{disk}}={\displaystyle {\int}_{{r}_{\mathrm{ISCO}}}^{{R}_{\mathrm{s}}}{r}^{2}{\mathrm{\Omega}}_{\mathrm{K}}\mathrm{d}m,}\end{array}$$(35)
with Ω_{K} the angular velocity given by Eq. (3) and dm the amount of mass in an infinitely thin annulus in the disk. Using mass conservation, we can express dm as a function of the mass accretion rate, Ṁ, and the accretion speed, v_{r}, and we get
$$\begin{array}{c}\hfill {L}_{\mathrm{disk}}=a\frac{\dot{M}c{r}_{\mathrm{g}}^{2}}{{v}_{\mathrm{r}}}{\displaystyle {\int}_{{x}_{\mathrm{ISCO}}}^{{x}_{\mathrm{s}}}\frac{{x}^{2}\mathrm{d}x}{1+{x}^{3/2}}\phantom{\rule{1em}{0ex}}\mathrm{with}\phantom{\rule{1em}{0ex}}x={(r/{r}_{\mathrm{g}})}^{3/2}/a.}\end{array}$$(36)
For convenience, we write this integral I. On the other hand, we measured ${\dot{L}}_{\mathrm{BH}\mathrm{disk}}$ that we write $\lambda {B}_{0}^{2}{r}_{\text{g}}^{3}$, with λ a number dependent on the spin and taken from the bottom plot in Fig. 5. The characteristic spinup timescale τ_{L} of the disk is thus
$$\begin{array}{c}\hfill {\tau}_{\mathrm{L}}=\frac{{L}_{\mathrm{disk}}}{{\dot{L}}_{\mathrm{BH}\mathrm{disk}}}=\frac{\mathit{ac}}{\lambda {v}_{\mathrm{r}}}\frac{{r}_{\mathrm{g}}}{c}\frac{\dot{M}{c}^{2}}{{B}_{0}^{2}{r}_{\mathrm{g}}^{2}c}I.\end{array}$$(37)
Assuming the BH accretes in the MAD regime, in the innermost regions, the magnetic field is maximal and given by Eq. (34). Then, the characteristic timescale τ_{L} is given by
$$\begin{array}{c}\hfill {\tau}_{\mathrm{L}}={\left(\frac{\pi}{25}\right)}^{2}\frac{\mathit{aI}}{\lambda}\frac{c}{{v}_{\mathrm{r}}}\frac{{r}_{\mathrm{g}}}{c}\xb7\end{array}$$(38)
We plotted this timescale as a function of the BH spin from 0.6 to 0.99 in Fig. 18 for a representative accretion speed v_{r} = c/200 (JacqueminIde et al. 2021). They are always longer than the duration of our simulations but it turns out that τ_{L} can be as short as a few hundred r_{g}/c for a = 0.99. In comparison, Chashkina et al. (2021) suggest that to activate the corona and trigger major reconnection events, the accretion of successive loops of opposite polarity and radial extent larger than 10 r_{g} is required. Such loops would need a couple thousand r_{g}/c to be fully accreted once the inner parts of the loop cross the horizon. During this time lapse, feedback through angular momentum is significant if the system is in the MAD regime and if the BH has a spin a ≳ 0.9. In this case, tremendous consequences on the structure of the innermost regions of the disk are to be expected. This mechanism would manifest on timescales of a few hundred r_{g}/c, which would be compatible with hourlong flares in Sgr A^{*}. Whether it could play a role in the flare dynamics deserves further investigation. Otherwise, if the system is not in the MAD regime or if the BH spin is lower than 0.9, the magnetic coupling between the disk and the BH is not sustained for an amount of time long enough to significantly modify the structure of the innermost regions of the disk via angular momentum deposit.
Fig. 18. Characteristic timescale of angular momentum deposited onto the coupled part of the disk, τ_{L}, as a function of the BH spin. 
5. Summary and perspectives
In this paper we have performed global PIC simulations of a collisionless and strongly magnetized electron–positron pair plasma embedded in a magnetic field sustained by a steady, aligned, and perfectly conducting disk in Keplerian rotation around a Kerr BH. We have shown that in this environment, often called the corona or the BH magnetosphere, magnetic field loops coupling the disk to the BH are maintained up to a maximal distance set by the BH spin, provided the disk rotation is prograde. Beyond the separatrix (i.e., the outermost closed magnetic field line), magnetic field lines open due to the strong shearing between the rotation speed at their footpoint on the disk and the shearing induced in the ergosphere by the LenseThirring dragging. For higher spins, the coupling region is smaller. The open magnetic field lines are twisted and are either anchored in the disk or thread the event horizon. The former ones are very inclined, which could favor the launching of a magnetocentrifugal disk outflow. The latter ones funnel an electromagnetic jet whose power dependence on the BH spin and magnetic flux is in agreement with the forcefree estimates of BlandfordZnajek jets by Tchekhovskoy et al. (2011), for spins ranging from 0.6 up to 0.99. Provided we rely on the magnetic field of a MAD disk rather than the one derived from onezone models, we reproduce the observed jet power in M 87^{*}. The two regions of open field lines are separated by a coneshaped current sheet, where the magnetic field reconnects and where plasmoids form and flow away from the Y point, the outermost point on the separatrix (which would be a Y ring in three dimensions). We find a quasiperiodic motion of the separatrix, which progressively stretches until a plasmoid detaches, with typical characteristic timescales on the order of 20 to 5 r_{g}/c for BH spins from 0.6 to 0.99. For the accreting stellarmass BH in Cyg X1, this quasiperiodicity corresponds to a few milliseconds. These episodes of sudden magnetic dissipation could thus be responsible for the millisecond flares observed by Gierliński & Zdziarski (2003) during the highsoft accretion states.
We identified a highly relativistic component in the energy distribution of particles from the current sheet. Its slope is close to −1, which indicates that these particles have been accelerated by magnetic reconnection. We retrieve the corresponding flat synchrotron component in the spectral energy distribution. This mechanism can qualitatively reproduce part of the nonthermal emission from accreting BHs, where the cooling timescale is long compared to the light crossing time (e.g., Sgr A^{*}). It could help to interpret their spectral properties, although the estimated dissipated powers and maximum Lorentz factor of the particles we derive might be underestimated by one to two orders of magnitude. This work may also provide an answer to the question of the origin of the coronal heating: in the current sheet, electromagnetic energy is tapped to continuously feed a relativistic and nonthermal population of electrons and positrons responsible for the upscattering of Xray photons from the disk up to γray energies. Thanks to the plausible magnetization we manage to reach in the corona, we are confident that these results provide reliable locations of the emission sites and photometric and spectral properties as a function of the BH spin: the higher the spin, the closer from the BH the Y ring is and the higher the dissipated electromagnetic power to fuel nonthermal emission. This physically motivated coronal emission could be used in Xray reflection spectroscopy models instead of semiempirical lamppost geometry, where the height of the point source above the BH is a degreeoffreedom. In our simulations, the current sheet coexists with the jetlaunching region that it engulfs, which leaves the door open for a contribution of a compact radio jet to the NIR excess observed in Sgr A^{*}, for instance (Moscibrodzka & Falcke 2013).
Finally, we computed synthetic images of synchrotron emission within ∼10 r_{g} around the BH and characterized the features due to the current sheet and the plasmoids it contains, for different viewing inclinations. If the mechanism we studied in this paper is at play around M 87^{*}, plasmoids could manifest in the images as transient and fainter sources of emission moving at mildly relativistic bulk speeds. These results also suggest that the hot spot reported by GRAVITY Collaboration (2018b) around Sgr A^{*} could be due to synchrotron emission from a large plasmoid propagating outward in the coneshaped current sheet along a helicoidal trajectory.
Coupling magnetic field line loops are anchored in the disk and close within the event horizon. We find that for the spins we explored (a > 0.6), they deposit electromagnetic energy onto the disk surface, between the ISCO and the footpoint of the separatrix, at a rate such that the associated power is comparable to the jet power. They also enable the extraction of angular momentum from the BH. For a > 0.9, this deposit of angular momentum in the innermost regions of the disk could dramatically alter its structure over a few hundred r_{g}/c. Relativistic particles accelerated at the Y ring are also a source of energy for the disk as they cascade onto the footpoint of the separatrix, whose location depends on the BH spin.
Our results remains qualitatively similar for thin (ϵ = 5%) and slim (ϵ = 30%) disks in prograde rotation, although magnetic reconnection is less vivid for thicker disks. In contrast, we find that for retrograde rotation (i.e., negative effective BH spin), the corona is not magnetically active. Reconnection is inhibited, and no coupling magnetic field line can be maintained. However, this conclusion strongly depends on the distribution of the magnetic flux provided by the disk and on the speed at which magnetic field lines are advected. It highlights the main limitation of our model: the disk is a fixed background whose dynamics we do not solve. Instead, we assume that over the few hundred r_{g}/c of our simulations, the disk and the magnetic field it brings are frozen. In order to capture the longer timescales of accreting BHs, the computationally expensive kinetic approach we followed will have to be coupled to GRMHD simulations, which are much more suitable to describe the dynamics of the disk (Ripperda et al. 2020).
Another aspect of our model to be questioned is the physical motivation for the magnetic topology we considered. Smallscale magnetic loops have been suggested to be a frequent outcome of the magnetorotational instability (MRI; Balbus & Hawley 1991). Indeed, the MRI fuels turbulence and leads to buoyantly emerging magnetic field loops that rise above the disk surface, as shown in shearing box and global MHD simulations (Davis et al. 2010; Zhu & Stone 2018). This confirms the semianalytic stochastic picture envisioned by Galeev et al. (1979) and Uzdensky & Goodman (2008). Depending on the extent of the loop and on its propensity to reconnect with neighboring loops, accretion could bring the innermost part of the loop into the event horizon while the outermost footpoint remains anchored in the disk. In our model and for the sake of simplicity, we only considered a single loop, by starting with an initially dipolar distribution of the magnetic field on the disk. The main motivation was to explore the dynamics of this configuration on timescales short enough to neglect the accretion of such a large loop, which would typically take a few thousand r_{g}/c, much longer than the timescale of our simulations, a few hundred r_{g}/c. In future PIC simulations, we intend to consider the effect of successive smallerscale loops on the corona but also on the formation of current sheets in a striped jet, analog to the general relativistic forcefree simulations of Parfrey et al. (2015) and Mahlmann et al. (2020).
We could not draw conclusions from synthetic light curves due to the axisymmetric assumption our twodimensional model relies on. This model artificially enhances the peaktopeak variability while lengthening peak duration due to the gravitational delay induced by the curved spacetime. Timing information is also flawed by the hypothesis we make about the symmetry of the corona with respect to the equatorial plane. For these reasons, we postpone an analysis of the synthetic light curves to future work, where the corona will be studied with fully threedimensional simulations. This will allow the possibility of studying the fragmentation of the current sheet into plasmoids of finite azimuthal extent and investigating whether magnetic reconnection can propagate along the Y ring and mimic an orbiting hot spot.
Acknowledgments
The authors thank the referee for the positive and enthusiastic feedback. This project has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No. 863412). Computing resources were provided by TGCC and CINES under the allocation A0090407669 made by GENCI. The authors thank Maïca Clavel, Antoine Strugarek, Guillaume Dubus and Geoffroy Lesur for fruitful discussions and constructive feedback.
References
 Algaba, J. C., Anczarski, J., Asada, K., et al. 2021, ApJ, 911, L11 [NASA ADS] [CrossRef] [Google Scholar]
 Arnowitt, R. L., Deser, S., & Misner, C. W. 1962, Recent Dev. Gen. Relativ., 127 [Google Scholar]
 Bacchini, F., Ripperda, B., Chen, A. Y., & Sironi, L. 2018, ApJS, 237, 6 [Google Scholar]
 Bacchini, F., Ripperda, B., Porth, O., & Sironi, L. 2019, ApJS, 240, 40 [Google Scholar]
 Baganoff, F. K., Maeda, Y., Morris, M., et al. 2003, ApJ, 591, 891 [Google Scholar]
 Balbus, S. A., & Hawley, J. F. 1991, ApJ, 376, 214 [Google Scholar]
 Ball, D., Özel, F., Christian, P., Chan, C.K., & Psaltis, D. 2021, ApJ, 917, 8 [NASA ADS] [CrossRef] [Google Scholar]
 Bambi, C., Brenneman, L. W., Dauser, T., et al. 2021, Space Sci. Rev., 217, 65 [CrossRef] [Google Scholar]
 Barkov, M. V., & Komissarov, S. S. 2016, MNRAS, 458, 1939 [NASA ADS] [CrossRef] [Google Scholar]
 Beloborodov, A. M. 2017, ApJ, 850, 141 [NASA ADS] [CrossRef] [Google Scholar]
 Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883 [CrossRef] [Google Scholar]
 Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433 [NASA ADS] [CrossRef] [Google Scholar]
 Bower, G. C., Wright, M. C. H., Falcke, H., & Backer, D. C. 2003, ApJ, 588, 331 [Google Scholar]
 Bransgrove, A., Ripperda, B., & Philippov, A. 2021, Phys. Rev. Lett., 127, 055101 [NASA ADS] [CrossRef] [Google Scholar]
 Cangemi, F., Beuchert, T., Siegert, T., et al. 2021, A&A, 650, A93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Cerutti, B., Werner, G. R., Uzdensky, D. A., & Begelman, M. C. 2013, ApJ, 770, 147 [Google Scholar]
 Cerutti, B., Philippov, A., Parfrey, K., & Spitkovsky, A. 2015, MNRAS, 448, 606 [Google Scholar]
 Chainakun, P., Watcharangkool, A., Young, A. J., & Hancock, S. 2019, MNRAS, 487, 667 [NASA ADS] [CrossRef] [Google Scholar]
 Chashkina, A., Bromberg, O., & Levinson, A. 2021, MNRAS, 508, 1241 [NASA ADS] [CrossRef] [Google Scholar]
 Chatterjee, K., Markoff, S., Neilsen, J., et al. 2021, MNRAS, 507, 5281 [NASA ADS] [CrossRef] [Google Scholar]
 Crinquand, B., Cerutti, B., Philippov, A., Parfrey, K., & Dubus, G. 2020, Phys. Rev. Lett., 124, 145101 [Google Scholar]
 Crinquand, B., Cerutti, B., Dubus, G., Parfrey, K., & Philippov, A. 2021, A&A, 650, A163 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Crinquand, B., Cerutti, B., Dubus, G., Parfrey, K., & Philippov, A. A. 2022, ArXiv eprints [arXiv:2202.04472] [Google Scholar]
 Davis, S. W., Stone, J. M., & Pessah, M. E. 2010, ApJ, 713, 52 [NASA ADS] [CrossRef] [Google Scholar]
 Del Santo, M., Malzac, J., Belmont, R., Bouchet, L., & De Cesare, G. 2013, MNRAS, 430, 209 [NASA ADS] [CrossRef] [Google Scholar]
 Dexter, J., & Agol, E. 2009, ApJ, 696, 1616 [Google Scholar]
 Dexter, J., Agol, E., Fragile, P. C., & McKinney, J. C. 2010, ApJ, 717, 1092 [NASA ADS] [CrossRef] [Google Scholar]
 Dexter, J., Tchekhovskoy, A., JiménezRosales, A., et al. 2020, MNRAS, 497, 4999 [Google Scholar]
 DoddsEden, K., Porquet, D., Trap, G., et al. 2009, ApJ, 698, 676 [NASA ADS] [CrossRef] [Google Scholar]
 Eckart, A., Baganoff, F. K., Schödel, R., et al. 2006, A&A, 450, 535 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 El Mellah, I., Sander, A. A. C., Sundqvist, J. O., & Keppens, R. 2019, A&A, 622, A189 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Event Horizon Telescope Collaboration (Akiyama, K., et al.) 2019, ApJ, 875, L1 [Google Scholar]
 Event Horizon Telescope Collaboration (Akiyama, K., et al.) 2021a, ApJ, 910, L13 [Google Scholar]
 Event Horizon Telescope Collaboration (Akiyama, K., et al.) 2021b, ApJ, 910, L12 [Google Scholar]
 Frank, J., King, A., & Raine, D. J. 1986, Phys. Today, 39, 124 [NASA ADS] [CrossRef] [Google Scholar]
 Fromang, S., Lyra, W., & Masset, F. 2011, A&A, 534, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Galeev, A. A., Rosner, R., & Vaiana, G. S. 1979, ApJ, 229, 318 [NASA ADS] [CrossRef] [Google Scholar]
 Gallo, E., Fender, R., Kaiser, C., et al. 2005, Nature, 436, 819 [Google Scholar]
 Gierliński, M., & Zdziarski, A. A. 2003, MNRAS, 343, 84 [Google Scholar]
 Goldreich, P., & Julian, W. 1969, ApJ, 157, 9 [NASA ADS] [Google Scholar]
 Gralla, S. E., & Jacobson, T. 2014, MNRAS, 445, 2500 [NASA ADS] [CrossRef] [Google Scholar]
 GRAVITY Collaboration (Abuter, R., et al.) 2018a, A&A, 615, L15 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 GRAVITY Collaboration (Abuter, R., et al.) 2018b, A&A, 618, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 GRAVITY Collaboration (Bauböck, M., et al.) 2020, A&A, 635, A143 [CrossRef] [EDP Sciences] [Google Scholar]
 GRAVITY Collaboration (Abuter, R., et al.) 2021, A&A, 654, A22 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Grinberg, V., Hell, N., Pottschmidt, K., et al. 2013, A&A, 554, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gutiérrez, E. M., Nemmen, R., & Cafardo, F. 2020, ApJ, 891, L36 [Google Scholar]
 JacqueminIde, J., Lesur, G., & Ferreira, J. 2021, A&A, 647, A192 [EDP Sciences] [Google Scholar]
 Kagan, D., Sironi, L., Cerutti, B., et al. 2015, Space Sci. Rev., 191, 545 [NASA ADS] [CrossRef] [Google Scholar]
 Kagan, D., Nakar, E., & Piran, T. 2018, MNRAS, 476, 3902 [NASA ADS] [CrossRef] [Google Scholar]
 Kara, E., Steiner, J. F., Fabian, A. C., et al. 2019, Nature, 565, 198 [Google Scholar]
 Kerr, R. P. 1963, Phys. Rev. Lett., 11, 237 [Google Scholar]
 King, A. R., & Pringle, J. E. 2021, ApJ, 918, L22 [NASA ADS] [CrossRef] [Google Scholar]
 Komissarov, S. S. 2004, MNRAS, 350, 427 [Google Scholar]
 Komissarov, S. S. 2022, MNRAS, 512, 2798 [NASA ADS] [CrossRef] [Google Scholar]
 Kuo, C. Y., Asada, K., Rao, R., et al. 2014, ApJ, 783, L33 [NASA ADS] [CrossRef] [Google Scholar]
 Laurent, P., Rodriguez, J., Wilms, J., et al. 2011, Science, 332, 438 [Google Scholar]
 Li, Z., Morris, M. R., & Baganoff, F. K. 2013, ApJ, 779, 154 [NASA ADS] [CrossRef] [Google Scholar]
 Loureiro, N. F., Schekochihin, A. A., College, K., & Cowley, S. C. 2007, Phys. Plasmas, 14, 100703 [NASA ADS] [CrossRef] [Google Scholar]
 Lucchini, M., Krauß, F., & Markoff, S. 2019, MNRAS, 489, 1633 [NASA ADS] [Google Scholar]
 MacDonald, D., & Thorne, K. S. 1982, MNRAS, 198, 345 [Google Scholar]
 Mahlmann, J. F., CerdáDurán, P., & Aloy, M. A. 2018, MNRAS, 477, 3927 [Google Scholar]
 Mahlmann, J. F., Levinson, A., & Aloy, M. A. 2020, MNRAS, 494, 4203 [Google Scholar]
 Malzac, J., & Belmont, R. 2009, MNRAS, 392, 570 [CrossRef] [Google Scholar]
 Markoff, S., Bower, G. C., & Falcke, H. 2007, MNRAS, 379, 1519 [NASA ADS] [CrossRef] [Google Scholar]
 Marrone, D. P., Moran, J. M., Zhao, J.H., & Rao, R. 2007, ApJ, 654, L57 [NASA ADS] [CrossRef] [Google Scholar]
 Matsumoto, T., Chan, C. H., & Piran, T. 2020, MNRAS, 497, 2385 [NASA ADS] [CrossRef] [Google Scholar]
 Mbarek, R., Haggerty, C., Sironi, L., Shay, M., & Caprioli, D. 2021, Phys. Rev. Lett., 128, 2022 [Google Scholar]
 Mehlhaff, J. M., Werner, G. R., Uzdensky, D. A., & Begelman, M. C. 2020, MNRAS, 498, 799 [NASA ADS] [CrossRef] [Google Scholar]
 Mehlhaff, J. M., Werner, G. R., Uzdensky, D. A., & Begelman, M. C. 2021, MNRAS, 508, 4532 [NASA ADS] [CrossRef] [Google Scholar]
 MillerJones, J. C. A., Bahramian, A., Orosz, J. A., et al. 2021, Science, 371, 1046 [Google Scholar]
 Moscibrodzka, M., & Falcke, H. 2013, A&A, 559, L3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mościbrodzka, M., Gammie, C. F., Dolence, J. C., Shiokawa, H., & Leung, P. K. 2009, ApJ, 706, 497 [Google Scholar]
 Mossoux, E., Grosso, N., Bushouse, H., et al. 2016, A&A, 589, A116 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Nathanail, A., & Contopoulos, I. 2014, ApJ, 788, 186 [Google Scholar]
 Nathanail, A., Fromm, C. M., Porth, O., et al. 2020, MNRAS, 495, 1549 [Google Scholar]
 Novikov, I. D., & Thorne, K. S. 1973, in Black Holes, eds. C. DeWitt, & B. DeWitt (New York: Gordon Breach), 343 [Google Scholar]
 Parfrey, K., Giannios, D., & Beloborodov, A. M. 2015, MNRAS, 446, L61 [Google Scholar]
 Parfrey, K., Spitkovsky, A., & Beloborodov, A. M. 2017, MNRAS, 469, 3656 [Google Scholar]
 Parfrey, K., Philippov, A., & Cerutti, B. 2019, Phys. Rev. Lett., 122, 35101 [NASA ADS] [CrossRef] [Google Scholar]
 Penrose, R., Królak, A., & Israel, W. 1969, Gen. Relativ. Gravit., 34, 1135 [Google Scholar]
 Ponti, G., George, E., Scaringi, S., et al. 2017, MNRAS, 468, 2447 [NASA ADS] [CrossRef] [Google Scholar]
 Porth, O., Chatterjee, K., Narayan, R., et al. 2019, ApJS, 243, 26 [Google Scholar]
 Porth, O., Mizuno, Y., Younsi, Z., & Fromm, C. M. 2021, MNRAS, 502, 2023 [NASA ADS] [CrossRef] [Google Scholar]
 Potter, W. J. 2021, MNRAS, 503, 5025 [NASA ADS] [CrossRef] [Google Scholar]
 Poutanen, J., & Vurm, I. 2009, ApJ, 690, L97 [NASA ADS] [CrossRef] [Google Scholar]
 Predehl, P., Sunyaev, R. A., Becker, W., et al. 2020, Nature, 588, 227 [CrossRef] [Google Scholar]
 Remillard, R. A., & McClintock, J. E. 2006, ARA&A, 44, 49 [Google Scholar]
 Ressler, S. M., Quataert, E., & Stone, J. M. 2018, MNRAS, 478, 3544 [NASA ADS] [CrossRef] [Google Scholar]
 Ressler, S. M., White, C. J., Quataert, E., & Stone, J. M. 2020, ApJ, 896, L6 [NASA ADS] [CrossRef] [Google Scholar]
 Ressler, S. M., Quataert, E., White, C. J., & Blaes, O. 2021, MNRAS, 504, 6076 [NASA ADS] [CrossRef] [Google Scholar]
 Reynolds, C. S. 2021, ARA&A, 59, 117 [NASA ADS] [CrossRef] [Google Scholar]
 Ripperda, B., Bacchini, F., & Philippov, A. A. 2020, ApJ, 900, 100 [NASA ADS] [CrossRef] [Google Scholar]
 Ripperda, B., Liska, M., Chatterjee, K., et al. 2022, ApJ, 924, L32 [NASA ADS] [CrossRef] [Google Scholar]
 Rodriguez, J., Grinberg, V., Laurent, P., et al. 2015, ApJ, 807, 17 [NASA ADS] [CrossRef] [Google Scholar]
 Ross, R. R., & Fabian, A. C. 2005, MNRAS, 358, 211 [CrossRef] [Google Scholar]
 Rushton, A., MillerJones, J. C., Campana, R., et al. 2012, MNRAS, 419, 3194 [NASA ADS] [CrossRef] [Google Scholar]
 Russell, D. M., Fender, R. P., Gallo, E., & Kaiser, C. R. 2007, MNRAS, 376, 1341 [Google Scholar]
 Russell, H. R., Fabian, A. C., McNamara, B. R., & Broderick, A. E. 2015, MNRAS, 451, 588 [NASA ADS] [CrossRef] [Google Scholar]
 Scepi, N., Dexter, J., & Begelman, M. C. 2022, MNRAS, 511, 3536 [NASA ADS] [CrossRef] [Google Scholar]
 Sen, K., Xu, X.T., Langer, N., et al. 2021, A&A, 652, A138 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337 [NASA ADS] [Google Scholar]
 Shapiro, S. L., & Teukolsky, S. A. 1983, Black Holes, White Dwarfs, and Neutron Stars: The Physics of Compact Objects (New York: Wiley), 645 [Google Scholar]
 Shreeram, S., & Ingram, A. 2020, MNRAS, 492, 405 [NASA ADS] [CrossRef] [Google Scholar]
 Sironi, L., & Spitkovsky, A. 2014, ApJ, 783, L21 [NASA ADS] [CrossRef] [Google Scholar]
 Sironi, L., Rowan, M. E., & Narayan, R. 2021, ApJ, 907, L44 [CrossRef] [Google Scholar]
 Sridhar, N., Sironi, L., & Beloborodov, A. M. 2021, MNRAS, 17, 1 [Google Scholar]
 Su, M., Slatyer, T. R., & Finkbeiner, D. P. 2010, ApJ, 724, 1044 [Google Scholar]
 Tchekhovskoy, A., Narayan, R., & Mckinney, J. C. 2011, MNRAS, 418, L79 [NASA ADS] [CrossRef] [Google Scholar]
 Uzdensky, D. A. 2002, ApJ, 574, 1011 [NASA ADS] [CrossRef] [Google Scholar]
 Uzdensky, D. A. 2004, ApJ, 603, 652 [NASA ADS] [CrossRef] [Google Scholar]
 Uzdensky, D. A. 2005, ApJ, 620, 889 [NASA ADS] [CrossRef] [Google Scholar]
 Uzdensky, D. A., & Goodman, J. 2008, ApJ, 682, 608 [NASA ADS] [CrossRef] [Google Scholar]
 Wald, R. M. 1974, Phys. Rev. D, 10, 1680 [Google Scholar]
 Wang, Q. D., Nowak, M. A., Markoff, S. B., et al. 2013, Science, 341, 981 [NASA ADS] [CrossRef] [Google Scholar]
 Werner, G. R., Uzdensky, D. A., Cerutti, B., Nalewajko, K., & Begelman, M. C. 2016, ApJ, 816, L8 [Google Scholar]
 Werner, G. R., Uzdensky, D. A., Begelman, M. C., Cerutti, B., & Nalewajko, K. 2018, MNRAS, 473, 4840 [Google Scholar]
 Witzel, G., Martinez, G., Hora, J., et al. 2018, ApJ, 863, 15 [NASA ADS] [CrossRef] [Google Scholar]
 Witzel, G., Martinez, G., Willner, S. P., et al. 2021, ApJ, 917, 73 [NASA ADS] [CrossRef] [Google Scholar]
 Yuan, F., Bu, D., & Wu, M. 2012, ApJ, 761, 130 [Google Scholar]
 Yuan, Y., Blandford, R. D., & Wilkins, D. R. 2019a, MNRAS, 484, 4920 [NASA ADS] [CrossRef] [Google Scholar]
 Yuan, Y., Spitkovsky, A., Blandford, R. D., & Wilkins, D. R. 2019b, MNRAS, 487, 4114 [NASA ADS] [CrossRef] [Google Scholar]
 Zdziarski, A. A., Shapopi, J. N. S., & Pooley, G. G. 2020, ApJ, 894, L18 [NASA ADS] [CrossRef] [Google Scholar]
 Zdziarski, A. A., Dziełak, M. A., Marco, B. D., Szanecki, M., & Niedźwiecki, A. 2021, ApJ, 909, L9 [NASA ADS] [CrossRef] [Google Scholar]
 Zenitani, S., & Hoshino, M. 2007, ApJ, 670, 702 [NASA ADS] [CrossRef] [Google Scholar]
 Zhao, X., Gou, L., Dong, Y., et al. 2021, ApJ, 908, 117 [NASA ADS] [CrossRef] [Google Scholar]
 Zhu, Z., & Stone, J. M. 2018, ApJ, 857, 34 [Google Scholar]
 Zoghbi, A., Miller, J. M., & Cackett, E. 2019, ApJ, 884, 26 [Google Scholar]
Appendix A: Numerical convergence
A.1. Numerical relaxation
We monitor the relaxation of the initial state by following reduced quantities representative of the dynamics of the fields and of the particles in the corona. For the former, we measure the BH net electric charge Q_{BH} defined based on the charges crossing the event horizon since the beginning of the simulation. Since the innermost regions within the event horizon are cropped (r_{min} > 0), charges entering the horizon can exit the simulation space. Instead of counting the charges, we thus applied the MaxwellGauss equation integrated over the event horizon:
$$\begin{array}{c}\hfill {Q}_{\mathit{BH}}=\frac{1}{4\pi}{\displaystyle {\u222f}_{\mathit{EH}}\mathit{D}\xb7\mathbf{d}\mathbf{\Sigma}.}\end{array}$$(A.1)
The electric field, D, is symmetric with respect to the equatorial plane, and as such we can integrate over the upper hemisphere and multiply by 2. We then compared Q_{BH} to the electric charge Q_{0} for which the time component of the 4potential cancels out (Wald 1974):
$$\begin{array}{c}\hfill {Q}_{0}=2a{B}_{0}.\end{array}$$(A.2)
In the top panel in Fig. A.1, we see that the net charge eventually plateaus, which indicates that numerical relaxation is achieved in the innermost regions of the simulation space after approximately 50r_{g}/c. The net charge of the BH is negative and for a = 0.99, its magnitude reaches almost a fourth of Q_{0}. It was recently alleged by King & Pringle (2021) that if the BH was to acquire an electric charge Q_{0}, the BlandfordZnajek mechanism would not be able to operate because the 4potential no longer contains a term function of the BH spin. However, the correspondence of the 4potential with the Wald solution for a nonspinning BH does not mean that the structure of the electromagnetic fields is the same due to the different underlying metrics. Komissarov (2022) showed that in contrast with a Schwarzschild BH, a maximally spinning BH with an electric charge Q_{0} has a magnetosphere where the component of the electric field parallel to the magnetic field does not cancel: particles are susceptible to be accelerated along the magnetic field lines.
Fig. A.1. Time evolution of monitored quantities. Top panel: Net amount of electric charge, Q_{BH}, within the event horizon as a function of time for different spin values and for a thin disk (ϵ = 5%) in prograde rotation. Bottom panel: Total number of particles within the simulation space as a function of time. 
Due to the injection procedure, the relaxation of the particle distribution takes more time. After a steep increase in the number of particles visible in the bottom panel in Fig. A.1, the number of particles keeps increasing slower, especially for simulations with higher BH spin values. A plateau is eventually achieved after ∼100r_{g}/c, once the rate at which particles exit the simulation space equals the injection rate. The relaxation of the outermost regions of the simulation space takes more time than the relaxation of the innermost parts, so we typically work between t = 100r_{g}/c and t = 150r_{g}/c to guarantee full time convergence. The total number of particles correspond to approximately 20 particles per grid cell, with more particles in the current sheet. It ensures a good signaltonoise ratio by lowering the particle shot noise.
A.2. Resolution and scale separation
We evaluated the impact of the resolution and of the scale separation on the dissipation by computing the radial profiles for simulations with a = 0.8, ϵ = 5% and a disk in prograde rotation. Figure A.2 shows the integral of J ⋅ E over a volume bounded by the pole and the disk in the longitudinal direction, and between the event horizon (vertical solid black line) and a given distance r to the BH in the radial direction. Apart from the innermost regions, energy is progressively transferred from the electromagnetic field to the particles until reaching a value of approximately 0.001${B}_{0}^{2}{r}_{g}^{2}c$ at 26r_{g}, consistent with Fig. 11. We compare a simulation with the standard resolution we used in this paper (NR in Fig. A.2, with 2,128 and 1,120 cells in the radial and longitudinal directions, respectively) with a simulation with a resolution twice higher (HR, 4,256×2,240 cells) and a highresolution simulation like HR but with a dimensionless magnetic field scale ${\stackrel{\sim}{B}}_{0}$ higher by a factor of 4 (HB). As emphasized by Eq. (11), the Larmor radius is resolved by the same number of cells in all 3 simulations. The profiles are averaged over 40r_{g}/c after the simulations have relaxed.
Fig. A.2. Timeaveraged dissipated electromagnetic energy between the event horizon and a given radius for simulations with normal resolution (NR), high resolution (HR), and high magnetic field (HB). The dashed and dotted black lines show the distance to the separatrix footpoint and to the Y point, respectively. 
The normal and highresolution cases NR and HR yield very similar profiles and prove that the number of grid cells we use is sufficient to resolve the magnetic reconnection process and the subsequent particle acceleration. The simulation HB at higher magnetic field gives a profile within 15% of the NR and HR ones. Due to our injection method (see Sect. 2.2.5), the higher magnetic field produces a denser corona. The magnetization is four times higher than in cases NR and HR but in the relativistic reconnection regime we explore, this does not significantly modify the reconnection rate. Instead, the slightly lower dissipation might be due to the higher plasma inertia.
Appendix B: Y point location
We provide fitting formulas for the coordinates (x_{Y}, y_{Y}) and the distance r_{Y} to the BH of the Y point as a function of the spin value. We fitted the measured median positions with the following power law:
$$\begin{array}{c}\hfill r={c}_{1}{a}^{{c}_{2}},\end{array}$$(B.1)
where c_{1} and c_{2} are the fitting parameters. The bestfit parameters for the range of spins we explore are given in Table B.1. Although reasonable in the range of BH spin values we sample, these values might lead to significant departure from the physical ones when the power law (Eq. (B.1)) is extrapolated to spins lower than 0.6.
Parameters of the bestfit power laws to locate the Y point as a function of the BH spin based on Eq. (B.1).
Appendix C: Synthetic images
In Fig. C.1 we show the synthetic images of synchrotron emission from the current sheet for different inclinations of the lineofsight with respect to the BH spin axis. The BH spin is a = 0.8 and the disk is thin and in prograde rotation around the BH. The inclination i = 55° corresponds to Fig. 13.
Fig. C.1. Intensity maps of synchrotron emission from the current sheet for different viewing angles. 
All Tables
Parameters of the bestfit power laws to locate the Y point as a function of the BH spin based on Eq. (B.1).
All Figures
Fig. 1. Fiducial radial profiles for the critical densities n_{δ}, n_{σ}, and n_{GJ} defined in the text. The numbers are the powerlaw exponents, and the arrows indicate the allowed regions. The green hatched region represents the region where the forcefree regime can be achieved while still resolving the skin depth. 

In the text 
Fig. 2. Poloidal magnetic field lines around the BH (solid black lines), with the event horizon (black disk) and the ergosphere (black dashed line), for a fiducial snapshot from a simulation with a = 0.8 and a thin disk (in chartreuse yellow in the equatorial plane). The thicker field lines delimit the region between the ISCO and the separatrix coupling the BH to the disk, visible in the zoomedin view in the bottom panel. The color maps represent the total plasma density, n, in units of n_{GJ} (top left), the mean Lorentz factor of the particles (top right), the toroidal component of the H field (bottom left), and the radial component of the current, J (bottom right). To compensate for spatial dilution, plasma density and current were multiplied by r^{2}. A smoothing Gaussian kernel a few cells wide was applied to the current map. In the plasma density map, the red circle locates the Y point. Distances to the BH on the x and y axes are given in units of r_{g}. The white line in the topright corner stands for the outer light surface. 

In the text 
Fig. 3. Angular speed of the magnetic field lines. Top panel: color map of the timeaveraged angular speed, Ω, with poloidal magnetic field lines overlaid, for a simulation with a = 0.6 and ϵ = 5%. The dashed white line is the outer light surface. Bottom panel: angular speed of magnetic field lines as a function of their magnetic flux function, A_{ϕ}, with the corresponding distance of the line footpoint on the disk at the top. The color map is the histogram of the (A_{ϕ}, Ω) pairs measured in each cell of the simulation space at a fiducial time. The left (respectively right) vertical line stands for the separatrix footpoint (respectively the ISCO). In orange are the Keplerian profiles (Eq. (2)) without (solid) and with (dashed) the smooth cutoff, and half of the BH angular speed (dashdotted). 

In the text 
Fig. 4. Threedimensional representation of the magnetic field lines around the BH (central black sphere), with the upper and lower disk surfaces in semitransparent black. Field lines are colorcoded from low magnetic potential (dark green) to high magnetic potential (magenta). Two additional field lines are represented in the zoomedin view in the bottom panel. 

In the text 
Fig. 5. Reduced quantities as a function of the BH spin. Top panel: ISCO (solid green line), event horizon (solid black line), and ergosphere at the equator (dashed black line) as a function of the BH spin. The blue line indicates the approximate location of the footpoint of the separatrix found by Uzdensky (2005), and the blue squares are the values measured in our simulations. Middle panel: magnetic flux Φ_{jet} of the open magnetic field lines threading the event horizon (green, left axis). In blue (right axis) are the jet electromagnetic powers measured (squares) and predicted by Tchekhovskoy et al. (2011) (crosses). Bottom panel: angular momentum (left axis) and energy (right axis) transferred from the BH to the disk per unit time, with and without the forcefree approximation (crosses and squares). 

In the text 
Fig. 6. Profiles along a line transverse to the current sheet and passing through an X point. The components of the drift velocity (solid black line) and of the positron and electron velocities (respectively blue circles and green squares) normal to the current sheet are represented. 

In the text 
Fig. 7. Plasmoid detachment from the Y point. Left panel: distance between the Y point and the footpoint of the separatrix, d_{Y}, as a function of time for different BH spin values. The origin of time is arbitrary, but over these time lapses the simulations have reached a numerically relaxed state. Right panel: estimates of the time periodicity of the motion of the Y point, τ_{Y}, as a function of the BH spin. 

In the text 
Fig. 8. Density maps of the locations of the Y point for different BH spin values. The inner solid lines correspond to the event horizons, and the dashed lines are the separatrix. Distances to the BH are given in r_{g}. 

In the text 
Fig. 9. Coordinates (x_{Y}, y_{Y}) and distance, r_{Y}, to the BH center of the Y point as a function of the BH spin. The dashed lines are bestfit power laws (see Appendix B). 

In the text 
Fig. 10. Dissipation of electromagnetic energy and particle acceleration. Left column: electromagnetic (top panel) and particle kinetic (bottom panel) energy flux through a surface extending up to a radius r (left axis) as a function of time and for a thin disk with a BH spin a = 0.99. Right column: slices at constant time t = 150 r_{g}/c (top panel) and constant radius r = 26 r_{g} (bottom panel) of the fluxes. The dashed lines in the upper panel are the timeaveraged profiles. The solid, dashed, and dotted yellow lines show the event horizon, the distance to the separatrix footpoint, and the distance to the Y point, respectively. 

In the text 
Fig. 11. Amount of electromagnetic energy dissipated per unit time in the current sheet as a function of the distance of the Y point to the BH (bottom x axis). The top x axis indicates the corresponding BH spin (analog to the black points in Fig. 9). The dashed line stands for the best semianalytic fit (see text). 

In the text 
Fig. 12. Map of synchrotron emission power for a BH spin a = 0.8, multiplied by r^{2} to highlight the contribution from the current sheet. Photons emitted from the region delimited with with a dotted line are those shown as originating from the current sheet in the synthetic images in Figs. 13 and C.1. 

In the text 
Fig. 13. Synthetic images of synchrotron emission from the current sheet for a simulation with a thin disk in prograde rotation around a Kerr BH with a spin a = 0.8. 

In the text 
Fig. 14. Map of the toroidal component, H_{ϕ}, for a simulation with a = −0.8 (retrograde) and ϵ = 5%. The magnetic field lines are those surrounding the nonreconnecting current sheet and the pair corresponding to the potential at the ISCO. The Keplerian disk is in light green and extends down to r_{ISCO} ∼ 8.4 r_{g}. 

In the text 
Fig. 15. Deposit of particle kinetic energy on the disk. Left panel: flux on the disk of particle kinetic energy per unit time per ring of infinitesimal width, dr, as a function of the distance, r, to the BH, for a = 0.8. The solid blue curve is the median profile, and the light blue shaded region is the 1σ variability. Vertical lines are the event horizon (solid) and the separatrix footpoint (dashed). Right panel: flux of particle kinetic energy per unit time near the separatrix footpoint compared to the rate of energy deposited onto the disk, Ė_{BH−disk}, per unit time via the coupling magnetic field lines. 

In the text 
Fig. 16. Distribution of the electrons and positrons in the current sheet as a function of their Lorentz factor. 

In the text 
Fig. 17. Spectral energy distribution of synchrotron emission from the current sheet. 

In the text 
Fig. 18. Characteristic timescale of angular momentum deposited onto the coupled part of the disk, τ_{L}, as a function of the BH spin. 

In the text 
Fig. A.1. Time evolution of monitored quantities. Top panel: Net amount of electric charge, Q_{BH}, within the event horizon as a function of time for different spin values and for a thin disk (ϵ = 5%) in prograde rotation. Bottom panel: Total number of particles within the simulation space as a function of time. 

In the text 
Fig. A.2. Timeaveraged dissipated electromagnetic energy between the event horizon and a given radius for simulations with normal resolution (NR), high resolution (HR), and high magnetic field (HB). The dashed and dotted black lines show the distance to the separatrix footpoint and to the Y point, respectively. 

In the text 
Fig. C.1. Intensity maps of synchrotron emission from the current sheet for different viewing angles. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.