Local stability and global instability in iron-opaque disks
aa r X i v : . [ a s t r o - ph . H E ] J un Astronomy & Astrophysicsmanuscript no. kotletymielone10 c (cid:13)
ESO 2018October 18, 2018
Local stability and global instability in iron-opaque disks
Mikołaj Grze¸dzielski ([email protected]) , Agnieszka Janiuk , and Bo˙zena Czerny Center for Theoretical Physics, Polish Academy of Sciences, Al. Lotnikow 32 /
46, 02-668 Warsaw, PolandReceived ...; accepted ...
ABSTRACT
The thermal stability of accretion disk and the possibility to see a limit-cycle behaviour strongly depends on the ability of the diskplasma to cool down. Various processes connected with radiation-matter interaction appearing in hot accretion disk plasma contributeto opacity. For the case of geometrically thin and optically thick accretion disk, we can estimate the influence of several di ff erentcomponents of function κ , given by the Roseland mean. In the case of high temperatures ( ∼ ) K, the electron Thomson scatteringis dominant. At lower temperatures atomic processes become important. The slope d log κ/ d log T can have locally stabilizing ordestabilizing e ff ect on the disk. Although the local MHD simulation postulate the stabilizing influence of the atomic processes, onlythe global time-dependent model can reveal the global disk stability range estimation. This is due to global di ff usive nature of thatprocesses. In this paper, using previously tested GLADIS code with modified prescription of the viscous dissipation, we examine thestabilizing e ff ect of the Iron Opacity Bump .
1. Introduction
The energy output observed in the Galactic X-ray Binaries(XRB), and Active Galactic Nuclei (AGN), suggests that thesource of emitted power in these sources must be connectedwith the gravitational potential energy of a compact object.In most of the former, and all of the latter, this object is ablack hole. Its mass can range from a few solar masses upto a few billion of M ⊙ . The material that is accreted onto ablack hole and emits radiation, may posses substantial angularmomentum. In this case, the accretion flow forms a geometri-cally thin disk, which is located in the equatorial plane of anXRB, or co-aligned with the plane perpendicular to the blackhole rotation axis in AGN, if the latter is powered by a Kerrblack hole. The classical solution of Shakura & Sunyaev (1973)with their α − viscosity prescription, describes a stationary ac-cretion disk where the dissipated heat is balanced by thermalradiation. As studied by Pringle, Rees & Pacholczyk (1973);Lightman & Eardley (1974) and Shakura & Sunyaev (1976),the viscous stress tensor scaling with a total (i.e., gas plus ra-diation) pressure leads to the runaway instability of the diskstructure. Alternatively, the viscosity prescription may be givenby a gas pressure only (in this case the thermal instability doesnot develop), or by some intermediate law, that invokes a com-bination of gas and radiation with di ff erent weights (see e.g.Szuszkiewicz (1990)). For instance, a general prescription thatwas recently discussed by Grze¸dzielski et al. (2017) reads: τ r φ = α P µ tot P − µ gas . (1)The above model leads to the unstable disk behaviour, whichmanifests in a limit-cycle type of oscillations of the emitted lu-minosity, characteristic for the periodically heated and cooled in-ner regions of the accretion disk. This kind of cycle is possible,if only the thermal runaway is captured by some stabilizing pro-cess. This might be advection of heat onto a black hole, as pro-posed for the so-called ’slim disk’ solution (Abramowicz et al.1988). The presence of radiation pressure instability in actionof cosmic sources has been a matter of debate (see e.g., re-view by Blaes (2014)). Nevertheless, there are strong observa-tional hints which support the limit-cycle type of behaviour in at least two well-known microquasars, GRS 1915 +
105 and IGRJ17091-324, in some of their spectral states (Belloni et al. 2000;Altamirano et al. 2011). The limit-cycle oscillations were de-tected also in the Ultraluminous X-ray source HLX-1, claimedto contain an intermediate-mass black hole (Farrell et al. 2009;Lasota et al. 2011; Servillat et al. 2011; Godet et al. 2012;Wu et. al 2016). Furthermore, some type of non-linear dynam-ics characteristic for an underlying unstable accretion flow wassuggested for a number of other XRBs (Sukova et al. 2015),while the statistical studies of a large sample of sources supportthe ’reactivation’ scenario in the case of compact radio sourceshosting supermassive black holes (Czerny et al. 2009). On theother hand, many of XRBs and AGN seem to be powered bya stable accretion, despite even large accretion rates. Therefore,some stabilizing mechanisms in the accretion process have beeninvoked, apart from the viscosity prescription itself. For instance,the propagating fluctuations in the flow (Janiuk & Misra 2012)may suppress the thermal instability, or at least delay the devel-opment of the instability (Ross, Latter & Tehranchi 2017). Asrecently discussed by Jiang, Davis & Stone (2016), a possiblestabilizing mechanism for the part of an accretion flow might bethe opacity changes connected with the ionization of heavy ele-ments. Using the shearing-box simulations of MRI-driven fluidin the gravitational potential of supermassive black hole (a par-ticular value of black hole mass, M = × M ⊙ , was used),Jiang, Davis & Stone (2016) have shown that the flow is stableagainst the thermal instability, if the opacity includes transitionsconnected with absorption and scattering on iron ions. This is be-cause the cooling rate, which includes now not only the Thom-son scattering (constant) term, but also the absorption and lineemission in the Roseland mean opacity, will depend strongly ondensity and temperature in some specific regions of the disk. Infact, as we discuss below in more detail, the dominant term fromopacity changes may completely stabilize the flow locally. Thesimulations of Jiang, Davis & Stone (2016) showed that e ff ectbut they did not describe the global evolution of the flow, whichis the subject of our present work. In our paper we consider thedisk stability for the range of black hole masses, characteristicfor either XRBs or AGN. We find that the influence of the opac- Article number, page 1 of 5 & Aproofs: manuscript no. kotletymielone10 ity changes on the global time evolution of the flow is essen-tial, although do not prevent the instability form developing. Weshow examples of lightcurves produced by our numerical simu-lation to illustrate this. The paper is organized as follows: in Sect.2 we present the analytical condition for local thermal instabil-ity, in Sect. 3 we present the values of κ opacities for the typicalaccretion disk densities and temperatures and range of the IronOpacity Bump, to reveal its influence on global disk behaviourin Sect. 4.
2. Local thermal stability in accretion disks
The domination of radiation pressure in the accretion diskleads to the thermal instability (Pringle, Rees & Pacholczyk1973; Lightman & Eardley 1974; Shakura & Sunyaev 1976;Janiuk, Czerny & Siemiginowska 2002). At the instability theheating rate Q + grows faster with temperature than cooling rate Q − . The appearance of local thermal instability is given by thecondition: d log Q − d log T < d log Q + d log T . (2)The analysis performed in Grze¸dzielski et al. (2017), underthe assumption on constant surface density during thermaltimescales leads to following formula on heating rate derivative: d log Q + d log T = + µ − β + β , (3)where β = P gas / P . For the case of this work, we assume α = . µ = .
56, being typical for the IMBH and AGN disk case.We have chosen these values since the dynamics of the outburstsof the disk matches the observed properties of the sources asshown in Grze¸dzielski et al. (2017). These values reproduce thecorrelation between the bolometric luminosity and outburst du-ration known for the observed sources, especially microquasarsand Intermediate Mass Black Holes. The radiative cooling ratedepends on disk surface density Σ and physical constants Stefan-Boltzmann σ b and speed of light. The radiative cooling rate isgiven by following formula (Janiuk, Czerny & Siemiginowska2002; Janiuk et al. 2015; Grze¸dzielski et al. 2017): Q − = σ b c Σ T κ . (4)For the radiation pressure dominated disk, combined withEqs.(2) and (3), from Eq. (4) for µ = .
56 we get: d log Q − d log T < . . (5)In our previous papers, only Thomson scattering was assumed,which resulted in the appearance of global radiation pressure in-stability among all scales of sub-Eddington accretion disks. Therecent results of Jiang, Davis & Stone (2016) are based on stabi-lizing influence of the iron opacity components. To confront theresults of their short time, local 3D MHD shearing-box simula-tion, we propose the global model, used previously for the caseof Intermediate Mass BH HLX-1 accretion disk (Wu et. al 2016;Grze¸dzielski et al. 2017). According to the lower temperaturesin radiation-pressure dominated areas of the accretion disks, weinclude in our model also the atomic opacity components. As-suming κ being function of ρ and T , we get the formula for logderivative of Q − : d log Q − d log T = + ∂ log κ∂ log T − − β + β ∂ log κ∂ log ρ . (6) The Eq. 6 shows possible stabilizing influence of the negativeslope of κ dependent on T and destabilizing influence of the pos-itive slope of κ . Also the dependence on ρ has inluence on diskstability. The value of κ itself is not important in the local stabil-ity analysis. In case of unstable disk, less e ffi cient cooling beingan e ff ect of the greater κ lowers the temperature of unstable equi-librium solution, and enlarge the temperature of stable solution,which can modify the duty cycle quantitatively but not qualita-tively. Similarly to Grze¸dzielski et al. (2017), we can derive thenecessary value for the thermal instability: β < µ + − ∂ log κ∂ log T + ∂ log κ∂ log ρ µ − + ∂ log κ∂ log T − ∂ log κ∂ log ρ (7)We can also define the thermal stability parameter s : s = d log Q + d log T − d log Q − d log T (8)Using the Eqs. (3) and (6), we can write Eq. (8) as follows: s = − + µ − β + β − ∂ log κ∂ log T + − β + β ∂ log κ∂ log ρ (9)The s parameter is connected with the Lyapunov exponent forthe system described by the energy equation in the accretion diskwith stress tensor given by Eq. (1). The value s > s ≤ µ = .
56 and β << ρ = − g cm − and T = K weget s = −
11 which corresponds to local thermal stability. Never-theless, for this density the parameter s gains positive values for T > × K. Below this temperature, for T > . × K,the disk is locally thermally stable because of the negative slopeof the bump. For temperatures in the range 1 . − . × K,the disk is locally thermally unstable.
3. The variable κ - Iron Opacity Bump The opacity κ is the local function describing interaction of pho-tons with matter from accretion disks. Under the assumptions oflocal thermal equilibrium, radiation and gas contribution to thetotal pressure, and local vertical hydrostatic equilibrium, boththe heating and cooling rates can be described as a function ofradius r , local density ρ and local temperature T . Although theradius r , a ff ecting the angular momentum transport is importantfor the heating rate in the α -disk model, it a ff ects the stabilityanalysis only indirectly, via the parameters of stationary solu-tions. In Fig. 1 we present the profiles of the total opacity forsolar metallicity, computed as a function of density and tempera-ture (Alexander et al. 1983; Seaton et al. 1994; Rozanska et al.1999). The combined conditions (3) and (4) to the opacity valuesresults in the significant local stabilization for the temperaturesof 1 − × K and densities about 10 − g cm − typical for theAGN accretion disks. We fitted the κ function with followingformulae: κ = κ Th + κ pl + κ bump , κ Th = . g − ,κ pl = . × ρ T − . , (10) Article number, page 2 of 5rzedzielski et al.: Local stability -2-1 0 1 2 3 4 4.5 5 5.5 6 6.5 7 7.5 8
Log κ [ c m g - ] Log T[K] ρ = 1 g cm -3 ρ = 0.1 g cm -3 ρ = 10 -2 g cm -3 ρ = 10 -3 g cm -3 ρ = 10 -4 g cm -3 ρ = 10 -5 g cm -3 ρ = 10 -6 g cm -3 ρ = 10 -7 g cm -3 ρ = 10 -8 g cm -3 ρ = 10 -9 g cm -3 ρ = 10 -10 g cm -3 Fig. 1.
The opacity κ functions including atomic components. Thedata are taken from (Alexander et al. 1983; Seaton et al. 1994;Rozanska et al. 1999). -0.4-0.2 0 0.2 0.4 0.6 0.8 5 5.5 6 6.5 7 7.5 Log κ [ c m g - ] ( t ab l e s ) Log T[K] ρ = 0.1 g cm -3 ρ = 10 -2 g cm -3 ρ = 10 -3 g cm -3 ρ = 10 -4 g cm -3 ρ = 10 -5 g cm -3 -0.4-0.2 0 0.2 0.4 0.6 0.8 5 5.5 6 6.5 7 7.5 Log κ [ c m g - ] ( ana l y t i c a l app r o x i m a t i on ) Log T[K] ρ = 10 -6 g cm -3 ρ = 10 -7 g cm -3 ρ = 10 -8 g cm -3 ρ = 10 -9 g cm -3 ρ = 10 -10 g cm -3 Fig. 2.
Top panel shows results from the tables (Alexander et al. 1983;Seaton et al. 1994; Rozanska et al. 1999), and bottom panel shows theanalytical approximation (Eq. 10). κ bump = . ρ . (0 . − ( T − . × K8 . × K ) + . − ( T − × K3 × K ) ) . The negative stabilizing slope of the Iron Opacity Bump isvisible in the Figure 1. In Fig. 2 in the upper panel we showthe opacities directly from the tables (Alexander et al. 1983; -10-8-6-4-2 0 2 4 5 6 7 8 9
Log ρ [ g c m - ] Log T [K] power lawTHIN DISKStoo low ρ IRON BUMPM = 1 MsM = 10 MsM = 10 MsM = 10 MsM = 10 MsM = 10 Ms-10-9-8-7-6-5 4.5 5 5.5 6 6.5
Log ρ [ g c m - ] Log T [K] M = 10 MsM = 10 MsM = 10 MsM = 10 MsM = 10 Ms Fig. 3.
Typical values of the ρ and T for Eddingtonian thin disks forcentral object mass between 1 and 10 M ⊙ . We set ˙ m = . Seaton et al. 1994; Rozanska et al. 1999), and in the lowerpanel we present the analytical approximation of opacity func-tion from Eq. ( ?? ). The detailed results of the dynamical modelare presented in Section 4.
4. Global model ρ and T In Fig. 3 we present typical values of ρ and T for the wide rangeof accretion disks, computed via the GLADIS code. For the val-ues of ρ and T in the upper left corner, the power law term of κ dominates, but matter with this parameter is too dense andtoo cold for central areas of sub-Eddington accretion disks. Theoblique belt below presents typical values of ρ and T for the ac-cretion disks with di ff erent masses. This area, for 5 . < log T < .
4, is covered by the Iron Opacity Bump, with stabilizing nega-tive slope (See Sect. 2). According to these results, the stabiliz-ing e ff ect of the negative slope of the Iron Opacity Bump can bevisible only for the Active Galactic Nuclei accretion disks with M ≈ − M ⊙ . We perform the simulations of the global disk behaviour usingthe time-dependent global code GLADIS (Grze¸dzielski et al.2017). The GLADIS code is a time-dependent code, whichsolves hydrodynamic equations, describing the long-time be-haviour of accretion disk under the assumption of the verticalhydrostatic equilibrium. The code models the time evolution ofthe flow in thermal and viscous timescales. Moreover, the codeassumes axial symmetry, and Keplerian angular velocity. In thispaper we set α = . µ = .
56, and M = × M ⊙ .The major change in comparison to Grze¸dzielski et al. (2017)is replacing constant Thomson κ with Eq. ( ?? ). Similarly toJiang, Davis & Stone (2016), we assumed the Eddington rate˙ m = .
03. The results of the time-dependent model are presented
Article number, page 3 of 5 & Aproofs: manuscript no. kotletymielone10 -1.2-1-0.8-0.6-0.4-0.2 0 0.2 0.4 0.6 0.8 1 20 40 60 80 100 120 S t ab ili t y pa r a m e t e r s r [R Schw ] instabilitys
Fig. 4.
The stability parameter, defined in Eq.(8) during the outburst for M = × M ⊙ , ˙ m = .
03 (same as in Jiang, Davis & Stone (2016)). in Fig. 4. The local shearing-box simulation resulted in the sig-nificant stabilization of the disk (Jiang, Davis & Stone 2016).However, the global model does not confirm these results. Thestabilization of the disk, which appears according to local pre-diction, is not found in global models considering a large rangeof radii. Fig. 4 presents the stability parameter profile s definedin Eq. 8. For the inner area of the disk, Thomson componentof opacity dominates and temperature is too large to expect anyform of stabilization. Outer area of the disk characterize thelarger value of total opacity temperature about 1 . − × K and negative values of s (the bump temperatures). The signifi-cant gradient of the s is correlated with the gradient of temper-ature and gradient of κ in the opposite direction. The stabilityparameter, presented in Fig. 4 can reach values between − − + µ (0 .
92 for our choice of µ ). In Fig. 4 the typical profileof the s parameter is presented. As the bump is approximatelyGaussian function, centered at 1 . × K, with standard de-viation σ = . × K, it is expected that the strongest e ff ectof the stabilizing slope would be visible for such temperatures.The combined outcome of the stabilizing influence of the neg-ative slope of the bump and destabilizing influence of the posi-tive slope of the bump for the dynamical model is presented inFig. 5. In contrast to the results of Jiang, Davis & Stone (2016),the iron bump does not stabilize the disk. However, it compli-cates the lightcurve pattern (many small short flares precedingmain outbursts instead of one simple flare), due to the complex-ity of the photon absorption process, but the inner regions ofthe disk remains hot enough to perform the limit-cycle oscilla-tions. In e ff ect, the bump partially stabilizes the disk - the am-plitude L max / L min decreases from 156 for model with constant κ from lower panel of 5 to 16.1 for model with bump (upperpanel of 5). The detailed analysis of the lightcurve shape is pre-sented in Table 1. Since the timescales presented in Fig. 5 aremuch longer than the duration of observation (lasting up to sev-eral decades), it is impossible to find such a lightcurve using thedirect method. For such a black hole mass, during the phase of Log ( L / e r g s - ) ( m ode l w i t h bu m p ) Time [yrs] with bump 44.5 45 45.5 46 46.5 47 0 50000 100000 150000 200000
Log ( L / e r g s - ) ( on l y T ho m s on opa c i t y ) Time [yrs] only Thomson opacity
Fig. 5.
Results of the time-dependent model for M = × M ⊙ , ˙ m = .
03 (same as in Jiang, Davis & Stone (2016)). model amplitude A period P [yrs] width ∆ κ with bump 16.1 12635 0.0023Thomson κ
156 70197 0.0044
Table 1.
The table describing the flare parameters for the model fromFig. 5 with κ described by formula ( ?? ) ( second row) and only Thom-son κ (third row). Those parameters were described more precisely inGrze¸dzielski et al. (2017). fastest growth of the luminosity, the luminosity change can reach1 per cent per year. The shape of the lighturve can be reflectedin the Eddington rate statistic - similar objects, being in the dif-ferent phase of the limit-cycle presented in Fig. 4 can emit theradiation with di ff erent luminosity and spectra. However, in caseof much smaller black hole masses new timescales are perhapsaccessible to observations. For example, digitalization of Har-vard plates (Grindlay et al. 2012) will bring lightcurves on theorder of a hundred years and perhaps the outbursts of AGN diskscan be discovered.
5. Conclusions
In our previous paper (Grze¸dzielski et al. 2017) we computedlarge grids of models confirming the universality of radiationpressure instability across the BH mass-scale. In this paper wechanged the opacity prescription to examine the heavy atoms in-fluence on the accretion disk instability. Comparison betweentwo models presented in Table 1 leads to the conclusions thatheavy atoms stabilize the disk partially, but do not imply thatthe variability vanishes. This stabilizing e ff ect manifests itselfrather in a significant period and amplitude decrease, without arelative broadening of the outbursts with respect to their separa-tion . Additionally, some mild precursors , being an outcome ofa non-monotonic profile of the s parameter distribution, are alsovisible. That partial stabilization, being an important e ff ect for http: // dasch.rc.fas.harvard.edu / Article number, page 4 of 5rzedzielski et al.: Local stability the Active Galactic Nuclei, has only weak influence on the radi-ation pressure among all BH mass-scale in accretion disk underthe assumption of solar metallicity. In case of sources with dif-ferent metallicity, this e ff ect can change its extent. Finally, weconclude that the radiation pressure driven limit cycle oscilla-tions, su ff ering some disturbances from the Iron Opacity Bump in case of the AGN disks are also expected, at least for moder-ately large supermassive black holes ( M = × M ⊙ ). Acknowledgments
This work was supported in part by the grants DEC-2012 / / E / ST9 / / / M / ST9 / References