Analytical Model of Disk Evaporation and State Transitions in Accreting Black Holes

State transitions in black hole X-ray binaries are likely caused by gas evaporation from a thin accretion disk into a hot corona. We present a height-integrated version of this process which is suitable for analytical and numerical studies. With radius $r$ scaled to Schwarzschild units and coronal mass accretion rate $\dot{m}_c$ to Eddington units, the results of the model are independent of black hole mass. State transitions should thus be similar in X-ray binaries and AGN. The corona solution consists of two power-law segments separated at a break radius $r_b \sim10^3 \,(\alpha/0.3)^{-2}$, where $\alpha$ is the viscosity parameter. Gas evaporates from the disk to the corona for $r>r_b$, and condenses back for $r<r_b$. At $r_b$, $\dot{m}_c$ reaches its maximum, $\dot{m}_{c,{\rm max}} \approx 0.02\, (\alpha/0.3)^3$. If at $r\gg r_b$ the thin disk accretes with $\dot{m}_d<\dot{m}_{c,{\rm max}} $, then the disk evaporates fully before reaching $r_b$, giving the hard state. Otherwise, the disk survives at all radii, giving the thermal state. While the basic model considers only bremsstrahlung cooling and viscous heating, we also discuss a more realistic model which includes Compton cooling and direct coronal heating by energy transport from the disk. Solutions are again independent of black hole mass, and $r_b$ remains unchanged. This model predicts strong coronal winds for $r>r_b$, and a $T\sim 5\times 10^8\,{\rm K}$ Compton-cooled corona for $r<r_b$. Two-temperature effects are ignored, but may be important at small radii.


INTRODUCTION
Black hole (BH) X-ray binaries have several distinct accretion states (see Tanaka & Shibazaki 1996;Done et al. 2007 for reviews). The most notable of these are the thermal (or high/soft) state and the hard (or low/hard) state. The thermal state is characterized by a dominant multitemperature blackbody spectrum in X-rays, consistent with a geometrically thin, optically thick, cool accretion disk, similar to the model described in Shakura & Sunyaev (1973) ;Novikov & Thorne (1973); Frank et al. (2002). The hard state, on the other hand, is characterized by a power-law X-ray spectrum extending beyond 100 keV, which is attributed to radiation from a geometrically thick, optically thin, hot accretion flow (also known as advectiondominated accretion flow, ADAF, or radiatively inefficient accretion flow, RIAF; Narayan & Yi 1994Abramowicz et al. 1995; see Yuan & Narayan 2014 for a review).
Several BH X-ray binaries (BHBs) have been observed during transitions from the thermal state to the hard state, and vice versa. From these observations, it is clear that the thermal state generally corresponds to a larger luminosity and mass accretion rate a few percent of Eddington, while the hard state corresponds to a lower luminosity and accretion rate (Tanaka & Shibazaki 1996;Done et al. 2007), though there is a notable hysteresis in the transitions Done et al. 2007). Esin et al. (1997) explained the thermal and hard spectral states as arising from different accretion disk configurations. According to their model, the thermal state corresponds to a configuration in which a thin accretion disk extends all the way down to the innermost stable circular orbit (ISCO), whereas in the hard state, the thin disk is truncated at some truncation radius outside the ISCO and the flow switches to an ADAF inside that radius. The Esin et al. (1997) model is broadly consistent with observations, though it does not readily explain the hysteresis phenomenon. Although some progress has been made on a theoretical explanation of hysteresis Liu et al. 2005), this topic still remains an open question. It is not the focus of the present paper.
Accreting supermassive BHs (SMBHs) in active galactic nuclei (AGN) again show different spectral states. The spectra of luminous AGN (quasars) are dominated by a big blue bump, consistent with a standard thermal thin accretion disk. Low luminosity AGN (LLAGN), on the other hand, have weak or absent blue bumps (Ho 1999), and appear to be consistent with a model in which the thin disk is restricted to large radii, with accretion at smaller radii occurring via an ADAF (see Yuan & Narayan 2014, and references therein). Compared to BHBs, whose M ∼ 10M black holes have short time scales, SMBHs with M ∼ 10 6 − 10 10 M have much longer time scales. Therefore we are less likely to witness a complete state transition within human time scales. Nevertheless, the class of "changing look" AGN (LaMassa et al. 2015) might well correspond to objects that are undergoing a thermal-to-hard (or vice versa) state transition (e.g., Noda & Done 2018). Interestingly, the state transition (if this is what one is seeing) in changing look AGN seems to occur at roughly the same Eddington-scaled luminosity as the corresponding transition in BHBs. Furthermore, in the hard state in both BHBs and AGN, the truncation radius between the outer thin disk and the inner hot ADAF, when scaled to Schwarzschild units, appears to show a similar dependence on the Eddington-scaled luminosity (Yuan & Narayan 2004). All this suggests that the state transition phenomenon is a robust feature of BH accretion, and that the basic physics is independent of BH mass.
While the model geometry described in Esin et al. (1997) is successful in explaining the spectra of the thermal and hard states (see, e.g., Poutanen et al. 2018), the question of how and why state transitions occur is not yet fully understood. It appears that a necessary ingredient is the presence of a hot magnetized corona above the cold thin disk. In seminal work, Galeev et al. (1979) suggested that amplified magnetic field in the disk causes loop-like magnetic structures to emerge due to buoyancy. In analogy with the solar corona, this will result in a very hot magnetically confined corona above the thin disk. Field amplification in the thin disk was later shown to be from the magneto-rotational instability (Balbus & Hawley 1991). Coronae have been widely invoked to explain the observation of intense X-rays from cold accretion disks (e.g., Haardt & Maraschi 1991;Svensson & Zdziarski 1994). Coronae also play an important role in studies of X-ray reflection spectra (Reynolds 2014). Coronae form spontaneously in MHD simulations of accretion disks, where they are generated by either the Parker instability (Machida et al. 2000) or dissipation of magnetic turbulence (e.g., Jiang et al. 2014). Coronae are also seen in general relativistic MHD simulations of black hole disks (e.g., De Villiers et al. 2003;Yuan & Narayan 2014).
Among several models for accretion state transitions, the best-studied is the evaporation model, which was first proposed by Meyer & Meyer-Hofmeister (1994, hereafter MMH94) for white dwarf accretion, and later adapted by Meyer et al. (2000, hereafter MLMH00) to explain state transitions in BHBs. In this model (hereafter "MM model"), mass and energy are exchanged in the vertical zdirection between the thin disk and the overlying corona. The hot corona thermally conducts energy downward to the thin disk and as a result matter in the disk "evaporates" and flows into the corona above. If the evaporation rate at some radius is large enough to remove all the gas in the thin disk, which can happen at sufficiently low mass accretion rates, the disk vanishes entirely and accretion at smaller radii occurs entirely via the corona. Liu et al. (1999) showed that the corona equations are equivalent to the ADAF model in the limit when there is no energy conduction to a cool disk. Therefore, once the thin disk is fully evaporated, the coronal flow automatically becomes a standard ADAF. The MM evaporation model and its subsequent extensions (e.g., Qiao & Liu 2009) have been successful in explaining X-ray observations of the state transition in BHBs. However, because of the complex vertical structure of the disk-corona system, as against the simpler height-integrated single layer structure of both the thin disk and ADAF models, studies of the evaporation model have been limited to numerically solving a set of complex differential equations in the vertical direction.
The goal of the present paper is two-fold. First, using simplifying approximations, we present a height-integrated version of the MM model which is suitable for analytical work as well as for simple numerical calculations. Using this height-integrated model, we analytically reproduce different regimes of the MM evaporation model and explain the relevant physics operating in each regime. We also study how the results depend on parameters such as the viscosity parameter α and the thermal conduction coefficient κ. Second, we show how to extend the height-integrated model to include additional physical effects beyond the original MM model. As an example, we consider the effect of direct magnetic heating of the corona by turbulent energy transport from the thin disk, and we present analytical results corresponding to this extension of the model. We also outline other extensions that could be explored in the future.
We provide an overview of the basic set-up of our model in Section 2, and we present analytical/numerical solutions in Section 3. We compare the analytical solutions with the original numerical approach of MMH94 and MLMH00 in Section 4, and describe the implications of our model for state transitions in Section 5. We include the effects of direct coronal heating and Compton cooling in Section 6 and present analytical and numerical results for this version of the model. We discuss applications to observations, comparisons to other models, and the limitations of our approach, in Section 7, and we summarize in Section 8.  The left half shows the black hole system in the hard state with a hot accretion flow (ADAF) at smaller radii and a cool thin disk with overlying hot corona at larger radii. The hot flow is shown as pink shaded circles and the cold disk is shown as a gray bar. The box in the right half is a zoomed-in cross-section of the thin accretion disk and corona at some radius. The disk mid-plane is at z = 0. The boundary between the thin disk and the transition region is at z d , the boundary between the transition region and the corona is at z t , and the scale-height corresponding to the top of the corona is at z = H. In our two-zone model, the equations are height-integrated separately for the transition region (z d to z t ) and the corona (z t to H). Quantities like Q + , Q − t , etc., are defined in Section 2.4 and Section 2.5. The total radial accretion rate in the corona isṀ c and the total thin disk accretion rate isṀ d . Each quantity is divided by 2 in the figure because the zoomed-in box shows only the top half of the system. The vertical mass flow rate per unit area from/to the corona is denoted asΣ z and is defined in Equation (7). Mass fluxes are shown with purple arrows and energy flows with gray arrows.
We begin by briefly reviewing the basics of the MM model. The left half of Figure 1 is a schematic of the accretion flow in the hard state. At large radii, accretion occurs via a sandwich structure in which a cool, geometrically thin equatorial disk (the gray bar) lies between two hot, geometrically thick coronal layers (shaded circles) above and below it. The cool disk evaporates fully and disappears at a certain radius (the inner edge of the gray bar), leaving the corona to continue inward as a hot accretion flow, labeled ADAF. In the thermal state, the sandwich structure extends down to the ISCO, and there is no separate ADAF zone.
The MM model deals with the vertical structure of the disk-corona region of the flow, indicated by the narrow vertical rectangle. Conservation laws are solved within this region, using appropriate boundary conditions at the bottom, where the corona meets the thin disk, and allowing for an outflowing wind at the top. The flow dynamics are simplified in a number of respects. Steady state (∂/∂t = 0) and axisymmetry (∂/∂φ = 0) are assumed. The angular velocity Ω and azimuthal velocity v φ in the corona are taken to be Keplerian, where M is the mass of the central BH and R is the cylindrical radius. The radial velocity v r is assumed to be given by where c s is the isothermal sound speed, P is the pressure, ρ is the density, and α is the dimensionless viscosity parameter. There is one additional key assumption in the model regarding radial derivatives, which is best explained by considering the mass conservation equation in cylindrical coordinates, The MM model retains the z-derivative in this equation but eliminates the R-derivative via an ansatz, where the particular coefficient −2 is motivated heuristically. A similar approximation is used to eliminate the R-derivative in the energy equation as well. The reason for this approximation is to reduce the problem to ordinary differential equations in z. The complete equations are given in Appendix A, where one more approximation, a convenient method of modeling the transition from cylindrical to spherical geometry for the outflowing wind, is also noted. With the above approximations, the MM model reduces to a set of four ordinary differential equations in z, which correspond to the conservation of mass, vertical momentum and energy, and an expression for the conductive energy flux in terms of the temperature gradient. Figure 2 shows a numerical solution of the equations for parameters similar to those used by MLMH00 in Figure 2 of their paper. Such numerical solutions may be computed at different cylindrical radii R, and then stitched together to explore the radial structure of the corona.

Vertically Integrated Two-Zone Model
Accretion disk models that use height-integrated equations are useful for many applications. Such "one-zone" models describe the disk at each radius in terms of a single average densityρ, single temperature T , etc. Thus, the sound speed in Equation (2) becomes c 2 s ≡ P/ρ in the one-zone model. Instead of solving for the detailed vertical structure, one computes an effective vertical scaleheight H using the condition of vertical hydrostatic equilibrium, and uses this to estimate the mass surface density Σ, Quantities likeρ, T , Σ are functions of the radius and satisfy ordinary differential equations in R. Often one can ignore the radial derivatives and solve directly for the quantities of interest at each radius, as is done successfully in the thin accretion disk model (Shakura & Sunyaev 1973). Alternatively, one can retain the derivatives and look for analytic (usually power-law) solutions, as is done in the ADAF model (Narayan & Yi 1994. Temperature T is in units of 10 9 K, vertical mass flux ρv z is in units of 10 −1 g cm −2 s −1 , vertical velocity v z is in units of sound speed c s , pressure P is in units of 10 8 dyn cm −2 , and conductive flux F c is in units of 10 16 erg cm −2 s −1 . A consistent solution satisfying all boundary conditions is obtained whenΣ z = (ρv z )(z = z d ) = 10 −1.30 g cm −2 s −1 and P (z = z d ) = 10 8.34 dyn cm −2 . Numerical calculations start from z = z d R and end at the sonic height z = z s , where v z ≈ c s and the conductive flux vanishes (F c ≈ 0). Note that the temperature rises steeply at low z, and becomes nearly constant (dT /dz ≈ 0) over a wide range of higher z. We identify the former as the transition region, and the latter as the corona. Also, note that the vertical mass flux ρv z at z = z s is very small compared to the flux at the bottom (z = z d ). This means that mass loss in a wind is fairly weak.
The MM disk evaporation model is not suitable for a one-zone approach since the numerical solution shown in Figure 2 has a complicated vertical structure. This is why previous efforts to utilize a onezone-like approach to describe re-condensation of a hot corona to an inner cool disk started with a one-zone ADAF solution and then modified it by considering the interaction with the underlying cool disk Meyer et al. 2007;Taam et al. 2008). However, a two-zone model appears promising. Consider the vertical temperature profile in Figure 2. Within a narrow range of z above the thin disk, the temperature T rises sharply until it is almost equal to the coronal temperature. Above that z, the temperature remains almost constant over a wide range of z. Thus we might be able to model the system in terms of two zones: a relatively thin transition zone lying just about the thin disk, and a geometrically thick corona above that. The two-zone structure is equally obvious in the downward conductive energy flux −F c . The flux starts at zero at the top of the corona and increases steadily with decreasing z until it reaches a peak value close to the thin disk, before plunging down to zero. In our work, we define the z coordinate corresponding to the peak of −F c as the boundary z t dividing the transition region and the corona.
Motivated by the above considerations, we consider below vertically integrated equations where we treat the transition region and the corona as two separate zones. This is shown schematically in the right-hand side of Figure 1. At a given R in our model, the thin disk extends vertically from z = 0 to z d , and the transition region from z d to z t , where both z d and z t are small compared to R. The corona zone extends from z = z t to H, where H is typically comparable to R. Mass and energy flow vertically between the different zones and also radially in the corona. Note in particular that there is a radial mass accretion rateṀ c in the corona, which is in addition to the mass accretion rateṀ d in the thin disk.

Vertically Integrated Continuity Equation
Using our simplifying approximations, we integrate the mass conservation Equation (3) vertically from z d to H. Contrary to the MM model, which eliminates the R-derivative via the ansatz Equation (4), we keep the R-derivative but eliminate the z-derivative by integrating vertically. In addition, since the numerical solution in Figure 2 shows that the vertical mass outflow rate (ρv z ) into a wind at the top of the corona is negligibly small, we assume in the height-integrated model that there is no wind: (ρv z )(z = H) ≈ 0. Therefore, the height integrated version of Equation (3) becomes is the mass evaporation rate (if the quantity is negative, its magnitude would represent the condensation rate) from the thin disk per unit area [g s −1 cm −2 ]. We neglect the radial mass flow in the transition region because of the small vertical thickness of this zone as well as the relatively low temperature, which suppresses the radial velocity (see Equation 2). Therefore, the entire evaporated mass fluxΣ z from the thin disk travels through the transition region and reaches the corona. Effectively, this means that the continuity equation in the transition region is trivial, and Equation (6) should be viewed as being applied only to the corona zone.
There are three unknowns at each radius R in our vertically integrated model: coronal temperature T (which determines v r ), surface density Σ =ρH, and mass evaporation rateΣ z . Equation (6) provides one ordinary differential equation to solve for these quantities. We need two more equations, which we obtain by vertically integrating the energy equation in the transition region and the corona, respectively, as we discuss next.

Vertically Integrated Energy Equation of the Transition Region
The full 2D partial differential equation for energy conservation in the entropy form is where (apart from an additive constant) is the entropy per unit mass, k is the Boltzmann constant, µ = 0.62m p is the mean molecular weight for a fully ionized gas of cosmic abundances (X=0.7 and Y=0.28), and γ = 5/3 is the ratio of the specific heats 1 . In the middle expression in Equation (8), the first term represents the radial advection of entropy and the second term describes the vertical transport of entropy. For the quantities in the right, q + and q − are the heating rate and radiative cooling rate per unit volume, and the last term is the energy loss rate via the divergence of the conductive heat flux F c . In the case of the transition region, we only consider the dominant energy terms. Radial advection and viscous heating are neglected because of the thinness of the transition region and the significantly cooler temperature compared to the corona. Vertical advection, the second term in the middle expression in Equation (8), is important. To compute the contribution from this term, we note that ρv z =Σ z is a constant in the transition region. Also, the gas is heated at constant pressure, and T (∂s/∂T ) P = c P 2 . Therefore, the integrated vertical advection term, which we call the evaporative cooling term, is given by Here, T represents the coronal temperature, T = T (z t ) = T (H), and we have assumed that the gas starts at zero temperature (the gas actually starts at the temperature T d of the thin disk, but T d T and is therefore negligible). Note that we use the symbol Q to represent heating/cooling rate per unit area, to distinguish it from the symbol q used earlier for the rate per unit volume. Also, note that we include a subscript 't' to indicate that the evaporative cooling rate computed here refers to the transition region. We later introduce a second evaporative (or vertical advection) contribution Q evap c for the corona. Considering next radiative cooling, we write the height-integrated term as Since the transition region spans a wide range of temperatures, cooling in this zone is very complicated and involves multiple atomic processes. Fortunately, we can use a simplifying approximation based on studies of the solar corona. Johnston et al. (2017) computed numerical solutions of the corona and transition region in the Sun and noted that the vertically integrated radiative loss in the transition region is approximately the same as the cooling loss in the corona (see Equation 13 and Figure 3 in their paper). Assuming their result to be valid for coronae in general, we set To verify the validity of this assumption, we compared thee magnitudes of Q − t and Q − c by vertically integrating the numerical solution of the MM model which uses a realistic cooling function. We confirmed that the two terms differ by only a factor of a few (∼ 2 − 3 at worst) over radii in the range r ∼ 10 3.5 − 10 4.5 where the MM model is most reliable. (As described in Section A, at radii r < 10 3 , the MM model does not self-consistently model the condensing corona, and at radii r > 10 5 , 1 The entropy formula in Equation (9) is derived from the first law of thermodynamics, T ds = du + P d(1/ρ), where u = P/((γ − 1)ρ) is the internal energy per unit mass for an ideal gas, and we treat γ as a constant.  describe a more careful approach to the entropy for the case when the pressure is not dominated by gas pressure alone. They present an effective adiabatic index γ in their Equation (17) which is a weighted sum of gas and magnetic adiabatic indices. 2 The constant pressure assumption in the transition region is in part motivated by the vertical structure from MMH94.
It is also an assumption used in the solar corona work of Shmeleva & Syrovatskii (1973) where MMH94's fourth boundary condition originates from (in Section A.2).
the coronal temperature approaches the lower temperature boundary condition of 10 6.5 K used by the MM model.) We note further that our assumption in Equation (12), which is inspired by numerical studies of the solar corona (Johnston et al. 2017), replaces a different approximation in the MM model, which also originates from solar corona research (Shmeleva & Syrovatskii 1973), where the fourth boundary condition in Section A.2 relates the conductive flux and pressure at the bottom of the transition region. We postpone discussion of the coronal cooling term Q − c to Section 2.5. The last term in the energy equation of the transition region is the conductive energy flux flowing into this region from the corona. As discussed earlier, the downward conductive flux −F c reaches its maximum value at the boundary between the corona and the transition region. Thus −F c (z t ) is the net conductive heating rate of the transition region per unit area. We therefore write where we have made the reasonable assumption that F c (z d ) = 0, i.e., the energy loss via conduction from the transition region into the thin disk at z = z d is negligible. Putting all the terms together, the height-integrated energy equation for the transition region is In words, the conductive heat flux flowing in from the corona Q cond provides energy to the transition region, and this energy is used partly to supply the evaporation energy Q evap t and the rest is lost via radiative cooling Q − t . Equation (14) is a second condition which the height-integrated model must satisfy and which we can use to solve for the three unknowns in the problem. Note that, compared to Equation (6), which is a differential equation, here we have an an algebraic equation with no radial derivatives.

Vertically Integrated Energy Equation of the Corona
When considering the corona, we keep all the terms in Equation (8) and integrate from z = z t to H. The first term is the radial advective energy term, which becomes after vertical integration, Note that quantities like Σ and H correspond to one side of the corona from the mid-plane. The temperature T is assumed to be independent of z in the corona, as described in Section 2.2. The second term in the middle of Equation (8) describes vertical energy advection. This is another component of evaporation energy, which we call Q evap c . As the gas rises in the corona, it moves at constant temperature T but gains entropy because its density declines steadily with increasing z. Using the entropy expression in Equation (9) and assuming a constant T , we estimate the evaporative energy per unit area to be It is hard to estimate the above integral precisely, but we expect its value to be similar to the corresponding term in the transition region, Q evap t . For convenience, and in the spirit of our toy model, we simply assume that the evaporation energy Q evap c in the corona is equal to Q evap t in the transition region multiplied by a factor of order of unity ζ. We thereby obtain the approximate expression given at the right of Equation (16). By vertically integrating the numerical solutions of the MM model, we find that ζ lies in the range 0.6 − 0.7 for radii in the range r ∼ 10 3.5 − 10 4.5 . Moreover, when we keep ζ and propagate it through the analytic solutions discussed later, the dependence on ζ is weak, e.g., T, Σ ∝ (2 + 5ζ) −1/3 for the direct magnetic heating model discussed in Section 6.3. Therefore, we assume ζ ∼ 1 for simplicity. With this approximation, we end up with Q evap c = Q evap t , hence we hereafter use Q evap without a subscript if referring to either term. Considering next the heating term in Equation (8), vertical integration gives In the present section, and in Sections 3 -5, we assume that the only source of heating in the corona is viscous dissipation (we consider an additional source of heating in Section 6). For this case, we have We write the pressure as where the coefficientp describes any additional contribution to the pressure over and above gas pressure P gas , andρ ≡ Σ/H is the average density of the corona. For example, Meyer & Meyer-Hofmeister (2002) consider the case when magnetic pressure is not negligible. In terms of the plasmaβ, they writep Considering next radiative cooling, the corresponding height-integrated term in the corona is where n e = (1 + X)ρ/(2m p ) and n = (X + Y /4)ρ/m p are the number densities of electrons and ions respectively (see below Equation (9) for the values of X and Y used), and Λ(T ) is the cooling function of the gas. In the present study we use where the coefficient corresponds to the bremsstrahlung curve in Figure 1 of Raymond et al. (1976). This estimate of the coefficient agrees with Equation (34.3) in Draine (2011) to within a factor of a few. We use bremsstrahlung cooling instead of a more detailed cooling curve. However, as we argue in Appendix C, this approximation is actually quite reasonable. In view of the approximation in Equation (12), the radiative cooling rate in the corona and the transition region are the same, and we refer to the corresponding rate as Q − Brem referring to height-integrated bremsstrahlung cooling rate.
We note that Compton cooling is likely to be important close to the black hole. We discuss this issue in Section 6.2 and we include an extra Compton cooling term in Section 6.3 when we discuss a generalized model where the corona is additionally heated by direct coronal heating.
The conduction term in Equation (8) is a divergence of the conductive heat flux. Therefore, when we vertically integrate, the result is the difference between the fluxes at the two boundaries. The conductive flux at the top of the corona is negligible (this is in fact one of the boundary conditions in the MM model). Therefore, we obtain the same formula for conduction as Q cond in Equation (13) In an unmagnetized plasma, the vertical heat flux is given by Allowing for a possible modest suppression of conduction, we thus write whereκ ≤ 1 is a numerical factor. An organized frozen-in field strongly suppresses cross-field conduction in a plasma. Since we expect the magnetic field in the corona to emerge from the underlying thin disk in a roughly vertical direction and to be swept back azimuthally, we do not expect much radial conduction, which is why we ignored the radial part of ∇ · F c . For the vertical direction of interest to us, we expect classical conduction to operate if the field is uniform, but the flux could be reduced if there are non-uniformities in the field. Narayan & Medvedev (2001) showed that a turbulent magnetic field can reduce conductivity from the classical value by a factor of a few, 1/5 κ 1. This is a reasonable range of values to consider in the present work. Combining all the terms, the vertically integrated energy equation for the corona is where the left-hand side corresponds to Q adv (R). Note that Q cond is an energy flux that flows out of the corona into the transition region. It is a cooling term for the corona, but becomes a heating term for the transition region. We now have three equations, viz., Equations (6), (14), (26), for the three unknowns at each radius, viz., T (R), Σ(R), andΣ z (R). Two of the equations are first-order differential equations in R, and the third (Equation 14) is an algebraic relation. The problem is thus reduced to a second-order differential equation in R.

Momentum equations
The three equations we derived above come from the continuity equation and the energy equation. While we have not explicitly mentioned the momentum equations, they are effectively used in the model. Following MLMH00, the azimuthal velocity is approximated to be Keplerian (Equation 1) using the R-component of the momentum equation, and the radial velocity formula Equation (2) is obtained from the φ-component of the momentum equation. Finally, the z-component of the equation is written as a condition for vertical hydrostatic equilibrium, whereby the formula for the scale-height H (Equation 5) is obtained.

Assumptions in the Current Model
Here we list the assumptions used for constructing vertically integrated equations in the model.
• In the transition region, we assume that the pressure is nearly constant in order to obtain the formula for the evaporative cooling term Q evap t in Equation (10).
• Due to the lower temperature of the transition region compared to the corona and the vertical thinness of this region, we neglect the radial advective term Q adv and the heating term Q + .
• We assume that the vertically integrated radiative cooling in the corona and the transition region are equal ( , except when we include Compton cooling as discussed in Section 6.3.
• For simplicity of the model, we also assume Q evap t ≈ Q evap c . From vertical solutions of MMH94, we confirm that this assumption is good.
• In the corona, we assume that the temperature is nearly independent of z, which is reasonable as can be seen in Figure 2.
• We neglect conductive flux at the bottom of the transition region F c (z d ) ≈ 0 to obtain Equation (13).
• We assume that ions and electrons are coupled such that the plasma is well described by a single temperature for both species.

RADIAL SOLUTIONS OF THE VERTICALLY-INTEGRATED CORONA MODEL
The three equations describing our model, Equations (6), (14), and (26), involve several terms and cannot be solved trivially. One approach is to focus only on the dominant terms, judiciously selected, and to look for analytic power-law solutions of the form Alternatively, we could retain all the terms and solve the equations numerically. We discuss both approaches in the following subsections and present some results. In the plots, we use M = 10M , α = 0.3,κ = 1,p = 1, unless otherwise stated. However, as we show, most of the results are independent of M when we scale quantities suitably.

Analytical Solutions: Evaporation-Dominated and Radiative Cooling-Dominated Regimes
We begin by considering the energy balance condition, Equation (14), in the transition region, which is an algebraic relation between a single heating term (conduction) and two cooling terms (evaporation and radiative cooling). It is natural to think that, depending on conditions, one of the cooling terms will dominate over the other. We can thus imagine two distinct regimes: (i) an evaporation-dominated regime where Q − Brem Q evap ≈ Q cond , and (ii) a radiative cooling-dominated regime where Q evap Q − Brem ≈ Q cond . We first consider the evaporation-dominated regime, where almost all the conductive heat flux from the corona is used to evaporate gas from the thin disk into the corona, and there is very little radiative cooling. Neglecting the cooling term in the equations, it is then possible to obtain a power-law solution with indices p = −1, q = −2, i.e., The analytic solution is written explicitly in Appendix Section B.1, where the radius R is scaled to the Schwarzschild radius R S as follows, From the analytical solution, we find that the two cooling terms scale as follows: We see that the evaporation term Q evap dominates at large R. Hence, the evaporation-dominated regime corresponds to radii far from the central black hole.
In the opposite case, when radiative cooling dominates over evaporative cooling, we can repeat the same analysis, but this time neglecting the evaporation term Q evap t . We find that we must also neglect the radial advection term Q adv in the corona (we discuss the validity of this approximation in Section 3.3). We then obtain a power-law solution with p = 0 and q = 0, corresponding to T ∝ R 0 = constant, Σ ∝ R 0 = constant. Once again, comparing the two competing cooling terms in this second regime, we find which is the same as before. It is clear that cooling dominates at small radii, i.e., closer to the black hole.
Using the analytic solutions for T and Σ in the two regimes, we can estimate the break radius r b ≡ R b /R S where the two zones meet. This can be done in two ways. First, we could look at the analytical solution for the evaporation-dominated regime, where Q − Brem is neglected, compute Q − Brem after the fact from the derived solution, and determine the radius at which Q − Brem becomes equal to Q evap . This is the radius at which the initial assumption Q − Brem Q evap breaks down. Similarly, we could do the converse, using the cooling-dominated analytical solution. These two independent calculations give where α 0.3 ≡ α/0.3. Note the identical scalings of the two results with respect to the parametersκ, α 0.3 andp, and the nearly identical coefficients. The close agreement indicates that the two radial zones we have identified are physically well-motivated. Hereafter, we define the break radius between the two regimes to be Note the interesting fact that r b is independent of the black hole mass M : the break between the cooling-dominated and evaporation-dominated regimes is located at the same Eddington-scaled radius for black holes of any mass. We now write the solutions for the coronal temperature T , scaling the radius by r b (note that the solutions written in Appendix Section B.1 do not scale by r b ). This gives for the two regimes We immediately see that, in Schwarzschild units, the temperature is mass-independent, meaning that the coronal states are the same for all BHs from stellar-mass BHs to SMBHs. In the left panel in Figure 3, the two analytic solutions for T are shown as red dashed lines. Also shown for comparison is the virial temperature, We see that the coronal temperature predicted by the model is subvirial over the entire range of radii, meaning that the corona is gravitationally well bound. This further supports our choice to ignore any wind escaping from the top of the corona. The analytic solutions for the surface density are which are again mass-independent. These solutions are shown in the right panel of Figure 3. However, since at a given r the scale-height H is proportional to BH mass, the average coronal densityρ = Σ/H does show a mass dependence, where m ≡ M/M is the dimensionless BH mass in solar units.

Radial numerical solutions
For comparison with the analytical solutions, we compute numerical solutions by solving the full set of equations without neglecting any terms. As a first step, we solve forΣ z from Equation (14) (energy equation for the transition region), which is the only algebraic equation. Then, we substitute this in the other two equations: continuity Equation (6) and the energy equation for the corona Equation (26). After substitution, we are left with two coupled differential equations, Σv r kT µ which we solve numerically 3 . Note that (RΣv r ) in the left-hand side of Equation (37) is a function of both T (via v r ) and Σ. Therefore, both Equations (37) and (38) contain radial derivatives of T and Σ on the left-hand side. These first-order differential equations are solved by starting with initial values of temperature T and surface density Σ at a large radius well into the outer evaporation-dominated regime. For the specific example shown in Figure 3, we chose the initial values to correspond to the analytic solutions Equations (33) and (36). However, the numerical solution is actually robust even when the initial values are changed by up to a factor of 10; in these cases, with decreasing r, the solution quickly adjusts and asymptotes to the solution shown in Figure 3. The numerical calculation is however sensitive to the direction of integration. Sometimes, we have found that integrating numerically from small to large radii fails, but integrating in the other direction is always stable. Figure 3 shows the numerical solution corresponding to the same parameters used for the analytical solutions. The numerical solution shows a break in slope at around r ∼ 10 3 , roughly where the analytical solution predicts the break r b between the evaporation-dominated and cooling-dominated regimes should be located (Equation 32). Also, the numerical solution agrees perfectly with the two analytical solutions, both in slope and normalization, at radii away from the break.  . Energy terms normalized by the conduction term Q cond from the radial numerical solution. Since we assume ζ ∼ 1 in Equation (16), we use Q evap without a subscript to refer to either Q evap c or Q evap t . In the inner radii, the dominant energy terms are viscous heating Q + , radiative cooling Q − Brem , and conduction Q cond . The advection term Q adv and evaporation term Q evap are negligible. In the outer radii, leaving aside Q − Brem , all the other terms, namely, Q + , Q cond , Q evap , and Q adv , are important.

Energetics
Here we investigate the contribution of each of the individual energy terms in the vertically integrated energy equations to obtain a better understanding of the physics of the system. Figure 4 shows how each energy term, normalized by the conduction energy Q cond , varies with radius in the numerical solution discussed in Section 3.2 (shown in Figure 3). It is evident that the system has two very distinct regimes as a function of radius. At radii r > 10 3 , all the terms are roughly comparable in magnitude, with the sole exception of the radiative cooling term Q − Brem , which is negligible. On the other hand, at radii r < 10 3 , the radiative cooling term Q − Brem quickly becomes one of the dominant terms, whereas both the radial advection term Q adv and the evaporative term Q evap become negligible. These trends are consistent with the approximations we used to obtain the two limiting analytical solutions described in Section 3.1.
Interestingly, with decreasing r, we see in Figure 4 that Q − Brem first overshoots before settling down to its asymptotic solution. The overshoot coincides with a region where Q evap changes sign, which signifies that there is a major change in the nature of the solution. Indeed, as we show in Section 4.2, the corona makes a transition at this point from evaporation to condensation, which is a major change in physics. This may explain in part why, with decreasing r in Figure 3, the numerical solution requires a factor of several in radius below the break before it settles down to the new power-law solution at small r.
The major differences between our model and the MM model are 1) the assumption in our model that there is no outflowing wind at the top of the corona, and 2) the manner in which we treat partial derivatives in the conservation laws. Each of these differences is discussed in the following subsections.

No-wind Approximation
The MM model considers disk-corona systems characterized by an outflowing wind at the top of the corona. This feature is reflected in their model as additional terms in the continuity Equation (A1) and energy Equation (A3), and also via boundary conditions at the top (see Section A.2). As described in Section 2, our model assumes that there is no wind escaping from the top of the corona, i.e., ρv z = 0 at z = H. We test the validity of this approximation by comparing numerical solutions as a function of z, similar to the work of MMH94, both with and without a wind. In Figure 5, the vertical flux of matter as a function of z at a fixed R is shown for both models. The "wind" solution is calculated using the same equations and boundary conditions as in MMH94, which are briefly summarized in Section A.1 and Section A.2. This solution passes through a sonic point at the top of the corona. For the "no wind" solution, we solve a similar set of equations, but modified so as to be consistent with the assumptions of our two-zone model. Specifically, we neglect the wind terms in the equations, and we modify the boundary condition such that the vertical speed goes to zero when z = H (see Section A.3 for details). Figure 5 shows that the two solutions are very similar at low values of z. They start to diverge at larger values of z near the top of the corona, but the density is low by this point. Note in particular that, at the sonic point v z ≈ c s , in the wind solution, the escaping vertical mass flux is only 2.0% of the initial mass fluxΣ z that flows into the corona at the bottom. Setting the mass outflow to zero in our model is unlikely to have a serious effect.

Mass Exchange between the Thin Disk and the Corona
As explained in Section 2.1, the MM model simplifies the problem by eliminating radial derivatives via the ansatz given in Equation (4). Here we examine what effect this approximation has on the solution.
Let us vertically integrate the mass conservation Equation (3). Ignoring mass loss in a wind (including it will not change anything), and applying ansatz Equation (4), we obtaiṅ where the final inequality follows from the fact that v r is negative (see Equation 2). As a direct result of using the approximation Equation (4), we find thatΣ z is always constrained to be positive, that is the numerical solution from the MM model always has mass evaporating from the disk into the corona.
Considering on the other hand our vertically integrated model, and using the power-law scalings Equation (27), the continuity Equation (6) becomeṡ The magnitude and sign of the mass fluxΣ z from the thin disk to the corona depend on the power-law indices p and q of T (R) and Σ(R). Using the solutions for these indices for the cooling-dominated (r < r b ) and evaporation-dominated (r > r b ) regimes, as discussed in Section 3.1, we finḋ We see that disk matter evaporates into the corona (Σ z > 0) in the outer regime (r > r b ). Moreover, the coefficient −3/2 is not very different from the coefficient −2 in the ansatz Equation (4). The MM model thus ought to give reasonable results in this regime.
However, for the cooling-dominated regime (r < r b ),Σ z is negative, i.e., gas from the corona condenses back onto the thin disk. This regime would appear to be inconsistent with the framework of the MM model. Nevertheless, the MM model does produce a change in the radial dependencies of quantities at about the correct break radius (see Figure 3 in MLMH00), and the logarithmic slopes at radii below the break appear to have similar values to those predicted by our vertically integrated model. Presumably, the MM model is able to achieve this, even though it assumes evaporation at all radii, because at radii below the break the cooling term Q − dominates over the evaporation term Q evap and so choosing the wrong sign for the latter term causes only a minor error.

TRANSITION BETWEEN DISK-CORONA SYSTEM AND HOT ACCRETION FLOW
Using our analytic solutions, we now consider the state transition problem. We begin by defining the mass accretion rate in the corona, where the extra factor of 2 is because we consider the total coronal accretion rate in both corona layers, whereas the Σ we defined earlier corresponds to only one side. When expressed in units of the Eddington mass accretion rate (defined with a fiducial radiative efficiency of 0.1), the dimensionless coronal accretion rates in the two analytical solutions arė As before, we see that the results are independent of the BH mass. The coefficients on the two sides differ by about 30%, suggesting that neither of the two analytical solutions is particularly accurate near r = r b . Equation (42) shows that the radial power-law index inṁ c is positive for r < r b and negative for r > r b . This means that the maximum coronal accretion rate occurs at around the break radius r b (consistent with MLMH00). Since our model ignores mass loss in a wind, the sum of the coronal mass accretion rateṀ c (R) and the thin disk accretion rateṀ d (R) should be independent of R. From Equation (42), the coronal accretion rate scales asṀ c ∝ R −3/2 at large R, which means that coronal accretion is negligible at very large distances. DefiningṀ d,0 ≡Ṁ d (∞) as the disk accretion rate at asymptotically large R, mass conservation givesṀ Since our model provides an estimate ofṀ c as a function of R, we can use Equation (43) to calculate how the accretion rate in the thin diskṀ d varies with radius. In Figure 6, the thin disk accretion rate,Ṁ d (R) =Ṁ d,0 −Ṁ c (R), is shown together with numerical and analytic coronal accretion rate estimatesṀ c for two values of the total mass accretion rateṀ d,0 . In panel (a), we consider the case whenṀ d,0 is 10% of the Eddington mass accretion rate. For this relatively high rate, the thin disk accretion rateṀ d decreases only slightly even at the peak of the coronal accretion rateṀ c . The thin disk loses a modest fraction of its mass with decreasing radius until r = r b and then it gains back the mass via condensation at smaller radii. The thin disk survives at all radii, so this model corresponds to the thermal state.
From Equation (42) and the numerical solution shown in Figure 6, we see that the maximum accretion rate in the corona is found at r ∼ r b , and is approximately equal tȯ The model shown in panel (a) in Figure 6 corresponds to a case where the Eddington-scaled total accretion rateṀ d,0 /Ṁ Edd = 0.1 is larger thanṁ c,max . Disk evaporation has only a minor effect on the system. Panel (b) shows a different example in which the total mass accretion rate is only 1% of the Eddington rate,Ṁ d,0 /Ṁ Edd = 0.01. In this case, the thin disk accretion rateṀ d decreases rapidly with decreasing r and goes to 0 at some radius outside r b ; in the example shown, this happens at r ≈ 2000. Interior to this radius (gray shaded region) there is no thin disk, and accretion occurs Here, the thin disk evaporates fully and truncates at r ≈ 2000 (which corresponds to R ≈ 10 10 cm for a 10M BH). At this radius,Ṁ c =Ṁ d,0 , andṀ d → 0. Radii below this are shaded in gray to indicate that accretion here is entirely in a hot accretion flow (there is no more thin disk). This model corresponds to the hard state.
entirely via a hot coronal flow. If the transition from the thermal state to the hard state is caused by a rapid decrease inṀ d,0 , one could imagine a scenario in which the inner thin disk does not have enough time to respond to the decreasing accretion rate and therefore survives temporarily. In this case, there is a possibility that some gas from the hot accretion flow might condense back on the the surviving inner disk. However, the inner disk will probably eventually disappear due to viscosity (Meyer et al. 2007;Liu et al. 2007), and once this happens, we will have the hard state.
In between the two choices ofṀ d,0 considered in Figure 6, the system under consideration will undergo a state transition if it has a time-varyingṀ d,0 (possibly caused by disk instability (Lasota 2001) or galaxy mergers). For the particular parameters chosen here, the transition would happen when the total mass accretion rate is around 2% of the Eddington rate.
We note an important implication of this model. Since matter condenses back onto the thin disk for r < r b , if the thin disk manages to survive down to r = r b with non-zeroṀ d , thenṀ d will only increase via condensation at smaller radii and so the disk-corona structure will survive all the way down to the BH. Therefore, the truncation radius where the thin disk transitions to a hot accretion flow should always be greater than r b . This prediction is in some tension with observations which suggest that the truncation radius is sometimes found at smaller radii (e.g., Yuan & Narayan 2004). We discuss this issue further in Section 7.1.

GENERALIZED MODEL WITH DIRECT MAGNETIC HEATING
The vertically-integrated two-zone description of the corona we have described so far is designed to match the basic version of the MM model. The main feature of the model is that it is simple and fast, while still giving similar results as the more detailed MM model (see Figure 12 in the Appendix). Taking advantage of this simplicity, we now generalize our model further.
It is believed that, in addition to viscous heating, coronae above thin disks are also directly heated by energy transported from the underlying thin accretion disk. The transport mechanism is likely via magnetic fields and the heating mechanism is probably through wave damping and magnetic reconnection of field loops, analogous to the situation in the solar corona (Haardt & Maraschi 1991Field & Rogers 1993). We now generalize our model such that, in addition to the viscous Q + visc heating term Equation (18) that we previously introduced, we include a second heating channel which we call "magnetic heating." For simplicity, we assume that a fixed fraction a of the viscously dissipated energy 3GMṀ d /8πR 3 [erg cm −2 s −1 ] in the thin disk (Frank et al. 2002) is transferred to the corona via magnetic field lines. The total coronal heating rate is then Liu et al. (2002bLiu et al. ( , 2016 investigated the value of the proportionality factor a (f in their notation). However, it is hard to obtain a first-principles estimate of a, so we treat it as a free parameter, constrained merely by the requirement 0 ≤ a ≤ 1. Even though we call the new direct heating term magnetic heating, there is no specific field strength that the model associates with this process; indeed, direct heating from the thin disk does not even have to be associated with magnetic fields (though, in analogy with the solar corona, it probably is). When obtaining solutions, once again we have the choice of either a numerical approach or an analytical approach. Numerically, the method is almost the same as in Section 3.2, the only difference being the addition of the new heating term Q + mag in the right-hand side of Equation (38). The thin disk accretion rateṀ d appearing in Q + mag is obtained as a function of radius through Equation (43): For the analytic solution, just as we previously solved a problem with only a single source of heating, viscous heating Q + visc (Section 3.1), we again solve the case where there is only one form of heating, now magnetic heating Q + mag (it is not possible to find a simple power-law solution when both heating terms are simultaneously present).

Analytical and Numerical Solutions when Heating is Dominated by Magnetic Heating
We begin by first identifying the radius range over which magnetic heating dominates. To do this, we use the analytical solutions obtained previously for the viscous-only heating problem (Section 3.1), and calculate from those solutions the ratio Q + mag (r)/Q + visc (r). For simplicity, we use the same fiducial parameter values as before, viz., α 0.3 =κ =p = 1, and add two new fiducial parameters, a = 0.1, M d = 0.1Ṁ Edd 4 . We then obtain We see that magnetic heating dominates over a wide range of radii, except in the vicinity of the break radius r b ≈ 10 3 . Our previous estimate of r b (Equation 32) is still valid since viscous heating dominates there, and our previous analytical solutions with pure viscous heating are roughly correct in the vicinity of r b (half an order of magnitude in r on either side). But everywhere else, magnetic heating dominates.
We now look for analytical power-law solutions of the model under the assumption that magnetic heating dominates. For this we substitute in Equation (26). As before, the solutions have two distinct regimes. There is a regime at smaller radii where radiative cooling Q − dominates over evaporation Q evap , and where the advection term Q adv is also negligible. Ignoring the less-important terms, we obtain an analytical solution. The second regime is at larger radii, where evaporation cooling is more important than radiative cooling, and where we obtain a different power-law solution. We do not write down the analytical solutions here because they suffer from an inconsistency which we discuss below.
The red dotted lines in Figure 7 show the analytical solutions for the temperature T and surface density Σ for the magnetic heating-dominated problem considered in the present subsection. We plot the analytical solutions only over the range r < 10 2.6 and r > 10 3.4 , omitting the region around r b ≈ 10 3 where viscous heating dominates.
We see an immediate problem with the large r solution shown in Figure 7, namely, the temperature exceeds the virial temperature T vir . Gas at supervirial temperature can easily escape a gravitationally bound system. Therefore, it seems that a magnetically heated corona with the fiducial parameters we have assumed will inevitably have a strong wind. The no-wind assumption that underlies our model thus breaks down. The model needs to be expanded to include wind physics before we can obtain trustworthy solutions at large r. We defer this to a future paper.
The small r solution does not suffer from a wind problem. However, the temperature increases with decreasing radius and so does the surface density (previously both were constant, see Equations 33,35). This means that, with decreasing radius, Compton cooling will become progressively more important and we will need to include this process in our radiative cooling term, which currently includes only bremsstrahlung cooling. We discuss the importance of Compton cooling in the following subsection. In addition, the temperature even reaches 10 9.5 K, where two-temperature plasma physics may need to be considered (e.g., Narayan & Yi 1995).

Compton Cooling
At radii r < r b , as kT approaches m e c 2 , Compton cooling could potentially be important. To check this, we estimate the Compton y-parameter (Rybicki & Lightman 1979), where τ es = κ es Σ is the electron scattering optical depth, with κ es = σ T n e /ρ ≈ 0.34 cm 2 g −1 (for our assumed composition).
The Compton y-parameter shows a large difference between the pure viscous heating model and the magnetic heating model. In the former case (Section 3.2), we find y < 1 at all radii (see Figure 8). Thus Compton cooling is less important 5 . However, with direct heating both T and Σ increase with decreasing r (see Figure 7), and so the Compton y-parameter becomes very large as we approach smaller radii (Figure 8). This is an indication that our neglect of Compton cooling is a serious limitation. In reality, as the temperature tries to increase with decreasing r, Compton cooling (and possibly pair production) will limit the temperature to lower values than the pure bremsstrahlung model predicts.
We therefore expand our model in this section to include a Compton cooling term in the equations. Since Equation (12) is based on the work of Johnston et al. (2017) which did not consider Compton cooling, the assumption of height-integrated bremsstrahlung in the corona being equal to heightintegrated bremsstrahlung in the transition region should still be valid. Thus, for the transition region we use In the corona, however, Compton cooling is important. Hence we write Comptonization is complex and non-linear if the scattering optical depth τ es is large or the Compton y-parameter exceeds unity (e.g., Rybicki & Lightman 1979;Sunyaev & Titarchuk 1980;Dermer et al. 1991). In the following, we simplify the problem by assuming that both τ es and y are reasonably small compared to unity (as we show later, this assumption is valid for our solutions). In this limit, Figure 8. Compton y-parameter as a function of cylindrical radius R as calculated from numerical solutions, for a = 0.1 andṀ d,0 = 0.1Ṁ Edd . The blue solid line represents the case when heating in the corona is purely from viscosity (Q + = Q + visc ), and the pink dashed line is when both viscous and magnetic heating are present (Q + = Q + visc + Q + mag ). The horizontal black dotted line shows y = 1, above which Compton cooling becomes dominant. Note that when magnetic heating is included with bremsstrahlung as the sole cooling mechanism, Compton cooling becomes important closer to the black hole. The new solution discussed in Sections 6.3 and 6.4, where the corona is heated by both Q + visc and Q + mag and cooled by both bremsstrahlung and Compton , is shown by the pink dotted line. After including Compton cooling consistently, the Compton y-parameter now remains less than 1 at all radii.
since the soft photon energy density, u ph , of the radiation above the photosphere of the thin disk is related to the total radiation flux, F d,tot , emerging from the disk as the Compton cooling term is We write F d,tot as the sum of the flux F d originally emitted by the thin disk and the backscattered flux F b which irradiates the disk from the corona, In the case of F d , we recall that a fraction a of the thin disk luminosity is assumed to be transported directly to the corona. The remainder emerges as radiative flux, hence The second equality in Equation (54) follows from standard thin disk theory: σT 4 eff,disk = 3GMṀ d /(8πR 3 ) (Frank et al. 2002).
In the case of F b , we first note that a fraction τ es (1 + τ es ) of the outgoing flux F d,tot from the disk is Compton-scattered by the corona (we keep the (1 + τ es ) correction factor for completeness, even though we have assumed that τ es is small). Each scattered photon gains in energy by a factor [1 + 4kT /m e c 2 (1 + 4kT /m e c 2 )], and half of these Comptonized photons move back toward the disk surface. These photons irradiate the disk and are in part reprocessed into blackbody radiation and are in part reflected (Haardt & Maraschi 1991). All of this irrradiation energy re-emerges from the disk with flux F b . We thus write We assume that the Compton-backscattered radiation remains soft and contributes fully to Compton cooling 6 . Substituting (55) in (53), we obtain and (52) then gives Having thus obtained an expression for Compton-cooling in the corona, the new set of equations describing the model consists of

Analytical Solutions when Heating is Dominated by Magnetic Heating and Compton Cooling is included
When Compton cooling is included and magnetic heating dominates, the analytical solution in the evaporation-dominated regime (outer radii) remains unchanged from the solution shown in Figure 7. This is because cooling is negligible at these radii, so it does not matter how we model the cooling.
However, Compton cooling strongly modifies the cooling-dominated regime (inner radii). Since Compton cooling dominates in the corona, we use the following approximate equations to calculate the analytical solution at smaller radii, In the transition region, the downward conductive energy flux balances bremsstrahlung radiative cooling, while in the corona, direct magnetic heating balances Compton radiative cooling. All other energy terms are negligible, as confirmed by the numerical results shown in Figure 11. Furthermore, since it turns out that the temperature is not relativistic (T ∼ few × 10 8 K (see below)) and the corona is optically thin τ es < 1, we use the simpler expression, and since y analytic 0.1, we use when deriving the analytical solution.
The full analytical solutions are given in the Appendix Section B.2. Here we show the solutions as a function of r/r b for the temperature, the surface density, In the above solutions, a 0.1 ≡ a/0.1 andṁ d,0.1 ≡Ṁ d,0 /0.1Ṁ Edd . Note that the solutions for both the temperature and the surface density are mass-independent, just as with the viscous heating dominated solutions (Section 3.1). The analytical solutions for the coronal temperature T and surface density Σ are shown as red doted lines in Figure 9. The results are consistent with observations of Seyfert 1 AGN coronae in that the mean temperature in these systems is measured to be kT = 65 ± 10 keV or T ∼ 7.5 × 10 8 K (Akylas & Georgantopoulos 2021) and no strong correlation has been observed between coronal temperature and black hole mass (Kamraj et al. 2022). From the above analytical solutions, we also obtain the coronal accretion rate in Eddington units, Correspondingly, the red dotted lines in Figure 10 show the analytical solutions forṁ c . Finally we note that 4kT /m e c 2 < 1 and τ es < 1 at all radii, hence the Compton y-parameter always stays below unity. This is confirmed by the numerical result shown in Figure 8. We note that, while the numerical model discussed in the next subsection is valid for any value of a over the range 0 ≤ a ≤ 1, the analytical solutions described here are consistent only when a is relatively small. This is because, as a approaches 1, the Compton cooling Q − Compt is suppressed by the (1 − a) factor (see Equations 54 and 64), while the direct heating Q + mag is enhanced by a factor a. The result is that the corona in the model heats up to a relativistic temperature and the approximation in Equation (63) breaks down. In addition, the Compton y-parameter approaches 1 so the simpler form of Q − Compt in Equation (64) is not valid (since (1 − y − τ (1 + τ )) −1 ≈ 1). For accurate results, the full expression for Q − Compt is required, and we need the numerical solution. However, the analytical solution presented here is within a factor of 2 of the numerical solution so long as a ≤ 0.5, which is a reasonable range for a.

Numerical Solution and Energetics
The blue solid lines in Figures 9 and 10 show the results we obtain by numerically integrating Equations (58)-(60), supplemented with Equation (43). In the corona energy equation we have included both the viscous and magnetic heating terms (Equation (45)), and we have set a = 0.1 anḋ M d,0 = 0.1Ṁ Edd , as for the analytical solutions. The numerical solutions agree fairly well with the analytical solutions, the agreement becoming better the farther we are from the break radius r b . At radii r > r b , the numerically computed temperature is even larger than the analytical solution, so the presence of a strong wind is even more certain. The small offset between the analytical (red dotted) and numerical (blue) lines at radii below r b is because the analytical model uses the simpler version of Compton cooling in Equation (64), whereas the numerical calculation uses the full expression in Equation (57) Figure 10 indicates that a model with an initial thin disk accretion rateṀ d,0 = 0.1Ṁ Edd does not evaporate fully when we include magnetic heating with a = 0.1 (just as in the previous case when only viscous heating was considered, see Figure 6). The maximum accretion rate in the corona for the case of Figure 10 isṁ c,max ≈ 0.064, but it is subject to change with different choices of a anḋ M d,0 . Magnetic heating causes the critical accretion rate below which the disk truncates to increase to aroundṁ crit ∼ 0.05 for the fiducial set of parameters,κ = α 0.3 = a 0.1 =p = 1 (compare with the critical rateṁ crit ∼ 0.02 for the pure viscous heating model in Equation 44). This revised estimate oḟ m crit should however be treated with caution since the regime we are considering likely has a strong wind in the outer radii, which we have not included in our model. Figure 11. Energy terms normalized by the direct magnetic heating term Q + mag for the numerical solution (blue lines) shown in Figures 9 and 10. The hatched region indicates radii with supervirial temperature. Figure 11 shows the variation of the different energy terms with radius, as determined from the numerical solution. The plot illustrates the complexity of the system when we include direct magnetic heating and Compton cooling. As predicted by the analytic work, viscous heating Q + visc dominates for radii near r b ≈ 10 3 , and magnetic heating dominates everywhere else. For r < r b , the only dominant terms are magnetic heating Q + mag and Compton cooling Q − Compt , which corresponds to the energy balance in the corona (Equation 62). Also note that in the inner radii bremsstrahlung cooling and conduction balance each other Q − Brem ≈ Q cond , which is the energy balance in the transition region (Equation 61). For r > r b , evaporation Q evap dominates over both radiative cooling terms Q − Brem and Q − Compt . All of these results are consistent with the assumptions we used in deriving the analytic solutions. Observational estimates of the thin disk truncation radius in hard state BHs (e.g., Yuan & Narayan 2004) show a wide range of radii from 10 − 10 4 R S . In the case of two low-luminosity active galactic nuclei (LLAGN),  argued that the truncation radius is at around 10 2 R S . According to our model, however, the thin disk can only truncate at a radius greater than the break radius r b (Section 5), which is ≈ 10 3 for our fiducial set of parameters, α 0.3 =κ =p = 1. As Equation (32) shows, r b would be smaller than 10 3 for other values of the parameters. For example, Rosner & Tucker (1989) and Narayan & Medvedev (2001) pointed out that the effective conduction coefficient in a collisionless turbulent magnetized plasma might be reduced by up to a factor of 5, i.e.,κ ≈ 1/5. This will cause r b to decrease by a similar factor. Additionally, our fiducial choicẽ p = 1 assumes that the pressure is entirely from the thermal gas. If there is a non-negligible magnetic pressure, then as Meyer & Meyer-Hofmeister (2002) noted, we will havep = (1 + 1/β) (Equation 20), where β is the plasma-β parameter. Shearing box 3D simulations of the magnetorotational instability show that the viscosity parameter α and the plasma-β are tightly correlated in differentially rotating disks such that αβ ∼ 0.5 (Blackman et al. 2008;Hawley et al. 2011). For our fiducial choice α = 0.3, this gives β ∼ 1.67 and a pressure correction factorp = 1.6. Since r b scales asp −3 , this reduces the break radius by a factor ∼ 4. Making α larger can produce even more extreme results, but our fiducial value of 0.3 is already rather large, and going higher would be unreasonable. Yuan & Narayan (2004); Cabanac et al. (2009) showed that for X-ray luminosities of BHBs L X 10 −2 L Edd , the truncation radius R tr 10R S , and as L X declines from 10 −2 L Edd to 10 −3 L Edd , R tr recedes farther and farther away from the ISCO. These observations are qualitatively consistent with our model (Section 5). The critical accretion rate above which the hard state disappears in our model when only viscous heating is present isṁ c,max ≈ 0.02 Eddington, and the limit is a factor of a few larger when we include magnetic heating (though predictions in this case are unreliable because the corona is likely to have a strong wind and we do not include this in the model). Furthermore, the analysis in Section 5 shows that the truncation radius should increase with decreasing mass accretion rate, and this is qualitatively in agreement with observations.

Strong Coronal Outflow
Strong accretion disk winds have been observed in many BHBs in the thermal state (Ponti et al. 2012). In the best-studied example, GRO J1655-40 (Miller et al. 2006(Miller et al. , 2008, the wind is launched at a radius of ≈ 10 3 R S . From the discussion of direct magnetic heating in Section 6.3, our model predicts for fiducial parameters that winds should be present for r > r b ≈ 10 3 , but not at smaller radii. Because the wind is not included self-consistently in our model, we cannot estimate the amount of mass loss in the wind. However, considering the highly supervirial temperature of our no-wind model, one suspects that the winds in these systems will be quite heavy. Similar winds are thought to be present also in AGN. Most AGN spectra show broad emission lines originating from what is called the broad line region (BLR) in the vicinity of the SMBH. A popular model is that the BLR lines are produced by outflowing gas in a disk wind (e.g. Lu et al. 2019;Matthews et al. 2020, and references therein). Kollatschny & Zetzl (2013a,b); Matthews et al. (2020) estimated the wind launching radius to be ∼ 10 3 R S , which is consistent with our model.
If the wind launching radius is around the break radius ∼ r b as our model predicts, then Equation (32) indicates an explicit dependence of the launching radius on the values of the parameters, α,κ,p, a, etc. In principle, it may be possible to use observations of disk winds to constrain these parameters. This might become practical if an extension of the current model that includes winds is developed.

Comparisons with Other Models
Various numerical models have been proposed based on the original work of MMH94 (see Liu & Qiao 2022 for a review). In later works, the same authors discussed the effect of changing the viscosity parameter (Meyer-Hofmeister & Meyer 2001;Qiao & Liu 2009) and reducing the conductivity of the plasma . A virtue of our analytic solutions is that they show explicitly how various quantities depend on the viscosity parameter and conduction coefficient. Our model is consistent with the previous numerical studies. An improved prescription for radial derivatives was introduced in Meyer-Hofmeister & Meyer (2003) in place of the approximation shown in Equation (4). However, they still limited themselves to evaporation, whereas our model shows that there is condensation at radii r < r b .
Building upon the MM evaporation model, Różańska & Czerny (2000a,b); Liu et al. (2002a); Qian et al. (2007) considered the effects of a two-temperature plasma and Compton cooling at small radii in the accretion system. A scenario where the Compton-cooled corona irradiates the underlying disk was studied by Liu et al. (2011);Qiao & Liu (2017), which is similar to our model set-up in Section 6.2. However, our model is different in that we additionally consider direct magnetic heating. Recently, Cheng et al. (2020) further improved the model by including radiation pressure and magnetic-reconnection heating in the model, and compared the predictions with observations of broadband spectra. The work presented here is different from these studies in that we provide simple and analytical solutions. Quantitative comparisons will be possible only after our model is expanded to include two-temperature effects.
Our disk-corona model is complementary to the work of Spruit & Deufel (2002) who investigated the evaporative process and transition to ADAF specifically at the truncation radius. Re-condensation of the ADAF back to an inner disk region closer to the black hole has been studied by Liu et al. (2006Liu et al. ( , 2007; Meyer et al. (2007); Taam et al. (2008) by starting from ADAF solutions and modifying them to consider disk-corona interactions.

Limitations of the Current Model
In §2.7, we listed a number of approximations we have made to derive the vertically-integrated diskcorona model presented in this paper. The resulting model is sufficiently simple that it permits us to obtain analytic solutions. These solutions reveal how the properties of the corona depend explicitly on parameters, and also how the physics varies between different regimes, e.g., evaporation-dominated vs. cooling-dominated regimes. While these are significant advantages, the model does have some limitations which demand more work.
At large radii r > r b , the version of our model with direct magnetic heating (Section 6), which we view as being more realistic than the simpler pure viscous heating version (Section 3), predicts supervirial temperatures. This is a sign that the model is not self-consistent. In practice, such systems will release the excess thermal energy by driving a wind to infinity, but our model explicitly assumes that there is no wind. We will need to modify the boundary conditions at the upper surface of the corona zone to allow mass and energy loss. If this could be accomplished without making the model overly complicated, it would allow a deeper understanding of disk winds in BHBs and AGN BLRs (see Section 7.2).
At radii r < r b , the ion and electron temperatures might deviate from each other, and twotemperature plasma physics effects may need to be considered. Energy transfer from ions to electrons could become inefficient, which might lead to a decrease in the conductive electron heat flux into the transition region and thus a decrease in the mass evaporation rate. A detailed investigation of this effect is left for future work. Additionally, if the magnetic field is strong enough, synchrotron cooling and synchrotron self-Compton (SSC) effects could become important as well. Simple one-zone prescriptions for synchrotron cooling have been developed for the ADAF model (Narayan & Yi 1995), and some of these may be of value for the corona problem.
In summary, the previous two paragraphs show that, when we include direct magnetic heating, there are serious problems at large radii (winds become important) and potential problems at small radii (two-temperature physics may need to be considered). However, the region around the break radius r b appears to be represented reasonably well by the present model, so our estimate of r b ≈ 10 3κ α −2 0.3p −3 (Equation 32) is probably safe.
Finally, we note that general relativistic gravity should replace our Newtonian gravity at sufficiently small radii, say r 10.

SUMMARY
The present work was inspired by the disk evaporation model proposed by Meyer & Meyer-Hofmeister (1994) and Meyer et al. (2000) to explain state transitions in black hole accretion disks. Their model, which we refer to as the MM model, considers a classical thin cold disk in the midplane, sandwiched between two hot coronae above and below it. Heat conduction from the the corona toward the thin disk causes gas in the thin disk to evaporate into the corona. The model calculates the structure of the corona and the mass evaporation rate by numerically solving a set of ordinary differential equations in the vertical direction. If the evaporation rate is large enough, the thin disk is truncated at a certain radius, and accretion proceeds via a hot accretion flow at smaller radii. This corresponds to the hard state (see also Esin et al. 1997). Alternatively, if evaporation is unable to eliminate the disk completely, the thin disk extends down to the black hole, and we obtain the thermal state.
In the present paper, we described a vertically-integrated version of the MM disk evaporation model, which allows both analytical and numerical exploration of the radial structure of the corona. Motivated by the distinct two-zone vertical structure of numerical solutions of the MM model (see Figure 2), we derived separate height-integrated equations for a vertically extended corona zone and a narrow transition zone between the corona and the thin disk. Assuming that cooling is dominated by bremsstrahlung emission, heating is purely by viscous dissipation (both as in the MM model), neglecting mass loss in a wind (the MM model includes winds), and applying conservation laws, the model reduces to three fundamental equations. These are Equations (6), (14), and (26), two of which are differential equations in radius, and one is algebraic. The main results of this pure viscous heating version of the model are as follows.
• For a fixed dimensionless radius r scaled to the Schwarzschild radius of the black hole, the solutions for the temperature T [K] of the corona, the vertically integrated surface density Σ [g cm −2 ], and the Eddington-scaled coronal mass accretion rateṁ c , are all independent of the black hole mass M (the volume average densityρ [g cm −3 ] scales inversely with M ). The model thus predicts that both corona formation and state transitions should be similar in stellar-mass and supermassive black hole systems, in rough agreement with observations.
• With suitable approximations, we obtain analytical solutions for the radial structure of the corona. As a function of r, we find two well-defined regions which are separated at a break radius, r b ≈ 10 3 (α/0.3) −2 , where α is the viscosity parameter in the corona. Equation (32) shows the scaling of r b with other model parameters.
• For r > r b , the analytical solution gives T ∝ r −1 , Σ ∝ r −2 ,ṁ c ∝ r −3/2 . We call this the evaporation-dominated region because the physics of the transition zone is dominated by the energetics of evaporation. For r < r b (cooling-dominated regime), we find that T and Σ are independent of r, andṁ c ∝ r 3/2 . In this region, gas condenses back from the corona onto the thin disk, and radiative cooling dominates the energetics of the transition zone. These analytical solutions agree with the results of detailed vertical structure solutions computed along the lines of the MM model ( Figure 12). The full analytical solutions are presented in Sections 3.1, 5, and Appendix B.1.
• Exact solutions as a function of radius, with no approximations, can be easily calculated numerically. They agree well with the analytical solutions in the respective regions ( Figure 3).
• The maximum mass accretion rate in the corona isṁ c,max ≈ 0.02 (α/0.3) 3 (Figure 6). If a thin accretion disk has an accretion rateṀ d,0 /Ṁ Edd <ṁ c,max at large r r b , the disk will be truncated and the system will be in the hard state. IfṀ d,0 /Ṁ Edd >ṁ c,max , the thin disk will extend down to the black hole and the system will be in the thermal state.
• The model predicts that in the hard state the disk truncation radius is always far from the black hole at r > r b . This disagrees with apparently smaller truncation radii observed in some systems (Section 7.1).
• The solutions of this version of the vertically-integrated model are consistent with our use of pure bremsstrahlung cooling and our neglect of mass loss in a wind. The model is thus fairly self-consistent at all radii.
Apart from providing analytical solutions, physical understanding, and an easy route for calculating numerical solutions, the vertically-integrated model has another advantage over the original MM model. It is simple to include additional physical effects in the model. In Section 6, we included a second source of heating (on top of viscous heating), namely, direct transfer of heat from the thin disk to the corona (presumably via magnetic fields). We find that with this additional heating, Compton cooling becomes important. We therefore include Compton cooling in addition to bremsstrahlung cooling. The resulting equations are nearly as simple as those of the pure viscous heating model and give the following results.
• The solutions for T , Σ andṁ c continue to be independent of the black hole mass, and the break radius remains roughly the same as before, r b ≈ 10 3 (α/0.3) −2 .
• For r > r b , in regions where direct heating dominates over viscous heating, we find an analytical solution in which T ∝ r −1/2 , Σ ∝ r −3/2 ,ṁ c ∝ r −1/2 . The full analytical solutions are presented in Section 6.3 and Appendix B.2. In this solution, T exceeds the virial temperature, which suggests that there should be strong mass loss in a wind (Section 7.2). The heaviest mass loss is likely to be at radii r ∼ r b , in agreement with observations of disk winds in stellar-mass black holes and the broad-line region in supermassive black holes. However, the model here is inconsistent since it does not include the effect of winds.
• The critical accretion rate which separates the hard state (at lowerṁ) and the thermal state (at higherṁ), increases toṁ crit ≈ 0.05 for α = 0.3,κ =p = 1, a = 0.1, but this result should be taken with caution.
The inconsistency identified in the second bullet point of the direct-heating model indicates that further improvements are needed; specifically, we need to include coronal winds at large radii. Given the simplicity of the framework, we anticipate that such an improvement will be tractable.

A.1. Equations
The numerical description of disk corona evaporation model first proposed by MMH94 involves solving differential equations in the vertical direction at a given radius. The full solution as a function of radius is obtained by repeating the vertical structure calculation at multiple radii and stitching the solutions together. In 3D cylindrical coordinates, they assume axisymmetry and introduce an approximation for radial derivatives to reduce the equations to a set of ordinary differential equations, only leaving z-derivatives. Specifically, the approximation for the radial derivatives of some quantity where one example for the mass conservation equation was shown in Equation (4). The mass conservation Equation (3) is then modified to The very last term on the right-hand side is added to self-consistently model the spherically expanding wind. Also, the z-component of the momentum equation simplifies to As described in Section 2.6, the radial and azimuthal momentum equations are used to define the angular velocity Ω and radial velocity v r respectively. The energy equation is modified with the similar approximation of the radial derivative to become where = ρ v 2 2 + γ γ − 1 P − ρ GM (R 2 + z 2 ) 1/2 is the energy density. Instead of bremsstrahlung, they used radiative cooling function Λ(T ) of an optically thin low-density plasma from Raymond et al. (1976). The last term in the Equation (A3) is again the due to the spherically expanding wind. Finally, the last equation is the conduction equation One can numerically solve for four unknowns P , T , v z , F c as functions of z with four differential Equations (A1)-(A4) and the following boundary conditions.

A.2. Boundary Conditions
In the four boundary conditions below, the subscript '0' indicates quantities at the bottom of the transition region (z = z d in our notation) and the subscript '1' indicates quantities at the top of the corona z = H. It is a free boundary problem where the height of the corona (H) is set to be at the sonic point (v z = c s ).

A.3. Vertical Numerical Solutions for No-wind Case
In Section 4.1, we compared our no-wind approximated model with the MM model. Here we describe the method used to obtain this no-wind solution, which is slightly modified from the MM model to be consistent with our vertically-integrated two zone model. We also confirm in this subsection that the obtained no-wind numerical solutions are in good agreement with our two-zone analytic solutions in Section 3.1.
In detail, the wind terms (the last term in Equations A1 and A3) are neglected and the first boundary condition in Section A.2 is adjusted such that the vertical speed is zero at scale-height v z (z = H) = 0. As we observed from Section 4.2, the radial derivative approximation of the MM model is not valid for r < r b . Therefore for the vertical numerical solutions here, we use new simplifications that are consistent with Equation (40 and the energy equation is also similarly tuned according to analytic scaling relations. The scaleheight H is set when the pressure is half of the pressure at the bottom of the transition region P (z = H) = P 0 /2. The lower boundary to start the numerical calculation is z d = R/300. The obtained no-wind numerical solutions (green solid line) are shown in Figure 12 along with the two-zone analytic solutions (red dashed lines) for the temperature and the average density. From the vertical numerical solutions, the surface density is approximately Σ ≈ ρ 1 kT 1 /µ/Ω K at each R so this quantity is plotted in panel (b). Because of the discontinuity at r b in the Equation (A5), there is a noticeable discontinuity around r b in numerical solutions. Note that similar equations are used for obtaining both solutions, but the applied simplifications are different. The no-wind numerical solutions are solved by making an approximation on radial derivatives (simplifying in R−direction similarly with the MM model) and the analytic solutions are obtained by vertically integrating the equations (simplifying in z−direction). The two solutions are in a very good agreement, meaning that our simple height-integration approach is reasonable and can reproduce the detailed calculation of the MM model. Note that at large R, the numerical solutions diverge from the analytic solutions. This is because as R increases, coronal temperature decreases and approaches the chromospheric temperature of T 0 = 10 6.5 K set at the bottom of the transition region (Section A.2) as a boundary condition.

B. EXACT FORMS OF ANALYTIC SOLUTIONS
The exact forms of analytic solutions are presented in this section. The solutions look slightly different here because in the main text, the dimensionless radius r is normalized by the dimensionless break radius r b . Each analytic solutions is presented in both regimes, with the first one for the inner cooling-dominated regime (r < r b ) and the second one for the outer evaporation-dominated regime (r > r b ).

B.1. Viscous heating only
The analytic solutions correspond to Section 3.1 with viscosity as the only source of heating Q + = Q + visc . In this section we discuss validity of using a bremsstrahlung cooling function in our analytic approach (Section 2) instead of a more realistic cooling function. In numerical studies of the MM model, the cooling function of a optically thin low-density plasma (Raymond et al. 1976) was used. In Figure 1 of Raymond et al. (1976), the cooling function Λ(T ) follows a bremsstrahlung curve above T 5 × 10 7 K.
As shown in Figure 3(a) or Figure 9(a), the coronal temperature interior to r b always have high temperature (T 10 8 K) regardless of which heating mechanism dominates. Thus, our assumption is reasonable in the inner regime (r < r b ) since cooling function asymptotes to bremsstrahlung curve at such high temperatures. However, exterior to r b , the temperature can get as low as 10 7 K very far from the black hole (for example, r 10 5 in Figure 3(a)). For the lower temperatures, atomic processes start to contribute to the cooling function so in principle realistic cooling curve should be used instead of a simple bremsstrahlung. However, analyses of energy (Figure 4 and Figure 11) showed that for r > r b , the radiative cooling (Q − Brem or Q − Compt ) becomes negligible compared to other energy terms so the choice of cooling functions is not as important. Therefore, we conclude it is reasonable to use bremsstrahlung as the only radiative cooling mechanism, except in the close vicinity of the BH in magnetically heated model where Compton cooling becomes important.

D. SATURATED CONDUCTION
The conductive heating flux formula in Equation (24) is taken from the classical theory of thermal conduction. However, the theory is based on an assumption that the mean free path λ is much smaller than the temperature scale length L T . Therefore, the classical conduction cannot be applied in certain cases where λ > L T and the conduction reaches an upper bound or saturates (Cowie & McKee 1977;Komarov et al. 2018). Cowie & McKee (1977) states that the upper bound of the conductive flux exists because the energy flux is limited to the thermal energy nkT times the thermal velocity. The mean free path of electrons mediating conduction is (Cowie & McKee 1977) where the electron-electron equipartition time t eq is (Spitzer 1962) t eq = 3m 1/2 e (kT e ) 3/2 4π 1/2 n e e 4 ln Λ .
Comparing the mean free path to the scale-height of the corona H, λ H ∼ 10 −2 for radial numerical solutions where the viscosity is the only source of coronal heating (Section 3.1) and λ H ∼ 10 −1 when magnetic direct heating is added. In principle λ should be compared with the temperature scale length L T which is approximately the height of the transition region z t . The scale length L T is around 0.1H, so we conclude here that the classical conduction theory is marginally applicable to our system. Różańska (1999) also confirmed that there is no need to consider saturation of conduction in the context of AGN or galactic black holes.