Active galactic nuclei feedback at parsec scale
aa r X i v : . [ a s t r o - ph . H E ] J u l Draft version July 19, 2019
Preprint typeset using L A TEX style emulateapj v. 12/16/11
ACTIVE GALACTIC NUCLEI FEEDBACK AT PARSEC SCALE
De-Fu Bu , Xiao-Hong Yang Draft version July 19, 2019
ABSTRACTWe perform simulations to study the effects of active galactic nuclei (AGNs) radiation and windfeedback on the properties of slowly rotating accretion flow at parsec scale. We find that whenonly radiative feedback is considered, outflows can be produced by the radiation pressure due toThomson scattering. The mass flux of outflow is comparable to that of inflow. Although strongoutflow is present, the luminosity of AGN can be easily super-Eddington. When wind feedback isalso taken into account, the mass flux of outflow does not change much. Consequently, the luminosityof the central AGN can still be super-Eddington. However, observations show that the luminosity ofmost AGNs is sub-Eddington. Some other mechanisms are needed to reduce the AGNs luminosity.Although the mass outflow rate does not change much by wind feedback, other properties of outflow(density, temperature, velocity and kinetic power) can be significantly changed by wind feedback.In the presence of wind feedback, the density of outflow becomes significantly lower, temperature ofoutflow becomes significantly higher, velocity of outflow is increased by one order of magnitude, thekinetic power of outflow is increased by a factor of 40 − Subject headings: accretion, accretion disks – black hole physics – galaxies: active – galaxies: nuclei INTRODUCTIONIt is widely believed that a super massive black hole(SMBH) is located at the center of almost each galaxy(see Heckman & Best 2014 for a review). The SMBHis believed to co-evolve with its host galaxy (Hopkins etal. 2005; Cattaneo et al. 2009). There are many obser-vations showing that the mass of the SMBH is tightlycorrelated with the properties of its host galaxy (Fabian2012). One of the example is that the mass of the SMBHis a constant fraction of the stellar bulge mass (H¨aring& Rix 2004). The even more remarkable evidence is thatthe mass of the SMBH is tightly correlated with the ve-locity dispersion of the host galaxy’s central bulge (e.g.,Magorrian et al. 1998; Ferrarese & Merritt 2000; Geb-hardt et al. 2000; Tremaine et al. 2002; G¨ultekin et al.2009; Kormendy & Ho 2013). AGNs feedback plays animportant role in the evolution of its host galaxy (Fabian2012; king & Pounds 2015).There are three kinds of output by an AGN, namely,radiation, wind and jet. The radiation from an AGNcan Compton heat/cool the interstellar medium (Ciotti& Ostriker 1997, 2001, 2007; Ciotti et al. 2009). Itis found that when the accretion flow is in hot mode,the gas in the region around the Bondi radius can beheated up to be above the virial temperature. Conse-quently, winds can be generated due to Compton heat-ing around Bondi radius (Bu & Yang 2018, 2019a; Yang& Bu 2018). For black hole X-ray binaries (BHBs), it isalso found that wind can be driven from accretion disk atlarge radii due to the Compton heating (Higginbottomet al. 2018, 2019). The line force due to the absorptionof UV photons by the low ionization gas can also launch Key Laboratory for Research in Galaxies and Cosmol-ogy, Shanghai Astronomical Observatory, Chinese Academyof Sciences, 80 Nandan Road, Shanghai 200030, China;[email protected] Department of physics, Chongqing University, Chongqing,400044; [email protected] strong wind (e.g., Proga et al. 2000; Murray et al. 1995;Liu et al. 2013; Nomura et al. 2016). AGN winds havelarge opening angle and can interact with the gas sur-rounding the AGN very efficiently. Gas surrounding anAGN can be directly blown away by winds, which willaffect the star formation rate and the AGN fueling (e.g.,Ostriker et al. 2010; Du et al. 2014; Weinberger et al.2017, 2018; Ciotti et al. 2017; Yuan et al. 2018; Bu &Yang 2019b). Due to the very small opening angle of jet,it is believed that jet can not effectively interact with gasin the region around an AGN. Jet should play an impor-tant role in galaxy cluster scale (Guo 2016; Guo et al.2018).Kollmeier et al. (2016) studied the luminosity distri-bution of a sample of 407 AGNs in the redshift range z ∼ . −
4. It is found that the distribution of theestimated Eddington ratios can be well described as log-normal, with a peak at L bol /L Edd ≈ / L bol and L Edd being bolometric and Eddington luminosity, re-spectively). In other words, almost all the AGNs inthis smaple have sub-Eddington luminosity. Later, Stein-hardt & Elvis (2010) used a even larger sample composedof 62185 quasars from the Sloan Digital Sky Survey Datato study the luminosity distributions of quasars. Theymade the same conclusion as that in Kollmeier et al.(2016). The fuelling gas at galactic scale has a very broadmass distribution (Liu et al. 2013). In other words, thegalactic scale inflow can provide sufficient gas to makethe AGN’s luminosity super-Eddington. Also, the blackhole slim accretion disk model predicts that the rotat-ing accretion disk can have well higher luminosity above L Edd (Abramowicz et al. 1988; Ohsuga et al. 2005; Yanget al. 2014). What is the reason for the sub-Eddingtonaccretion of the AGNs?AGN radiative and wind feedback may be possibleways to solve the ‘sub-Eddington puzzle’. In the inner-most region very close the central black hole, the numer-ical simulations of slim disk (Ohsuga et al. 2005) with Bu & Yangradiative transfer show that the accretion disk can radi-ate well above L Edd . Therefore, the radiative feedbackin the innermost region can not be effective in reduc-ing the mass accretion rate (or equivalently luminosity).On the parsec scale, the effects of AGN radiative feed-back on the accretion flow have been studied (kurosawa& Proga 2009; Liu et al. 2013). It is shown that eventhough the line force can drive strong wind, the accre-tion rate (or luminosity) of AGN is still super-Eddington.The sub-Eddington puzzle can not be solved by radiativefeedback.The motivations of the present paper are twofold. Thefirst one is to study whether AGN wind feedback at par-sec scale can solve the sub-Eddington puzzle. The secondmotivation is to study the effects of wind feedback on theproperties of the accretion flow at parsec scale. The ac-cretion flow at parsec scale connects galactic inflow andthe accretion flow onto the central black hole. It is vitallyimportant to study the accretion flow in this region, be-cause it is in this region that the inflow rate to the centralblack hole is determined.In Bu & Yang (2019b), we studied the effects of feed-back of wind from low-luminosity AGN (LLAGN) onproperties of accretion flow at parsec scale. In the presentpaper, we study the high accretion rate luminous AGNcase by setting much higher gas density at the outerboundary. The properties of wind from LLAGN are quitedifferent from those of wind from luminous AGNs (seethe details in Sections 2.2 and 2.3). Therefore, we expectthat the effects of wind feedback by LLAGN should bequite different from those of wind feedback by luminousAGN.The structure of this paper is as follows. In section2, we present the numerical settings of simulations andthe physical assumptions. The results are presented insection 3; We summarize our results in Section 4. NUMERICAL METHODThe mass of the central black hole is set to be M BH =10 M ⊙ with M ⊙ being solar mass. We study the lowangular momentum accretion flow in the region 1500 r s ≤ r ≤ . × r s ( r s is Schwarzschild radius). In this paper,we only calculate the region above the midplane, we have0 ≤ θ ≤ π/
2. We set the specific angular momentum ofgas to equal to the Keplerian angular momentum at r c =900 r s , which is smaller than our inner radial boundary.By using the ZEUS-MP code (Hayes et al. 2006), werun two-dimensional hydrodynamic simulations in spher-ical coordinates ( r, θ, φ ). The equations describing theaccretion flow are as follows, dρdt + ρ ∇ · v = 0 , (1) ρ d v dt = −∇ p − ρ ∇ Φ + ρ F rad (2) ρ d ( e/ρ ) dt = − p ∇ · v + ˙ E (3)where, ρ is density, v is velocity and e is internal energy.Ideal gas equation is adopted p = ( γ − e with γ = 5 / F rad is the radiation pressure due to Thomson scattering.˙ E is the net heating rate per unit volume. For the radiation force term F rad , we only considerthe radiation pressure due to Thomson scattering. Thecentral AGN emits most of its energy in the region veryclose to the black hole. The inner radial boundary of oursimulations is significantly larger than the AGN emittingregion, therefore, we can treat the AGN emitting regionas a point source. The radiation force only has the radialcomponent, F rad , r = κ es c [ F X + F UV ] (4) F X and F UV are the X-ray flux and UV flux from thecentral AGN, respectively. We set κ es = 0 . g − .We define the ionization parameter ξ as, ξ = 4 πF X n = f X L bol exp( − τ X ) /nr (5) f X is the fraction of the total luminosity ( L bol ) in X-ray. τ X = R r ρκ es dr is the X-ray scattering optical depth.The net heating rate ˙ E in Equation (3) includes Comp-ton heating/cooling rate ( S c), photoionization heating-recombination cooling rate ( G x), bremsstrahlung coolingrate ( B r) and line cooling rate ( L line ). The Comptonheating/cooling rate is: S c = 8 . × − n (T X − T) ξ erg cm − s − (6) T X is the Compton temperature of the X-ray photonsemitted by the central AGN, T is the temperature of theaccreting gas. n = ρ/ ( µm p ) is number density of gas,with µ = 0 . m p being mean molecular weight andproton mass, respectively. The sum of photoionizationheating-recombination cooling rate is: G x = 1 . × − n ξ / T − / (1 − T / T X ) erg cm − s − (7)The bremsstrahlung cooling rate is: B r = 3 . × − n √ T erg cm − s − (8)The line cooling rate is: L line = 1 . × − exp( − . × /T ) /ξ/ √ T +10 − erg cm − s − (9)In this paper, in all of the models, we take into accountthe radiative feedback by the central AGN. The radiativefeedback includes radiation force F rad , r (Equation (4)),Compton heating/cooling (Equation (6)) and photoion-ization heating/recombination cooling (Equation (7)). Insome models, in addition to radiative feedback, we alsotake into account wind feedback. The detailed propertiesof wind are introduced in Sections 2.2 and 2.3.2.1. Critical accretion rate
Observations of BHBs show that when accretion flowluminosity is lower than 2% L Edd , the BHBs is in hardstate and the accretion flow is hot accretion flow (Mc-Clintock & Remillard 2006). When luminosity is higherthan 2% L Edd , the BHBs is in soft state and the accretiondisk model is standard thin disk (Shakura & Sunyaev1973). We use 2% L Edd as the boundary between LLAGNand luminous AGN. Therefore, we defined LLAGN asAGN which has luminosity lower than 2% L Edd . The ac-cretion disk model for LLAGN is hot accretion flow (seeGN feedback 3Yuan & Narayan 2014 for reviews). Correspondingly,we define luminous AGN as AGN which has luminos-ity higher than 2% L Edd . The accretion disk model forluminous AGN is standard thin disk.The mass accretion rate corresponding to the criticalluminosity is, ˙ M c = 2% L Edd ǫ cold c (10) ǫ cold is the radiative efficiency of a cold thin disk. In thispaper, we set ǫ cold = 0 . M in )is just 2 times the accretion rate at the black hole hori-zon. Therefore, in order to determine the accretion modeof the central AGN, we directly compare ˙ M in with thecritical accretion rate ˙ M c . When ˙ M in ≥ ˙ M c , the stan-dard thin disk is operating in the central AGN. When˙ M in < ˙ M c , the central AGN is a LLAGN and hot accre-tion flow operates.2.2. Cold disk feedback
Wind
Strong winds can be generated from a cold accretiondisk. There are many observational evidences showingthat strong winds are present in both luminous AGNs(e.g., Crenshaw et al. 2013; Tombesi et al. 2010; Goffordet al. 2015; King & Pounds 2015; Liu et al. 2014; Liuet al. 2015; Wang et al. 2016; Sun et al. 2017; Park etal. 2018; He et al. 2019) and the soft state of black holeX-ray binaries (e.g., Neilsen & Homan 2012; Homan etal. 2016; D´ıaz Trigo & Boirin 2016; You et al. 2016).In this paper, we set the mass flux and velocity of windgenerated by a cold thin disk according to the observa-tions introduced in Gofford et al. (2015). The reasonsfor adopting the observational results in Gofford et al.(2015) are as follows. First, the winds are found to begenerated at a distance ∼ − r s from the blackhole (Gofford et al. 2015). The inner boundary of oursimulation domain is 1500 r s . Therefore, it is reasonableto directly adopt the observational result in our simu-lations. Second, Gofford et al. (2015) give formulas todescribe the properties of wind. Therefore, it is mucheasier to set the wind properties in simulations.Gofford et al. (2015) analyzed the properties of windin a sample composed of 51 Suzaku-observed AGNs. Itis found that the wind velocity is a function of the bolo-metric luminosity of the AGNs, v wind = 2 . × (cid:18) L bol erg s − (cid:19) km s − (11)Observations indicate that the wind velocity saturates ata value of 10 km s − . Accordingly, in our simulations,we set a upper limit of wind velocity at this value. Themass flux of wind is found to be approximately equalto the accretion rate onto the black hole (Gofford et al.2015). According to the observation results, we set themass flux of wind to be equal to the black hole accretion rate. Once we get the mass accretion rate at the innerboundary ( ˙ M in ) of the simulation domain, we directly setthat the both the wind mass flux and black hole accretionrate as, ˙ M BH = ˙ M wind = 12 ˙ M in (12)The final question to describe wind is its opening an-gle. In the sample of AGNs analyzed by Gofford et al.(2015), there are some radio-loud AGNs. The propertiesof the wind in these radio-loud AGNs are also analyzedby Tombesi et al. (2014). Because the inclination angleof jet in these radio-loud AGNs is known, it is easy toknow the angle at which wind is present. It is reportedin Tombesi et al. (2014) that the winds are detected ata wide range of jet inclination angles ( ∼ ◦ to 70 ◦ , seealso Mehdipour & Costantini 2019). In this paper, wesimply inject wind at the inner boundary in the region10 ◦ ≤ θ ≤ ◦ . We also assume that the density of windis uniformly distributed in this region.2.2.2. Radiation
For the luminous AGN ( ˙ M in ≥ ˙ M c ), we can ob-tain the black hole accretion rate according to Equation(12). The bolometric luminosity is obtained by using L bol = ǫ cold ˙ M BH c . We study low-angular momentumgas accretion in this paper. The circularization radiusof the gas r c is smaller than the inner radial boundaryof the simulation. We assume that once the gas fallsinto the inner boundary. A cold thin disk will form at r c . Observations show that a hot corona is present atapproximately 10-20 gravitational radius (Reis & Miller2013; Uttley et al. 2014; Chainakun et al. 2019), whichemits X-ray photons. In this paper, we also assume thata spherical X-ray emitting hot corona is present in theinner most region above and below the thin disk. Weassume the X-ray flux is spherically distributed. The X-ray luminosity is L X = f X L bol . The cold thin disk isassumed to emit UV photons. As done in Proga et al.(2000) and Liu et al. (2013), we assume that the flux ofUV photons is F UV = 2 f UV cos( θ ) L bol exp( − τ UV ) / πr ,with f UV being the ratio of UV luminosity to the bolo-metric luminosity. We set τ UV = τ X . We neglect theemissions in radio, optical and infrared bands, and wehave f X + f UV = 1. We assume f X = 0 .
05 in this pa-per. The Compton temperature of the X-ray photons isassumed to be T X = 8 × K (Sazonov et al. 2004).2.3. Hot accretion flow feedback
Wind
It is very hard to directly observe wind from a hot ac-cretion flow through blue-shifted absorption lines. Thereason may be that hot accretion flow is fully ionized andno absorption line is present. There are some indirect ob-servational evidences for the presence of wind producedby hot accretion flow (e.g., Crenshaw & Kramemer 2012;Wang et al. 2013; Weng & Zhang 2015; Cheung et al.2016; Homan et al. 2016; Almeida et al. 2018; Parket al. 2019; Ma et al. 2019; Riffel et al. 2019). Thedetailed properties of wind can not be given by observa-tions. However, we have much better theoretical under-standing of the properties of wind from hot accretion flowdue to the numerical simulations (e.g., Tchekhovskoy et Bu & Yangal. 2011; Yuan et al. 2012, 2015; Narayan et al. 2012;Li et al. 2013; Almeida & Nemmen 2019) and analyticalworks (e.g., Li & Cao 2009; Gu et al. 2009; Cao 2011; Wuet al. 2013; Gu 2015; Nemmen & Tchekhovskoy 2015; Yiet al. 2018; Kumar & Gu 2018, 2019).We set the wind properties from hot accretion flowaccording to the numerical simulation work of Yuan etal. (2015). In Bu & Yang (2019b), we have introducedthe details of how to set physical properties of wind. Forconvenience, we briefly introduced here. Simulations ofviscous hot accretion flow (e.g., Yuan et al. 2015) findthat wind is present outside 10 r s . Because wind can takeaway mass, mass accretion rate decreases inwards (e.g.,Yuan et al. 2012, 2015),˙ M inflow ( r ) = ˙ M R H (cid:18) rR H (cid:19) . , for 10 r s < r < R H (13)with R H being the outer boundary of the viscous hotaccretion flow; ˙ M R H is the mass inflow rate at R H . Bu etal. (2013) studied low angular momentum hot accretionflow. In Bu et al. (2013), we find that when gas fromlarge radii has very small angular momentum, the gaswill free fall from large radii to the circularization radius( r c ). In the region outside r c , no outflows can be formed(Bu et al. 2013). When gas arrives at r c , a viscous hotaccretion flow will form and winds can be generated bythe viscous hot accretion flow. Therefore, in the presentpaper, we set R H = r c and ˙ M R H = ˙ M in . The accretionrate at black hole horizon is,˙ M BH = ˙ M R H (cid:18) r s R H (cid:19) . (14)The mass flux of wind is,˙ M wind = ˙ M in − ˙ M BH (15)The inner boundary of the simulation (1500 r s ) is veryclose to r c . Therefore, even we set R H = 1500 r s , accord-ing to equations (14) and (15), the mass flux of windcan not differ much from that in the case R H = r c . Wehave done calculations and find that if R H = 1500 r s ,the wind mass flux will be larger than that used in thepresent paper by a very small factor of 0.02. Therefore,we expect that the results will not be differ much if weset R H = 1500 r s .Comparing Equations (12) and (14), it is clear thatthere is a gap in the accretion rate (or luminosity) ofthe AGN. The maximum accretion rate for hot accretionflow is ∼ . L Edd / . . According to equation (22),when accretion rate is ∼ . L Edd / . , the radiativeefficiency is ∼ .
05. Correspondingly, the maximum lu-minosity of hot accretion flow is 0 . L Edd . The gap is0 . − L Edd .Observationally, there should be no gap in the lumi-nosity of AGNs. The reason may be as follows. Whenthe AGN is in the low luminosity phase, a thin disk istruncated at a radius
Rtr . Inside
Rtr , hot accretion flowoperates. With the increase of accretion rate,
Rtr de-creases (Liu et al. 1999; Gu & Lu 2000; Yuan & Narayan2014). The accretion rate onto the black hole is equal to˙ M in (10 r s /Rtr ) . . With the decrease of Rtr , the accre-tion rate of the black hole will gradually reach the value of ˙ M in . Therefore, there should be no gap in reality. Inthis paper, we do not consider the details of this point.However, we note that according to Figures (1), (4) and(5), the AGN spends most of their time at the luminousphase (with luminosity higher than 2% L Edd ). There areonly several snapshots at which AGN is in hot accretionflow mode. Therefore, the gap in luminosity of the AGNcan not affect the results much.The radial velocity of wind depends on the locationwhere it is produced (Yuan et al. 2015), v wind , r = 0 . v k ( R wind ) (16) v k ( R wind ) is the Keplerian velocity at wind launchingradius ( R wind ). According to Equation (15), we have, d ˙ M wind /dr = d ˙ M inflow /dr (17)By using Equations (13), (16) and (17), we can derive themass flux weighted radial velocity of wind. Given that R H = r c = 900 r s , the mass flux weighted wind radialvelocity is, v wind = 0 . v k (900 r s ) (18)with v k (900 r s ) being the Keplerian velocity at 900 r s . v θ of wind is found to be much smaller than both radial androtational velocities. Therefore, we set v θ of wind to bezero. Wind rotational velocity is also given by Yuan etal. (2015), V φ, wind = 0 . v k ( R H ) (19)Wind internal energy is found to be about 0.6 times thegravitational energy (Yuan et al. 2015), e wind = 0 . GM BH R H (20)The final issue is the opening angle of wind. Wind isfound to be distributed in the regions θ ∼ ◦ − ◦ and θ ∼ ◦ − ◦ (Yuan et al. 2015). In the present work,only the region above midplane is calculated. Hence,wind is injected into the computational domain within30 ◦ < θ < ◦ . We assume that density of wind is inde-pendent of θ . 2.3.2. Radiation
For a LLAGN ( ˙ M in < ˙ M c ), hot accretion flow is oper-ating. We assume that the hot accretion flow only emitsX-ray photons. In this case, f X = 1 and L X = L bol .We also assume that the X-ray photons are sphericallyemitted. Xie & Yuan (2012) find that the radiative ef-ficiency ( ǫ hot ) depends on black hole accretion rate andthe parameter δ which describes the fraction of the directviscous heating to electrons. The radiative efficiency isas follows (Xie & Yuan 2012): ǫ hot ( ˙ M BH ) = ǫ ( 100 ˙ M BH ˙ M Edd ) a , (21)˙ M Edd = L Edd / . c is the Eddington accretion rate; ǫ and a for the case of δ = 0 . ǫ , a ) = (1 . , .
65) if ˙ M BH ˙ M Edd . . × − ;(0 . , . . × − < ˙ M BH ˙ M Edd . . × − ;(0 . , .
12) if 3 . × − < ˙ M BH ˙ M Edd . . × − . (22)GN feedback 5When ˙ M BH ˙ M Edd > . × − , the radiative efficiency issimply set to be 0.1. The bolometric luminosity is L bol = ǫ hot ˙ M BH c . The flux of X-ray photons is F X = L bol exp( − τ X ) / πr .The Compton temperature ( T X ) of the X-ray photonsfrom hot accretion flow is calculated by Xie & Yuan(2017). It is found T X ∼ K. We set T X = 10 K,when the central engine is a LLAGN.2.4.
Initial and boundary conditions
Initially, we put gas with uniform density ( ρ ) and tem-perature ( T ) in the whole computational domain. Ourresolution is 140 ×
80. In r direction, in order to well re-solve the inner region, we employ logarithmically spacedgrid. Grids in θ direction are uniformly spaced.In the models with wind feedback, at the inner radialboundary, we inject wind into the computational domain.The injection θ angle and properties of the injected windare introduced above. For the θ angle at which no windis injected, we use outflow boundary conditions. Out-flow boundary conditions do not permit gas flows intothe computational domain. The physical variables in theghost zones are just copied from the first active zone. Inthe models without wind feedback, at the inner radialboundary, in the whole θ angle (0 ◦ ≤ θ ≤ ◦ ), we useoutflow boundary conditions.In order to continuously supply gas at the outer bound-ary, we set the outer radial boundary conditions as fol-lows. When we find that gas falls inward at the last ac-tive zone ( v r <
0) at a given θ angle, then at this angle,we inject gas into the computational domain. Density,temperature, specific angular momentum of the injectedgas are equal to those for the initial condition. If we find v r > θ angle, then at this angle, we employoutflow boundary conditions.At θ = 0, we use axis-of-symmetry boundary condi-tions. At θ = π/
2, reflecting boundary conditions areemployed. RESULTSAll of the models in this paper are summarized in Table1. In order to discriminate the outflow from the centralAGN (not resolved in this paper) and outflow found inthe simulation domain, we use “wind” to represent themass outflow from the central AGN, which is injectedinto the computational domain at inner boundary in sim-ulations with wind feedback. “Outflow” is used to denotethe outward moving gas found in our simulation domain.3.1.
Effects of wind feedback
We take the two models (D19T6rad andD19T6windrad) with ρ = 10 − g / cm and T = 2 × K as our fiducial models.Figure 1 shows the time evolution of luminosity of thecentral AGN in unit of L Edd . The top and bottom pan-els are for models D19T6rad and D19T6windrad, respec-tively. The horizontal solid line in each panel is the time-averaged value. The luminosity in both models fluctuatessignificantly. The highest luminosity in model D19T6radis ∼ L Edd . The lowest luminosity can be lower than1% L Edd . In model D19T6rad, the fluctuation is due tothe radiative feedback (radiation pressure) from the cen-tral AGN. The reason for the fluctuation is as follows.
Fig. 1.—
Time evolution of the luminosity of the central AGNnormalized by Eddington luminosity for models D19T6rad (top-panel) and D19T6windrad (bottom-panel). The horizontal axis istime in unit of free-fall timescale measured at the outer boundary.The horizontal solid line in each panel corresponds to the time-averaged (from t = 3 . When the luminosity is significantly higher than L Edd ,radiation pressure due to Thomson scattering can be sig-nificantly exceeding gravity. Outflows can be launchedand gas around the central AGN is taken away. Thenthe accretion rate (or luminosity) decreases and radiativefeedback becomes weak. The gas around the AGN willfall back and the AGN will be triggered again. The time-averaged luminosity for model D19T6rad is 28 . L Edd .The luminosity evolution patten in modelD19T6windrad is very similar as that in modelD19T6rad. The fluctuation of luminosity in modelD19T6windrad is due to both wind and radiative (pres-sure) feedback. The wind feedback does decrease theluminosity of the central AGN. The luminosity of thecentral AGN in model D19T6windrad has never beinghigher than 100 L Edd . The time-averaged luminosity inmodel D19T6windrad is 21 . L Edd . After taking intoaccount wind feedback, the decrease of luminosity is notsignificant.In order to study the detailed properties of the flow,we plot Figure 2. In this figure, we show the time-averaged (from 3.5 to 25 free-fall timescales measured atthe outer boundary) two-dimensional structure for mod-els D19T6rad (top-panels) and D19T6windrad (bottompanels). The left column panels plot logarithm density(color) and velocity vector (arrows). The right columnpanels plot logarithm temperature (color) and unit vec-tor of velocity (arrows). In model D19T6rad, inflow ex-ists in the region around midplane (65 ◦ < θ ≤ ◦ ). Bu & Yang TABLE 1Models and results
Model wind feedback ρ T L/L
Edd ˙ M ( r out ) / ˙ M Edd P K ( r out )(10 − g cm − ) 10 K ( L Edd )(1) (2) (3) (4) (5) (6) (7)D19T6rad OFF 1 2 28.5 190 0 . . M Edd = L Edd / . c ( c is speed of light). Col. 7: time-averaged mechanical energy flux (in unit of Eddington luminosity) of outflow measured at the outerboundary. The time-average are done from 3.5 to 25 free fall timescales measured at the outer boundary. Fig. 2.—
Two-dimensional time-averaged (from t = 3 . In the region θ < ◦ , outflow is present. The outflowvelocity increases with decreasing of θ . In this model,outflows are launched due to radiation (mainly UV pho-tons) pressure due to Thomson scattering. Because theUV photon flux from the central AGN has a sin θ de-pendence (see Section 2.2.2), the radiation force per unitmass increases with decrease of θ . Therefore, the outflowvelocity increases with decrease of θ . The maximum ve-locity of outflow is ∼ . c .In model D19T6windrad, we take into account wind feedback. From the bottom panel of Figure 1, we seethat the central AGN spends most of its time in luminousAGN phase. For the luminous AGN, wind is injected intothe computational domain in the region 10 ◦ < θ < ◦ (see Section 2.2.1). In addition to radiation pressure,outflows can also be driven by wind feedback. There areseveral differences between the properties of outflows inmodels D19T6rad and D19T6windrad. First, in modelD19T6windrad, the outflows have much larger openingangle (0 ◦ < θ < ◦ ). Second, the density in the outflowGN feedback 7 Fig. 3.—
Radial profiles of time-averaged (from 3.5 to 25 free-falltimescales measured at the outer boundary) mass inflow rate (solidlines) and outflow rate (dashed lines) for models D19T6rad (blacklines) and D19T6windrad (red lines). The mass fluxes are in unitof Eddington accretion rate ( ˙ M Edd = L Edd / . c ). Fig. 4.—
Same as Figure 1, but for models D19T7rad (top-panel)and D19T7windrad (bottom-panel). region in model D19T6windrad is much lower. The rea-son is that wind from the central AGN can directly blowaway gas, which results in a much lower density in theoutflow region. Third, the gas temperature in the outflowregion in model D19T6windrad is significantly higherthan that in the outflow region in model D19T6rad. Theare two reasons. The first one is that some portion of theenergy of the wind from the central AGN is transferred to the gas and becoming gas internal energy. The gasin the outflow region is heated by the compression workdone by wind. The second one is that lower gas density inoutflow region in model D19T6windrad results in lowerradiative cooling rate. Finally, the maximum velocity ofoutflow in model D19T6windrad is ∼ . − . c , whichis one order of magnitude higher than that of outflow inmodel D19T6rad.In Figure 3, we plot the radial profiles of time-averaged (from 3.5 to 25 free-fall timescales measuredat outer boundary) mass inflow (solid lines) and outflow(dashed lines) rates for models D19T6rad (black lines)and D19T6windrad (red lines). In both of the models,due to the presence of outflow, the mass inflow rate de-creases inwards. The mass outflow rate in both mod-els increases with increasing of radius. This means thatwhen outflow moves outwards, more and more inflow gasjoins outflows. The mass inflow rate in model D19T6radis larger. The reasons are as follows. First, in modelD19T6rad, the inflow occupies larger solid angle at anyradius (see Figure 2). Second, the density of the inflowgas in model D19T6rad is higher (see Figure 2). Finally,the infall velocity of inflowing gas in the two models iscomparable. In the region r > r s , the mass outflowrate in model D19T6windrad is lower. The reason is asfollows. Even though the velocity of outflow is higher inmodel D19T6windrad, the significant lower outflow den-sity (see Figure 2) makes the mass outflow rate in thismodel much lower.From Figure 3, we can see that at the inner boundary,the mass inflow rate in model D19T6rad does not differmuch from that in model D19T6windrad. The accretionrate onto the black hole is half of the mass inflow rate atthe inner boundary (see Section 2.2.1). The mass accre-tion rates onto the black hole in these two models do notdiffer much. Consequently, the luminosity of the centralAGN does not change much by wind feedback (see Fig-ure 1). Liu et al. (2013) studied the accretion flow withexactly same initial and boundary conditions as thosein models D19T6rad and D19T6windrad. In Liu et al.(2013), wind feedback has not been taken into account.They find that the radiative feedback (including radia-tion pressure due to Thomson scattering and line force)can not make the luminosity of the central AGN sub-Eddington. The results in the present paper demonstratethat the presence of wind and radiative (without lineforce) feedbacks can also not solve the ‘sub-Eddington’puzzle. Therefore, in order to make the AGN’s luminos-ity sub-Eddington, we may simultaneously consider windfeedback, radiation pressure due to Thomson scatteringand line force. In addition, in this paper, we do notconsider dust. In the region at parsec scale, dust mayexist (Wang 2008; Wang et al. 2017; Czerny et al. 2016,2019). Dust can absorb the radiation from the centralAGN and then emits infrared photons. The emitted in-frared photons can exert a strong radiation pressure ondust (Dorodnitsyn & Kallman 2012; Dorodnitsyn et al.2016). This may also help produce strong wind and de-crease the luminosity of the central AGN. In future, it isvery necessary to perform simulations by taking into ac-count wind feedback, radiation pressure due to Tomsonscattering and line force, and radiation pressure on dustsimultaneously to study this problem again.One of the most important properties of outflows in Bu & Yang Fig. 5.—
Same as Figure 1, but for models D20T6rad (top-panel)and D20T6windrad (bottom-panel).
AGN feedback study is the kinetic power. The kineticpower of outflow is as follows, P K ( r ) = 2 πr Z ◦ ◦ ρ max( v r ,
0) sin θdθ (23)We list the kinetic power (in unit of Eddington Lumi-nosity) of outflows calculated at outer boundary for allof our models in column 7 of Table 1. The kinetic power P K ∝ ˙ M outflow v , with ˙ M outflow and v outflow beingmass flux and velocity of outflow, respectively. FromFigure 3, we see at the outer boundary, the outflowmass flux in model D19T6rad is 2 ∼ Dependence on gas temperature at the outerboundary
In models D19T7rad and D19T7windrad, if gas flowsinwards at the outer boundary, then we set the gastemperature at outer boundary to be one order ofmagnitude higher than that in models D19T6rad andD19T6windrad. In Figure 4, we plot the time evolutionof luminosity of the central AGN. As in the models intro-duced above, due to the AGN feedback, the luminosityin models D19T7rad and D19T7windrad fluctuates sig-nificantly. The time-averaged luminosity of the centralAGN is just slightly decreased by AGN wind feedback.As in models with other values of density and tempera- ture, the mass flux of outflow at the outer boundary doesnot be affected much by wind feedback (see the sixth col-umn of Table 1). The kinetic power of outflow measuredat the outer boundary is significantly enhanced by windfeedback (see the seventh column of Table 1). We alsostudied the properties of outflow (density, temperatureand velocity). The effects of wind feedback are also verysimilar as those in models with ρ = 10 − g cm − and T = 2 × K. For example, compared to model withoutwind feedback (D19T7rad), in the model with wind feed-back (D19T7windrad), density in the outflow region issignificantly lower, outflow temperature is much higher,velocity of outflows is significantly higher.3.3.
Dependence on gas density at the outer boundary
In models D20T6rad and D20T6windrad, if gas flowsinwards at the outer boundary, then we set the gas den-sity at outer boundary to be one order of magnitude lowerthan that in models D19T6rad and D19T6windrad. InFigure 5, we plot the time evolution of luminosity of thecentral AGN. As in models with other different valuesof density and temperature, the wind feedback has verysmall effects on the time-averaged luminosity of the cen-tral AGN. Also, the mass fluxes of outflows at the outerboundary do not differ much in models with and with-out wind feedback (see the sixth column of Table 1). Thekinetic power of outflow measured at the outer bound-ary in model with wind feedback is roughly 2 orders ofmagnitude higher than that of outflow in model with-out wind feedback (see the seventh column of Table 1).The effects of wind feedback on properties of outflow arealso similar as those in the models introduced above. Inthe model with wind feedback (D20T6windrad), temper-ature in the outflow region is much higher, density in theoutflow region is much lower, the velocity of outflow issignificantly higher. SUMMARYTwo-dimensional hydrodynamic simulations are per-formed to study the effects of AGN radiation and windfeedback on the properties of slowly rotating accretionflow at sub-parsec to parsec scale. We find that whenonly radiation feedback is considered, outflows can belaunched by the radiation pressure force due to Thom-son scattering. The maximum velocity of the outflow is ∼ . c . Due to the presence of outflow, the mass inflowrate decrease inwards. The mass flux and kinetic powerof outflows depend on the density and temperature ofgas at parsec scale.We find that when wind feedback is taken into account.the mass flux of outflow can not be significantly affected.Accordingly, the luminosity of the central AGN does notchange much after taking into account wind feedback.The ‘sub-Eddington’ puzzle can not be solved by consid-ering wind feedback. However, wind feedback can affectother properties of outflow significantly. For example,maximum velocity of outflow can be increased by a fac-tor of ∼
100 by the wind feedback. Consequently, thekinetic power of outflow can be increased by a factor of40 −100. Other properties of outflow such as density andtemperature are also significantly affected by wind feed-back. In the model with wind feedback, the density ofoutflow is significantly lower, the temperature of outflowis significantly higher.GN feedback 9ACKNOWLEDGMENTSD.-F. Bu is supported in part the Natural Sci-ence Foundation of China (grants 11773053, 11573051,11633006 and 11661161012), and the Key Research Pro-gram of Frontier Sciences of CAS (No. QYZDJSSW-SYS008). X.-H. Yang is supported in part by the Funda- mental Research Funds for the Central Universities (No.106112015CDJXY300005 and No. 2019CDJDWL0005)and the Natural Science Foundation of China (grant11847301). This work made use of the High PerformanceComputing Resource in the Core Facility for AdvancedResearch Computing at Shanghai Astronomical Obser-vatory.