Numerical simulation of hot accretion flows (IV): effects of black hole spin and magnetic field strength on the wind and the comparison between wind and jet properties
DDraft version February 8, 2021
Preprint typeset using L A TEX style AASTeX6 v. 1.0
NUMERICAL SIMULATION OF HOT ACCRETION FLOWS (IV): EFFECTS OF BLACK HOLE SPIN ANDMAGNETIC FIELD STRENGTH ON THE WIND AND THE COMPARISON BETWEEN WIND AND JETPROPERTIES
Hai Yang , Feng Yuan , Ye-fei Yuan , and Christopher J. White Shanghai Astronomical Observatory, Chinese Academy of Sciences, 80 Nandan Road, Shanghai 200030, [email protected] (HY), [email protected] (FY) University of Chinese Academy of Sciences, 19A Yuquan Road, Beijing 100049, China Department of Astronomy University of Science and Technology of China Hefei, Anhui, China; [email protected] Department of Astrophysical Sciences, Princeton University, Peyton Hall, Princeton, NJ, USA; [email protected]
ABSTRACTThis is the fourth paper of our series of works studying winds from hot accretion flows around blackholes. In the first two papers, we have shown the existence of strong winds in hot accretion flows usinghydrodynamical and magnetohydrodynamical (MHD) simulations. In the third paper, by using threedimensional general relativity MHD numerical simulation data of hot accretion flows and adopting a“virtual particle trajectory” data analysis approach, we have calculated the properties of wind, such asits mass flux and velocity. However, that paper focuses only on a non-spinning black hole and SANE(standard and normal accretion). In the present paper, we extend the third paper by including casesof a rapidly rotating black hole and MAD (magnetically arrested disk). We focus on investigatingthe effect of spin and magnetic field on the properties of wind and jet. It is found that a larger spinand stronger magnetic field usually enhance the wind and jet. The formulae describing the massflux, poloidal velocity, and fluxes of momentum, kinetic energy, and total energy of wind and jet arepresented. One interesting finding, among others, is that even in the case of very rapidly spinningblack hole where the jet is supposed to be the strongest, the momentum flux of jet is smaller than thatof wind, while the total energy flux of jet is larger than that of wind by at most a factor of 10. Theimplications of this result to the importance of wind relative to jet in active galactic nuclei feedbackare discussed.
Keywords: accretion, accretion disks — black hole physical : jet — hydrodynamics INTRODUCTIONDepending on the value of mass accretion rate, blackhole accretion is divided into two modes, namely coldmode above roughly
2% ˙ M Edd and hot mode below thisrate. Here ˙ M ≡ L Edd /c is defined as the Eddingtonaccretion rate, with L Edd being the Eddington luminos-ity. The cold mode includes standard thin disk (Shakura& Sunyaev 1973; Pringle 1981) and super-Eddington ac-cretion (Abramowicz et al. 1988; Sądowski et al. 2014;Jiang et al. 2014), bounded by the Eddington accretionrate. It corresponds to the cold mode (or quasar modeor radiative mode) in the study of active galactic nu-clei (AGN) feedback. The hot mode includes advection-dominated accretion flow (Narayan & Yi 1994, 1995;Yuan & Narayan 2014) and luminous hot accretion flows(Yuan 2001), bounded by ∼ . α ˙ M Edd . It correspondsto the hot mode (or radio mode or jet mode or mainte-nance mode) of AGN feedback.In the present paper we continue our study of windfrom hot accretion flows. The motivation of such a study is two-fold. First, wind is a fundamental ingredient ofblack hole accretion flow, affecting both the dynamicsand radiation of the accretion flow. Second, active galac-tic nuclei feedback is now believed to play a crucial rolein galaxy formation and evolution (Fabian 2012; Kor-mendy & Ho 2013; Naab & Ostriker 2017), while windfrom hot accretion flow is potentially one of the mostimportant mediums of feedback (Weinberger et al. 2017;Yuan et al. 2018; Yoon et al. 2019). For example, recentIllustrisTNG cosmological simulations (Weinberger et al.2017) find that to overcome some serious problems ingalaxy formation, e.g., reducing star formation efficiencyin the most massive halos, winds launched from hot ac-cretion flows must be invoked to interact with the inter-stellar medium on the galaxy scale. Yuan et al. (2018) re-cently comprehensively include feedback by wind and ra-diation from AGNs in both cold and hot feedback modesand find that wind plays a dominant role in controllingthe star formation and black hole growth, although ra-diative feedback cannot be neglected. a r X i v : . [ a s t r o - ph . H E ] F e b The study of wind was started by Stone et al. (1999).In this milestone paper they find for the first time thatthe mass accretion rate of the hot accretion flow de-creases inward. Two models have been proposed toexplain this surprising result. In the adiabatic inflow-outflow solution (ADIOS; Blandford & Begelman 1999,2004; Begelman 2012; see also an earlier speculationfor the existence of strong wind in hot accretion flowsin Narayan & Yi 1994), it is assumed that the inwarddecrease of accretion rate is because of mass lost inthe wind. In the convection-dominated accretion flow(CDAF) model, it is thought to be because of convectivemotion in the hot accretion flow (Narayan et al. 2000;Quataert & Gruzinov 2000; Abramowicz et al. 2002).To solve this debate, in the first paper of this seriesof works, Yuan et al. (2012b) first significantly extendthe radial dynamical range of Stone et al. (1999) basedon a “two-zone” approach and confirm the Stone et al.(1999) result. In the second paper of this series (Yuanet al. 2012a), by systematically comparing the statisti-cal value of some physical quantities such as radial ve-locity and angular momentum of “inflow” and “outflow”and analyzing the convective stability of the MHD accre-tion flow, they show that strong wind must exist and isthe reason for the inward decrease of mass accretion rate(see also Narayan et al. (2012)). In addition to numericalsimulations, some analytical solutions of accretion withwind have also been obtained (Gu 2015; Mosallanezhadet al. 2016), with the results broadly consistent with thenumerical simulations. Bu et al. (2016a,b) address thequestion of where the wind can be produced, they findthat winds can only be produced within the Bondi radius.This is likely because of the change of the gravitationalpotential beyond the Bondi radius due to the contribu-tion of stellar cluster. The theoretical prediction of theexistence of strong wind from hot accretion flow has beenconfirmed by some observations in recent years, includ-ing the supermassive black hole in the Galactic center(Wang et al. 2013; Ma et al. 2019), low-luminosity AGNs(Tombesi et al. 2014; Cheung et al. 2016; Park et al.2019), and black hole X-ray binaries in the hard state(Homan et al. 2016).The remaining question is what are the main proper-ties of the wind such as mass flux and velocity. Thisquestion is difficult to answer because the accretion flowis strongly turbulent; thus it is hard to discriminate theturbulent outflow and the real wind. The widely adoptedsolution in literature is to time-average the simulationdata first to filter out turbulence. However, since windis always instantaneous, such a procedure will also in-evitably eliminate the real wind; thus usually wind willbe significantly underestimated.This difficult is solved in Yuan et al. (2015), the thirdpaper of this series. In this work, they propose a “virtual particle trajectory” approach, which can faithfully reflectthe motion of fluid elements and thus clearly distinguishturbulent outflow and wind. Using this approach, andbased on the three dimensional general relativity MHDsimulation data of hot accretion flows, they successfullyobtain the main properties of wind. They have alsocompared the mass flux of wind obtained by the tra-jectory approach with that obtained by the usual time-average streamline approach. The exact difference be-tween the two results depends on the radius and timeinterval adopted in the time average calculation. Thelonger the interval, the larger the difference. Taking r g as an example, it is found that the flux of wind obtainedby the time-average streamline approach is smaller by afactor of 10 compared to that obtained by the trajectoryapproach (refer to Fig. 5 in Yuan et al. 2015).The current work is a direct extension of Yuan et al.(2015). Yuan et al. (2015) only deals with SANE (stan-dard and normal evolution; Narayan et al. 2012) aroundnon-spinning black holes. In reality, it is possible thatthe magnetic field in the flow may be much stronger, i.e.,the accretion is in the MAD (magnetically arrested disk;Narayan et al. 2003) state, and the black hole spin ismore likely nonzero. In the present work, we investigatehow the properties of wind will change with the strengthof magnetic field and black hole spin . Since a rapidlyspinning black hole will power a relativistic jet, it is in-teresting to compare the properties of winds with the jetsuch as their momentum and energy fluxes. Such infor-mation is also valuable for us to evaluate the respectiverole of wind and jet in AGN feedback.The main structure of the paper is as follows. In Sec-tion 2, we will describe our equations and simulationmethod (§2.1), the three models we study (§2.2), anda brief overview to the trajectory approach we use toanalyze the simulation data (§2.3). The results will bepresented in Section 3. We then summarize our work inSection 4. The implications of our results in the contextof AGN feedback are discussed in Section 5. MODEL2.1.
Equations and numerical method
The equations of ideal MHD describing the evolutionof the accretion flow read (e.g., Gammie et al. 2003), (cid:53) µ ( ρu µ ) = 0 , (1) (cid:53) µ T µν = 0 , (2) (cid:53) µ ∗ F µν = 0 . (3) In fact, Sadowski et al. (2013) have studied the wind in thecases of MAD and rapidly spinning black hole. However, as in mostworks, the time-average streamline approach is adopted. This iswhy they find the wind very weak. with T µν = ( ρh + b λ b λ ) u µ u ν + ( p gas + 12 b λ b λ ) g µν − b µ b ν , (4) ∗ F µν = b µ u ν − b ν u µ . (5)Here T µν is the stress-energy tensor of the MHD, h isthe specific enthalpy of the fluid. We use the Athena++(White et al. 2016; Stone et al. 2020) code to solve theabove GRMHD equations in the Kerr metric. This codeuses a finite-volume Godunov scheme to ensure totalenergy conservation, with the flux of conserved quanti-ties obtained by solving the Riemann problem at eachinterface. In our simulation, we use the HLLE Rie-mann solvers (Einfeldt 1988), which is commonly used inGRMHD simulation. In order to satisfy the divergence-free constraint to prevent spurious production of mag-netic monopoles, the staggered-mesh constrained trans-port (CT) method is applied. For the spatial reconstruc-tion, we use the piecewise linear method (PLM; van Leer1974). For the mesh grid, we use the static mesh refine-ment (SMR). It allows us to easily use higher resolutionin areas of interest, and can give a good balance betweenaccuracy and performance, which is very useful for large-scale 3D simulations.All our simulations are performed in the Kerr–Schild(horizon penetrating) coordinates ( t , r , θ , ϕ ). The radiusof the black hole horizon r H is r H = (1 + √ − a ) r g ,with r g ≡ GM BH /c being the gravitational radius and a the black hole spin parameter. The determinant ofthe metric g ≡ det( g µν ) = − Σ sin θ , where Σ ≡ r + a cos θ (McKinney & Gammie 2004). The comovingrest-mass density is denoted as ρ , u µ is the componentof the coordinate-frame 4-velocity, the equation of thestate of the gas is u = p gas / (Γ − , with p gas being thegas pressure of the comoving mass, u being the internalenergy of the gas, and Γ being the adiabatic index. In ourwork, Γ is taken to be 4/3, and time is measured in unitof r g /c . The units we use is Heaviside-Lorentz, both thelight speed and gravity constant are set to be unity, andthe metric sign convection is ( − , + , + , +) . The metric ofour simulation is stationary, the self-gravity is ignored.2.2. Three models and their setup
In this paper, we consider three models, i.e., MAD00,SANE98, and MAD98. They denote the MAD accretionflow around a black hole of spin a = 0 , SANE accretionflow around a black hole of a = 0 . , and MAD accretionflow around a black hole of a = 0 . , respectively. Wewill combine these three models with the SANE arounda black hole of a = 0 presented in (Yuan et al. 2015) toanalyze the dependence of wind properties on black holespin and magnetic field strength.Our simulations starts with a torus that rotates around Table 1 . The static mesh refinementgrid of MAD00.
Level r / r g θ / π ϕ / π [1 . , [0.1305, 0.8695] [0, 2] [8 . , [0, 0.1210] [0, 2] [8 . , [0.8790, 1] [0, 2]2 [1 . , [0.1942, 0.8058] [0, 2] [44 . , [0, 0.0573] [0, 2] [44 . , [0.9427, 1] [0, 2] Table 2 . The static mesh refinementgrid of MAD98 and SANE98.
Level r / r g θ / π ϕ / π [1 . , [0.1305, 0.8695] [0, 2]1 [6 . , [0, 0.1210] [0, 2] [6 . , [0.8790, 1] [0, 2] [1 . , [0.1942, 0.8058] [0, 2]2 [36 . , [0, 0.0573] [0, 2] [36 . , [0.9427, 1] [0, 2] Table 3 . The setup of three models.
Model a N r N θ N ϕ DurationMAD00 0 288 128 64 40 000SANE98 0.98 352 128 64 80 000MAD98 0.98 352 128 64 40 000 a black hole. The torus is initially in hydrostatic equi-librium, described by Fishbone & Moncrief (1976). Theinner edge of the torus is at r = 40 . r g , and the radiusof pressure maximum is at r = 80 r g . The gas torus isthread by a poloidal magnetic field (Penna et al. 2013).In our simulation, we use two different initial magneticfield configurations for MAD and SANE. For the MAD,we set one poloidal loop threading the whole torus. Thiswill result in rapid accretion. The magnetic flux wouldbe accumulated and finally impede the accretion of mass.The magnetic field flux will quickly saturate at a max-imum value on the black hole for a given mass accre-tion rate, reaching the MAD state (Narayan et al. 2003;McKinney et al. 2012; Sadowski et al. 2013). For theSANE, we initially set a seed field that consists of multi-ple poloidal loops of magnetic field with changing polar-ity. Such a configuration makes magnetic reconnectioneasy to occur, thus the magnetic field will always stayweak, preventing the accumulation of magnetic flux sothe accretion flow is in the SANE state (Narayan et al.2012). Figure 1 shows the initial magnetic field configu-ration of MAD and SANE.Three parameters are required to determine the spe-cific magnetic field of the torus, namely r start , r end , and λ B . The first two denote the inner and outer edges ofthe magnetized region, while the last parameter controlsthe size of the poloidal loops, or equivalently the num-ber of the loops. The normalization of the magneticfield is determined by the gas-to-magnetic pressure ratio β ≡ p gas /p mag . It has the minimum value in the equato-rial plane for each loop. It peaks at loop edges and dropsto the loop center. For MAD, we set r start = 25 r g , r end =810 r g , λ B = 25 (it just has one loop), and β min = 0 . .For SANE, we set r start = 25 r g , r end = 550 r g , λ B = 3 . (it has five loops), and β min = 0 . .Our grid uses static mesh refinement with a root gridplus a more refined grid. The root grid is × × cells in radial, polar, and azimuthal direction. In the r direction, the inner and outer edges are located at 1.1 r g and 1200 r g respectively. Logarithmic spacing is adopted,with the ratio r i +1 /r i being 1.0827. The grid in the polarand azimuthal directions are uniform. Because of the useof Kerr–Schild coordinates, the inner edge of our simu-lations are within the black hole horizon. For differentmodels we use different static mesh refinement. Table 1and Table 2 give the details of the three models. For arefinement region, each additional level means that thegrid of this area is refined by a factor of 2 for all direc-tions based on the previous level grid. Table 3 is thefinal effective grid of the different model. Figure 2 showsthe different zoom levels of the mesh grid in the poloidalplane of MAD98.For the inner and outer boundaries we use the outflowboundary conditions. We use the polar axis boundarycondition in the θ direction and periodic boundary con-ditions in the ϕ direction. Since numerical MHD cannot deal with vacuum, like all GRMHD numerical sim-ulations, we need to impose the density and gas pres-sure floors: ρ min = max(10 − r − . , − ) , p gas , min =max(10 − r − . , − ) . Following White & Chrystal(2020), we also enforce σ < and γ < . Here σ = 2 p mag /ρ is the magnetization parameter, γ ≡ αu t with α ≡ ( − g tt ) − / the lapse. The former is an extralimitation for the density and gas pressure floors. Thelatter is to limit the velocity to keep the normal-frameLorentz factor from becoming too large.2.3. Virtual particle trajectory approach
We briefly overview the “virtual particle trajectory” ap-proach proposed in Yuan et al. (2015). Trajectory isrelated to the Lagrangian description of fluid, obtainedby following the motion of fluid elements at consecutivetimes. Since the motion of fluid is not strictly steady butturbulent, trajectory is very different from streamlines.To get the trajectory, we first need to choose a set of “test particles” in the simulation domain. Note that theyare not real particles, but a collection of spatial coordi-nates as the starting points for the trajectory calculation.Usually we put some “test particles” in different θ and ϕ at a given initial radius r and a given initial time t .Their velocity can be obtained from the interpolation ofthe simulation data. With this information, we can ob-tain their location at time t + δt . We use the “VISIT”software to perform the calculation. Repeating this pro-cess, we can obtain the trajectory of these particles. It isimportant to choose an appropriate time step δt to ob-tain a convergent trajectory. After some tests, we choose δt = 100 r g /c in our work; it corresponds to the Keple-rian timescale at r ≈ r g . In this way, we can obtain thetrajectories of fluid elements or “virtual test particles”based on our high time resolution numerical simulationdata. The readers can refer to Figures 1 & 2 in Yuanet al. (2015) for examples of various trajectories we haveobtained in the case of SANE00. Once we obtain thetrajectories, we can judge whether the motion is turbu-lence or real wind, and calculate the wind properties suchas mass flux and poloidal velocity, as we will discuss in§3. Obviously, obtaining the trajectory of fluid elementsis much more time consuming than the streamlines ap-proach, but it can faithfully reflect the motion of fluidcompared to the latter. RESULT3.1.
Definitions of wind and jet
Let us first define some terminology before describingour results. “Outflow” means the flow with a positive ra-dial velocity v r , including both “turbulent outflow” and“real outflow.” In the former case the test particle willfirst move outward but eventually will return after mov-ing outward for some distance. In the latter case thetest particle continues to flow outward and eventuallyescapes the outer boundary of the simulation domain.“Real outflow” consists of two components, i.e., jet andwind. The jet region is defined as the region occupiedby magnetic field lines connected to the ergosphere (de-scribed by r erg ≡ r g + √ − a cos θ r g (Visser 2007))of the black hole. Thus the jet region is bounded by themagnetic field line whose foot point is rooted at the blackhole ergosphere with θ = 90 ◦ , i.e., the boundary betweenthe black hole ergosphere and the accretion flow (refer tothe red lines in Figure 5). In this case all magnetic fieldlines in the jet region are anchored to the black hole ergo-sphere and thus can extract the spin energy of the blackhole via the Blandford & Znajek mechanism (Blandford& Znajek 1977) to power the jet (“BZ-jet”). Real outflowsoutside of this boundary are powered by the rotation en-ergy of the accretion flow and we call them wind. Notethat our definition of wind adopted here is different from Figure 1 . The initial magnetic field configuration of MAD (left) and SANE (right). The dashed line in the SANEmeans that it has a different polarity compared with the solid line. The color is the logarithm of the gas-to-magneticpressure ratio β .that adopted in some literature, where they require thatthe Bernoulli parameter of wind must satisfy Be > .We do not add this requirement, because we find thatfor non-steady accretion flow Be is not constant alongtrajectories, but usually increases outward at least un-til the radius within which turbulence is well developed(Yuan et al. 2015). We find that even though Be is neg-ative at a small radius, it can become positive when itpropagates outward.In Yuan et al. (2015), although the black hole a = 0 ,we still find jet-like outflow in terms of high velocity.This is confirmed by the results of the present paper, aswe will show in the left panels of Figure 8 for MAD00.We can see in this figure that the poloidal speed of theplasma within the red line is significantly larger than thatof the wind. Of course, these outflows must be poweredby the rotation energy of the underlying accretion diskrather than by the black hole, so we call them “disk-jet”(Yuan & Narayan 2014; Yuan et al. 2015; see also Ghosh& Abramowicz 1997; Livio et al. 1999). These resultssuggest that in AGNs and black hole X-ray binaries jetsmay still be present even though the black hole is non-spinning. In the present work, we focus on the differencebetween outflows powered by the spinning black hole andby the rotating accretion disk, so we simply call the for-mer “BZ-jet” or “jet” and the latter “wind”.3.2. Overview of the simulation results
Both MAD00 and MAD98 reach steady state after t = 5000 ; for SANE98, it reaches the steady state after t = 50 , . Figure 3 shows the t - and ϕ -averaged two di-mensional distribution of density and velocity in the r − θ plane for the three models at the chosen time chunks of t = 10 , – , , , – , and – , , re-spectively. Figure 4 shows in a more quantitative wayvarious quantities of the flow averaged over ϕ and timeas a function of θ at r g . From the top to bottom panels we have density, the ratio of the gas pressure to magneticpressure β , temperature, and B . From the figure wecan see that the accretion flows of the two MAD modelsare geometrically thinner than that of the SANE model,consistent with previous works (e.g., Tchekhovskoy et al.2011; McKinney et al. 2012). The reason is because in theMAD model the magnetic field is so strong that it com-presses the accretion flow vertically. We also find thatthe velocity field of the two MAD models are much moreordered than that of the SANE model. This is becausethe magnetorotational instability (MRI) is suppressed inthe MAD model so there is no turbulence, while MRI isstill present in the SANE (e.g., McKinney et al. 2012).Figure 5 shows the ϕ - and t -averaged magnetic fieldlines of the three models. The field configuration of thetwo MAD models are very similar and both are very or-dered. This is again because MRI is suppressed so thereis no turbulence. For the SANE model, the field linesare less ordered due to the existence of turbulence. Theplasma β is also much smaller in the two MAD modelsthan in the SANE model, as expected. The red line isthe magnetic field line anchored at the ergosphere with θ = 90 ◦ , i.e., the boundary between the BZ-jet and windregions. From the figure, we can see that the BZ-jetregion is strongly dominated by magnetic pressure, espe-cially for the two MAD models. The shape of the red lineis parabolic, i.e., z ∝ R s in the cylindrical coordinates ( R, z ) with s > . There are many discussions about theshape of this magnetic field line or streamline. Readerscan refer to, e.g., Nakamura et al. (2018) and Chen &Zhang (2021) for details.3.3. Mass fluxes of wind and jet
To calculate the mass flux of the wind at a given time t , we first put test particles at a given radius r withdifferent θ and ϕ and obtain their trajectories. The massflux of wind is calculated by summing up the correspond- Figure 2 . Different zoom levels of the mesh grid of the poloidal plane of the MAD98.
Figure 3 . Two dimensional distribution of density (by color) and velocity vector (by arrows) of the accretion flowaveraged in ϕ direction and in time interval of t = 55 , – , , , – , and – , for SANE98,MAD00 and MAD98 respectively. The bottom panels are the zoom-in pictures of the top panels.ing mass flux carried by test particles whose trajectoriesbelong to the real outflow (i.e., wind or jet) using thefollowing formula: ˙ M wind(jet) ( r ) = (cid:88) i ρ i ( r ) u ri ( r ) √− gδθ i δϕ i . (6)Here the subscript represent the different test particles, ρ i and u ri are the mass density and four-velocity at thelocation where the test particle “ i ” originates, and δθ i and δϕ i are the ranges of θ and ϕ the particles occupy.In our calculations, we usually choose 10 different ini-tial times t to obtain the trajectories of these test par-ticles and the mass flux corresponding to each choice of t . We then do time-average of these 10 groups to ob-tained the averaged mass flux. These 10 initial timesare uniformly selected after the simulation has reachedsteady state. Because the calculation of trajectory re-quires simulation data after time t , and it needs “testparticles” to move away from their initial location farenough to distinguish the real outflow from turbulence,the initial time t should not be too late. For MAD00,MAD98 and SANE98, the time intervals we choose forthe time-average are t = 10 , to , , t = 8000 to23000, and t = 55 , to , , respectively. To calcu-late the wind properties (such as mass flux, velocity, andenergy and momentum fluxes) as a function of radius,we have chosen a series of initial radii, namely r = 20 r g , r g , r g , r g , r g , r g , r g for MAD98 andMAD00, and r = 10 r g , r g , r g , r g , r g , r g , r g , r g for SANE98. The initial position of the“test particle” in θ and ϕ directions are consistent withour simulation grid, i.e., 128 grid points in the θ direction [0 , π ] and 64 grid points in the ϕ direction [0 , π ] . Wesum the mass flux of all real outflow as wind mass fluxfor MAD00 since there is no BZ-jet in this case; but forthe a = 0 . model, i.e., SANE98 and MAD98, becauseof the existence of the BZ-jet, we calculate the mass fluxof the wind and jet separately.Figure 6 shows the results. The black, red, and bluesolid lines show the inflow, outflow, and net rates, re-spectively. The inflow and outflow rates are calculatedusing the following equation (e.g., Tchekhovskoy et al.2011; White 2019), ˙ M ( r ) = − (cid:90) π (cid:90) π ρu r √− gdθdϕ, (7)where u r is the radial contravariant component of the4-velocity, and the integral is over all θ, ϕ at fixed r . Weadd the negative sign so that the flux is positive whenflow moves inward. The net rate is the difference betweeninflow and outflow rates. Note that here the outflow rate ˙ M out includes both real outflow (jet and wind) and theturbulent outflow. The red dashed and red dot-dashedlines denote the mass fluxes of wind and jet respectively, obtained using our trajectory approach. Combined withthe result of SANE00 from Yuan et al. (2015), the radialprofiles of the mass flux of wind for SANE00, SANE98,MAD00 and MAD98 can be described by ˙ M w ( r ) = (cid:18) r r s (cid:19) . ˙ M BH , (8) ˙ M w ( r ) = (cid:18) r r s (cid:19) . ˙ M BH , (9) ˙ M w ( r ) = (cid:18) r r s (cid:19) . ˙ M BH , (10) ˙ M w ( r ) = (cid:18) r r s (cid:19) . ˙ M BH , (11)respectively. Here ˙ M BH is the mass accretion rate at theblack hole horizon. From these, we can find the followingresults:• By comparing SANE00 and MAD00, we can seethat close to the black hole, the wind flux becomesweaker when the magnetic field becomes stronger.This is likely because turbulence helps the produc-tion of wind. The turbulence in SANE is caused byMRI; but in the case of MAD, the magnetic fieldis too strong close to the black hole thus MRI issuppressed and there is no longer turbulence (e.g.,Narayan et al. 2003; McKinney et al. 2012).• However further away from the black hole, wherethe magnetic field is not too strong so that MRIis present in both SANE and MAD, the magneticfield helps the production of wind. This is whythe power-law index for the MAD00 model is 1.54,larger than that of the SANE00, which is 1. Similarresults can be found when we compare SANE98and MAD98.• By comparing SANE00 and SANE98, we can findthat the wind flux become stronger with the in-crease of black hole spin if the radius is not toolarge. This is confirmed by comparing MAD00 andMAD98. This indicates that the black hole spinhelps the production of wind.For comparison purpose, we have also used the time-averaged approach calculating the mass flux of wind forthe three models. As in Yuan et al. (2015), we find thatthe mass flux obtained by this approach is smaller thanthat obtained by the trajectory approach, and the dis-crepancy depends on the radius and time interval of av-erage. For example, at 50 r g and 200 r g , for the time in-terval we use, the former is found to be smaller than thelatter by a factor of − . The difference for SANE98 is Figure 4 . Different quantities of the flow averaged over ϕ and time as a function of θ for MAD00 (left), SANE98(middle), and MAD98 (right) at r = 40 r g . The red dashed line represents the boundary between the BZ jet and wind. Figure 5 . Magnetic field distribution of MAD00, SANE98 and MAD98. The white lines denote magnetic field lines;the colors denote the logarithm of the β of the plasma. The red line is the magnetic field line that is rooted at the blackhole ergosphere with θ = 90 ◦ , i.e., the boundary between the ergosphere and the accretion flow. It is the boundarybetween the BZ-jet and wind.slightly larger than MAD00 and MAD98, which is likelybecause turbulence in MAD is weaker. Such a discrep-ancy is smaller than that found in Yuan et al. (2015). Wespeculate that this is because different simulation codesand different time interval are used.Now let us discuss the mass flux of jet. From Figure 6,we can see that for SANE98, the jet is produced at smallradii and gradually increases until it saturates at r ≈ r g . For MAD98, the jet mass flux increases faster, andsaturates at a larger radius of r ≈ r g . The saturatedmass fluxes of the BZ-jet in SANE98 and MAD98 are ˙ M jet ≈ . M BH , (12) ˙ M jet ≈ . M BH , (13)respectively. Again, ˙ M BH denotes the mass accretionrate at the black hole horizon. The mass flux of the jetin MAD98 is larger than that of the SANE98. implyingthat magnetic field helps the formation of jet.Figure 7 shows the mass flux of wind per unit θ inte-grated over all ϕ and averaged from t = 55 , – , , , – , and – , for SANE98, MAD00and MAD98, respectively. The positive value denotesreal outflow while negative value denotes the inflow. Thefigure shows that for most cases, the mass flux of inflow isconcentrated on the equatorial plane. But at r = 160 r g of the SANE98 model, the inflow flux has some bimodaldistribution, the accretion is mainly via the coronal re-gion. Such a result has been discussed and termed “coro-nal accretion” in Zhu & Stone (2018), and is also similarto the accretion model of Quadrupole Topology QDPa inBeckwith et al. (2008).3.4. The poloidal speed of wind and jet
In this section, we discuss the poloidal speed of windand jet. This is the dominant component of the outflow velocity at large radii where the magnetic field becomesweak. Figure 8 shows the poloidal speed of wind and jetaveraged over ϕ and time from , to , , , to , and to , for SANE98, MAD00, andMAD98 respectively. We can see the same with theSANE00 model presented in Yuan et al. (2015), for thethree models the velocity is higher in the polar directionbut smaller closer to the equatorial plane. By comparingMAD98 and SANE98, we can see that the largest veloc-ity of the jet of the two models are similar. However, thevelocity shown in this figure is time-averaged. We findthat the instantaneous velocity of the two models arequiet different when examining the snapshot result. Thevelocity of jet in SANE98 varies little with time, whilethat in MAD98 varies greatly. This means that the veloc-ity fluctuation in MAD98 is stronger than SANE98. Forillustration purposes, Figure 9 shows the poloidal speedof wind and jet at a certain time and radius as functionsof θ and ϕ . For MAD00, SANE98, and MAD98, the cho-sen time is t = 19 , , , , and , respectively.We can see that the peak of the poloidal speed of the realoutflow is about . c, . c, . c for MAD00, SANE98,and MAD98, respectively.Figure 10 shows the mass-flux-weighted poloidal veloc-ity of wind and jet as a function of radius. The blue andred dots are for jet and wind respectively, while the redline is the fitting curve for the wind. The poloidal veloc-ity of wind, combined with the result for SANE00 takenfrom Yuan et al. (2015), are described by: v p ( r ) = 0 . v k ( r ) , (14) v p ( r ) = 0 . v k ( r ) , (15) v p ( r ) = 0 . v k ( r ) , (16)0 Figure 6 . The mass fluxes of inflow, outflow, wind, andjet for MAD00 (top), SANE98 (middle) and MAD98(bottom) respectively. The values have been normalizedby the net rate ˙ M BH . v p ( r ) = 0 . v k ( r ) (17)for SANE00, MAD00, SANE98 and MAD98, respec-tively. Here v k ( r ) is the Keplerian speed at radius r .We can find that the poloidal speed of wind of MADmodels is larger than that of SANE models. This is alsothe case for jet. As we can see from the figure, the mass- flux-weighted velocities of jet at r g are v p , jet ≈ . c, (18) v p , jet ≈ . c, (19)for SANE98 and MAD98, respectively. This is becauseLorentz force is one of the main acceleration forces ofwind (Yuan et al. 2015) . In other words, the Poyntingflux of outflow is converted into the kinetic energy flux,as we will discuss in §3.6. Thus when the magnetic field isstronger, as in the MAD model, the poloidal speed will belarger. Another result we can find is that the black holespin has little effect on the velocity of wind, although theblack hole spin can strongly strengthen the BZ-jet. Thisis because, unlike the jet, wind is produced at relativelylarge radii where the effect of general relativity is weak.The mass-flux-weighted poloidal speed of wind de-creases with increasing radius. But we note that this doesnot mean that if we follow the trajectory of a given “windparticle”, its polodal speed will decrease when propagat-ing outward. The study in Yuan et al. (2015) has indi-cated that the poloidal speed actually increases outwardor at least keeps constant. This is because of the accelera-tion due to the gradient of magnetic and thermal pressureof the wind . They find that the speed of wind launchedfrom radius r can be described by v p ( r ) ≈ (0 . − . v k ( r ) for SANE00 (Yuan et al. 2015). The outward decrease ofthe mass-flux-weighted wind speed shown in eqs. 15-17is because, when the wind propagates outward, more andmore wind material whose velocity is becoming smallerwill join in, which makes the speed decreases with in-creasing radius. In fact, we can see in Figure 6 that thewind mass flux increases rapidly with radius. Therefore,the mass-flux-weighted wind velocity v ( r ) mainly reflectswind launched close to radius r .Unlike with the wind, the blue dots of Figure 10 showthat the mass-flux-weighted poloidal velocity of the jetincreases with radius. On the one hand, this is because ofthe strong acceleration of jet matter during its outwardmotion. On the other hand, unlike the case of wind,there is little low-speed matter joining into the jet whenit propagates outward. In this sense, the acceleration mechanism is similar to theBlandford & Payne (1982). Cui et al. (2020) and Cui & Yuan (2020) have studied thelarge-scale dynamics of wind after they are launched from the ac-cretion disk scale, with and without taking into account the ef-fect of magnetic field. The boundary conditions of wind are takenfrom the numerical simulations of wind launching from accretionflows combined with observational results. Both the black hole andgalaxy potentials are included. During the outward propagation,the enthalpy and rotational energy compensate for the increaseof gravitational potential. As a result, the wind can travel for along distance with roughly constant speed. This is true for windslaunched from both hot accretion flow and cold thin disks, althoughthe travel distance of cold wind is shorter than that of hot wind. Figure 7 . The mass flux of wind per unit θ as a function of θ for MAD00 (left), SANE98 (middle), and MAD98(right), respectively. The values are integrated over all ϕ and averaged from time , – , , , – , and – , for SANE98, MAD00 and MAD98 respectively, and have been normalized by ˙ M in ( r = 2 r g ) of thecorresponding models. The top and bottom panels are for r = 40 r g and r = 160 r g respectively. The positive valuedenotes the wind flux, negative value denotes the total inflow ˙ M in . The red dashed line represents the boundarybetween the BZ-jet and wind.3.5. The kinetic energy and momentum fluxes of windand jet
The kinetic energy and momentum fluxes of wind andjet are calculated by the following equations, ˙ E jet(wind) ( r ) = 12 (cid:90) γρ ( r, θ, ϕ ) v r ( r, θ, ϕ )( r + a cos θ ) sin θdθdϕ, (20) ˙ P jet(wind) ( r ) = (cid:90) γρ ( r, θ, ϕ ) v r ( r, θ, ϕ )( r + a cos θ ) sin θdθdϕ, (21)where γ = 1 / (cid:112) − v r . Figures 11 & 12 show the ra-dial profiles of the momentum and kinetic energy fluxesof wind and jet. From these results, combined withthe results of SANE00 shown in Yuan et al. (2015), wefind that both momentum and kinetic energy fluxes ofwind increase with the black hole spin. This is becausethe mass flux of outflow increases with spin (§3.3) al-though the poloidal velocity of wind is almost indepen-dent of spin (§3.4). As for the effect of magnetic field,we find that, at r g , the momentum flux of outflow(i.e., wind and jet) of MAD00/MAD98 is slightly largerthan that of SANE00/SANE98, while the kinetic en- ergy flux of MAD00/MAD98 is much larger than thatof SANE00/SANE98. This is because, at r g themass flux of wind of MAD00/MAD98 is less than that ofSANE00/SANE98 (eqs. 8-11 in §3.3) by about a factorof 1.5, while the velocity of the wind of MAD00/MAD98is greater than SANE00/SANE98 by a factor of about 3(eqs. 14-17 in §3.4). For jet, we can easily understandthe results by combining eqs. 12, 13, 18, & 19.When the black hole spin is large, both BZ-jet andwind will be present. It is interesting in this case to com-pare the momentum and energy fluxes of jet and wind.Such information is needed when we study the respec-tive roles of wind and jet in AGN feedback. The resultscan be obtained from Figures 11 & 12. For SANE98,the comparison of momentum and kinetic energy fluxesbetween jet and wind at r g are, ˙ P wind ≈ . P jet , (22) ˙ E jet ≈ . E wind . (23)For MAD98, the results are, ˙ P wind ≈ ˙ P jet , (24)2 Figure 8 . The poloidal speed of wind and jet in unit of speed of light averaged over ϕ and time from , – , , , – , and – , for MAD00 (left), SANE98 (middle) and MAD98 (right), respectively. The top andbottom panels are for r = 40 r g and r = 160 r g . The red dashed line represents the boundary between jet and wind. ˙ E jet ≈ E wind . (25)From these results, we can see that even in the case ofrapidly spinning black hole where the jet is supposed tobe the strongest, no matter the accretion flow is SANE orMAD, the power of jet is only a factor of ∼ larger thanthat of wind, while the momentum of jet is even smalleror at most similar to that of wind. In this section, weonly consider the kinetic energy flux. When the totalenergy is considered, i.e., the Poynting flux and enthalpyare also included, as we will discuss in §3.6, the jet powerwill be larger than that of wind by a factor of ∼ and10 for SANE98 and MAD98, respectively (refer to eqs.40 and 41).In the case of a slowly rotating black hole, the windbecomes relatively more important compared to jet. Inthis case, the jet is mainly powered by the rotation ofthe accretion flow. Yuan et al. (2015) show that themomentum and energy fluxes of wind are much largerthan that of jet in the case of SANE00 (refer to theirFigures 11 & 12): ˙ P wind ≈
15 ˙ P jet , (26) ˙ E wind ≈ E jet . (27)These results have important implications in the study of AGN feedback, indicating the importance of windcompared to jet, which are often neglected in many AGNfeedback studies. We will discuss this issue in §5.3.6. The Poynting flux of BZ-jet
Following Tchekhovskoy et al. (2011) and White(2019), the total energy fluxes (as measured at infinity)of the BZ-jet is calculated as follows, ˙ E tot ( r ) = (cid:90) π (cid:90) π T rt ( r + a cos θ ) sin θdθdϕ (28)here T rt is a component of the stress-energy tensor T µν describing the radial flux of energy: T rt = ( b + u + p + ρ ) u r u t − b r b t . (29)Note that T rt represents the total energy transported bythe fluid and the magnetic field, including the rest massenergy of the gas. The positive sign of ˙ E tot means the en-ergy flux is inward. Given that rest mass energy doesn’tplay a role in AGN feedback, following Sadowski et al.(2013), we define the following “energy flux” by eliminat-ing the rest mass energy, ˙ e ( r ) = − T rt − ρu r . (30)Positive values of ˙ e ( r ) mean energy is lost from the sys-tem. The magnetic flux that threads the hemisphere of3 Figure 9 . The poloidal speed of wind and jet as a function of θ at a given time of for MAD00 (left; t = 19 , ),SANE98 (middle; t = 61 , ), and MAD98 (right; t = 15 , ) respectively. Different color denotes different ϕ . Thetop and bottom panels are for 40 r g and 160 r g .black hole horizon, Φ BH , is calculated as follows, Φ( r ) = 12 (cid:90) π (cid:90) π √ π | B r | ( r + a cos θ ) sin θdθdϕ. (31)A dimensionless magnetic flux that is normalized withthe mass accretion rate is φ ( r ) = Φ (cid:112) ˙ M cr g . (32)The energy extraction efficiency, or the energy outflow ef-ficiency, η , is defined as the energy return rate to infinitydivided by the time-average rest-mass accretion power: η = ˙ M ( r H ) c − ˙ E tot ( r H ) (cid:104) ˙ M ( r H ) c (cid:105) ×
100 per cent . (33)Here ˙ M ( r H ) ≡ ˙ M BH is the mass flux at the black holehorizon r H calculated by eq. 7. Figure 13 shows the timeevolution of mass flux ˙ M ( r H ) (top panel), energy outflowefficiency (middle panel), and dimensionless magneticflux (bottom panel) for MAD00, SANE98, and MAD98.They all strongly fluctuate with time, especially the twoMAD models. For example, the energy outflow efficiencyof MAD98 ranges from 3 to 0. The time-averaged valuesof some physical quantities for SANE98 and MAD98 are given in Table 4. The dimensionless magnetic fluxes φ are ∼ and 40. The energy outflow efficiency of MAD00is zero, as expected. The time-averaged energy outflowefficiency η of SANE98 and MAD98 models are and . We can see that for MAD98, the variability of η and φ are synchronous, and both of them are larger com-pared to SANE98, which is consistent with our expecta-tion (e.g., Narayan et al. 2003; Tchekhovskoy et al. 2011).The values of η and φ for MAD98 are also roughly quanti-tatively consistent with the result of Tchekhovskoy et al.(2011). While SANE98 has a larger value of φ than isoften seen in SANE simulations in the literature (usuallyless than , see Porth et al. 2019), the lack of short-termvariability in this quantity clearly distinguishes it froma true MAD state. Such large values of φ occur wheninitial field configurations like ours are run beyond thecommon stopping point of t = 10 , or , (Whiteet al. 2020).Tchekhovskoy et al. (2010, see also BZ77) give the BZ-jet power as, P BZ = κ π Ω Φ f (Ω H ) , (34)where κ is a numerical constant that depends on thegeometry of the magnetic field, we adopt κ = 0 . f l u x - w e i g h t p o l o i d a l v e l o c i t y a=0 MAD0 25 50 75 100 125 150 175 200R / Rg0.00.10.20.30.40.5 f l u x - w e i g h t p o l o i d a l v e l o c i t y a=0.98 SANE0 25 50 75 100 125 150 175 200R / Rg0.00.10.20.30.40.5 f l u x - w e i g h t p o l o i d a l v e l o c i t y a=0.98 MAD Figure 10 . Radial profiles of mass-flux-weighted and timeaveraged poloidal velocity of wind (red dots) and jet(blue dots) for MAD00, SANE98, and MAD98.
Table 4 . Time-averaged quantities of BZ-jet
Model ˙ M φ P P P BZ η SANE98 0.91 23.81 0.275 0.29 35.3%MAD98 47.27 40.53 42.60 44.22 110% in this paper. Ω H = ac/ r H is the angular velocity ofthe black hole horizon, Φ BH is the magnetic flux thread-ing the hemisphere of black hole horizon, and f (Ω H ) is Figure 11 . The radial profiles of the momentum fluxes ofwind and jet. The values have been normalized by ˙ M ( r =2 r g ) c of MAD00, SANE98 and MAD98 respectively.a modifying factor for high spin a , which is f (Ω H ) ≈ . H r g /c ) − . H r g /c ) (Tchekhovskoy et al.2010). Using this formula, combined with the magneticflux from our simulations, we obtain P BZ = 0 . , . for SANE98 and MAD98 respectively. To compare thesevalues with our simulation results, we have also calcu-lated the Poynting flux of the BZ jet based on the fol-5 Figure 12 . The radial profiles of the kinetic energy fluxesof wind and jet. Those values have been normalized by ˙ M ( r = 2 r g ) c of the three models respectively.lowing equation, P P = (cid:90) π (cid:90) θ BZ − T rt (EM)( r + a cos θ ) sin θdθdϕ, (35)where T rt (EM) is the electromagnetic component ofstress-energy tensor describing the radial flux of thePoynting flux: T rt (EM) = b u r u t − b r b t . (36) −1 ̇ M SANE98 MAD98 MAD00 η g ̇ c]050100 ϕ Figure 13 . Time evolution of mass flux, energy extractionefficiency, and dimensionless magnetic flux for MAD00(green), SANE98 (blue), and MAD98 (red), respectively.
Figure 14 . The radial profiles of various components ofthe energy flux of the BZ-jet for SANE98 (upper panel)and MAD98 (bottom panel). They are the total energyflux F tot , the Poynting flux F p , the rest mass energysubtracted energy flux F tot − F m , the flux of enthalpy F u , the sum of the enthalpy and Poynting fluxes F u + F p ,and the rest mass energy flux F m = ρu r . Those valueshave been normalized by ˙ M c at the black hole horizonof the respective models.6We obtain P P = 0 . , . for SANE98 and MAD98respectively. These values are in good agreement withthat predicted by eq. (34), only slightly smaller. Wehave also calculated the “Poynting flux jet efficiency” η P defined by η P = (cid:104) P P (cid:105)(cid:104) ˙ M BH c (cid:105) (37)For MAD98, we obtain η P = 90 . , while the total en-ergy efficiency η = 99 . , so the Poynting flux of jet ac-counts for ∼ of the total energy flux. For SANE98, η P = 32 . while η = 35 . , so the Poynting flux of jetaccounts for ∼ of the total energy flux. Thus thetime-averaged Poynting powers of jet for SANE98 andMAD98 are, ˙ E jet = 0 .
32 ˙ M BH c , (38) ˙ E jet = 0 . M BH c , (39)respectively.To understand the conversion among various compo-nents of the energy fluxes in the jet, we have calculatedthe radial profiles of these components. The results areshown in Figure 14. They are the total energy flux(blue line) F tot ( ≡ − (cid:2) ( b + u + p + ρ ) u r u t − b r b t (cid:3) ) , thePoynting flux (red line) F p ( ≡ − ( b u r u t − b r b t )) , the restmass energy-subtracted energy flux (blue dotted line) F tot − F m ( ≡ − (cid:2) ( b + u + p + ρ ) u r u t − b r b t (cid:3) − ρu r ), theenthalpy flux (green line) F u ( ≡ − ( u + p ) u r u t ) , the sumof the enthalpy and Poynting flux (yellow line) F u + F p ,and the rest-mass energy flux (black line) F m ( ≡ ρu r ) .A positive value of flux means that the direction of theflux is outward. The following results can be found fromFigure 14:• The general “shape” of the radial profile of the totalenergy flux F tot is very similar to that of the restmass energy flux F m , implying that the change ofthe total energy is mainly due to the change of therest mass energy.• The rest mass energy-subtracted total energy flux F tot − F m roughly remains constant within r g .• Both the enthalpy and rest-mass energy fluxes ( F u and F m ) are pointing inward within the stagnationradius (Broderick & Tchekhovskoy 2015; Nakamuraet al. 2018) of the black hole, while the Poyntingflux F p is always pointing outward.• When the jet is initially launched from the blackhole, its energy is dominated by the Poynting flux F p . With the increase of radius, F p decreases while This value is slightly different from that given in Table 4. Thisis because, for simplicity, here the time interval used when we dothe time-average is shorter than that used in Table 4.
Figure 15 . The ratio of wind to jet energy flux forSANE98 and MAD98.the enthalpy flux F u increases. Their sum, F u + F p ,decreases with radius, and such a decrease is fasterthan that of F tot − F m . This means that the reducedelectromagnetic energy is not only converted intoenthalpy of the jet, but also converted into kineticenergy. This is consistent with the increase of jetkinetic energy with radius shown in Figure 12. Thedifference between the blue dashed line and yellowline is the magnitude of the increased kinetic en-ergy.• The increase of enthalpy flux with radius is largerthan that of the kinetic energy flux.To compare the energy flux between wind and jet, wehave calculated the ratio of the wind to jet total energyfluxes. The results for SANE98 and MAD98 are shownin Figure 15. We find that in general, the ratio increaseswith increasing radius, implying that the wind energyflux increases faster than that of jet. This is because,unlike the jet, which is produced only at small radii, windis produced in a large range of radius from small to largeradii. Thus the total energy flux of jet roughly remainsconstant but the total energy flux of wind increases withradius. At ∼ r g , the ratio of wind to jet total energyfluxes reaches F tot , wind /F tot , jet ≈ (40)and F tot , wind /F tot , jet ≈ (41)for SANE98 and MAD98, respectively. Referring to eqs.(23) & (25), we can see that the ratio of the wind tojet total energy flux is smaller compared to the ratio oftheir kinetic energy fluxes. Such a result implies that therelative contribution of kinetic energy to the total energyin jet is smaller than that in wind. SUMMARYIn a previous work (Yuan et al. 2015), based on threedimensional general relativity MHD (3D GRMHD) sim-7ulations of black hole hot accretion flows, we have inves-tigated the wind launched from the accretion flow. A“virtual particle trajectory” approach has been adoptedin that work, which can faithfully reflect the motion offluid elements and thus discriminate turbulence and realoutflows. This approach is superior to the time-averagedapproach of dealing with wind often adopted in the lit-erature. Due to the instantaneous nature of wind, thetime-average approach strongly underestimates the massflux of wind by an order of magnitude or so. The Yuanet al. (2015) paper only focuses on the SANE and non-spinning black holes (i.e., SANE00). In the present work,we extend the Yuan et al. (2015) work by consideringwind and jet in both SANE and MAD and rapidly spin-ning black holes. Since the magnetic field lines are quiteordered close to the rotation axis, we define jet as the re-gion occupied by field lines connected to the ergosphereof the rotating black hole; while outflow beyond the jetregion is called wind. The boundary between jet andwind is shown by, e.g., the red line in Figure 5. Threemodels have been considered, namely SANE98, MAD00,and MAD98, which denotes SANE accretion flow arounda black hole with spin a = 0 . , MAD accretion flowaround a black hole with a = 0 and a = 0 . , respec-tively. Our main results can be summarized as follows.• The radial profiles of mass flux of wind and jet invarious models are presented in Figure 6 and in eqs.(8-11). Close to the black hole, the wind becomesweaker in MAD compared to SANE, which is dueto the suppression of turbulence when the magneticfield becomes too strong. Further away from theblack hole, the wind mass flux becomes strongerin MAD compared to SANE, which is because atlarge radii turbulence is present and in general themagnetic field helps the formation of wind.• When the radius is not too large, the wind becomesstronger when the black hole spin becomes larger.• For rapidly spinning black holes, the mass flux ofjet for SANE98 and MAD98 can be described inunits of the black hole accretion rate by eqs. 12 &13.• The time and ϕ − averaged poloidal velocity of windand jet as a function of θ are presented in Fig-ure 8. It is much larger in the jet region than in thewind region, as expected. The mass-flux-weightedpoloidal velocity of wind and jet are shown in Fig-ure 10. The quantitative results are given by eqs.15-17 and eqs. 18-19 for wind and jet, respectively.Magnetic field significantly enhances the poloidalvelocity of wind and jet, as expected, but the effectof black hole spin on the wind velocity is weak. • The radial profiles of the momentum and kineticenergy fluxes of wind and jet are shown in Fig-ures 11 & 12. Both the momentum and kineticenergy fluxes of wind increase with the black holespin and the strength of magnetic field. The val-ues of momentum and kinetic energy fluxes of windcan be easily calculated by combining the above-mentioned mass flux and poloidal velocity results.• The comparisons of the fluxes of momentum andkinetic energy between wind and jet for SANE98and MAD98 are described by eqs. 22–25. Evenin these two cases with extremely rapidly rotat-ing black hole where the jet is supposed to thestrongest, the momentum flux of wind is larger orsimilar to that of jet; the kinetic energy flux of jetis only stronger than that of wind by a factor of ∼ . In the case of a non-spinning black hole suchas SANE00, both the momentum and kinetic fluxesof wind will be significantly larger than that of jet(eqs. 26–27). When other energy components suchas enthalpy and Poynting flux are included, the jetpower will be larger than that of wind by a factorof ∼ and 10 for SANE98 and MAD98. Theseresults have important implications for the studyof AGN feedback, when we consider the respectiverole of wind and jet. We will discuss this issue in§5.• The time-averaged Poynting power of jet forSANE98 and MAD98 is described by eqs. 38 and39. As shown by Figure 14, the power of jet is dom-inated by Ponyting flux close to the black hole andis gradually converted into enthalpy and kinetic en-ergy when the jet propagates outward. DISCUSSION: THE RELATIVE IMPORTANCEOF WIND AND JET IN AGN FEEDBACKIn many papers studying the effect of hot mode (alsoradio or jet or maintenance mode) AGN feedback ingalaxy evolution, the authors only take into account jetsbut neglect wind. Our results suggest that wind shouldalso be included, and they may even be more importantthan jet.First, the momentum of wind is in general larger thanthat of jet. In the study of AGN feedback, some impor-tant aspects, including the growth of the black hole mass,the determination of the AGN mass accretion rate, andthe degree of suppression of star formation in the centralregion of the galaxy, are mainly controlled by the mo-mentum rather than energy feedback (e.g., Ostriker et al.2010). Even in the case that energy feedback is more rel-evant than momentum feedback, since the power of jetis larger than wind only by a factor of – even in theextreme case of a very rapidly spinning black hole, and8since the opening angle of wind is much larger than thatof jet, which makes the energy deposition efficiency ofwind to the interstellar medium potentially much largerthan that of jet, the role of wind must be included.On the other hand, in the extreme case of rapidly spin-ning black holes, since the jet power is much larger thanwind and the jet is well collimated, it can easily penetratethrough the galaxy and propagate into much further dis-tance than wind. In this way, the jet is able to heat thecircum-galactic or intergalactic medium, which is hardfor wind to do so.In the above discussions, one caveat is that the com-parisons of momentum and power of wind and jet areconducted at r = 200 r g due to limitation of our GRMHDsimulations. It is not clear whether the relation roughlyremains correct or not at much larger radii when the jetpropagates deeply into the galaxy. In the case of mo-mentum, we can see from Figure 11 that the flux of windseems to keep increasing with increasing radius while theflux of jet seems to saturate. This will make the role ofwind relatively more important.ACKNOWLEDGEMENTSHY thanks Dr. Z. Gan for his help in the initial setupof the simulations. HY and FY are supported in partby the National Key Research and Development Pro-gram of China (grant 2016YFA0400704), the NaturalScience Foundation of China (grant 11633006), and theKey Research Program of Frontier Sciences of CAS (No.QYZDJSSW-SYS008). YFY is supported by the Na-tional Natural Science Foundation of China (Grant No.11725312, 11421303). This work has used the High Per-formance Computing Resource in the Core Facility forAdvanced Research Computing at Shanghai Astronomi-cal Observatory. APPENDIX A. KERR–SCHILD COORDINATES
Since our simulation output data is in Kerr–Schild co-ordinates, we need to know the tetrad carried by the lo-cally non-rotating frame (LNRF) observer in Kerr–Schildcoordinates. We can then use the tetrad to convert thesimulation data into physical quantities in LNRF coor-dinates. We know that in Boyer–Lindquist coordinatesthe Kerr metric form is ds = − (cid:18) − r Σ (cid:19) dt + Σ∆ dr + Σ dθ + A sin θ Σ dϕ − ar sin θ Σ dϕdt (A1)with the definitions Σ = r + a cos θ , ∆ = r − r + a , A = ( r + a ) − a ∆ sin θ Both r and a have units of the black hole mass M . In contravariant form we have (cid:18) ∂∂s (cid:19) = − A Σ∆ (cid:18) ∂∂t (cid:19) − ar Σ∆ (cid:18) ∂∂t (cid:19) (cid:18) ∂∂ϕ (cid:19) + ∆Σ (cid:18) ∂∂r (cid:19) + 1Σ (cid:18) ∂∂θ (cid:19) + ∆ − a sin θ Σ∆sin θ (cid:18) ∂∂ϕ (cid:19) (A2)The frame carried by the LNRF observer has basis vec-tors: e µ ( t ) = (cid:114) A Σ∆ (cid:18) , , , arA (cid:19) (A3) e µ ( r ) = (cid:114) ∆Σ (0 , , , (A4) e µ ( θ ) = (cid:114)
1Σ (0 , , , (A5) e µ ( ϕ ) = (cid:114) Σ A θ (0 , , , (A6)in contravariant form and e ( t ) µ = (cid:114) Σ∆ A (1 , , , (A7) e ( r ) µ = (cid:114) Σ∆ (0 , , , (A8) e ( θ ) µ = √ Σ(0 , , , (A9) e ( ϕ ) µ = (cid:114) A Σ sin θ (cid:18) − arA , , , (cid:19) (A10)in covariant form.In Kerr–Schild coordinates, the line element of the Kerrspace time is ds = − (cid:18) − r Σ (cid:19) dt + (cid:18) r Σ (cid:19) drdt + (cid:18) r Σ (cid:19) dr + Σ dθ + sin θ (cid:20) Σ + a (cid:18) r Σ (cid:19) sin θ (cid:21) dϕ − (cid:18) ar sin θ Σ (cid:19) dϕdt − a (cid:18) r Σ (cid:19) sin θdϕdr (A11)In contravariant form, it is (cid:18) ∂∂s (cid:19) = − (cid:18) − r Σ (cid:19) (cid:18) ∂∂t (cid:19) + 4 r Σ (cid:18) ∂∂t (cid:19) (cid:18) ∂∂r (cid:19) + ∆Σ (cid:18) ∂∂r (cid:19) + 2 a Σ (cid:18) ∂∂t (cid:19) (cid:18) ∂∂ϕ (cid:19) + 1Σ (cid:18) ∂∂θ (cid:19) + 1Σ sin θ (cid:18) ∂∂ϕ (cid:19) (A12)9The transformation matrix from Boyer–Lindquist co-ordinates to Kerr–Schild coordinates is: ∂x µ (KS) ∂x ν (BL) = r ∆ a ∆ (A13)The transformation matrix from Kerr–Schild coordi-nates coordinates to Boyer–Lindquist coordinates is: ∂x µ (BL) ∂x ν (KS) = − r ∆ − a ∆ (A14)Through the vector transformation formula: e µ ( α ) (KS) = ∂ x µ (KS) ∂ x ν (BL) e ν ( α ) (BL) (A15)and e ( α ) µ (KS) = ∂ x ν (BL) ∂ x µ (KS) e ( α ) ν (BL) (A16)It can be obtained that in Kerr–Schild coordinates, theLocally Non-Rotating Frame observers carry the follow-ing frames respectively: In contravariant form, e µ ( t ) = (cid:114) A Σ∆ (cid:18) , , , arA (cid:19) (A17) e µ ( r ) = (cid:114) ∆Σ (cid:18) r ∆ , , , a ∆ (cid:19) (A18) e µ ( θ ) = (cid:114)
1Σ (0 , , , (A19) e µ ( ϕ ) = (cid:114) Σ A θ (0 , , , (A20)in covariant form, e ( t ) µ = (cid:114) Σ∆ A (cid:18) , − r ∆ , , (cid:19) (A21) e ( r ) µ = (cid:114) Σ∆ (0 , , , (A22) e ( θ ) µ = √ Σ(0 , , , (A23) e ( ϕ ) µ = (cid:114) A Σ sin θ (cid:18) − arA , ar A ∆ − a ∆ , , (cid:19) . (A24)We can use this frame to convert the simulation datainto LNRF coordinates. For instance, we can use the be-low formulas to have the velocity vector in the LNRF co-ordinates. Athena++ outputs velocity components ν i = ( ν , ν , ν ) in Kerr–Schild-like coordinates where the timebasis vector is orthogonal to hypersurfaces of constanttime and has unit length. Then ν j = g ij ν i for i, j =1 , , ; γ = √ ν i ν i , α = (cid:112) − /g ; u t (KS) = γ/α ; u r (KS) = ν − γαg (KS) ; u θ (KS) = ν − γαg (KS) ; u ϕ (KS) = ν − γαg (KS) ; u a (LNRF) = e ( a ) µ (KS) u µ (KS); (A25) v r = u r (LNRF) u t (LNRF) ; v θ = u θ (LNRF) u t (LNRF) ; and v ϕ = u ϕ (LNRF) u t (LNRF) .REFERENCES Abramowicz, M. A., Czerny, B., Lasota, J. P., & Szuszkiewicz, E.1988, ApJ, 332, 646Abramowicz, M. A., Igumenshchev, I. V., Quataert, E., &Narayan, R. 2002, ApJ, 565, 1101Beckwith, K., Hawley, J. F., & Krolik, J. H. 2008, ApJ, 678, 1180Begelman, M. C. 2012, MNRAS, 420, 2912Blandford, R. D., & Begelman, M. C. 1999, MNRAS, 303, L1—. 2004, MNRAS, 349, 68Blandford, R. D., & Payne, D. G. 1982, MNRAS, 199, 883Blandford, R. D., & Znajek, R. L. 1977, MNRAS, 179, 433Broderick, A. E., & Tchekhovskoy, A. 2015, ApJ, 809, 97Bu, D.-F., Yuan, F., Gan, Z.-M., & Yang, X.-H. 2016a, ApJ, 818,83—. 2016b, ApJ, 823, 90Chen, L., & Zhang, B. 2021, ApJ, 906, 105Cheung, E., Bundy, K., Cappellari, M., et al. 2016, Nature, 533,504Cui, C., & Yuan, F. 2020, ApJ, 890, 81Cui, C., Yuan, F., & Li, B. 2020, ApJ, 890, 80Einfeldt, B. 1988, SIAM Journal on Numerical Analysis, 25, 294Fabian, A. C. 2012, ARA&A, 50, 455Fishbone, L. G., & Moncrief, V. 1976, ApJ, 207, 962Gammie, C. F., McKinney, J. C., & Tóth, G. 2003, ApJ, 589, 444Ghosh, P., & Abramowicz, M. A. 1997, MNRAS, 292, 887Gu, W.-M. 2015, ApJ, 799, 71Homan, J., Neilsen, J., Allen, J. L., et al. 2016, ApJL, 830, L5Jiang, Y.-F., Stone, J. M., & Davis, S. W. 2014, ApJ, 796, 106Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511Livio, M., Ogilvie, G. I., & Pringle, J. E. 1999, ApJ, 512, 100Ma, R.-Y., Roberts, S. R., Li, Y.-P., & Wang, Q. D. 2019,MNRAS, 483, 5614McKinney, J. C., & Gammie, C. F. 2004, ApJ, 611, 977McKinney, J. C., Tchekhovskoy, A., & Bland ford, R. D. 2012,MNRAS, 423, 3083Mosallanezhad, A., Bu, D., & Yuan, F. 2016, MNRAS, 456, 2877Naab, T., & Ostriker, J. P. 2017, ARA&A, 55, 59Nakamura, M., Asada, K., Hada, K., et al. 2018, ApJ, 868, 146Narayan, R., Igumenshchev, I. V., & Abramowicz, M. A. 2000,ApJ, 539, 798—. 2003, PASJ, 55, L69Narayan, R., SÄ dowski, A., Penna, R. F., & Kulkarni, A. K.2012, MNRAS, 426, 3241Narayan, R., & Yi, I. 1994, ApJL, 428, L13—. 1995, ApJ, 452, 710Ostriker, J. P., Choi, E., Ciotti, L., Novak, G. S., & Proga, D.2010, ApJ, 722, 642Park, J., Hada, K., Kino, M., et al. 2019, ApJ, 871, 257Penna, R. F., Kulkarni, A., & Narayan, R. 2013, A&A, 559, A116Porth, O., Chatterjee, K., Narayan, R., et al. 2019, ApJS, 243, 26Pringle, J. E. 1981, ARA&A, 19, 137Quataert, E., & Gruzinov, A. 2000, ApJ, 539, 8090