An axisymmetric, hydrodynamical model for the torus wind in AGN

An axisymmetric,hydrodynamical model for the torus wind inAGN.A.Dorodnitsyn 1,2,T.Kallman 1,and D.Proga 3ABSTRACTWe report on time-dependent axisymmetric simulations of an X-ray excited flow from a parsec-scale,rotating,cold torus around an active galactic nucleus.Our simulations account for radiative heating and cooling and radiation pressure force.The simulations follow the development of a broad bi-conical outflow induced mainly by X-ray heating.We compute synthetic spectra predicted by our simulations.The wind characteristics and the spectra support the hypothesis that a rotationally supported torus can serve as the source of a wind which is responsible for the warm absorber gas observed in the X-ray spectra of many Seyfert galaxies.Subject headings:acceleration of particles –galaxies:active –hydrodynamics –methods:numerical –quasars:absorption lines –X-rays:galaxies 1.Introduction Approximately 50%of low-red-shift AGNs exhibit spectra rich in atomic absorption features in the X-ray band.These spectra,referred to as warm absorbers,contain numerous lines and bound-free edges from ions of intermediate-Z elements,and indicate blue-shifts of the absorbing material of 102−103km s −1.The term refers to the fact that the X-ray absorbing gas has an electron temperature and ionization degree which are interme-diate between that expected for neutral and fully ionized material.Warm absorber gas has been known for more than two decades,beginning with observations using the Einstein Ob-servatory satellite (Halpern 1984),although the limited spectral resolution of early X-ray experiments hampered a detailed study of its properties.1Laboratory for High Energy Astrophysics,NASA Goddard Space Flight Center,Code 662,Greenbelt,MD,20771,USA2Space Research Institute,Profsoyuznaya st.,84/32,117997,Moscow,Russia 3Department of Physics and Astronomy,University of Nevada,Las Vegas,NV 89154,USAa r X i v :0801.1470v 1 [a s t r o -p h ] 9 J a n 2008Recent X-ray telescopes,Chandra and XMM-Newton,provide unprecedented spectral resolution up to∼10keV,and have revealed much of the complexity of warm absorber spectra.The ubiquity of X-ray warm absorbers and the implication of outflow in objects which are likely powered by accretion,combine to make warm absorbers of great interest for the understanding of active galactic nuclei.While observationally the reality of the X-ray absorbing gas isfirmly established,the source of this gas remains uncertain.If the wind originates from an accretion disk(Konigl&Kartje1994;Murray et al.1995),then warm absorber gas can occupy a vast spa-tial region from the inner accretion disk up to and beyond the broad line region.In addition, it is likely to have a speed which is comparable to the escape velocity in the region where the wind originates.If so,the line widths favor a wind originating far outside the inner accretion disk,at∼1pc,for an AGN containing a106M black hole.The natural source of gas at this distance is the molecular torus which is responsible for obscuring the broad line region in Seyfert2galaxies,and which is thought to exist in some form in many low and intermediate luminosity AGN(Antonucci&Miller1986). The existence of an outflow from the torus has been suggested by(Krolik&Begelman1986; Krolik&Begelman1988),and as the source of warm absorberflows by(Krolik&Kriss1995; Krolik&Kriss2001).Direct interferometric observations of this dusty structure have recently become avail-able:mid-infrared high spatial resolution studies of the active nucleus of NGC1068,revealed a dust torus-like structure that is2.1pc thick and3.4pc in diameter(Jaffe et al.2004). Evidence of at least two temperature components has been established:thefirst one is ex-tended and warm(320K)and conceals the second,which is compact and hot(>800K).It is natural to assume that the latter component can be attributed to the inner part of the torus,heated by the radiation from the compact nucleus.Although such evidence is most solid for nearby active galaxies it is likely that the same obscuring torus paradigm may ap-ply in quasars whose central regions are heavily obscured by gas and dust(Type II quasars) (Zakamska et al.2006).The work done so far in modeling torus warm absorberflows employ either one-dimensional models e.g.(Chelouche&Netzer2005)or two dimensional models which omit rotational forces and also employ an artificial shape for the torus,which is adopted as a boundary for theflow(Balsara&Krolik1992).These calculations support the idea that the outflow originates at r∼1pc-the estimate supported both from virial arguments and also deduced from the warm absorber line widths.The content of the warm absorber gas is highly ionized cosmically abundant heavy elements.This gas can occupy a region from0.01pc to∼10pc from the nucleus.The column densities are in the range1021−1024cm−2.These estimatesare consistent with constraints on the position of the warm absorber gas,which are deduced from the absence of correlated variability of the line equivalent widths to changes of the continuum (Netzter et al.2003;Behar et al.2003).The goal of the current work is to show that a state-of-the-art hydrodynamic model of the evaporative torus wind possesses key characteristics of the the warm absorber.To approach this goal,in this paper we present results of models containing the following in-gredients:1)Time-dependent axisymmentric numerical simulations that take into account rotational forces ,X-ray heating and cooling and radiation pressure both in spectral lines and continuum;2)Exploration of the resulting distributions of density,velocity etc.;and3)A calculation of the synthetic X-ray absorption spectra using the density and velocity distributions as input.Although calculating the spectra is crucial to check the relevance of the particular hydrodynamic solution to the warm absorber phenomenon,the primary goal of the work reported here is the outflow dynamics .2.Basic equationsThe equations to be solved include:∂ρ∂t +∇·(ρv )=0,(1)ρ ∂v ∂t+(v ·∇)v =−∇p −ρ∇Φ+ρg rad ,(2)∂ ∂t +∇·(v ( +p +ρΦ))=H ,(3)these are the conservation equations for:mass,momentum and energy.Heating and cooling processes are described by the function H (erg cm −3s −1); -is the sum of the kinetic and internal energy densities: =ρv 2/2+e .Three components of the velocity are taken into account.However,the adopted azimuthal symmetry implies that all ∂/∂φ=0.The gas is assumed to have a constant adiabatic index,so that e =P/(γ−1)with γ=5/3.Equations (1)-(3)are cast in a non -dimensional form with the characteristic scales set by the properties of the plasma orbiting at a distance of R 0=1pc from a black hole with massM 6in units of 106M :time in measured units of t 0(s)=4.7·1011r 3/2pc M −1/26,velocity inunits of V 0(cm s −1)=6.6·106M 1/26r −1/2pc .We anticipate that the conditions at the interface of a torus and an outflow driven by X-ray heating will resemble the X-ray self-excited wind in X-ray binaries.That is,an X-rayheated skin is created at the outer layers of the companion star atmosphere as a result of X-ray heating by the compact object.The temperature there rises almost discontinuously to some value that allows the gas to escape the gravity of the star.Theflow below the jump is subsonic and is determined by the details of the heating,cooling,and transfer processes in this region(Basko et al.1977).In such a situation,the distribution of temperature depends on the heating and cooling,including both adiabatic effects and radiation,in the equations of hydrodynamics.This can only be calculated accurately if the interior of the star is included as part of the computational domain,and the same applies to the problem of interest to us: the torus X-ray excited wind.We,therefore include in our computational domain the entire volume containing the torus and we start our simulations from a rotating torus configuration in equilibrium in the external gravitationalfield of the black hole.A stationary,polytropic,toroidal structure with constant angular momentum was sug-gested by(Papaloizou&Pringle1984)(PP-torus),and was shown to be stable in2.5di-mensions against radial or axial perturbations,an found to be unstable to non-axisymmetric perturbations(this effect cannot be numerically investigated in axisymmetry).In our dimen-sionless units the locus of the maximum torus densityρmax occurs at R0,so thatρmax=1. The distortion of the torus is controlled by the parameter d=( ++ −)/(2R0),where ±refers to the inner or outer edge of the torus.No replenishing of the torus gas is assumed, i.e.gas can only leave computational domain.The drawback of this approach is that it is not possible to obtain a true stationary solution of the problem.Nonetheless,as we will show,the torusflow can be regarded as a quasi-equlibrium on timescales shorter than the depletion timescale,and can thus provide representative snapshots of the character of the flow.2.1.X-ray exited wind:evolution of the wind and the torusThe thermodynamic properties of X-ray heated gas depend on the spectrum of the incident radiation as well as on the local atomic physics.Under the assumption of photo-ionization equilibrium the thermodynamic state of the gas can be parameterized in terms of ratio of radiation energy density to baryon density(Tarter,Tucker,&Salpeter1969):ξ=4πF x/n,where F x=L x e−τ/(4πr2)is the local X-rayflux,L x is the X-ray luminosity ofthe nuclei,τ=r0κeρdr-is the optical depth.We assume that the attenuation is dominatedby Thomson scatteringκ=0.2(1+X H) 0.4cm2g−1,where X H is the mass fraction of hy-drogen,and the factor e−τ-describes approximately the attenuation of the radiation on the way from the source toward afiducial point.Assuming that the there is a fractionηx of the total luminosity L BH available in X-rays and that the disk radiates a fractionηedd of its Ed-dington luminosity L edd =1.25·1044M 6we estimate:ξ 4·102·ηx ηedd M 6/(N 23r pc ),where r pc -is the distance in parcecs and N 23is the column density in 1023cm −2.If the wind is driven by thermal evaporative mechanism,then its pressure and temperature are determined by the condition of photo-ionization equilibrium,if the dynamical time within the flow is much larger than the characteristic time of the heating-cooling.In a stationary flow the rate of the Comp-ton and photo-ionization heating should be balanced by the Compton,radiative recombina-tion,bremsstrahlung and line cooling respectively.In the calculations we adopt approximate formulas of (Blondin 1994)for these processes.The radiative acceleration,g rad ,accounts for both the acceleration due to lines g rad =(F UV κ/c )k t −α (1+τmax )(1−α)−1 /τ1−αmax ,where L UV is the UV flux,t =σe v th ρ/|dv/dr |,v th is the thermal velocity,σe is the electron-scattering opacity,and α=0.5(Castor,Abbott,&Klein 1975;Stevens &Kallman 1990),and for radiation pressure on electrons g cont =F UV κ/c .In our approach,the radiation force which is due to lines and also k and τmax depend on the ionization parameter ξ,and we adopt the description made by (Stevens &Kallman 1990).The ionization parameter ξis calculated self-consistently at each grid zone,at each time step taking into account the attenuation of the X-ray flux due to finite radial optical depth.We adopt spherical-polar (r,θ,φ)coordinates for calculation of a wind heated by a central source of radiation.We extend the computational domain r i from r in =0.01pc to r out =50pc and the polar domain θi from 0to πwith no assumption about the equatorial symmetry.Models presented here have been calculated with a N r ×N θ=100×100grid.We have also performed test calculations on the 300×300grid finding no noticeable difference with the result presented here.We postpone the discussion of this and related issues to a separate publication.Both angle and radial grids are non-uniformly spaced:the angle gradually increases towards θ=π/4with fewer angular points towards the pole (θ=0)and equator (θ=π/2),and analogously spaced in the southern hemisphere.The use of a logarithmically spaced radial grid allows for the ”throat”of the torus to be resolved.We set outflowing boundary conditions at both ends of the radial grid and axially symmetric at θ=0,π.This initial torus is imbedded in a low density gas with ρ=10−4ρmax ,where ρmax is the maximum density inside the initial torus.To solve the hydrodynamical part of the system of equations (1)-(3)we use the ZEUS-2D code described by (Stone &Norman 1992).As a test we have evolved a torus distribution of matter for t =2(two rotational periods in our units)and found the whole configuration stable.We include the heating term into the hydrodynamical equations using a fully implicit approach because the characteristic heating time scale can be much shorter than the dynamical time.We assume M 6=1,ηedd =0.4,ηx =0.5and also that 0.5L BH is radiated in UV and contributes to the flow energetics through the radiation force g rad .The initial torus config-uration chosen here is only mildly distorted,the distortion parameter:d =1.25.Following(Papaloizou&Pringle1984),we assume that the distribution of specific angular momentum inside the initial torus,is constant.To ensure that our initial torus is in equilibrium,we follow(Papaloizou&Pringle1984)with the exception that we reduce the effective gravity GM BH/r2(1−ηedd)due to the effect of the radiation pressure.However the real continuum radiation pressure,which is adopted in the calculations,is a factor of e−τr less because of the optical depth effect.Our initial torus is Thomson thickτr( +,θ=π/2) 4.35.Although the torus itself is not the primary topic of this study,the X-ray illumination does affect the interior and average properties of the torus.During the evolution the locus of the maximum torus density(located initially at =1)shifts towards larger .This happens primarily because of the over-pressured region created by the evaporating gas and also because of radiation pressure.This allows for the energy to be deposited directly to the torus interior. In our case,direct radiative heating produces gradients of gas pressure in the torus interior. Thus,by making the initial torus marginally optically thick,we crudely mimic the model in which the torus is sustained against vertical collapse by the gradients of infrared radiation energy density(Krolik2007).In what follows,we describe the torus outflow behavior at three different times after the initial turn-on of X-rays.After one orbital period,an over-pressured region extends to r 4pc throughout the area that is not shadowed by the high density torus.At this time the non-shadowed area lies atθ<56◦.In this space,an axisymmetric region exists between0.7pc< <2.5pc and z<7where the temperature is everywhere,of the order of the local virial temperature (except in the torus interior):T 5T vir(r).Closer to the black hole,the heating is much stronger,leading to T∼107K.High temperature,low density gasfills the space in the ”throat”of the torus.The outer edge of the torus has extends to∼5.3pc in temperature, and to∼4.5pc in density contours.The ionization parameter in the same region isξ 600.A wind is well-developed:the velocityfield in the vicinity of the high density torus,follows the equal pressure contours;the maximum radial velocity is observed close to the axes at v max(θ=3◦)=1180km s−1,the second local maximum is atθ=12◦:v max=850km s−1. The maximum velocity gradually decreases towards higher inclinations v max(θ=27◦)= 343km s−1;and atθ=60◦it is v max=110km s−1.Closer to black hole,at about0.2pc,the vθcomponent of the velocity is more pronounced,until heights z 1.4pc are reached forming a base of the evaporative wind.After this time,theflow is approximately symmetrical in both hemispheres.Because we are solving equations of ideal hydrodynamics(with only a small numerical viscosity),accretion through inner boundary(at r=0)is almost negligible:˙M in (M yr−1)<10−8.The mass-loss rate˙M maxΩ(M yr−1sterrad−1)peaks atθ 10◦,i.e.atmuch higher inclinations than v max.The total mass-loss rate at this stage is˙M(M yr−1) 3·10−4.Figure1shows theflow streamlines and density contours at three orbital periods.Atthis time the high pressure region expands vertically,to z 10pc from the inner torus.The torus inner edge can be inferred from the temperature and density maps: − 1.3pc.Inside a throat of this torus,T 10T vir.The contours of T become moreflattened:z max 4pc, + 8pc.The conservation of angular momentum helps to keep the outer edge of the torus at∼5pc.However,there is an envelope around this dense torus:z max 4pc, + 8.1pc. The values of the minimum ionization parameter and the column density vary significantly with the inclination angle:ξmin(θ=8◦)=2000,N22=0.3;atθ 25◦,the ionization parameter decreases significantly:ξmin(θ=26◦)=63,N22=2.3;and at higher inclination, the ionization parameter is quite low:ξmin(θ=53◦)=6.6,and the obscuration is quite high:N22=75.The mass loss rate changes slowly d˙M/dt((M yr−1)/yr) 1.53·10−6.The mass-loss rate peaks at˙M maxΩ=0.05atθ 20◦,where v max(θ=20◦)=616km s−1;the maximum value of v max(θ=5◦)=1000km s−1;The total mass-loss rate is˙M(M yr−1) 0.08.It is this evolved model that gives mass-loss rate in accord with the estimates of (Krolik&Begelman1986).Figure2shows theflow atfive orbital periods.At this time the density maximum is at =3.6.Remarkably,the inner edge of the inner,dense torus is at the same distance from the black hole: − 1.3pc,as it was at three orbital periods,while the outer edge is at + 6.2pc.The position of the outer edge of the torus,as infered from the temperature maps,extends to11pc.Outside this region theflow temperature is T 2.5T vir.The significant drop of ionization parameterξfrom∼103to∼10,occurs again atθ=25◦,where the column also rises to a∼1022.Aspect ratio of the torus:∆r/r∼1in accord with what is inferred from observations(Krolik&Begelman1986;Jaffe et al.2004).Almost everywhere in the wind the poloidal component of the velocity is determmined by v r.However,closer to the black hole,the vθcomponent is important,vθ∼v r.The mass-loss rate is approximatelyconstant at15◦<θ<60◦:˙M maxΩ=0.07,and has a peak atθ=74◦:˙M maxΩ=0.16.The totalmass-loss rate is˙M(M yr−1) 0.25and d˙M/dt((M yr−1)/yr) 10−6.The maximum value of v max(θ=3◦)=1032km s−1;another local peak is atθ=15◦,where the maximum velocity is v max=606km s−1,and the minimum v max(θ=90◦)=180km s−1.The torus is losing mass in all directions,although with very different speed at different inclinations.We have performed multidimensional simulations of a wind excited by X-ray heating of a rotationally supported torus.Our dynamical calculations predict the velocityfield and the spatial distribution of opacity and emissivity of the gas in both continuum and lines as a function of position.(Synthesizing the spectrum provides a key test for our model).We have calculated synthetic spectra by integrating over the wind and using the formal solution of the transfer equation.This is done by passing the results of the hydrodynamical calculations into a code which calls the subroutines from the xstar code(Kallman&Bautista2001) to calculate the transmittedflux.In this paper we present only absorption spectra.Thisapproximately represents the case when warm absorber gas lies on the line of sight to the observer and we see only transmitted spectra.These are shown in Figure3,for various anglesθ=26◦,35◦,44◦at three orbital periods.These show that our simulations provide gas with intermediate ionization and sufficient column to produce observable warm absorber features.This supports the idea that it is this wind that is responsible for the warm absorber gas observed in Seyfert galaxies.More detailed results,including dependence of these results on different parameters as well as limitations of the adopted model and comparisons with observations,will be discussed in a separate paper.This research was supported by an appointment to the NASA Postdoctoral Program at the NASA Goddard Space Flight Center,administered by Oak Ridge Associated Universities through a contract with NASA,and by grants from the NASA Astrophysics Theory Program 05-ATP05-18.REFERENCESAntonucci,R.R.J.,Miller,J.S.ApJ,297,621Balsara,D.S.,Krolik,J.1992,ApJ,402,109Blondin,J.M.1994,ApJ,435,756Basko,M.M.,Hatchett,S.,McCray,R.,Sunyaev,R.A.1977,ApJ,215,276Behar,E.,Rasmussen,A.P.,Blustin,A.J.,et al2003,ApJ,598,232Castor,J.I.,Abbott,D.C,Klein,R.1975,ApJ,195,157Chelouche,D.,Netzer,H.2005,ApJ,625,95Crenshaw,D.M.,et al.1999,ApJ,516,750Chartas,G.,Brandt,W.N.,Gallagher,S.C.,Garmire,G.P.2002,ApJ,279,169 Halpern,J.P.1984,ApJ,132Jaffe et al.2004,Nature,429,47Kallman,T.,Bautista,M.2001,ApJ,Supp.,133,221Kinkhabwala,A.et al.,2002,ApJ,575,732Konigl,A.,Kartje,J.1994,ApJ,434,466Krolik,J.H.2007,ApJ,661,52KKrolik,J.H.,Kriss,G.A.1995,ApJ,447,512Krolik,J.H.,Begelman,M.C.1986,ApJ,308,L55Krolik,J.H.,Begelman,M.C.1988,ApJ,329,702Krolik,J.H.,Kriss,G.A.2001,ApJ,561,684Murray,N.,Chiang,J.,Grossman,S.A.,Voit,G.M.1995,ApJ,451,498 Netzter,H.,Kaspi,S.,Behar,E.,et al.2003,ApJ,599,933 Papaloizou,J.C.B.,Pringle,J.E.1984,MNRAS,208,721Pounds,K.A.,Reeves,J.N.,King,A.R.,et al.2003b,MNRAS,346,1025 Reynolds,C.S.1997,MNRAS,286,513Stone,J.M.,Norman,M.L.1992,ApJS,80,753Stone,J.M.,Pringle,J.E.,Begelman,M.1999,MNRAS,310,1002 Stevens,I.,Kallman,T.1990,ApJ,365,321Tarter,C.B.,Tucker,W.,Salpeter,E.E.1969,ApJ,156,943 Zakamska,N.,et al.2006,ApJ,132This preprint was prepared with the AAS L A T E X macros v5.2.Fig. 1.—Velocity streamlines superimposed on contours of log(n/n0)after three orbital periods.Note,that streamlines carry information about the direction only i.e.the magnitude of the velocity at low inclination is much higher than near the equatorial plane.orbital periods;at northern hemisphere superimposed with velocity vectors.Fig. 3.—The transmitted spectra observed at three orbital periods from2.5D model of the evaporative torus wind.Each panel corresponds to a different viewing angle,θ=26◦,35◦,44◦;x-axis:energy;y-axis:normalized,transmittedflux.The initial spec-trum is a power law with energy indexα=1.。