MODELING THE X-RAYS RESULTING FROM HIGH-VELOCITY CLOUDS

We use similar hydrodynamic and magnetohydrodynamic algorithms as in Kwak et al. (2009) and Kwak & Shelton (2010). To wit, we use the FLASH computer code, version 2.5 (Fryxell et al. 2000) with adaptive mesh refinement (AMR) to model the hydrodynamics and magnetohydrodynamics of fast-moving clouds and ambient gas. For our suite of cloud shock simulations (Section 3.2), we test the effects of radiative cooling by modeling it in some simulations but not others. The radiative cooling calculations are done with the FLASH module, which uses cooling rates for CIE plasmas. We do not simulate gravity in our cloud shock simulations, but instead start the cloud with a large initial velocity. We do the same in our CIE testing (Section 3.1) and resolution experiments (Section 3.5). In our turbulent mixing simulations (Section 3.3), we use the gravity module in FLASH and the expression for gravitational acceleration presented in Ferriére (1998) to establish hydrostatic equilibrium in the background gas and to simulate the Milky Way's gravitational pull on the cloud. In order to maintain hydrostatic balance in the background gas, radiative cooling must be disallowed and is so in the turbulent mixing simulations. It is not needed in the CIE testing and so is disabled in them as well. In two of our simulations, we model a magnetic field oriented perpendicular to the cloud's motion. The others have no magnetic field. For the purposes of the magnetohydrodynamics, we treat all of the gas as if it is fully ionized in all of our simulations. Thermal conduction is not modeled.

The predictions of the X-ray count rates and very high ion column densities require predictions of the fractions of ions in any given state in the gas. We calculate the ion fractions using CIE and/or partially non-equilibrium ionization models. In the CIE models, the ionization levels are calculated from the temperature of the simulated gas, using the Raymond and Smith code (Raymond & Smith 1977, with updates). For the NEI models, we enable FLASH to track the ionization levels of select elements as a function of time by accounting for their collisional ionization and recombination during each time step. Ionization and recombination rates from Summers (1974) are used. We limit our NEI calculations to only two elements, silicon and oxygen, in order to economize on CPU and memory. Silicon is a major contributor to X-ray spectra in the 1/4 keV band, while oxygen is an important contributor in the 3/4 keV band. In Section 3.1, we show that CIE models adequately approximate the NEI models.

In our CIE tests (Section 3.1), shock-heating analysis (Section 3.2) and turbulent-mixing analysis (Sections 3.3), we use a three-dimensional Cartesian coordinate system with width (span in the x-direction), and depth (span in the y-direction) of 1.5 kpc and 1.5 kpc, respectively. The footprint is centered on x = 0 kpc, y = 0 kpc. Most of these simulations use a domain that is 6 kpc tall. The exceptions are the high-resolution, radiatively cooling, shock scenario simulations, for which we reduced the domain height to 1.5 kpc. In the turbulent mixing simulations, hydrostatic equilibrium is achieved by decreasing the thermal pressure with height above the midplane. We set the physical conditions at the base of the domain to simulate those in the ISM at z = 8 kpc, thus in these simulations the domain runs from z = 8 to 14 kpc. In the remaining simulations, there are no gravitational, density, pressure, or magnetic field gradients in the modeled ambient gas. Thus, an arbitrary z offset can be added to the height along the z-axis. For these simulations, we choose to call the lower boundary z = 0 kpc, only for convenience. In the shock models that had no radiative cooling, each domain is initially segmented into four blocks. Each block is subdivided by the AMR routine an additional three to five times, with the number of subdivisions depending upon the density gradient. With five levels of refinement, the resulting blocks are 1500/(2^5–1) pc, on a side, while with three levels of refinement they are 1500/(2^3–1) on a side. Each of these blocks is segmented into 8³ zones, which, consequently, range in size from (12 pc)³ to (47 pc)³. Like the shock models that had no cooling, some of our shock models that had cooling also had moderate resolution and 6 kpc tall domains initially made from four blocks, while others used higher resolution and a single 1.5 kpc tall block that was subdivided an additional three to seven times depending upon the density gradient. These blocks were then subdivided into 8³ zones that range in size from (2.9 pc)³ to (47 pc)³. For our resolution tests (Section 3.5), we use both three-dimensional Cartesian domains and two-dimensional cylindrically symmetric domains and test maximum refinement levels between four and eight.

We calculate each zone's emission spectrum using the Raymond and Smith spectral code (Raymond & Smith 1977, with updates), the gas density calculated by FLASH, the ionization levels from either the CIE or the NEI calculations, and Allen (1973) cosmic abundances. The Allen abundances are greater than the gas-phase abundances in the halo and clouds. For example, the gas-phase metal abundances in the warm halo are ∼4/10 of the Allen solar abundances (Savage & Sembach 1996), while the gas-phase oxygen abundance in Complex C is about ∼1/10 to ∼1/4 of the Allen solar values (Fox et al. 2004), and the gas-phase metal abundances in Magellanic Stream clouds are ∼1/14 of the Allen solar values (Fox et al. 2010). The strength of the emitted spectra should be scaled accordingly, bearing in mind that in the shock scenario (the most X-ray productive scenario), the X-rays originate in the ambient gas and that over time, collisions in hot gas break down dust.

In order to obtain the spectrum pertaining to any given line of sight, we sum the spectra from the intersected zones. If the zones have differing sizes, then we weight the terms accordingly. The X-ray spectra are convolved with the ROSAT response matrix in order to determine the count rate in the ROSAT R12 band, the 1/4 keV energy band. For the brighter models, we also report intensities at the photon energies of the O vii triplet (∼570 eV) and O viii Lα line (653 eV). During our line intensity calculations, we do not subtract the continuum or the pseudo-continuum composed of faint, unresolved lines, because they account for ∼5% of the emission in these energy bins during the brightest phases. They account for a greater fraction of the intensity when the gas is dimmer, but at these times, the reported intensity is too dim to be observed. In addition, for comparison with the observed X-ray excess associated with a Magellanic Stream cloud (Bregman et al. 2009; 0.64 ± 0.10 counts ks⁻¹ arcmin⁻² seen by the XMM-Newton pn detector in 0.4–1.0 keV X-rays), we calculate the XMM pn count rate for four of our shock models.

We calculate the density of O vii ions for select models by combining the fraction of oxygen in the O vii ionization state with the oxygen abundance and the gas density in each zone. We integrate the density along sightlines through the domain in order to obtain O vii column densities. Again, the results can be scaled if elemental abundances other than those of Allen (1973) are preferred.

We initialize the 6 kpc tall domains such that a spherical cloud is located in the upper portion of the grid and is surrounded by stationary ambient gas. In the 1.5 kpc tall domains, the center of the cloud is initially located at 2/3 of the height of the domain. At the beginning of each simulation, the cloud has a radius of 0.2 kpc, which corresponds to 11 if the cloud is seen from the current distance to Complex C (10 kpc; Thom et al. 2008). In most of our simulations, the cloud's initial volume density of hydrogen nuclei, n_{cl, H}, is 0.0645 cm⁻³. Note that the simulations assume a 10 to 1 ratio of hydrogen to helium, making the mass density 0.0903 amu cm⁻³, and note that the number density shown in the figures and quoted in the text is the number of hydrogens per unit volume. We choose an initial hydrogen volume density of 0.0645 cm⁻³ because clouds with this density and a radius of 0.2 kpc have hydrogen column densities along sightlines through their centers of 8 × 10¹⁹ H cm⁻², which is similar to the observed H i column density in small features within Complex C. Initial cloud temperatures range from 100 K to ∼11, 000 K and are chosen from the constraint that the cloud's initial thermal pressure must balance the initial thermal pressure of the gas around it. In our cases with magnetic fields, the initial magnetic field pressures also balance.

Table 1 lists the cloud and ISM parameters for our simulations. They are grouped into three general cases, Case A (clouds falling under the influence of gravity through a hot ambient medium that is in hydrostatic equilibrium; all of these test turbulent mixing, and some also test shock heating), Case B (fast clouds moving through hot, rarefied gas before colliding with warm, moderate-density ambient gas; these test shock heating and there are two categories, those with radiative cooling (Br) and those without (Ba)), and Case C (fast clouds moving through warm, moderate-density ambient gas; these test shock heating).

In our suite of Case A simulations, the ambient medium is either T_ISM = 1 × 10⁶ or 3 × 10⁶ K and, as mentioned in Section 2.1, radiative cooling is disabled in order to maintain hydrostatic equilibrium in the ambient gas. We determine the initial density gradient from the equation for hydrostatic equilibrium with the constraints that the magnetic field strength and temperature are constant. We constrain the midplane density by requiring the midplane thermal pressure to be reasonable, i.e., roughly 6000 K cm⁻³, and the 1/4 keV X-ray surface brightness of the ambient medium above z = 2 kpc to be comparable to observationally derived values. As an example, we note the X-ray surface brightness in the 1/4 keV ROSAT band produced by the gas above 2 kpc in our Model A reference run, A1. Its intrinsic surface brightness is 1510 × 10⁻⁶ counts s⁻¹ arcmin⁻². If the column density of intervening H i is 1 × 10²⁰ cm⁻², then the absorbed surface brightness is 583 × 10⁻⁶ counts s⁻¹ arcmin⁻², which is consistent with measurements presented in Snowden et al. (1998, 2000). Our constraints yield the ambient density functions plotted in Figure 2; thus, the density of ambient gas, expressed in units of hydrogens per unit volume, at the cloud's initial height in most of our Case A simulations is between ∼5 × 10⁻⁵ and 2.5 × 10⁻⁴ cm⁻³. The values are consistent with those used or determined elsewhere (e.g., density of ≳ 2.4 × 10⁻⁵ cm⁻³ (Blitz & Robishaw 2000, ram pressure stripping of nearby dwarf galaxies); ∼5 × 10⁻⁵ cm⁻³ (Moore & Davis 1994, ram pressure stripping of Magellanic Clouds); 2 × 10⁻⁴ cm⁻³ in model by Bland-Hawthorn et al. 2007; 1 and 3 × 10⁻⁴ cm⁻³ in models by Heitsch & Putman 2009). In order to balance the ambient and cloud pressures at the beginning of each simulation, we allow the cloud temperature to vary slightly with height. We make a suite of simulations in order to sample various initial cloud heights and velocities, ambient densities and temperatures, and strengths of the magnetic field lying parallel to the Galactic midplane. All of the Case A models use CIE calculations in order to determine the fraction of atoms in any given ionization state. The hydrodynamical results from Model A are stored at 2 Myr intervals.

Zoom In Zoom Out Reset image size

In our suite of Case B simulations, the clouds pass through a hot ambient medium before hitting a warm extended medium. The hot medium may represent halo or intergalactic gas, while the warm medium may represent material shed from a preceding HVC, a high-z cloud, or, given the proximity of some HVCs to the Galactic disk, it may represent warm gas in the Galaxy's thick disk. In our setup, the hot and warm sectors are in pressure balance with each other. In most of our Model B simulations, the hot medium has a temperature of T_ISM = 1 × 10⁶ K, a density of hydrogen of n_{ISM, H} = 6.45 × 10⁻⁵ cm⁻³, i.e., n_{cl, H}/1000, and occupies the uppermost 1.5 kpc of the grid. The density of the hot portion of the domain is equal to that at z  ∼  12 kpc in Models A1, 2, 3, 6, and 7. In most of our Model B simulations, the warm part of the domain has a temperature of T_ISM = 10⁴ K, a density of hydrogen of n_{ISM, H} = 6.45 × 10⁻³ cm⁻³, i.e., n_{cl, H}/10, and occupies the lower 4.5 kpc of the grid. At the beginning of each Case B simulation, the environmental gas is in pressure balance with the cloud. There is no magnetic field and therefore no magnetic pressure. It is not practical to model hydrostatic equilibrium in Case B because the weight of the warm medium cannot be supported unless the pressure has a strong gradient. By giving up hydrostatic equilibrium, we also give up gravity. Because both are absent, an arbitrary offset can be added to our domain's z without affecting the results. On this arbitrary scale, the initial location of the center of the cloud is z = 5 kpc in simulations having 6 kpc tall domains (the moderate-resolution simulations) and 0.5 kpc from the top in the simulations having 1.5 kpc tall domains (the high-resolution simulations). Lacking gravity to accelerate them, the clouds must be given an initial velocity. We sample several velocities; see Table 1. Cases with speeds of 200 and 300 km s⁻¹ can be used for some Galactic HVCs. Cases with speeds of 300 and 400 km s⁻¹ can be compared with the Magellanic Stream, whose orbital speed is 378  ±  18 km s⁻¹ (Bland-Hawthorn et al. 2007). The simulations with faster clouds are useful for identifying trends. We test the effects of NEI by creating a small suite of models for which both NEI and CIE ionization fractions are calculated (these models have names that begin with Br and include radiative cooling). We test the effects of radiative cooling at the CIE rate by comparing models with and without radiative cooling. Those without radiative cooling have names that begin with Ba. The computational results for Case Ba models are stored at 2 Myr intervals. We performed several Model Br simulations with the same resolution and archival time intervals as the Model Ba simulations, but in order to better model the radiative cooling in this scenario, simulations having smaller zone sizes and stored time intervals are also needed. Thus, we also present simulations of Br models in which up to seven levels of refinement are allowed, the spacing between archived epochs is 40,000 years, and the simulations are stopped at 2 Myr. It is this set of simulations that has 1.5 kpc domain heights.

Case C is a variation on Case B, with the primary difference being that all of the ambient gas is warm (T_ISM = 1 × 10⁴ K) in Case C. In Case C, we also sample two ambient densities, a rarefied case having n_{ISM, H} = 6.45 × 10⁻⁴ H cm⁻³ and a moderate case having the same density as the warm gas in Case B, i.e., n_{ISM, H} = 6.45 × 10⁻³ H cm⁻³. In the former case, we also reduce the clouds' temperature in order to achieve pressure balance at the beginning of the simulations. All of the Case C models use CIE calculations, omit radiative cooling, and are archived every 2 Myr.

In order to evaluate the CIE approximation, we compare spectra that we calculate from the NEI and CIE ionization fractions of silicon and oxygen. The NEI ionization fractions, i.e., the fractions of the atoms at any given ionization level, were obtained by having FLASH track the ionization levels in a time-dependent fashion in moderate-resolution versions of Models Br3, Br4, and Br5. The CIE ionization fractions were obtained from the hydrodynamic information for these same models; specifically, they were calculated by the Raymond and Smith code using the temperature of the gas as an input. For each of these three models, we calculate the spectra produced by NEI silicon, CIE silicon, NEI oxygen, and CIE oxygen for various vertical sightlines through the domain at each epoch, 2 Myr, 4 Myr, etc. (Although the shocked gas has already radiated way much of its energy by 2 Myr, it is still sufficiently emissive for these experiments.) In order to reduce the quantity of information, we then convolve each silicon spectrum with the ROSAT 1/4 keV band response function and each oxygen spectrum with the ROSAT 3/4 keV band response function, yielding count rates in these bands. In each case, the CIE count rates in the 1/4 keV band are fairly similar to the NEI count rates. For example, the ROSAT 1/4 keV count rate calculated from the CIE and NEI silicon spectra produced along the central sightline of the Br4 model at 2 Myr are within 5% of each other. The relationship between the CIE and NEI oxygen spectra were not as consistent from model to model and epoch to epoch. For Model Br4 at 2 Myr, the CIE and NEI oxygen spectra are within 14% of each other. The emission line spectra for this case are displayed in Figure 3. Sometimes, the O vii intensities from the models are less consistently aligned with each other than are the silicon intensities and frequently the O viii quantities (which are always small) are poorly aligned, presumably due to delayed ionization and recombination. Similarly, the high-resolution Model Br simulations also show poor alignment between the CIE and NEI O viii column densities.

Zoom In Zoom Out Reset image size

Theoretically, if the clouds collide with neighboring gas at fast enough speeds, then shocks should develop in both media. The temperature (T₂) of the shocked plasma can be estimated from the equation for plane parallel shocks (Shu 1992):

$$\begin{eqnarray} &&T_2 = T_1 \frac{\big[(\gamma +1)+2\gamma \big(M_1^2-1\big)\big] \big[(\gamma +1)+(\gamma -1)\big(M_1^2-1\big)\big]}{(\gamma +1)^2 M_1^2},\nonumber \\ \end{eqnarray} \tag{ 1 }$$

where T₁ is the temperature of the unshocked gas, M₁ is the Mach number calculated from the speed at which the unshocked gas is overrun by the shock, and γ is the adiabatic constant for monotonic gas (5/3). Here, we compare the theoretical post-shock temperature with those of Models Ba3 and C3 at the first epoch (t = 2 Myr). These models begin with cloud velocities of v_z = −300 km s⁻¹ and do not model radiative cooling (although, in post-processing, we can estimate the X-ray spectra that would be emitted by gas at the simulated temperatures). By the first archived epoch at 2 Myr of simulated time, the cloud has shock heated the environmental gas beneath it while a reverse shock has been driven into the lower part of the cloud. In both Models Ba3 and C3, the shocked environmental gas was initially warm (T = 1 × 10⁴ K), ionized ISM. Although the unshocked portion of the cloud still travels with approximately its initial velocity, v_z = −300 km s⁻¹, at this epoch, the shocked portion of the cloud moves downward at speeds between 240 and 300 km s⁻¹ and the swept up ISM moves downward at speeds up to 240 km s⁻¹ with an average of about 200 km s⁻¹ and moves sideways with speeds exceeding 100 km s⁻¹. Thus, the impacted material does not simply pile up in front of the cloud, as in a shock-tube simulation, but instead partially skirts around the cloud. The shockfront travels downward at 4/3 of the shocked ISM's downward speed, which is confirmed by the displacement of the shockfront between t = 2 and 4 Myr. Thus, M₁ = 18. From Equation (1), we find that T₂ should be 1.0 × 10⁶ K, which is consistent with the post-shock temperatures in Models Ba3 and C3: 0.95 and 1.0 × 10⁶ K, respectively; see Figure 4. Since the clouds decelerate in our gravityless simulations, the post-shock velocity decreases with time and thus the post-shock temperature decreases with time as well. This can be seen by comparing the temperatures at 10 Myr with those at 2 Myr in Figure 4.

Zoom In Zoom Out Reset image size

A reverse shock propagates through the cloud, but due to the density contrast, it propagates slower than the forward shock in the ISM. In our simulations of Model C3 at 2 Myr, the reverse shock moves into the cloud at a speed of ∼90 km s⁻¹ from the cloud's reference frame. Given an initial cloud temperature of T_cl = 1000 K, the reverse shock should be strong and, according to Equation (1), it should heat the cloud to ∼1.2 × 10⁵ K. Our simulational results find the temperature to range from the cloud temperature to this value; see Figure 4. While the reverse shock-heated material is far hotter than the unshocked portion of the cloud, it is still not hot enough to produce X-rays. The X-rays that result from this geometry originate in the shock-heated ambient medium. In the following subsections, we explain that their intensity ranges from unobservably dim to significantly bright, depending upon the time since interaction, collision speed, and strength of radiative cooling.

3.2.1. Effect of Radiative Cooling

The simulations and theory discussed above did not include the effects of radiative cooling. When it is included in the calculations, the shocked ambient gas quickly cools. This can be seen by comparing the temperature and density distributions in Models Ba3 and C3 with those from Model Br3 in Figure 4. These are fair comparisons because each model has the same parameters for the cloud and lower ISM and, although the ionization fractions of some elements were calculated in a time-dependent manner in Model Br3 (but not in the others), the NEI ionization levels did not affect the hydrodynamics and we replace them with the CIE ionization levels for the following spectral calculations. Furthermore, we note that the ambient medium in Models Ba3 and Br3 include a layer of hot gas, but this layer contributes only 1 × 10⁻⁶ counts s⁻¹ arcmin⁻² which is trivially small and is subtracted before we report the results.

Without radiative cooling, the shock-heated gas is very bright for millions of years. As can be seen in Figure 5, even the moderately fast clouds (v_z = 200 km s⁻¹) induce surface brightnesses across their footprints of ∼1000 × 10⁻⁶ counts s⁻¹ arcmin⁻² in the ROSAT R12 band at their peak. Faster clouds produce ≳ 10, 000 × 10⁻⁶ counts s⁻¹ arcmin⁻² in this band at their peaks. These models are much brighter than the typical intrinsic intensity of 1/4 keV photons emitted above the Galactic disk, which is ∼1000 to ∼2000 × 10⁻⁶ counts s⁻¹ arcmin⁻² with significant directional variation (Snowden et al. 2000). All of these simulated clouds continue to make >100 × 10⁻⁶ counts s⁻¹ arcmin⁻² until at least 14 Myr after the collision.

Zoom In Zoom Out Reset image size

When radiative cooling is included in the simulation, the hot gas in front of the cloud sheds most of its energy via radiation, only some of which is in the X-ray band. Figures 5 and 6 include plots of the resulting 1/4 keV surface brightnesses from the high-resolution Model Br simulations. The average surface brightnesses across footprints that are 200 pc in radius peak around 3, 500, and 2000 × 10⁻⁶ counts s⁻¹ arcmin⁻² in the ROSAT 1/4 keV band within 1 to ∼1.5 Myr of the collision for the v_z = 200, 300, and 400 km s⁻¹ collisions, respectively.

Zoom In Zoom Out Reset image size

The surface brightness profiles of Models Ba3 and C3 in Figure 6 show that the surface brightness peaks near each cloud's axis and drops non-monotonically with radius. Model Br3's profile is more complicated. The peaks and valleys in the profiles are associated with the density and temperature structure of the off-axis gas. When the gas is very bright, the profile extends beyond the footprint of the cloud. Sensitive X-ray observations of the clouds would see extended disks.

The shock-heated gas is less conspicuous in the 3/4 keV band than in the 1/4 keV band, in both the radiative cooling and the non-cooling scenarios. For both Models Ba3 and C3, the O vii triplet (photon energy ∼570 eV) intensity across a disk of radius = 200 pc averages to ∼1 photon s⁻¹ cm⁻² sr⁻¹ for the first several million years (see Figure 7). It averages to far less than 1 photon s⁻¹ cm⁻² sr⁻¹ for Model Br3. The former intensity is about ∼1/5 of the typical observed intensity on randomly chosen sightlines away from the Galactic center (Henley & Shelton 2010). Because the collisions weaken with time, the O vii intensities fade sooner than the 1/4 keV surface brightnesses. Furthermore, none of the v_z = −300 km s⁻¹ models produce appreciable O viii intensities at any time.

Zoom In Zoom Out Reset image size

The average O vii column density within a footprint of radius 200 pc is 3 × 10¹⁵ cm⁻², for both Models Ba3 and C3 at 2 Myr. It remains near this level for 10 Myr, thus outliving the duration of the bright O vii emission. In Model Br3, the average O vii column density is ∼4 × 10¹⁴ cm⁻² from about 1.3 to about 1.5 Myr. The O vii and O viii column densities for these and other speeds simulated for Model Ba and the high-resolution versions of Model Br are plotted in Figure 8.

Zoom In Zoom Out Reset image size

The gas in these models is also UV emissive. When we divide the predicted energy spectra into 0.25 dex wide bins, we find that the 10–18 eV photon energy range is the most luminous part of the ∼5 to ∼2000 eV spectrum for Models Ba3, Br3, and C3. Model Br3, for example, radiates more than 100 times more power in 10–18 eV photons than in 100–180 eV photons when it is at its peak 1/4 keV X-ray brightness (age ∼1.3 Myr). At earlier and later times, the ratio is even larger. With its comparatively large radiative power, the UV plays a more important role in cooling the shock-heated gas than does the X-ray.

We have presented results for both radiatively cooling and adiabatic models. In the following subsections, we continue to discuss both cases. While, on the face of it, it appears contradictory to consider the X-ray count rates expected from hot gas whose hydrodynamics were calculated in models that excluded radiative cooling (Models Ba and C), we note that (1) the true cooling rate is unknown, (2) if the gas-phase abundances are between zero and the solar values, then the cooling rate should be between that used in the Br suite of models and those used in the Ba and C suites, and (3) X-ray emission accounts for very little cooling. Because X-ray emission accounts for such a small fraction of the cooling and because the relative contribution to the spectrum made by any given element varies from one energy band to the next, it is possible for variations in gas-phase metallicities to preferentially lower the cooling rate compared with those in Case Ba or raise the X-ray count rate compared with those from Case Br. Iron provides a good example of the effect of the gas-phase metallicity. Almost 60% of the energy radiated in 5 to 100 eV photons by a T = 10⁶ K plasma having Allen (1973) abundances are emitted by iron ions, while slightly less than 20% of the ROSAT R12 counts are due to photons emitted by these ions. Thus, if the true iron abundance in the shock-heated gas were half that expected by Allen (1973), then the overall radiative cooling rate would be noticeably affected (it would be about 70% of the CIE rate), while the ROSAT R12 count rate would be insignificantly affected (it would be about 90% of CIE rate). Reducing the iron abundance to zero decreases the radiative cooling rate by 3/5, while reducing the 1/4 keV luminosity by only 1/5. Similar phenomena occur at higher temperatures. Silicon provides another example, because it accounts for 13% of the 5–100 eV power, but almost 30% of the ROSAT R12 count rate. Thus, doubling its relative abundance would increase the radiative loss rate by only 1/8, but increase the 1/4 keV count rate by almost third. Here, we are not suggesting that the halo has supersolar abundances of silicon, but are pointing out that if the ratio of silicon to other elements is larger than in solar abundance gas, the effect would preferentially benefit the X-ray luminosity. The interested reader is referred to Sutherland & Dopita (1993, Figure 18) for relative cooling rates from various elements. Other reasons why the rate of radiative cooling might differ from the CIE rate for solar abundances include disequilibrium between electron and ion temperatures and delayed ionization and recombination.

3.2.2. Effect of HVC Velocity

3.2.3. Effect of Density

Here, we consider mixing between cool and warm clouds as they travel through hot, rarefied ambient gas. As the cloud passes through the hot gas, Kelvin–Helmholtz instabilities develop and mixing between the cloud and ambient material creates a zone of intermediate-temperature, intermediate-density gas (Esquivel et al. 2006). Analytical calculations and computational simulations have already shown that mixed layers are rich in high ions, such as O vi (Slavin et al. 1993; Esquivel et al. 2006; Kwak & Shelton 2010; Kwak et al. 2011). In order to determine whether or not mixing zones between HVCs and very hot gas can also be X-ray productive, we simulate the hydrodynamics of three-dimensional clouds passing through hot, rarefied gas. This approach differs from the more common approach, which treats the two regions as blocks that slide past each other, but better replicates the effects of the cloud's rounded shape and motion. The names of these simulations begin with A.

In order to model the interaction adequately, the simulations must resolve the mixing length scale. Interstellar turbulence may exist across a wide range of length scales, but the largest length scale dominates the mixing (Esquivel et al. 2006). In our case, the largest length scale is the height of the distorted cloud. This height and width are adequately resolved in our Model A simulations, which, when fully resolved, use ∼30 zones to model a 400 pc span (the nominal height of the interaction zone in Figure 9).

Zoom In Zoom Out Reset image size

Our simulations do not include thermal conduction, but de Avillez & Breitschwerdt (2007) found that turbulent diffusion is much more effective than thermal conduction at leveling the temperature gradient. Thus, our omission of thermal conduction is not problematic. Note, also, that radiative cooling has not been allowed in the Model A simulations. It was disallowed in order to maintain hydrostatic balance in the background gas. The ramifications of cooling are too complicated for us to simply say that adding cooling would lower the radiation rates of the turbulently mixed gas. This is because radiative cooling would not only allow the ∼1 × 10⁶ K gas to cool to lower temperatures, but it would allow hotter gas to cool to ∼1 × 10⁶ K.

Model A1 serves as an example of Case A simulations and as our reference simulation for comparison with the others. In this simulation, we track the cloud as it falls from its initial height of 12 kpc to near the bottom of our grid (z = 8 kpc), which takes 32 Myr of simulated time. We record the simulation results at 2 Myr intervals and here describe the results at two of the sampled epochs. By 16 Myr, the cloud has accelerated to v_z = −120 km s⁻¹ and fallen to a height of z  ∼  11 kpc. By 30 Myr, it has accelerated to v_z = −220 km s⁻¹ and fallen to z  ∼  8.2 kpc. We obtain the cloud-induced excess intrinsic count rate in the 1/4 keV band by subtracting the count rate of the undisturbed hot halo from that toward the disturbed gas. Vertical sightlines throughout the domain are sampled. At 16 Myr, the brightest X-ray excess is 0.1 × 10⁻⁶ counts s⁻¹ arcmin⁻². At 30 Myr, it is 0.4 × 10⁻⁶ counts s⁻¹ arcmin⁻². The latter is the brightest rate for Model A1 at any time. Note that these are intrinsic count rates; absorption by N_H  ∼  10²⁰ cm⁻² would decrease them by factors of >2. Even without absorption, the small count rate enhancements would be inconspicuous against the diffuse X-ray background.

The timescale for the growth of the Kelvin–Helmholtz instabilities decreases as the cloud falls and the density contrast between the cloud and the ambient gas decreases. For the last ∼14 Myr of the simulation, we estimate that the timescale for the growth of instabilities is ≲ 14 Myr, from Chandrasekhar (1961, Section 101), and assuming that instabilities grow on similar timescales in plane parallel geometries as in our geometry. Hence, the lack of X-rays due to turbulent mixing in our simulations is unlikely to be due to the instabilities having insufficient time to grow.

Models that simulate a larger ambient temperature (Models A4 and A11), a greater ambient pressure (Models A5 and A11), faster initial speeds (Models A2, A8, A9, and A10), denser clouds (Model A11), less dense clouds (Models A3 and A10), and nonzero magnetic fields (Models A6 and A7) were also simulated. Of the slow models (initial |v_z| ⩽ 50 km s⁻¹) the greatest X-ray enhancement seen anywhere in the domain at any time during the simulations is 35 × 10⁻⁶ counts s⁻¹ arcmin⁻². This occurred in Model A11 at 30 Myr, just before the cloud crossed the domain's lower boundary. The next brightest model was Model A4, with a peak enhancement of 14 × 10⁻⁶ counts s⁻¹ arcmin⁻² at 30 Myr, also just before the cloud crossed the domain's lower boundary. Models A4 and 11 have the largest ambient gas temperature (T_ISM = 3 × 10⁶ K) and, through most of the domain, have the largest ambient density of all of the Case A models. These attributes prime Models A4 and 11 to be bright; the hotter ambient medium produces hotter mixed gas, while the denser ambient and cloud gas produce denser mixed gas than the other models. Although Model A11 is within a factor of a few of being comparable with X-ray observations, suggesting that further adjustment of the model parameters could result in a sufficiently bright model, the parameters in Model A11 are already near reasonable limits.

Some of the fast models (initial |v_z| = 300 or 400 km s⁻¹) create brighter disturbances. Model A9, for example, has a peak enhancement of 78 × 10⁻⁶ counts s⁻¹ arcmin⁻² at 8 Myr, which is the last epoch before the cloud leaves the domain. The bright region is fairly extended, with a radius of 740 pc (with the cut off being defined where the count rate exceeds the background by 5%) over which the count rate averages 75 × 10⁻⁶ counts s⁻¹ arcmin⁻². While an excess of this magnitude and spatial extent may be observable, it is not due to turbulent mixing. It is due to the shock. Likewise, the shock in Model A8 created a bright region, whose maximum average brightness was 20 × 10⁻⁶ counts s⁻¹ arcmin⁻² over a 200 pc radius footprint at 10 Myr, the last epoch before the cloud left the grid. Again, the X-rays are due to the shock, not the mixed material.

Here, we consider the dynamics of the mixed gas and how they lead to X-ray dimness. Although our models develop a mixed zone, as expected, and although the mixed zone contains some hot gas, in several models (A1–A3 and A5–A10) too little of the mixed gas was sufficiently hot (i.e., T nearing 10⁶ K) to be X-ray emissive and even when the temperature was sufficiently high, as it was in Models A4 and A11, the mixed gas fell behind the clouds into the semi-vacuum created by the cloud's passage. Here, the density was not sufficient for great emissivity. Model A11 provides an example of the varying conditions and low pressure in the ablated material (see Figure 9). Within the trailing stream, the temperature and density vary from hot (T = 2.6 × 10⁶ K) but diffuse (stream density ∼90% of the density in the undisturbed halo gas at this height) near the surface of the stream to only mildly hot (T = 1.2 × 10⁶ K) but denser (stream density ∼180% of the density in the undisturbed gas) in the core. The pressures in these example locations are ∼70% of those of the ambient gas.

The enhancement due to turbulent mixing by a single cloud is modest when compared with the typical X-ray count rate for high-latitude sightlines. Scattered clouds may contribute unidentifiably to the soft X-ray background and multiple, aligned clouds could create a non-negligible X-ray surface brightness. However, obtaining bright X-rays from individual clouds requires fast speeds as discussed in the preceding section.

Like the suites of Ba, Br, and C models, the suite of A models radiates more prolifically in the UV than in the X-ray. The broadband spectra are several orders of magnitude brighter in the far UV than in the soft X-ray. The ultimate source of this energy is the reservoir of thermal energy in the hot ambient gas. Mixing with the ablated cloud material has lowered the temperature and ionization level of the neighboring hot gas such that it has become highly emissive, especially in the UV.

To the observer who looks straight upward at an incoming HVC, the X-ray bright region would extend for at least 200 pc from the center of the cloud, thus 400 pc in diameter. If the observer is not located directly beneath the cloud, then the bright region would be somewhat displaced from the cloud itself and would subtend a smaller angle in the direction along the cloud's motion. If the clouds were as far away as Complex C, which is located ∼12 kpc from Earth, then the 400 pc diameter X-ray bright footprint due to a single cloudlet would subtend only a 2° angle. Such a small feature cannot be easily examined with ROSAT All Sky Survey (RASS) data. However, the overlapping footprints of many bright clouds in an ensemble could create an extended X-ray bright region that should be compared with the RASS maps.

If Complex C is composed of such an ensemble, then we can estimate the average 1/4 keV count rate across the complex. Here, we examine the cases in which the X-rays result from shock heating by assuming that the individual clouds within Complex C are like Model Ba3 or Model Br3 clouds; faster models would result in greater X-ray count rates while slower models or more cooling would result in lower count rates. Any gas that was ablated from the clouds and mixed with the ambient medium would also be dim.

The expected count rate is a product of the count rate of a single cloud when viewed from below (r), the number of such clouds ($\mathcal {N}$), and a scaling factor (f) that accounts for the dilution of the surface brightness over the larger area of Complex C. For Model Ba3, at its brightest epoch, the average 1/4 keV surface brightness, r, within a circular extraction region of radius 400 pc (this is twice as large as the previously mentioned extraction region, in order to capture more of the photons) is 3600 × 10⁻⁶ counts s⁻¹ arcmin⁻² in the ROSAT R12 band. In Model Br3, the emission is dimmer, but more concentrated. The surface brightness in the ROSAT R12 band for a 200 pc radius footprint is 560 × 10⁻⁶ counts s⁻¹ arcmin⁻². The number of individual clouds is the ratio of Complex C's mass (M_CC = 8.2^+4.6_{− 2.6} × 10⁶ M_☉; Thom et al. 2008) to the initial mass of one Model Ba3 or Model Br3 cloud (M = 7.5 × 10⁵ M_☉), thus $\mathcal {N}$ is 10.9^+6.1_{− 3.5}. If the area of each model extraction region, A = π × (400 pc)² for Model Ba3 and A = π × (200 pc)² for Model Br3, were to be diluted so as to encompass Complex C, whose area is A_CC = 3 kpc ×15 kpc (Thom et al. 2008), then conservation of luminosity would require the average count rate of each cloud to be reduced by a factor of f = A/A_CC. Not only does this factor account for dilution, but it also accounts for the concentration of the X-rays when the viewing angle results in a foreshortened cross section. Combining r, $\mathcal {N}$, and f yields a theoretical intrinsic R12 count rate of 440⁺²⁴⁰_{− 140} counts s⁻¹ arcmin⁻² for the case in which the clouds are like those in Model Ba3 and 17.0^+9.6_{− 5.4} × 10⁻⁶ counts s⁻¹ arcmin⁻² when they are like Model Br3.

Attenuation by intervening material will reduce the count rate. Assuming that N_H  ∼  10²⁰ cm⁻², which is the typical column density of galactic gas in the directions toward the brighter parts of Complex C, ∼2/3 to ∼3/4 of the original photons will be absorbed or scattered before reaching the observer. Thus, the observer would see an average cloud-induced X-ray surface brightness of ∼130 × 10⁻⁶ counts s⁻¹ arcmin⁻² from a complex of Model Ba3 clouds. This is a significant enhancement, though not as large as that seen along some directions in the region of Complex C (see Figure 1). For a complex of Model Br3 clouds, the predicted count rate is 4 × 10⁻⁶ counts s⁻¹ arcmin⁻², which is negligible. Irrespective of the cooling rate, faster clouds would result in brighter X-rays.

The Magellanic Stream provides another point of comparison. Bregman et al. (2009) report an enhancement of 0.4–1.0 keV X-rays on the leading side of the MS30.7-81.4-118 cloud within the Magellanic Stream. With the XMM-Newton pn detector, they found an excess of 0.64 ± 0.10 counts ks⁻¹ arcmin⁻². This is a 6.4σ effect. (Bregman et al. 2009 also found an excess in their Chandra data, but at the 1.6σ level.) The bright region is roughly 6 arcmin across, equating to roughly 90 pc in width. Here, we compare with the count rates for 100 pc wide circular footprints from our Models Ba3, Ba4, Br3, and Br4, which, with velocities of 300 and 400 km s⁻¹, bracket the stream velocity (roughly 380 km s⁻¹; Bland-Hawthorn et al. 2007). Our non-radiative models, Models Ba3 and Ba4, produce 0.25 and 3.2 counts ks⁻¹ arcmin⁻² at their brightest epochs (6 and 4 Myr, respectively), thus bracketing the observed value. Meanwhile, our radiative models, Models Br3 and Br4, produce only 7.0 × 10⁻³ and 0.42 counts ks⁻¹ arcmin⁻² at their brightest epochs (1.16 and 0.92 Myr, respectively). Only the faster non-radiative model is within range of the Magellanic Stream observations. Thus, collisions between fast HVC gas and relatively dense warm gas can account for the observed X-rays, but do so more easily if the radiation rate is quenched and/or the cloud is moving very fast.

As shown in Figure 8, our faster Model Ba simulations yield O vii column densities of ≳ 4 × 10¹⁵ cm⁻². This is similar to the median column density for sightlines through the Galactic halo (Lei 2010, excluding sightlines through the Galactic center soft X-ray enhancement), but is also of the same order of magnitude as the sightline-to-sightline variation in observed values. Thus, in the absence of radiative cooling, HVC-shocked gas might be observable with future, high-resolution instrumentation that could distinguish fast-moving from slow-moving ions. Again, if the shock-heated gas cooled at the CIE rate, then far fewer O vii ions would result, making the region unobservable.

In order to examine the effects of computational resolution, we calculate additional versions of Models Ba3, 4, 5, and 6 and Model Br3 using lesser and greater numbers of refinement levels than in the foregoing simulations, i.e., in FLASH, we use lrefine_max = 4 and 6, for the additional Model Ba3, 4, 5, and 6 simulations, rather than the value of 5 used in the primary simulations discussed in previous sections of the paper. We follow a similar pattern for comparison with the moderate-resolution Model Br3 simulations, but also use lrefine_max = 5 and 6 when making simulations for comparison with the high-resolution Model Br3 simulations, which used lrefine_max = 7. We find that within this range, the refinement level does not affect the timescale on which the clouds fragment, although it does affect the shapes of the clouds and of the X-ray bright regions. In order to compare the X-ray productivities of the various simulations, we extract the 1/4 keV X-ray count rates within circular regions of radius equal to 400 pc for the Case Ba simulations and 200 pc for the comparisons with the Case Br simulations. (We set the footprints for the Case Ba simulations to be greater than those used in earlier parts of this paper in order to capture all or nearly all of the downward directed flux.) We extract these count rates from every epoch in the Model Ba4-like simulations, every epoch from the Model Br3-like simulations that were made for comparison with the primary moderate-resolution Model Br3 simulation, every epoch from the simulations that were made for comparison with the high-resolution Model Br2, 3, and 4 simulations, and the t = 10 Myr epoch from the Model Ba3-, Ba4-, Ba5-, and Ba6-like simulations. We then compare our resolution experiment simulations with the control simulations having the same initial cloud velocity. The X-ray count rates vary somewhat between our test simulations, but in almost all cases are within 45% of that of the relevant control simulation. Frequently, they are much closer to those of the control simulations.

We also examined two-dimensional cylindrically symmetric simulations using maximum refinement levels of 4, 5, 6, 7, and 8. The principle advantage of this configuration is that it allows much smaller zone sizes. The principle disadvantage is that two-dimensional simulations are not able to track azimuthal instabilities (Korycansky et al. 2002). Lacking azimuthal modulation, the cloud material concentrates along the symmetry axis and resists fragmentation until later times than in the three-dimensional simulations. In general, the two-dimensional simulations predict similar X-ray count rates as the three-dimensional simulations, but with greater variation between runs having different refinement levels.