2 The evaporation process

Any mechanism heating the protons in the disk, even if only small, can potentially give rise to a thermal instability since the disk protons do not loose energy efficiently. Their energy loss by radiation is negligible on account of their high mass, and the transfer of energy to the disk electrons by Coulomb interaction is inefficient in a sufficiently hot plasma. In addition, by the nature of the Coulomb interaction, the time scale for the disk protons to equilibrate with the disk electrons increases with temperature (cf. Sect. 1.3). Depending on the strength of the mechanism heating the protons and its temperature dependence, a runaway may occur in which the disk protons continue to heat up as their cooling rate continues to decrease.

One heating process to be taken into account is internal viscous energy release (due to friction). At the low surface densities in the region of interest, the amount of energy released by friction is small compared with the heating flux by the external ISAF protons. If we may assume that a reasonable fraction of the viscous energy release goes into the protons, this small amount can still be quite important for the energy balance of the disk protons because of their low energy loss rate, and must therefore be taken into account. A second heating process is energy transfer from the incident ISAF protons to protons in the disk by Coulomb interaction. Though we have seen (Sect. 1.3) that the incident protons loose their energy mainly to the electrons, the small amount transferred to the disk protons may conceivably be relevant for their energy balance. In the following we investigate these heating mechanisms in turn, an conclude that viscous heating can give rise to a runaway, while heating of the disk protons by Coulomb interaction with the incident protons does not.

   
2.1 Viscous heating of the disk protons, instability

The energy balance of the protons in the warm disk is conveniently described in terms of heating and cooling time scales due to the processes involved. First we consider the cooling time scale by transfer of energy from protons to electrons in the disk.

The energy exchange timescale by Coulomb interactions between nonrelativistic protons and electrons with Maxwellian distributions at temperatures $T_{\rm p}$ and $T_{\rm e}$, respectively, is given by Spitzer's (1962) classical result:

$$t_{\rm ep}=\frac{3~m_{\rm p}~k^{3/2} (T_{\rm e} + \frac{m_{\... ...3~ \frac{T_{\rm e}^{3/2}} {n_{\rm e} \ln\Lambda}\;\;{\rm s},$$ (8)

where $m_{\rm p}$, $m_{\rm e}$ are the mass of the proton and electron, respectively, k is the Boltzmann constant, e the charge of the electron, $n_{\rm e}$ the electron density and $\ln \Lambda \approx 20$the Coulomb logarithm. The approximate equality is used since we consider cases where initially the electron temperature is of the order of the proton temperature and therefore the contribution of $\frac{m_{\rm e}}{m_{\rm p}} T_{\rm p}$ can safely be neglected.

The relativistic correction to this result is small for the conditions of interest here. In Deufel et al. (2002) we have shown that Spitzer's formalism is accurate to better than 5% compared to the relativistic treatment by Stepney & Guilbert (1983), for proton energies below 100 MeV and electron temperatures $kT_{\rm e}\la 50$ keV.

Next, the viscous heating timescale in an accretion disk is given by

$$t_{\rm th}=\frac{1}{\alpha \Omega}\;,$$ (9)

where $\Omega = (G M / R^3)^{1/2}$ is the Kepler angular velocity and $\alpha$ is the viscosity parameter (Shakura & Sunyaev 1973).

Comparing Eqs. (8), (9) yields a critical electron density

$$n_{\rm e}^{*}= \frac{3~m_{\rm p}k^{3/2}}{8\sqrt{2\pi~m_{\rm e}}~e^4~\ln\Lambda}~T_{\rm e}^{3/2} ~\alpha~\Omega\;.$$ (10)

For electron densities $n_{\rm e}<n_{\rm e}^{*}$ viscous heating in an accretion disk works faster than proton-electron coupling. In those regions the protons can not loose their energy fast enough, and heat up. This increases $n_{\rm e}^*$. At the same time, the higher proton temperature causes the layer to expand, by hydrostatic equilibrium. On both accounts, the Coulomb coupling between protons and electrons decreases, and the heating of the protons accelerates. Thus we expect viscous heating of the protons to lead to thermal instability wherever the density is less than given by Eq. (10). Since the density decreases steeply with height in the atmosphere, there is always a level above which it is unstable. A new equilibrium is reached only when the protons reach the virial temperature. The unstable part of the atmosphere has then expanded to $H/R\sim 1$, and it is then, effectively, part of the ISAF in which the disk is embedded. The unstable part of the atmosphere has evaporated, feeding the ISAF. Whether such an instability has an effect on the global accretion properties depends on the mass in the unstable region.

The amount of mass involved in the instability is the mass in the atmosphere above the level where the density has dropped to $n_{\rm e}^*$. Let this mass be ${\mathcal M}^*$ (per unit surface area of the disk). The time scale on which the protons heat up in the unstable region is just the viscous heating time scale $t_{\rm th}$. The rate at which the atmosphere evaporates is thus approximately

$$\dot\mathcal M \approx \alpha\Omega \mathcal M^* \;.$$ (11)

The unstable mass ${\mathcal M}^*$ is found by integrating the density in the atmosphere upward from the level z^* where n=n^*:

 

$$\mathcal{M}^* = \int_{z^*}^\infty m_{\rm p}~n_{\rm e}(z)~{\rm d}z = m_{\rm p}\;n_{\rm e}^*\; ~H \sqrt{2}~ f(u^*) \;\;,$$ (12)

where H is the scale height of the disk atmosphere (assumed isothermal) with density stratification, $n_{\rm e} =n_0 \exp[-u^{2}]$ with u² = z²/2H². The function $f(u^\star)$ is

 

$$f(u^\star)\! = \!\int_{u^\star}^{\infty}\exp[u^{\star2}-u^2]~{\rm... ... \frac{\sqrt{\pi}}{2}{\rm e}^{u^{\star2}}\!\left(1-{\rm erf}[u^{\star}]\right).$$ (13)

To evaluate Eq. (13) we need to know the dimensionless critical height above the midplane, $u^\star$. The level at which $n_{\rm e}=n_{\rm e}^*$ depends on the surface density of the disk, being higher at large surface density. We measure the surface density by the Thomson scattering optical depth $\tau_{1/2}$, measured from the midplane to the surface of the disk. In terms of the density at the midplane n₀, this is:

$$\tau_{1/2} = n_0 ~ \sigma_{\rm T} H~\sqrt{\pi/2}\;.$$ (14)

With Eq. (10) we can then calculate $u^\star$:

$$u^{\star}=\left[-\ln\left(\frac{\pi~m_{\rm p}^{1/2}~k^2} {\sq... ...rac{\alpha ~T_{\rm e}^2}{\tau_{1/2}}\right)\right]^{1/2}\cdot$$ (15)

Above a certain temperature (for a given $\alpha$) u^* is not defined, as the value within the square brackets in Eq. (15) drops below zero. This means that the whole atmosphere is subject to the instability, and we can set $u^\star=0$. In this cases $f(u^\star)=\sqrt{\pi}/2$ (the maximum value for f).

Using Eq. (12) to replace $\mathcal{M}^\star$ in Eq. (11) this yields, with Eq. (10):

 

$$\dot\mathcal{M} = \frac{3}{4\sqrt{2\pi}}\frac{m_{\rm p}^{3/2}~ k^2}{m_{\rm e}^{1/2}~e^4~\ln\Lambda} ~\alpha^2~\Omega~T_{\rm e}^2~f(u^\star).$$ (16)

Here $H=c_{\rm s}/\Omega$ and $c_{\rm s}=\sqrt{2~k~T/m_{\rm p}}$ is the isothermal sound speed.

An evaporating part exists in every disk atmosphere, but at low temperatures only the highest layers evaporate. The mass loss is unimportant in such a case and can not change the properties of the disk. Equation (16) shows that the evaporation rate is a relatively strong function of R, $\alpha$ and $T_{\rm e}$. Whether situations exist where this evaporation rate becomes important is investigated in Sect. 3.

   
2.2 Heating by proton-proton interactions

In this subsection we show that Coulomb interaction of the incident protons with the resident disk protons is relatively unimportant for the conditions encountered. It is included in the numerical results reported in the next section, however.

The main energy channel is from the hot protons to the electrons in the warm surface layer (see Sect. 1.3 and Paper II). The properties of the Coulomb interactions at the temperature of the incoming protons are such that only a small fraction of their energy is transferred to the protons in the layer (which we will call "field protons'' here following the terminology in Spitzer 1962). The energy budget of the field protons is also small, however (cf. Fig. 1), so the effect of this channel on the proton temperature in the warm layer has to be considered.

The interaction of a fast proton with the protons in a much cooler plasma is given by Spitzer's result discussed in Sect. 1. In terms of the rate of change of the energy of the incoming proton, it can be expressed as

$$\frac{{\rm d}E_{\rm p'}}{{\rm d}t}=-\frac{8\pi~e^4}{m_{\rm p} v_{\rm p'}}~n_{\rm p}~\ln\Lambda~[\psi(x)-x\psi'(x)]\;,$$ (17)

where $n_{\rm p}$ is the density of the (cool) field protons. The subscript ${\rm p}'$ indicates the hot penetrating protons and ${\rm p}$ the field protons. Here $\psi(x)$ and $\psi'(x)$ are the error function and its derivative and $x^2=m_{\rm p}v_{\rm p}^2/2kT_{\rm þ}$ is the ratio of the velocity of the incoming protons to the thermal velocity of the field protons. In the following we set $\psi-x\psi'=1$. This approximation is valid if the velocity of the incoming protons much exceeds the velocity of the field protons. In our case we consider virialized protons penetrating into the considerably cooler disk plasma, so the approximation is valid.

The heating rate $W_{\rm p'\!p}$ per unit volume due to the interaction of hot protons with the field protons is then

$$W_{\rm p'\! p} = n_{\rm p'}\cdot\frac{{\rm d}E_{\rm p'}}{{\rm... ...p} v_{\rm p'}^2}~\ln\Lambda~n_{\rm p} ~n_{\rm p'}~v_{\rm p'}.$$ (18)

This can be compared with the rate of energy transfer from the field protons to the electrons in the warm layer as given by Eq. (5-30) from Spitzer (1962),

 

$$W_{\rm pe} -~n_{\rm e}~k~\frac{T_{\rm e}-T_{\rm p}}{t_{\rm ep}}$$ =

$$\simeq - \frac{8(2\pi~m_{\rm e})^{1/2}~e^4~ n_{\rm e}^2 \ln\Lambda}{3~m_{\rm p}} \frac{kT_{\rm e}-kT_{\rm p}}{(kT_{\rm e})^{3/2}}\;,$$ (19)

where we have used Eq. (8) to replace $t_{\rm ep}$. The approximate equality indicates that the contribution from ${m_{\rm e}}/{m_{\rm p}}T_{\rm p}$ can again be neglected. Using $v_{\rm p'}=\sqrt{2kT_{\rm p'}/m_{\rm p}}$, we compare Eqs. (18), (19) and obtain an estimate for the equilibrium temperature of the disk protons exposed to penetrating hot protons:

 

$$T_{\rm p} \approx T_{\rm e} + 114 \left(\frac{n_{\rm p'}}{n_{\rm e}} \right) \left(\frac{T_{\rm e}}{T_{\rm p'}} \right)^{\!\frac{1}{2}} ~T_{\rm e}.$$ (20)

If the mass flux in the ISAF and its viscosity parameter $\alpha$are of the same order as in the cool disk, a thin disk approximation for both flows yields $n_{\rm p'}/n_{\rm e}\sim (T_{\rm e}/T_{\rm p'})^{3/2}$. In the inner region of the accretion flow, where $T_{\rm e}\sim 100$ keV and $T_{\rm p'}>10$ MeV, the second term in (20) is then negligible.

We conclude that heating of the ions in the warm disk by the incident ions alone is not sufficient to drive the ion temperature away from equilibrium with the electrons. For this to happen, viscous heating has to be included.

2.3 Numerical simulation of the evaporation

We test the inferred heating instability numerically, starting with a warm disk of moderate optical depth and investigating its temporal evolution due to internal viscous heating, and p'-p and p-e energy exchange as described in the previous subsection.

We use a plane-parallel, one-dimensional model geometry. The vertical density distribution through the atmosphere is found from the equation of hydrostatic equilibrium as in Paper II, except that here we do not account for the pressure exerted by the decelerating protons.

The temporal evolution of the temperature profile is computed with a fourth order Runge-Kutta method. The timestep $\Delta t$ of each integration is set to the shortest timescale of the heating and cooling processes involved. The change of temperature per timestep within a certain volume is due to the rate of change of the enthalpy there, $\Delta_{\rm t} w = \rho c_{\rm p} \Delta_{\rm t} T$.

The temperature increase per timestep as a function of optical depth due to viscous heating of the protons is then given by

$$\frac{\Delta T_{\rm visc}(\tau)}{\Delta t}= \frac{9}{4}~\alph... ...\tau)~ \Omega_{\rm K} \cdot \frac{1}{\rho(\tau)~c_{\rm p}}\;,$$ (21)

where $c_{\rm p}$ is the specific heat at constant pressure.

The temperature change of the field protons due to the impinging hot protons is computed similar to the proton illumination method described in detail in Papers I and II. Thus we follow protons from a Maxwellian distribution with virial temperature through the warm atmosphere of the disk and record their energy losses. The energy flux in the illuminating protons is related to the mass condensation rate (see Eq. (26) and Sect. 3 for details). The kinetic energy decrease of the fast protons (p') due to their interactions with the field protons (p) as a function of optical depth can be derived from Eq. (17). This gives the heating rate $\Lambda_{\rm p'p}$ and the contribution of p'p exchange to the temperature change per time step is

$$\frac{\Delta T_{\rm p'p}(\tau)}{\Delta t}= \frac{\Lambda_{\rm p'p}}{\rho(\tau)~c_{\rm p}}\cdot$$ (22)

Finally the rate of change of the temperature due to the electron-proton [ep] coupling is (according to Spitzer (1962))

$$\frac{\Delta T_{\rm ep}(\tau)}{\Delta t}= \frac{T_{\rm e}(\tau)-T_{\rm p}(\tau)}{t_{\rm ep}(\tau)}\cdot$$ (23)

To complete the model the electron temperature in the disk needs to be specified, which is regulated by the radiation processes. We have considered these in some detail in Paper II, where we found that for the present conditions electron temperatures of 100-500 keV result. The radiation losses by inverse Compton scattering increase quite rapidly with $T_{\rm e}$ and with the increasing optical depth when pair production sets in. This makes the electron temperature relatively insensitive. For the present purpose, it is sufficient to keep the electrons at a constant temperature, which we treat as a parameter of the model.

Figure 2: Evolution of the temperature profile due to the viscous instability of the protons. The model layer has optical depth $\tau _{1/2}=1.5$ and initial temperature $T_{\rm0,p} = T_{\rm0,e} = 500$ keV (dashed line). Dashed-dotted lines show the profile after 10 thermal time scales, solid lines after t=20. Values of the viscosity parameter $\alpha$ are 0.1, 0.3 and 0.5 (numbers at the lines). For comparison the temperature profile after $t=20~t_{\rm th}$ due to ${\rm p'p}$ interactions alone is also shown (dotted line).

When the disk protons in our simulation have reached their local virial temperature, we do not allow a further temperature increase. In a more realistic calculation, this limit on the temperature would come about through the advection of internal energy. Such a calculation requires a more detailed model of the ISAF and its sources and sinks of mass and energy, which is beyond the scope of the present calculation.

As an example we consider a warm disk with an optical depth (from the surface to the midplane) $\tau _{\rm 1/2}=1.5$, an initially isothermal temperature profile $T_{\rm e}=T_{\rm p}=500$ keV, illuminated by an ISAF with accretion rate $\dot M_{\rm I}=0.1 \dot M_{\rm Edd}$, where $\dot M_{\rm Edd} = 4 \pi c R_\star m_{\rm p} /\sigma_{\rm Th}$ is the Eddington accretion rate.

Figure 3: Local evaporation and condensation rates from a warm (100-500 keV) layer of optical depth $\tau _{\rm 1/2}=1.5$, formed by proton illumination of an initially cool disk. Solid line shows the condensation rate onto the layer from an ISAF accreting at $0.1~\dot{M}_{\rm E}$, as a function of distance from the hole. Dashed lines show the evaporation rates from the layer, for three values of its temperature. In the innermost region, the net effect is evaporation, feeding the ISAF, further out the ISAF condenses onto the layer. At lower ISAF accretion rates ( $0.01~\dot{M}_{\rm E}$, dotted) the evaporating region becomes wider (see text for details).

Figure 2 shows the temperature profiles for various values of the viscosity parameter $\alpha$ after $t=10~t_{\rm th}$(dashed-dotted lines) and $t=20~t_{\rm th}$ (solid lines), where $t_{\rm th}=1/(\alpha\Omega)$ is the thermal time scale of the disk. For $\alpha=0.5$ all disk protons heat up to the virial temperature, whereas for $\alpha=0.3$ and $\alpha =0.1$ only a part of the disk is subject to the instability. The extent of the unstable part is in excellent agreement with our estimate in Sect. 2.1. In test calculations in which only ${\rm p'p}$ heating and ${\rm ep}$energy coupling are included, the temperature increase of the disk protons is considerably smaller than with the viscous heating included. This confirms our conclusion from Sect. 2.2 that ${\rm p'p}$interactions alone can not heat the bulk of the internal protons considerably above the ambient electron temperature (cf. Eq. (20)).

Summarizing our investigation of the energy channels in a warm disk we conclude that viscous energy release in a hot plasma causes a runaway temperature increase in the unstable upper part of the disk atmosphere, until the protons have reached their local virial temperature. The process affects a significant part of the stratification if the electron temperature is above $\sim$100 keV.