The effect of rotation on the thermal instability of stratified galactic atmospheres - II. The formation of High Velocity Clouds

TL;DR

This study investigates how rotation influences thermal instabilities in galactic coronae and suggests that rotation can lead to the formation of High Velocity Clouds (HVCs) in the galaxy's inner regions, aligning with observations.

Contribution

It demonstrates that rotation is a key factor in the formation of HVCs via thermal instabilities, supported by stability analysis and hydrodynamical simulations.

Findings

Condensations occur in rotating regions of the corona.

Locations of condensations match observed HVC complexes.

Rotation supports HVC formation in inner galactic regions.

Abstract

Whether High Velocity Clouds (HVCs) can form by condensation of the hot ($T \sim 10^6 \, {\rm K}$) Galactic corona as a consequence of thermal instabilities has been controversial. Here we re-examine this problem and we suggest that rotation of the corona might be a missing key ingredient. We do this by studying the evolution of the models of rotating galactic coronae presented in Sormani et al. (2018) under the presence of cooling and thermal conduction. We combine a linear stability analysis with the results of local and global hydrodynamical simulations. We find that condensations are likely to occur in regions where the corona has substantial rotational support. Under reasonably general assumptions on the rotation profile of the corona, the locations where condensations are expected are in remarkable agreement with the observed location of the major non-magellanic HVCs complexes in our Galaxy (namely, at distances $< 15 \, \rm kpc$ from the Sun and within $30^\circ$ from the disc plane). We conclude that HVCs can form by thermal instabilities provided that (i) the corona is rotating substantially in the inner ($R < 50 \, \rm kpc$) parts, as suggested by current observational data and predicted by cosmological simulations of galaxy formation; (ii) close to the disc the corona is well-represented by a nearly-equilibrium stratified rotating structure (as opposed to a fast cooling flow). Our results also suggest that a better understanding of the disc-halo interface, including supernova feedback, is critical to understand the origin of HVCs.

Equations

$$\frac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho \mathbf{v}\right)=0;,$$

$$\left[\frac{\partial }{\partial t}+\left(\mathbf{v}\cdot \nabla \right)\right]\mathbf{v}=-\frac{\nabla P}{\rho }-\nabla \Phi ,$$

$$\frac{P}{\gamma -1}\left[\frac{\partial }{\partial t}+\left(\mathbf{v}\cdot \nabla \right)\right]\sigma =\nabla \cdot \left(f{\kappa }_{\mathrm{S}}{T}^{5/2}\nabla T\right)-{\left(\frac{\rho }{\mu m_{\mathrm{p}}}\right)}^2\Lambda \left(T\right),$$

$$t_{\mathrm{dyn}}=\sqrt{2}\frac{c_{\mathrm{s}}}{|{\mathbf{g}}_{\mathrm{eff}}|}$$

$$t_{\mathrm{rot}}=\frac{2\pi }{\Omega },$$

$$t_{\mathrm{cool}}=\frac{3k_{\mathrm{b}}T}{2n\Lambda (T)},$$

$$t_{\mathrm{cond}}=\frac{P{L}^2}{f{\kappa }_{\mathrm{S}}{T}^{7/2}}.$$

$${\tilde{\omega }}^3+i{\omega }_{\mathrm{d}}{\tilde{\omega }}^2-\left({\omega }_{\mathrm{BV}}^2+{\omega }_{\mathrm{rot}}^2\right)\tilde{\omega }-i{\omega }_{\mathrm{rot}}^2{\omega }_{\mathrm{d}}=0,$$

$${\omega }_{\mathrm{BV}}^2=-\frac{k_z^2}{k^2}\frac{\mathcal{D}\left(P_0\right)}{\gamma {\rho }_0}\mathcal{D}\left({\sigma }_0\right),$$

$${\omega }_{\mathrm{rot}}^2=-\frac{k_z^2}{k^2}\frac{1}{R^3}\mathcal{D}\left(R^4{\Omega }^2\right),$$

$${\omega }_{\mathrm{d}}={\omega }_{\mathrm{c}}+{\omega }_{\mathrm{th}}.$$

$$\mathcal{D}=\frac{k_R}{k_z}\frac{\partial }{\partial z}-\frac{\partial }{\partial R},$$

$${\omega }_{\mathrm{c}}=\frac{\gamma -1}{\gamma }\frac{k^{2}f{\kappa }_{\mathrm{S}}{T}_0^{7/2}}{P_0},$$

$${\omega }_{\mathrm{th}}=-\frac{\gamma -1}{\gamma }{\left(\frac{{\rho }_0}{\mu m_{\mathrm{p}}}\right)}^2\frac{\Lambda \left(T_0\right)}{P_0}\Delta \left(T_0\right),$$

$$\Delta \left(T\right)=2-\frac{\mathrm{d}\log \Lambda \left(T\right)}{\mathrm{d}\log T},$$

$$\tilde{\omega }=-i{\omega }_{\mathrm{th}}.$$

$$\tilde{\omega }=-i\frac{{\omega }_{\mathrm{rot}}^2}{{\omega }_{\mathrm{BV}}^2+{\omega }_{\mathrm{rot}}^2}{\omega }_{\mathrm{th}}.$$

$${\lambda }_{\mathrm{crit}}=\frac{2\pi \mu m_{\mathrm{p}}}{{\rho }_0}{\left[\frac{f{\kappa }_{\mathrm{S}}{T}_0^{7/2}}{\Lambda (T_0)\Delta (T_0)}\right]}^{1/2}.$$

$$\rho =-\frac{\partial P/\partial z}{\partial \Phi /\partial z},$$

$$v_{\varphi }^2=R\left(\frac{\partial \Phi }{\partial R}-\frac{\partial P/\partial R}{\partial P/\partial z}\frac{\partial \Phi }{\partial z}\right).$$

$$\frac{\delta \rho }{{\rho }_0}\sim \frac{\partial \delta P/\partial z}{\partial P_0/\partial z}\sim \frac{\delta P}{P_0}\frac{\partial \log \delta P/\partial z}{\partial \log P_0/\partial z}\gg \frac{\delta P}{P_0},$$

$$\frac{2v_{0\varphi }\delta v_{\varphi }}{R}=\frac{1}{{\rho }_0}\frac{\partial \delta P}{\partial R}-\frac{\delta \rho }{{\rho }_0^2}\frac{\partial P_0}{\partial R}.$$

$$\frac{1}{\delta \rho /{\rho }_0}\frac{\partial \delta P/\partial R}{\partial P_0/\partial R}\sim \frac{\partial P_0/\partial z}{\partial \delta P/\partial z}\frac{\partial \delta P/\partial R}{\partial P_0/\partial R}\sim 1.$$

$$\frac{\delta v_{\varphi }}{v_{0\varphi }}\sim \frac{c_{\mathrm{s}}^2}{v_{0\varphi }^2}\frac{\delta \rho }{{\rho }_0}.$$

$$\delta P/P_0\ll \delta \rho /{\rho }_0\sim (v_{0\varphi }^2/c_{\mathrm{s}}^2)\delta v_{\varphi }/v_{0\varphi }.$$

Full text

1 Introduction

It is now well established that the Milky Way (MW) is surrounded by a hot gaseous atmosphere, which is known as the Galactic corona (e.g. Putman et al., 2012; Tumlinson et al., 2017; Bregman et al., 2018, and references therein). The Galactic corona has a typical temperature of $T\simeq 2\times 10^{6}\rm\,K$ and is presumed to extend up to the virial radius of the MW dark matter halo, $r_{\rm vir}\simeq 250\,{\rm kpc}$, although due to the difficulty of observations most direct detections only probe the innermost $r\lesssim 50\,{\rm kpc}$ (e.g. Yoshino et al., 2009; Gupta et al., 2012; Miller & Bregman, 2013, 2015).

High Velocity Clouds (HVC) are a number of cold discrete gas clouds which are found throughout the inner portions of the Galactic corona (e.g. Putman et al., 2012; Westmeier, 2018).111In this paper we consider HVC to be a population distinct from that of the Intermediate Velocity Clouds (IVCs), which are found closer to the Galactic disc and constitute the so-called extraplanar HI gas (e.g. Fraternali, 2017, for a recent review). Most of the gas that constitutes IVCs has been probably pushed off the Galactic disc by stellar feedback associated with star formation (e.g. Fraternali & Binney, 2008; Marasco & Fraternali, 2011; Marasco et al., 2012), as opposed to HVCs that are probably of external origin as implied by their low metallicities (see below). They are primarily detected through the HI 21cm line and are characterised by radial velocities that are inconsistent with the circular rotation of the Galactic disc. Analogous populations of cold gas clouds are observed in M31 and other nearby spiral galaxies (e.g. Thilker et al., 2004; Westmeier et al., 2008; Lockman, 2017).

The HVCs have the following properties:

They are much colder than the surrounding hot medium which envelopes them; 2. 2.

They have pressures which are roughly consistent with pressure equilibrium with the surrounding hot medium (e.g. Stanimirović et al., 2002; Fox et al., 2005); 3. 3.

They typically have low metallicities, in the range of $0.1\mathchar 45\relax 0.3\,Z_{\odot}$ (e.g. van Woerden & Wakker, 2004); 4. 4.

Current observational constraint place most of the HVCs in proximity of the Galactic disc, at heights $z\lesssim 10\,{\rm kpc}$, distances from the Sun between $2<d<15\,{\rm kpc}$, and within $30^{\circ}$ from the disc plane as viewed from the Galactic centre (e.g. Putman et al., 2012). The exception is the Magellanic stream, which is clearly associated with the Magellanic clouds and is found at distance of $\sim 55\,{\rm kpc}$;222It should be noted that while all non-Magellanic HVCs are indeed found in proximity of the Galactic plane, a significant amount of low velocity cold gas in the CGM (even close to the viral radius) may be missed by current observations (e.g. Zheng et al., 2015). 5. 5.

They seem to be infalling with typical velocities of the order of $v\sim\text{tens of }\,{\rm km\,s^{-1}}$ (Fox et al., 2004; Putman et al., 2012).

The origin of HVCs is not fully understood. Several hypotheses have been put forward including:333Interestingly, most of the hypotheses still discussed today were already considered by Oort (1966) in the very early days (see his section 2).

1. (a)

They are made of gas ejected from the Galactic disc by supernova explosions through the galactic fountain mechanism (e.g. Bregman, 1980). This possibility is currently disfavoured by the low metallicities measured in the main HVC complexes (see item iii above), which strongly point to an extragalactic origin; 2. (b)

They are gas that has been stripped from satellites, and/or from gas that is coming onto the Milky Way along cold cosmological filamentary streams (e.g. Mayer et al., 2006; Olano, 2008; Fernández et al., 2012); 3. (c)

They are formed through condensation of the Galactic corona via thermal instabilities (e.g. Maller & Bullock, 2004); 4. (d)

They are formed through condensation of the Galactic corona, which however do not occur spontaneously but need to be ‘seeded’ by the blowout of a particularly powerful superbubble that reaches high altitudes above the disc (Fraternali et al., 2015).

The interesting possibility (c) was supported by early cosmological simulations (Kaufmann et al., 2006; Peek et al., 2008). This possibility was examined critically by Binney et al. (2009) which, by performing a linear stability analysis of a non-rotating, stratified, pressure-supported atmosphere, found that thermal instability is strongly suppressed by a combination of buoyancy and thermal conduction and therefore concluded that it is an unlikely mechanism for the formation of the HVCs. In support of this it has been discovered that the cold clumps seen in the earlier simulations were numerical artefacts (Hobbs et al., 2013). However, newer state of the art cosmological simulations at high resolution show the widespread presence of neutral hydrogen cloudlets in the circumgalactic medium (e.g. van de Voort et al., 2019).

One of the critical assumptions in the work by Binney et al. (2009), as mentioned in their conclusions, is that the corona does not rotate. However, there is robust evidence that the Galactic corona rotates. For example, by measuring the Doppler shifts of the OVII absorption lines, Hodges-Kluck et al. (2016) found that the innermost parts of the Galactic corona rotate at $\sim 180\,{\rm km\,s^{-1}}$ and contain an amount of angular momentum comparable to that in the stellar disc. Marinacci et al. (2011) used wind-tunnel simulations of cold clouds travelling through a hot medium to show that if the inner corona of a Milky Way-like galaxy did not rotate, then it would be rapidly accelerated by interaction with the extra-planar HI gas. They concluded that the corona is expected to rotate and to lag in the inner regions by $\sim 80\,{\rm km\,s^{-1}}$ with respect to the cold disc. By analysing a sample of simulated galaxies from the Eagle cosmological simulations, Oppenheimer (2018) concluded that galactic coronae have significant rotational support in the inner parts (at radii $\ll$ than the virial radius).

Nipoti (2010) extended the work of Binney et al. (2009) and investigated whether rotation might change the stability properties of the Galactic corona. Although this author formally derived the relevant dispersion relations, his analysis was inconclusive because of a lack of explicit models of rotating Galactic coronae to which he could apply his analysis to. Moreover, it is now clear that the results of the linear theory are of limited applicability to quantify the condensation of non-linear structures (see for example McCourt et al., 2012; Sharma et al., 2012b).

This motivates a re-examination of the question of whether condensation can form via thermal instabilities from a rotating Galactic corona. In this paper we use the models of Sormani et al. (2018) to study this question. In the absence of cooling and thermal conduction these models are stable equilibria which rotate in the inner parts. We study how they evolve when cooling and thermal conduction are added to the equations of motion. This study follows a companion paper (Sobacchi & Sormani, 2019, hereafter Paper I), in which we investigated the effects of rotation on the thermal instability using an idealised plane-parallel setup, and we found that rotation enhances condensation by increasing the threshold value of the ratio between the cooling and the dynamical time below which thermal condensations can be formed.

The paper is structured as follows. In Sect. 2 we describe our physical model by writing down the equations that govern the evolution of our system and the important timescales of the problem. In Sect. 3 we apply the linear stability analysis of Nipoti (2010) and the results of Paper I to our models. To confirm our findings and study the non-linear behaviour, in Sect. 4 we perform 2.5D axisymmetric simulations which include cooling and thermal conduction using the models of Sormani et al. (2018) as initial conditions. We discuss our results in Sect. 5, and we sum up in Sect. 6.

2 Physical model

The Galactic corona is expected to be approximately in equilibrium, stratified in the gravitational potential of the Galaxy, and rotating in the inner parts (see introduction). The question is whether condensations due to the thermal instability can form in such a system. The question is not whether the system is formally unstable according to a linear analysis, i.e. whether there exist at least one unstable mode that grows exponentially in time. Linear perturbation analysis shows that the system is always formally unstable but it has been previously argued that, despite the presence of unstable modes, the system will not develop condensations in the non-rotating case because buoyancy strongly suppresses their formation during the non-linear evolution of the system (Malagoli et al., 1987; Binney et al., 2009). Numerical simulations have later confirmed that indeed buoyancy can suppress the production of condensations and that non linear effects play a crucial role in the development of the instability (e.g. McCourt et al., 2012). Nipoti (2010) suggested that conclusions similar to those of Binney et al. (2009) (i.e. the corona should not develop condensations) also hold in the rotating case. However (i) the presence of additional unstable modes makes the final result more uncertain; (ii) we have found in Paper I that condensations are more likely to form in the rotating case; (iii) the presence of rotation changes the global structure of the corona, which was not taken into account in the discussion of Nipoti (2010); (iv) Nipoti’s analysis was limited to the linear regime, but we have demonstrated in Paper I that non-linear effects play a key role also in the rotating case.

It is therefore important to (i) revisit the analysis of whether condensation can form in the Galactic corona using realistic rotating models and (ii) follow the evolution of the system in the non-linear regime. We address these items in Sect. 3 and 4. In the remainder of this section, we lay out the fundamental equations and the important timescales that govern the evolution of the system.

2.1 Fundamental equations

We assume that the system is governed by the equations of fluid dynamics:

[TABLE]

where $\rho$ is the density, ${\bf v}$ is the velocity, $P$ is the pressure, $\Phi$ is the external gravitational potential, $T$ is the temperature and $\sigma~{}=~{}\log\left(P\rho^{-\gamma}\right)$ is the entropy. The energy equation has been modified with respect to the ideal case by the addition of cooling and thermal conduction to the right hand side of Eq. (3), where $\Lambda$ is the cooling function, $\kappa_{\rm S}=(5/2)k_{\rm b}^{7/2}/(\pi m_{\rm e}^{1/2}{\rm e}^{4})\simeq 1.53\times 10^{-5}\;{\rm erg}\;{\rm s}^{-1}\;{\rm cm}^{-1}\;{\rm K}^{-7/2}$ is the coefficient of thermal conductivity (Spitzer, 1962), $m_{\rm e}$ is the electron mass, ${\rm e}$ is the electron charge, $k_{\rm b}$ the Boltzmann constant, and $f\leq 1$ is the fraction by which thermal conduction is suppressed by a tangled magnetic field (e.g. Binney & Cowie, 1981).444Of course, the suppression factor $f$ is a crude effective prescription, while a self-consistent treatment of the problem should include the effect of anisotropic thermal conduction. However, this requires to include magnetic fields in the analysis, which lies beyond the scope of the paper. We take a fiducial $f=0.01$ as indicated by studies of galaxy clusters (e.g. Nipoti & Binney, 2004) and to provide a direct comparison with the work by Binney et al. (2009) which adopted the same value. In this paper we adopt an adiabatic index $\gamma=5/3$, the value for monoatomic ideal gases, and an ideal equation of state $P=nkT$, where $n=\rho/(\mu m_{\rm p})$ is the number density, $\mu=0.58$ is the mean molecular weight and $m_{\rm p}$ is the proton mass.

Throughout this paper, we adopt the collisional ionisation equilibrium cooling function from Sutherland & Dopita (1993) for ${\rm[Fe/H]}=-0.5$, which accounts for the sub-solar metallicities found in the HVCs and which are believed to characterise the Galactic corona (e.g. van Woerden & Wakker, 2004; Bregman et al., 2018).

We use standard cylindrical coordinates $\left(R,\phi,z\right)$, and we use $r$ to denote the spherical radius.

2.2 Important timescales

The evolution of the corona is determined by the interplay of the following timescales (see also Paper I). The first is the dynamical timescale:

[TABLE]

where $c_{\rm s}=\sqrt{P/\rho}$ is the sound speed (apart from an unimportant factor $\gamma=5/3$) and $\mathbf{g}_{\rm eff}=-\nabla\Phi+\Omega^{2}R\;\hat{\mathbf{e}}_{R}$ is the effective gravity. In a nearly-equilibrium corona we have $\mathbf{g}_{\rm eff}\simeq\nabla{P}/\rho\sim P/(r\rho)\sim c_{\rm s}^{2}/r$ so $t_{\rm dyn}\sim\sqrt{2r/|\mathbf{g}_{\rm eff}|}$. Therefore, in the case of a non-rotating corona, $t_{\rm dyn}$ is essentially the free-fall time. In the case of a rotating corona, $t_{\rm dyn}$ it is a useful generalisation of the free-fall time that is based only on local quantities and therefore allows for a direct comparison with the local analysis of Paper I.

The second is the rotational timescale

[TABLE]

where $\Omega=v_{\rm\phi}/R$ is the local angular velocity. The third is the cooling timescale

[TABLE]

where $k_{\rm b}$ is the Boltzmann constant. This is ratio between the internal energy per unit volume, $(3/2)nk_{\rm b}T$, and the amount of energy lost due to cooling per unit time and volume, $n^{2}\Lambda(T)$. It therefore represents the time a blob of gas would take to emit all its energy if the cooling continued at the current rate.

Finally, we have the thermal conduction timescale:

[TABLE]

This is the time it takes to thermal conduction to be effective over a region of size $L$.

3 Thermal instability in rotating coronae

3.1 Linear analysis

The linear evolution of axisymmetric perturbations of a rotating atmosphere described by equations (1)-(3) has been studied by Nipoti (2010). Here we follow his treatment and refer to his paper for a detailed derivation of the dispersion relation.

The unperturbed corona is assumed to be axisymmetric and in a steady state over the time-scales of interest, even in the presence of radiative cooling and thermal conduction. Thus, the unperturbed state can be described by the time-independent axisymmetric pressure $P_{0}$, density $\rho_{0}$, temperature $T_{0}$ and velocity field $\mathbf{v}_{0}=(v_{0R},v_{0z},v_{0\phi})$, which satisfy Eqs. (1)-(3) with vanishing partial derivatives with respect to $t$.

We linearise Eqs. (1)-(3) with Eulerian perturbations of the form $F_{0}+\delta F\exp\left(ik_{R}R+ik_{z}z-i\omega t\right)$, where $F_{0}$ is the equilibrium value and $\left|\delta F\right|\ll\left|F_{0}\right|$. We work in the short-wavelength (WKB) limit, i.e. we assume that $kH\gg 1$ where $H=c_{\rm s}^{2}/|\mathbf{g}_{\rm eff}|$ is the typical lengthscale of the system. We also use the Boussinesq approximation, which allows to filter out sound wave modes from the dispersion relation so that only the relevant low-frequency modes with $\omega\ll|\mathbf{k}|c_{\rm s}$ remain (see for example Nipoti 2010; Balbus & Potter 2016).

Under these assumptions one finds the following dispersion relation:

[TABLE]

where $\tilde{\omega}=\omega-k_{R}v_{0R}-k_{z}v_{0z}$, so that a mode is exponentially growing and thus formally unstable if and only if $\operatorname{Im}(\tilde{\omega})>0$. Equation (8) is formally identical to the dispersion relation obtained in the local analysis of Paper I and therefore confirms that we can use the results obtained in that paper for the problem considered here. Indeed, the problem studied in Paper I can be considered a simplified version of the problem studied by Nipoti (2010) that makes the following additional assumptions (not present in Nipoti’s paper): (i) no differential rotation; (ii) the surfaces of constant pressure and density coincide.

In the above we have defined the following characteristic frequencies:

[TABLE]

where

[TABLE]

and $\Omega=v_{0\phi}/R$.

The quantity $\omega_{\rm BV}$ is the characteristic frequency of buoyant oscillations, while $\omega_{\rm rot}$ is the characteristic frequency related to the rotation of the system. Note that $\omega_{\rm BV}$ is related to the dynamical time (4) by $\omega_{\rm BV}\sim t_{\rm dyn}^{-1}$ and similarly $\omega_{\rm rot}\sim t_{\rm rot}^{-1}$ (the exact relations depend on the direction of ${\bf k}$). The combination $\omega_{\rm BV}^{2}+\omega_{\rm rot}^{2}$ is positive when the equilibrium is convectively stable, as we assume in this paper,555The condition that $\omega_{\rm BV}^{2}+\omega_{\rm rot}^{2}>0$ for all the possible $\left(k_{R},k_{z}\right)$ gives the classical Solberg-Høiland stability criteria (e.g. Tassoul, 2000). while the individual terms $\omega_{\rm BV}^{2}$ and $\omega_{\rm rot}^{2}$ are not necessarily positive. The quantity $\omega_{\rm d}$ is the characteristic frequency of dissipative processes, and is composed of the sum of the cooling ($\omega_{\rm th}$) and thermal conduction ($\omega_{\rm c}$) frequencies. We note $\omega_{\rm th}$ is negative if the plasma is thermally unstable in the classical sense of Field (1965) (as we assume throughout this paper) while $\omega_{\rm c}$ is always positive, so the effect of thermal conduction always goes in the direction of stabilising the system (see below).

3.1.1 Analysis of the dispersion relation

Since the dispersion relation (8) is formally identical to the one derived in Paper I, we will not repeat the analysis here but simply state some points that are relevant for this paper.

• •

In the fast cooling regime ($t_{\rm cool}\ll$ all other timescales so that we can set $\omega_{\rm BV}=\omega_{\rm rot}=\omega_{\rm c}=0$), we recover the classical Field (1965) instability, and the dispersion relation reduces to

[TABLE]

Hence, condensations are expected in this regime.

• •

As $t_{\rm cool}$ (and in particular the ratio $t_{\rm cool}/t_{\rm dyn}$) is increased, the number of unstable modes decreases (i.e. the modes corresponding to more and more $\mathbf{k}$ turn into overstable modes). For example, in the limit $t_{\rm cool}/t_{\rm dyn}\gg 1$ in the non-rotating case only the modes with ${\bf k}$ nearly parallel to $\nabla P_{0}$ remain unstable ((in the non-rotating case this means that ${\bf k}$ is radial, so that $\omega_{\rm BV}=0$ for these modes, see also Binney et al. 2009). This suggests that the instability becomes weaker and weaker as $t_{\rm cool}/t_{\rm dyn}$ is increased and leads to the expectation that there may exist a threshold value of $t_{\rm cool}/t_{\rm dyn}$ above which condensation cannot occur.

• •

A new unstable mode appears due to the presence of rotation. This was already pointed out by Nipoti (2010). In the slow cooling regime ($t_{\rm cool}/t_{\rm dyn}\gg 1$), neglecting the effect of thermal conduction ($f=0$), this mode grows with a frequency of (see Paper I):

[TABLE]

Since $\omega_{\rm BV}/\omega_{\rm rot}\sim t_{\rm rot}/t_{\rm dyn}$, this mode becomes important when $t_{\rm dyn}/t_{\rm rot}$ is large (note that in this limit $\tilde{\omega}\sim-i\omega_{\rm th}$), i.e. when the corona has substantial rotational support. We discuss the physical interpretation of this new mode in Appendix A.

• •

When thermal conduction is added to the system, the dispersion relation is formally identical with the substitution $\omega_{\rm th}\to\omega_{\rm d}$. Thus the effect of thermal conduction is to stabilise all modes with wavelengths shorter than:

[TABLE]

This condition is obtained by requiring that $\omega_{\rm d}=0$.

3.2 When do we expect condensations?

From the linear analysis it follows that (in the absence of thermal conduction) condensation will certainly occur in the fast cooling regime $t_{\rm cool}/t_{\rm dyn}\ll 1$ (because we recover the classical Field 1965 instability, see previous section). Considerations that use the linear analysis as a guide to infer what happens in the non-linear regime lead to the expectation that condensations should not occur in the opposite regime $t_{\rm cool}/t_{\rm dyn}\gg 1$ (Malagoli et al., 1987; Binney et al., 2009, see also Sect 3.1.2 and 4.3 of Paper I).666Intuitively, the physical reason why condensations do not form if $t_{\rm cool}/t_{\rm dyn}\gg 1$ essentially boils down to the fact that an overdense clump will sink and mix with the hot gas before having time to cool significantly. Hence there should exist a threshold value of $t_{\rm cool}/t_{\rm dyn}$ above which condensation cannot occur. While the linear analysis constitutes a useful guide to what we can expect in the non-linear regime, the physics involved is inherently non linear, so the exact value of $t_{\rm cool}/t_{\rm dyn}$ above which condensations cannot occur cannot be inferred from the linear analysis but must be determined using numerical simulations (McCourt et al. 2012, Paper I).

In the absence of rotation ($t_{\rm dyn}/t_{\rm rot}=0$), McCourt et al. (2012) and Sharma et al. (2012b) suggested that condensation will occur when $t_{\rm cool}/t_{\rm dyn}\lesssim 10$, the threshold being quite sharp.777McCourt et al. (2012) found a somewhat more stringent criterion, $t_{\rm cool}/t_{\rm dyn}\lesssim 1$, than Sharma et al. (2012b), who found $t_{\rm cool}/t_{\rm dyn}\lesssim 10$. Previous claims attributed the difference to the different geometry (plane parallel vs spherical), but Choudhury & Sharma (2016) has shown that the difference is due to the fact that McCourt et al. (2012) looked at the condensation at the location of $\min(t_{\rm cool}/t_{\rm ff})$, while Sharma et al. (2012b) considered condensations anywhere in the box, which happens further inside the location of the initial $\min(t_{\rm cool}/t_{\rm ff})$. Also note that the threshold becomes even less stringent if large density fluctuations are assumed to be present in the initial conditions (e.g. Choudhury et al., 2019).

In Paper I we have found that the presence of rotation can increase the threshold value of $t_{\rm cool}/t_{\rm dyn}$ by a factor up to $f_{\rm rot}\sim 10$. The effect becomes largest when (i) rotation is dynamically important (namely, in the regime $\Omega c_{\rm s}/g_{\rm eff}\gtrsim 1$, or equivalently $t_{\rm dyn}/t_{\rm rot}\gtrsim 0.2$); (ii) $\theta$ is small, where $\theta$ is defined as the angle between ${\bf g}_{\rm eff}$ and the direction perpendicular to ${\boldsymbol{\Omega}}$ (see fig. 1 of Paper I), which occurs preferentially near the Galactic plane. The strength of this effect as a function of the parameters is quantified in fig. 4 of Paper I.

While item (i) could be interpreted as a result of the new unstable mode that appears as a result of rotation (Eq. 17), the fact that condensation depends on $\theta$ (item ii) must be an entirely non linear effect because the dispersion relation (8) does not depend on $\theta$. In Paper I, we have interpreted this as a result of the fact that the motion in the direction perpendicular to ${\boldsymbol{\Omega}}$ is hampered by the conservation of the angular momentum, so contrary to the non-rotating case an overdense region will not sink and mix if ${\bf g}_{\rm eff}$ and ${\boldsymbol{\Omega}}$ are nearly perpendicular (see Sect. 4.3 in Paper I).

In conclusion, we may expect the thermal instability in a rotating corona in the absence of thermal conduction to occur when $t_{\rm cool}/t_{\rm dyn}\lesssim 10f_{\rm rot}$ where $f_{\rm rot}\geq 1$ is the factor by which the threshold is increased by the presence of rotation, which depends on $t_{\rm dyn}/t_{\rm rot}$ and $\theta$ (see fig. 1 of Paper I for the definition of $\theta$). The factor $f_{\rm rot}$ can be read off the right panel in fig. 4 of Paper I. The presence of thermal conduction has a stabilising effect up to a certain maximum spatial scale given by Eq. (18).

3.3 Application to the equilibrium models of Sormani et al. (2018)

We now apply the considerations of the previous sections to the models of Sormani et al. (2018), to see where whether we may expect them to produce condensations.

Fig. 1 shows the fiducial model of Sormani et al. (2018). According to Sect. 3.2, the parameters that control the thermal stability of the system are $t_{\rm cool}/t_{\rm dyn}$, $t_{\rm dyn}/t_{\rm rot}$ and $\theta$ (the angle between the effective gravity and the direction perpendicular to the rotation axis, see fig. 1 of Paper I). These are shown in Fig. 2. The regions where the condensations are most likely to happen are those where $t_{\rm cool}/t_{\rm dyn}$ becomes small, $t_{\rm dyn}/t_{\rm rot}$ becomes large and $\theta$ becomes small. Hence these are the regions we should look for. Excluding the region $z\lesssim 5\,{\rm kpc}$, which is the transition zone where the corona is interacting with the Galactic disc (of which our model cannot be taken as a faithful representation), we see that the most favourable area is the region $5\lesssim z\lesssim 10\,{\rm kpc}$, $10\lesssim R\lesssim 30\,{\rm kpc}$. Interestingly, this roughly coincides with the locations where HVCs are found in the MW: as mentioned in the introduction, most of the HVCs are found in proximity of the Galactic disc, at heights $z\lesssim 10\,{\rm kpc}$, distances from the Sun between $2<d<15\,{\rm kpc}$ and within $30^{\circ}$ from the disc plane as viewed from the Galactic centre. Note that, in contrast, in a non-rotating spherical model the most favourable region for the formation of the HVCs (for example, model 1 or 4 of Sormani et al. 2018 or the models discussed in Binney et al. 2009) is a spherical region within few kpc of the Galactic centre, in disagreement with observations.

We see from Fig. 1 that the presence of rotation is expected to help the production of condensation in the favourable area mentioned above by lowering the threshold $t_{\rm cool}/t_{\rm dyn}$ by a factor of $f_{\rm rot}\sim$ a few according to the right panel of fig. 4 of Paper I. This is a direct effect by which condensations are enhanced by rotation.

There is also an indirect effect that makes condensation more likely in a rotating corona compared to a non-rotating corona in the same gravitational potential. Fig. 3 shows the ratio $|\mathbf{g}_{\rm eff}|/|\mathbf{g}|$ for the same model as in Fig. 1. From this figure we see that compared to a non-rotating corona in the same external gravitational potential (i.e. with the same value of $\mathbf{g}$) and with the same density and temperature, the presence of rotation increases the value of $t_{\rm dyn}$ (Eq. 4) by roughly a factor of 2 in the relevant regions, and therefore reduces the ratio $t_{\rm cool}/t_{\rm dyn}$ (which is the parameter that matters for the production of condensations) by the same factor. This makes condensations more likely.

The fact that the formation of HVCs is enhanced in a region in proximity of the Galactic plane ($z\lesssim 10\,{\rm kpc}$) at radial distances of few tens of $\,{\rm kpc}$ is not a particular feature of our fiducial model, but is a quite a general feature of any plausible rotating model. Any model which rotates at a substantial fraction of the circular velocity of the disc will generate (i) a similar enhancement in density and therefore of $t_{\rm dyn}/t_{\rm cool}$ in the relevant regions, and (ii) a relative enhancement of $T$ in a cone above the Galactic centre (i.e. whose axis is $R=0$) compared to a spherical model, which strongly suppresses the formation of HVCs in this region. Indeed, other models such as for example the models of Pezzulli et al. (2017) and all the rotating models in Sormani et al. (2018) have similar characteristics.

Fig. 4 shows the critical wavelength below which modes are stabilised by thermal conduction for the same model, assuming that the suppression factor of thermal conduction due to a tangled magnetic field is $f=0.01$ (the same value used in Binney et al. 2009). This shows that thermal conduction is likely to be very effective in further suppressing the formation of condensations in a cone above the Galactic centre. In the relevant region where HVCs may be expected to form ($5\lesssim z\lesssim 10\,{\rm kpc}$, $10\lesssim R\lesssim 30\,{\rm kpc}$), the situation is less clear. The figure suggests that suppression due to thermal conduction may be substantial but, given that the values of $\lambda_{\rm crit}$ are too large for the WKB approximation which underlies the linear analysis to be valid, it is unclear whether this will actually stop condensations from forming.

In conclusion, the analysis in this section suggests that condensation might form in a rotating corona due to a combination of the effects discussed in Sect. 3.2. It is unclear whether thermal conduction is able to suppress the formation of condensations. To determine this it is necessary to run some numerical experiments.

4 Numerical experiments

To better understand the non-linear behaviour of the thermal instability, we run axisymmetric 2.5D hydrodynamical simulations in the presence of cooling and thermal conduction using the models of Sormani et al. (2018) as initial conditions. We do not impose any initial perturbation. In these simulations, the evolution of the system is governed by Eqs. (1)-(3). Note that our initial conditions are not equilibria with respect to the energy equation (Eq. 3). These simulations are 2.5D in the sense that we assume axisymmetry, so the grid in the $(R,z)$ plane is two-dimensional, but the velocity vector has three components $(v_{R},v_{z},v_{\phi})$.

4.1 Numerical setup

We have used the public Eulerian code PLUTO version 4.3 (Mignone et al., 2007). This is a free software for the numerical solution of systems of conservation laws targeting high Mach number flows in astrophysical fluid dynamics. We use a two-dimensional static cartesian grid in the region $R\times z=[0,200]\,{\rm kpc}\times[0,200]\,{\rm kpc}$. The grid is uniformly spaced in both $R$ and $z$ with $800\times 800$ points, which corresponds to a spatial resolution of $0.25\,{\rm kpc}$. We use the following parameters: rk2 time-stepping, no dimensional splitting, roe Riemann solver and the mc flux limiter. Boundary conditions are reflective on the plane $z=0$ and on the outer boundaries at $z=200\,{\rm kpc}$ and $R=200\,{\rm kpc}$ and axisymmetric (i.e. $v_{\phi}$ is reversed on the ghost cells just beyond the axis $R=0$) on the axis $R=0$. Thermal conduction is treated with the Explicit scheme.

We adopt the collisional ionisation equilibrium cooling function from Sutherland & Dopita (1993) for ${\rm[Fe/H]}=-0.5$. This accounts for the sub-solar metallicities found in the HVCs, which are believed to characterise the Galactic corona (e.g. van Woerden & Wakker, 2004; Bregman et al., 2018). The cooling is switched off for $T\leq 10^{4}\,\rm K$, and for $z\leq z_{0}=5\,{\rm kpc}$ since in this region the corona is likely interacting with the Galactic disc and our model cannot be taken as a faithful representation of the real physical conditions expected in this region. Indeed the cooling timescale of our model close to the plane becomes very short and the presence of cooling would lead to a cooling catastrophe. To make the transition smoother, for $z>z_{0}$ we add a position-dependent and time-independent heating function $H(R,z)=\exp[-(z-z_{0})/z_{0}]\,H_{0}(R,z)$ where the function $H_{0}(R,z)$ is such that at $t=0$ the heating exactly balances the cooling. We have experimented with different choices for the way in which the cooling is switched off near to the plane, and all lead to similar results (see discussion in Sect. 5).

Finally, we have verified by running simulations in which both cooling and thermal conduction are switched off that all the models of Sormani et al. (2018) are indeed stable equilibrium models. In this setup, the initial state does not change noticeably during the course of the simulations, which lasts for several Gyr.

4.2 Results

Fig. 5 shows the evolution of our fiducial model under the presence of cooling but in the absence of thermal conduction. The model slowly contracts as the consequence of cooling until at $t\simeq 1\,{\rm Gyr}$ a condensation forms at $z\sim 6\,{\rm kpc}$, $R\sim 16\,{\rm kpc}$. This leads to a big drop of condensed gas which falls onto the Galactic disc. This is quickly followed a sequence of similar drops which fall on the Galactic disc in a way reminiscent of a leaky faucet. The location where the condensation forms is in agreement with the expectations from the linear analysis (Sect. 3.3) and is roughly consistent with the observationally determined locations of HVCs. This suggests that indeed thermal instability can produce condensation resembling HVCs. We discuss in more detail in Section 5 the dependence of these results on our assumptions, and in particular on our parametrisation of the heating and cooling functions at the disc-halo interface.

Fig. 6 shows the evolution of the same model under the presence of both cooling and thermal conduction. The latter has the effect of slightly delaying the formation of the condensation and to slightly smooth out its density and temperature contrast (for a detailed study of the dynamics of cold clouds in the presence of thermal conduction, see e.g. Armillotta et al., 2017). However, it seems that thermal conduction is not able to halt the formation of the condensation. It must be mentioned however that this relies on a crude model for the suppression of thermal conduction due to the presence of a tangled magnetic fields and on the value of $f=0.01$. A proper treatment of the problem would require magneto-hydrodynamics rather than hydrodynamic simulations. We have tested that in a simulation with $f=1$ thermal conduction is so strong that not only suppresses the condensation, but also changes the entire structure of the corona making it almost isothermal within a few hundred of Myr.

We have tested what happens to all the other models presented in Sormani et al. (2018). All the rotating models are unstable and produce condensations. The instability is much stronger in the isothermal rotating models (models 2 and 3) on account of their higher densities in proximity of the Galactic plane, while it is similar to our fiducial model in the adiabatic rotating model 5. The non-rotating spherical adiabatic model 4 is stable, and produces no condensations, as expected since its cooling times are very long everywhere ($t_{\rm cool}/t_{\rm dyn}\gtrsim 100$ - see also the similar model in the central column of fig. 5 of Binney et al., 2009). In the non-rotating spherical isothermal model 1 cooling times are short within a few kpc of the Galactic centre compared to the duration of the simulations (similarly to the model in the central column of fig. 4 of Binney et al., 2009), thus the models cools very quickly near the centre, producing condensations. However, these condensation are all within a roughly spherical regions of $r\lesssim 8\,{\rm kpc}$ from the Galactic centre, so they are simply a central accumulation of cold gas and do not agree with the observationally determined locations of HVCs.

5 Discussion

The results presented in Sections 3.3 and 4.2 suggest that condensation may form in a rotating stratified corona. However, there are some aspects that deserve further discussion.

The first aspect is the turn off of the cooling near the Galactic plane for $z\leq z_{0}$. We have experimented with different values of $z_{0}$ and with different ways of turning off the cooling close to the Galactic plane. The latter always produce similar results for comparable values of $z_{0}$. For example, we can (i) turn off the cooling gradually following a function $F(z)=\tanh[(z-z_{0})/z_{1}]$, where $z_{0}$ and $z_{1}$ are two parameters, instead of using a heating function as in Sect. 4; (ii) assume that heating balances the cooling on large scales (the so called “global thermal balance”) if $z<z_{0}$, and turn off the cooling at $z<z_{1}$. This produces results that are very similar to those described in Sect. 4.2 for comparable values of $z_{0}$ and, say, $z_{1}=2\,{\rm kpc}$ (results are not very dependent on this value, provided that $z_{1}\ll z_{0}$). The effects of changing the value of $z_{0}$ are more tricky. If $z_{0}$ is increased to, say, $z_{0}=7.5\,{\rm kpc}$, then condensations will form at higher $z$ and at later times than shown in Fig. 5. If we raise $z_{0}$ to even higher values, say $z_{0}=20\,{\rm kpc}$, condensations will stop altogether forming because the cooling time becomes comparable to the Hubble time at those heights. The trickiest point is that, for example, while condensations form at $z\gtrsim 7.5\,{\rm kpc}$ if $z_{0}=7.5\,{\rm kpc}$, they will form at lower heights if $z_{0}$ is lower. This is because they form first where the cooling time is shortest, and the upper parts of the atmosphere then collapse before having time to produce condensations. If $z_{0}\sim 0$ our setup reduces to a pure cooling flow and no condensations are observed outside of the disc: the gas is falling too fast and collapses onto the disc before having time to cool significantly. We may argue that near the Galactic plane some form of heating should be present, and that excessive cooling in the corona at $z\lesssim 5\,{\rm kpc}$ must be somehow prevented, and thus condensations should form at those heights. We conclude that condensations resembling HVCs may form only if the corona is well represented by a nearly-equilibrium stratified rotating structure near the plane, while they cannot form if it is well-represented by a fast cooling flow. The detailed study of this question requires better understanding of the conditions in the transition zone between the disc and the Galactic atmosphere, which is a complicate and poorly understood problem at the moment.

The second aspect is that it is difficult to disentangle which of the two effects discussed in Sect. 3.3 is actually responsible for the production of condensations. The presence of rotation enhances the condensations (i) indirectly, by affecting the global density and temperature structure of the corona and thus changing the global profile of $t_{\rm cool}/t_{\rm dyn}$, which is the parameter that controls the formation of condensations; (ii) directly, by increasing the threshold value of $t_{\rm cool}/t_{\rm dyn}$ below which condensations can form. However, the locations where these two effects become relevant are almost identical in the models. We conclude that it is probably a combination of both effects that is at work here and makes the corona thermally unstable. Note, however, that effect (ii) relies on the results of Paper I which does not include the effect of differential rotation, whose impact on condensations is yet to be explored in the non-linear regime.

It should also be mentioned that throughout this paper we assume axisymmetry, but we cannot exclude a priori that the formation of condensation is suppressed when this assumption is relaxed (Nipoti, 2010; Nipoti & Posti, 2014). To check this would require 3D simulations, which are however much more computationally expensive than the ones presented in this paper.

The pressure within the condensation in Fig. 5 and 6 is slightly lower than in the surrounding hot medium. This is in agreement with the fact that the measured pressures of HVCs are slightly below what is expected from pressure equilibrium with the hot medium, where the pressure of the hot medium is based on its density and temperature determinations (see for example fig. 1 and remarks in Sect. 4.2 in Sormani et al. 2018).

We have crudely modelled the thermal conduction by taking into account its suppression due to a tangled magnetic field using a factor $f$. However, we note that this problem should be properly studied in the context of magnetohydrodynamics. This lies outside of the scope of the current paper, but is a worthwhile direction for future investigations. Indeed, magnetic field are known to have non trivial effects on the dynamics of cold gas clouds travelling through the hot coronal medium (e.g. Grønnow et al., 2017, 2018).

Finally, since density determinations of the Galactic corona have large uncertainties, we have experimented by running simulations of our fiducial model in which $T_{0}(R,z)$ and $v_{0\phi}(R,z)$ are unchanged but the density is scaled down by a factor of 2. This lowers the cooling time and thus has the effect of delaying the formation of condensation by approximately the same factor. Hence condensations may be suppressed only if the density is lower by a factor of $\sim 10$, which would lead to condensation forming on times longer than the Hubble time. However, based on current constraints on the density of the Galactic corona this seems unlikely albeit possible. In relation to this it is worth mentioning that, while in this paper we have assumed that the densities and pressure of the corona are a ‘given background’, they may be set by the requirement that the gas the value of $t_{\rm cool}/t_{\rm dyn}$ is close to the threshold value for the production of condensations (e.g. Sharma et al., 2012a; Voit et al., 2015, 2019). The consequences of this requirement are yet to be explored in the context of rotating coronae.

6 Conclusions

We have examined the effects of rotation on the formation of HVCs in the Galactic corona. We have used both a linear stability analysis, and the results of local and global hydrodynamical simulations applied to the models of Sormani et al. (2018).

Rotation enhances the formation of condensation through two effects:

An indirect effect. Rotation changes the global density and temperature structure of the corona, thus changing the global profile of $t_{\rm cool}/t_{\rm dyn}$, the parameter that controls the formation of condensations (see Sect. 3.3 and Figs. 1, 2 and 3). 2. 2.

A direct effect. The presence of dynamically significant rotation ($\Omega c_{\rm s}/g_{\rm eff}\gtrsim 1$ or equivalently $t_{\rm dyn}/t_{\rm rot}\gtrsim 0.2$) enhances the condensation by directly changing the threshold value of $t_{\rm cool}/t_{\rm dyn}$ below which condensations can form. The effect is largest in the proximity of the Galactic plane, where the effective gravity is nearly perpendicular to the rotation axis (see Sect. 3.2 and 3.3 and fig. 4 in Paper I). This is the effect we discovered in Paper I and assumes that we can neglect the impact of differential rotation, whose importance for the condensation remains to be tested in the non-linear regime.

It is difficult to disentangle the two effects because the locations where they are important are essentially the same. In any reasonable rotating model of the Galactic corona this makes the formation of condensation more favourable precisely in the location where HVCs are found in the observations (see Sect. 3.3). Although this may be a coincidence and it may be that current observations are missing some neutral gas further out, given that several independent lines of research point to a corona that is significantly rotating in the inner parts (see Introduction) at the moment it seems reasonable to interpret this accordance as evidence that HVCs are formed by the thermal instability.

Thermal conduction does not suppress the formation of condensations, assuming that it is lower by a factor $\lesssim 10^{-2}$ with respect to the Spitzer’s standard value due to the presence of a tangled magnetic field. However this treatment is rather crude, and a proper study requires explicitly taking explicitly into account the often subtle effects of the magnetic fields which is outside the scope of this work.

Condensations at $z\gtrsim 5\,{\rm kpc}$ are absent in our simulations if we completely turn off any form of heating near the Galactic plane. In this case a fast cooling flow is established near the disc and the gas does not have time to cool significantly before falling onto the disc. Hence, the formation of HVCs requires some amount of heating (e.g. supernova feedback) that prevents excessive cooling near the disc. In this sense, our results indicate that the disc-corona interaction is important, similarly to what has been suggested by Fraternali et al. (2015). The main difference between Fraternali et al. (2015) and the scenario presented here is that we do not require the condensations to be ‘seeded’ by a powerful superbubble that reaches the corona, but can happen spontaneously. We suggest that the formation of HVCs occurs under much more general circumstances, and only require that the Galactic corona is approximately in global thermal balance relatively close to the Galactic disc.

We conclude that if the corona rotates:

1. (a)

HVCs can occur as a consequence of the thermal instability in the Galactic corona. Our conclusion differs from that of Nipoti (2010) because his analysis was limited to the linear regime, while the direct effect (ii) mentioned above is largely a consequence of non-linear effects (Paper I), and because he did not apply his analysis to any concrete global model of rotating corona so he also partly missed the importance of the indirect effect (i). 2. (b)

The predicted locations for the formation of the condensations agrees with the observationally determined locations of most HVCs, which are observed at $z$-heights $\lesssim 10\,{\rm kpc}$, at distances between $2\mathchar 45\relax 15\,{\rm kpc}$ from the Sun and within $30^{\circ}$ from the disc plane as viewed from the Galactic Center (e.g. Putman et al., 2012). This does not require much fine tuning, and is a general property of models of the MW corona with a substantial rotational support. 3. (c)

The formation of HVCs requires the corona to be in global thermal balance (namely, heating roughly balances cooling on average) close to the disc. Hence, further understanding of the properties of HVCs requires a detailed modelling of the disc-halo interface, including the effect of supernova feedback.

Acknowledgements

The authors thank Lucia Armillotta, Filippo Fraternali, Uri Keshet, Ralf Klessen, Federico Marinacci, Carlo Nipoti, Gabriele Pezzulli and Steve Shore for useful comments and discussions. We are grateful to the referee, Prateek Sharma, for constructive suggestions that improved the paper. MCS acknowledges support from the Deutsche Forschungsgemeinschaft via the Collaborative Research Centre (SFB 881) “The Milky Way System” (sub-projects B1, B2, and B8). ES acknowledges support from the Israeli Science Foundation (grant 719/14) and from the German Israeli Foundation for Scientific Research and Development (grant I-1362-303.7/2016).

Appendix A Physical interpretation of the unstable mode that arises due to rotation

In the absence of rotation, the modes described by the dispersion relation (8) essentially describe blobs undergoing buoyant oscillations while cooling is acting. When rotation is added to the system, as noted by Nipoti (2010), a new unstable mode appears, which in the slow cooling regime has a frequency given by Eq. (17). What is the physical interpretation of this new mode?

In the limit $t_{\rm cool}\to\infty$ the frequency of the new unstable mode $\omega\to 0$ (Eq. 17). Hence the perturbed state is just a new equilibrium state in this limit and fluid elements associated with this type of perturbations maintain their equilibrium values in the $(R,z)$ plane at all times. Hence, these modes are clearly not buoyant oscillations. Indeed, it can be shown from the linearised equations of motion that the perturbations associated with this new unstable mode, in the slow cooling regime ($t_{\rm cool}\to\infty$), have $\delta v_{R}\to 0$ and $\delta v_{z}\to 0$. We now show how such modes can be constructed explicitly in this limit.

Imagine taking an equilibrium model such as those described in Sormani et al. (2018). Can we find another equilibrium model which is exactly the same as the previous one everywhere except within a region $A$ (see Fig. 7)?

Since as discussed in Sormani et al. (2018) an equilibrium is completely characterised by the pressure distribution $P_{0}(R,z)$, this is in principle easy to do: just take a function $P(R,z)=P_{0}(R,z)+\delta P(R,z)$ such that $P(R,z)=P_{0}(R,z)$ outside $A$ and $P(R,z)\neq P_{0}(R,z)$ within $A$. This will result in an equilibrium corona that is identical to the original one except within the region $A$. It is clear that if the associated density perturbation $\delta\rho/\rho_{0}$ is non-negligible then in the presence of cooling this region will cool at a substantially different rate than its surroundings (because cooling depends on the density as $\rho^{2}$), leading to a thermal instability. However, the perturbation analysis which underlies the derivation of the dispersion relation (8) requires that (i) $\delta F/F_{0}\ll 1$ where $F$ is any of $P,\rho,v_{\rm\phi}$ (ii) perturbations to be quasi-isobaric so that $\delta P/P_{0}\ll\delta\rho/\rho$ as a consequence of the Boussinesq approximation which filters out sound waves (iii) $\delta P$ varies on short spatial scales, so that $\left|\nabla\log\delta P\right|\gg\left|\nabla\log P_{0}\right|$, as a consequence of the WKB approximation. Thus the question is, given a pressure perturbation $\delta P$ that satisfies these requirements, can we achieve a non-negligible density perturbation within the region $A$? Only if this can happen will the system be unstable according to the linear perturbation analysis that leads to Eq. (8).

Given an equilibrium characterised by a pressure distribution $P\left(R,z\right)$ and an external potential $\Phi\left(R,z\right)$, the density $\rho(R,z)$ and azimuthal velocity $v_{\phi}(R,z)$ can be readily calculated using the following formulas (Sormani et al., 2018):

[TABLE]

Writing $P(R,z)=P_{0}\left(R,z\right)+\delta P\left(R,z\right)$, substituting into (19) and using assumptions (i)-(ii) above we find:

[TABLE]

which means that we may indeed regard the rapidly varying component of the solution as being isobaric, as required by the Boussinesq approximation. Substituting $P(R,z)=P_{0}\left(R,z\right)+\delta P\left(R,z\right)$ into Eq. (20) and linearising we find

[TABLE]

The ratio of the first to the second term on the right hand side is of order

[TABLE]

Therefore, using the fact that $P_{0}\sim\rho_{0}c_{\rm s}^{2}$ we see that

[TABLE]

Putting everything together, we have that

[TABLE]

Therefore, we can achieve a significant perturbation of $\delta\rho/\rho_{0}$ only if $(v_{0\phi}^{2}/c_{\rm s}^{2})\gtrsim 1$, i.e. only if the corona rotates significantly. This explains why the unstable mode only becomes important when the corona is rotating significantly. If the corona does not rotate significantly instead, $(v_{0\phi}^{2}/c_{\rm s}^{2})\ll 1$, the density perturbation $\delta\rho/\rho_{0}$ is suppressed by a factor $(v_{0\phi}^{2}/c_{\rm s}^{2})$ compared to the velocity perturbation. Therefore, if the velocity perturbation is small as required by the linear regime, the density perturbation is negligible, which also explains why the frequency of the unstable mode given by Eq. (17) goes to zero as the rotation velocity $\Omega\to 0$.

In conclusion, the new unstable mode that appears due to the fact that only in the presence of substantial rotational support it is possible to find a new equilibrium through a nearly-isobaric perturbation which has a non-negligible density perturbation. In the limit in which the corona is not rotating instead, a nearly-isobaric perturbation cannot achieve a substantial density perturbation unless the velocity perturbation $\delta v_{\phi}/v_{0\phi}$ is large, which is unlikely to happen spontaneously due to random fluctuations.