The density variance – Mach number relation

The density variance – Mach number relation in isothermal and non-isothermal adiabatic turbulence

C. A. Nolan, C. Federrath & R. S. Sutherland
Research School of Astronomy Astrophysics, The Australian National University, Canberra, ACT 2611, Australia
E-mail: (ANU)E-mail: (ANU)E-mail: (ANU)

The density variance – Mach number relation of the turbulent interstellar medium is relevant for theoretical models of the star formation rate, efficiency, and the initial mass function of stars. Here we use high-resolution hydrodynamical simulations with grid resolutions of up to cells to model compressible turbulence in a regime similar to the observed interstellar medium. We use Fyris Alpha, a shock-capturing code employing a high-order Godunov scheme to track large density variations induced by shocks. We investigate the robustness of the standard relation between the logarithmic density variance () and the sonic Mach number () of isothermal interstellar turbulence, in the non-isothermal regime. Specifically, we test ideal gases with diatomic molecular () and monatomic () adiabatic indices. A periodic cube of gas is stirred with purely solenoidal forcing at low wavenumbers, leading to a fully-developed turbulent medium. We find that as the gas heats in adiabatic compressions, it evolves along the relationship in the density variance – Mach number plane, but deviates significantly from the standard expression for isothermal gases. Our main result is a new density variance – Mach number relation that takes the adiabatic index into account: and provides good fits for . A theoretical model based on the Rankine-Hugoniot shock jump conditions is derived, , and provides good fits also for . We conclude that this new relation for adiabatic turbulence may introduce important corrections to the standard relation, if the gas is not isothermal ().

equation of state – galaxies: ISM – hydrodynamics – ISM: clouds – ISM: structure – turbulence

1 Introduction

The interstellar medium (ISM) is a complex, turbulent, multi-phase gaseous medium, which permeates the space between stars in the galactic plane (Ferrière, 2001). It is an essential part of the evolutionary cycle in stars, recycling the products of nucleosynthesis from dying stars and creating the stellar nurseries for a new generation of star formation (Mac Low & Klessen, 2004; Elmegreen & Scalo, 2004; McKee & Ostriker, 2007; Krumholz, 2014; Padoan et al., 2014). The ISM interacts with supernova explosions, protostellar jets, winds and outflows, which shape its structure and drive the turbulence we observe via atomic and molecular line observations of the ISM.

In many simulations that include an ISM (e.g., Bland-Hawthorn et al., 2007; Wagner & Bicknell, 2011; Cooper et al., 2009; Fischera & Dopita, 2005), a statistical construction with isotropic properties related to turbulent statistics have been used, as a proxy for the turbulent ISM. Causal models of the turbulent ISM will drastically increase the accuracy of these models, but this first involves an in-depth study of the statistics and evolution of fully-developed turbulence.

In purely isothermal gas, the probability density function (PDF) of the gas densities may be approximated by a lognormal distribution (Vázquez-Semadeni, 1994), which in log-space has the form


where , and are the mean and variance of the logarithm of the density , scaled by the mean density, . The logarithmic density variance is a function of the root-mean-squared (rms) sonic Mach number (), and is given by


The coefficient is known as the turbulence driving parameter and depends on the mode mixture induced by the turbulent forcing mechanism (Federrath et al., 2008). Purely solenoidal (divergence-free) driving leads to , while purely compressive (curl-free) driving corresponds to .

Equation (2) has been studied extensively for isothermal gases (Padoan et al., 1997; Passot & Vazquez-Semadeni, 1998; Kritsuk et al., 2007; Beetz et al., 2008; Federrath et al., 2008; Price et al., 2011; Burkhart & Lazarian, 2012; Seon, 2012; Konstandin et al., 2012), with investigation into different simulation techniques (Price & Federrath, 2010) and stirring methods (Federrath et al., 2008, 2010). It has also been studied by employing a heating and cooling curve (Wada & Norman, 2001; Kritsuk & Norman, 2002; Audit & Hennebelle, 2005, 2010; Hennebelle & Audit, 2007; Seifried et al., 2011; Gazol & Kim, 2013). Recently an investigation has been done in the magnetohydrodynamic (MHD) regime (Molina et al., 2012), and on polytropic gases (Federrath & Banerjee, 2015).

In our work, we investigate the robustness of this well-established density variance – Mach number relation, Equation (2), in the non-isothermal regime, specifically in ideal gases with diatomic molecular () and monatomic () adiabatic indices.

This is relevant because theoretical models of the star formation rate (Krumholz & McKee, 2005; Padoan & Nordlund, 2011; Hennebelle & Chabrier, 2011; Federrath & Klessen, 2012), the star formation law (Federrath, 2013b), the star formation efficiency (Elmegreen, 2008), and the initial mass function of stars (Hennebelle & Chabrier, 2008; Hopkins, 2013a; Chabrier et al., 2014) heavily rely on Equation (2).

Section 2 summaries our simulation and analysis methods, Section 3 first presents results for the isothermal case, in order to make contact with previous studies and to verify our analysis techniques. Then we present a numerical resolution study to determine the minimum resolution required in order to measure the density variance – Mach number relation in simulations with and present our main results for adiabatic indices and . We provide a theoretical model for the relation in Section 4 and discuss the discrepancies that we find compared to the standard Equation (2). Section 5 summarises our conclusions.

2 Simulation and analysis methods

Simulation name
(01) AD-TURB-256-A200-G1 1.0001 200
(02) AD-TURB-256-A400-G1 1.0001 400
(03) AD-TURB-256-A800-G1 1.0001 800
(04) AD-TURB-256-A1600-G1 1.0001 1600
(05) AD-TURB-256-A100-G7/5 7/5 100
(06) AD-TURB-256-A200-G7/5 7/5 200
(07) AD-TURB-256-A400-G7/5 7/5 400
(08) AD-TURB-1024-A200-G5/3 5/3 200
(09) AD-TURB-512-A200-G5/3 5/3 200
(10) AD-TURB-256-A100-G5/3 5/3 100
(11) AD-TURB-256-A200-G5/3 5/3 200
(12) AD-TURB-256-A400-G5/3 5/3 400
(13) AD-TURB-128-A200-G5/3 5/3 200
(14) AD-TURB-64-A200-G5/3 5/3 200
Table 1: Simulation parameters. Notes. Column 1: simulation name. Column 2: grid resolution. Column 3: adiabatic exponent in Equation (5). Column 4: dimensionless driving amplitude of the turbulence. Column 5: Resulting time-averaged velocity dispersion in code units in the regime of fully developed turbulence. Column 6: turbulent box crossing time: in code units.

To simulate the turbulent ISM we use the high-resolution, shock-capturing code Fyris Alpha (Sutherland, 2010) to solve the equations of compressible hydrodynamics across a three-dimensional, periodic domain with side length , initial uniform density , pressure of 1/2 (), and zero initial velocities. Unlike previous studies, the goal here is to test the density and velocity statistics of purely adiabatic turbulence with an ideal gas EOS, rather than a purely isothermal or polytropic EOS, and simpler than employing a cooling curve or running chemo-hydrodynamical simulations (Glover et al., 2010). Ultimately, simulations with multiple species including all relevant chemical reactions, as well as radiative heating and cooling would be the most realistic, but their complexity might not allow us to reduce the results to some simple rules of thumb that can be used in practical applications. We thus simplify the problem significantly by studying the turbulence in purely adiabatic, ideal gases with the aim of extracting results that might be applicable to a wider range of cases, including terrestrial experiments and atmospheric turbulence, in addition to the ISM. Table LABEL:tab:sim lists the key parameters of all our adiabatic turbulence simulations.

2.1 Ideal gas equation of state

The ideal gas equation of state relates the pressure , density and temperature , and is given by


with the volume , the Boltzmann constant , the particle mass , the total number of particles , and the particle number density . The ratio of the specific heat capacities at constant pressure and constant volume defines the adiabatic index,


where denotes the number of degrees of freedom. For monatomic gas, and , while for diatomic molecular gas, and , because diatomic molecules have two rotational degrees of freedom in addition to the three translational degrees of freedom. Note that at typical molecular cloud temperatures (about ), oscillatory degrees of freedom cannot be excited by collisions, which is why—although theoretically present—they do not contribute to increase in such cases. The specific internal energy of an ideal gas, , is only determined by its temperature. Inserting this equation into Equations (3) and (4), leads to


expressed via the adiabatic index , which serves as the equation of state. In order to determine the statistics of turbulence in this adiabatic regime, we use isothermal, diatomic molecular and monatomic equations of state ( and respectively).

Note that in order to model isothermal gases, is often set close to unity (e.g., ), as if the gas had an extremely large number of degrees of freedom . This trick produces a gas that approximately stays at constant temperature, because any excess heat from dissipation (e.g., by shocks) is absorbed in such a big internal energy reservoir that the temperature of the gas does not notably change.

2.2 Driving of turbulence, time evolution, and definition of the Mach number

The driving of turbulence in the gas is performed by stirring with random purely solenoidal (divergence-free) forcing at low wavenumbers for the duration of the simulation. All wavenumbers in the range were driven. The driving pattern is evolved with an Ornstein-Uhlenbeck process similar to the methods explained in Eswaran & Pope (1988), Schmidt et al. (2009), and Federrath et al. (2010).

From the work done by Federrath et al. (2008, 2010), we expect a proportionality constant of in Equation (2) for solenoidally driven isothermal gas and we will test that in both isothermal and adiabatic gases. The rms Mach number of the gas is modified by varying the stirring amplitude of the driving force, allowing each to be tested at a range of Mach numbers.

Figure 1: The velocity dispersion , as a function of simulation time for driving amplitudes , and , and for fixed adiabatic . The times for the onset of turbulence in each case are shown as vertical dashed lines and are approximated with the box crossing time , where is the linear size of the computational domain. The box crossing time for each simulation is listed in Table LABEL:tab:sim.

All simulations are run for several turbulent box crossing times to test the density variance – Mach number relation in the regimes of transient as well as fully-developed turbulence. The turbulent crossing time is defined as , where is the asymptotic velocity dispersion. The velocity dispersion and hence the crossing time depend on the driving amplitude . We show an example of this dependence in Fig. 1. We assume that the turbulence becomes fully developed after one crossing time in each simulation, indicated by the dashed vertical lines in Fig. 1. At this point the gas properties no longer vary drastically but change smoothly, indicating a statistically stable configuration. This allows us to distinguish regimes of transient () and fully-developed () turbulence.

Given the velocity dispersion and sound speed , we have different Mach numbers , depending on the driving amplitude and depending on the value of adiabatic . This is because the sound speed depends on the derivative of the EOS, Equation (5). It furthermore depends on the simulation time, because the internal energy and thus the mean temperature of the gas keeps increasing during the course of the adiabatic simulations with and . This is in stark contrast to the isothermal and polytropic simulations performed in previous studies (Padoan et al., 1997; Passot & Vazquez-Semadeni, 1998; Federrath et al., 2008; Price et al., 2011; Konstandin et al., 2012; Federrath & Banerjee, 2015), where the sound speed did not change systematically, after the turbulence was fully developed. Here, however, the total energy is conserved, which means that all dissipated energy is conservatively added to the internal energy. Thus, the injected energy from the driving is converted into internal energy and heats the gas continuously, leading to an ever increasing average sound speed and to a continuously decreasing rms Mach number. We thus use instantaneous measurements of and in the following to determine the density variance – Mach number relation in adiabatic gases.

2.3 Measuring the density variance

The density variance of the gas is calculated using method 4 in Section 2.3 of Price et al. (2011), but instead of fitting a lognormal distribution, we fit the more appropriate Hopkins (2013b) distribution. The advantage of the Hopkins fit is that it takes turbulent intermittency effects into account and provides excellent fits to the density PDFs over a wide range of physical parameters, including different Mach numbers, driving amplitudes and mixtures (Federrath, 2013a), as well as magnetic field strengths and variations in the polytropic exponent for simulations that employ a polytropic EOS (Federrath & Banerjee, 2015). The Hopkins (2013b) density PDF is defined as


where is the modified Bessel function of the first kind. Equation (2.3) is motivated and explained in detail in Hopkins (2013b). It contains two parameters: 1) the volume-weighted standard deviation of logarithmic density fluctuations , and 2) the intermittency parameter . In the zero-intermittency limit (), Equation (2.3) becomes the lognormal distribution from Equation (1), .

In order to measure the density variance , we fit our simulation density PDFs in a restricted range around the mean (from to , where is the second moment of the density distribution, directly computed by summation over all simulation data points) with Equation (2.3) and determine the best-fit parameter . In agreement with the conclusions drawn in Price et al. (2011) and Hopkins (2013b), we find that the fitted is the same within a few percent as (computed by summation over all simulation grid cells).

3 Results

3.1 Isothermal comparison

Figure 2: The density variance – Mach number relation for isothermal turbulence (approximated by setting ). In order to reach a wide range of Mach numbers to test the relation, we use stirring amplitudes between and (models 1–4 in Table LABEL:tab:sim), leading to Mach numbers between 3 and 12. We fit Equation (2) to the four points and find the proportionality constant with a goodness of fit of R. The fit is shown as a solid line, while cases with and are shown as dashed lines for comparison. Our best fit is consistent with the expectation value for purely solenoidal driving (Federrath et al., 2008, 2010).

The isothermal gas relation has been studied extensively and is therefore a good comparison test for our hydrodynamics code, setup and post-processing methods to determine (see §2). For purely solenoidal forcing of the turbulence, we expect a proportionality value of in Equation (2) (Federrath et al., 2008, 2010). Four simulations were performed at a grid resolution of with to prevent non-isothermal effects brought on by high stirring amplitudes. The four time-averaged points lie between and 0.4, with a fit of and goodness-of-fit parameter R (Fig. 2). Our measurement of spans the expected value, thus our methods produce reasonable results for isothermal turbulence. Now that we have established that our simulation and analysis techniques reproduce previous results for isothermal turbulence, we can now move on to study non-isothermal turbulence in the adiabatic regime.

3.2 Resolution study

We now test the resolution requirements of density variance and Mach number measurements by performing a series of identical simulations at increasing resolutions of , , , and grid cells for . Fig. 3 shows the density variance – Mach number relation for each simulation and simulation time. We see that first, the density variance and Mach number increase and reach a maximum after about (the lower part of the correlation). In this first part of the evolution, the kinetic energy of the gas increases due to the driving until the kinetic energy power spectrum is established (Schmidt et al., 2009). After the kinetic energy and rms velocity have reached a saturated state (see Fig. 1), only the sound speed keeps increasing monotonically, because the dissipated energy heats the gas. This leads to a continuously decreasing Mach number and density variance in the regime of fully developed turbulence (the upper part of the correlation).

Figure 3: Left panel: Evolution of simulations with increasing resolution in density variance – rms Mach number space. The dashed lines represent functions of Equation (2) with and , for comparison. Initially the rms Mach number increases substantially, then decreases smoothly once the turbulence becomes fully developed (at about ), as a result of the continuously increasing internal energy, temperature and sound speed for our ideal gas EOS, Equation (5). Right panel: A magnified region of the left panel, displaying the convergence of statistics. Simulations with grid cells are representative of the converged system.

The higher the numerical grid resolution, the greater the maximum of for , which is seen to converge for . A zoom within the focus area of the density variance – Mach number relation is shown in the right-hand panel of Fig. 3. These points are independent of the initial jump, but are still seen to converge at increased resolution. For resolutions equal to and above , the values are almost identical. Therefore the statistics of density variance and rms Mach number may be approximated well by numerical resolutions cells. This is consistent with the resolution requirements established in Kitsionas et al. (2009), Federrath et al. (2010), Kritsuk et al. (2011), and Federrath (2013a).

Looking closely at Fig. 3, we see that the gas evolves along a curve in the density variance – Mach number plane. This is due to the continuous heating of the gas via the turbulent driving, which lowers the Mach number continuously. The evolution of this curve seems to correlate somewhat with Equation (2), but is steeper than the theoretical prediction for isothermal turbulence. The behaviour of this curve is quantified in detail in the following subsections, §3.3 and §3.4 for and , respectively.

3.3 Diatomic molecular gas:

Figure 4: Evolution of simulations with and (red, small arrows), 200 (purple, normal arrows) and 400 (blue, large arrows). The arrows indicate the direction of time evolution.

In the first case we look at a diatomic equation of state, i.e., . We perform three separate simulations with different stirring amplitudes , and plot their evolutionary curves in Fig. 4. The best fit to Equation (2) from the combined points of all three simulations is with a goodness of fit parameter of R. Note that the increase in sound speed due to gas heating is slow compared to the time it takes to establish a statistically steady state, which is shown in Fig. 1, where we see that the velocity dispersion is fully developed after one crossing time. However, the continuous heating in the fully developed regime of turbulence leads to a continuously increasing sound speed, the consequences of which are discussed in more detail in Section 4.

As in Fig. 3, we see in Fig. 4 that simulations with produce a somewhat steeper relation compared to the isothermal relation (cf. Fig. 2) as time progresses and both and decrease due to the continuous heating of the gas. The different times of each simulation are connected by a line with arrows in Fig. 4, indicating increasing time. Therefore a new functional fit might be more appropriate to describe the behaviour in our non-isothermal, adiabatic turbulence simulations.

The simplest modification to the existing model function, Equation (2), is to allow for variations in the exponent on the Mach number. Our data in Fig. 4 indicate that the exponent is somewhat higher than the standard dependence from Equation (2). Thus we use the following new fit function to determine the exponent :


We do not necessarily expect that the coefficient in this new relation is the same as in Equation (2), but we will test that below. The new function is fitted to the data from each of the three simulations with different driving amplitude and to the combined set of data points. We determine the goodness of fit parameter for the fits to Equations (2) and (7) and compare them. The results are summarised in Table LABEL:y140_fits.

Fit function
0.90 0.96
0.90 0.99
0.92 1.0
All data 0.90 0.99
Table 2: Statistical fit parameters for different functions , for . Notes. denotes the goodness-of-fit parameter. A value of indicates a perfect fit to the given data.

Table LABEL:y140_fits shows that for the individual simulations as well as for the combined data set, the modified power-law function from Equation (7) provides the better fit to the data as quantified by the goodness of fit . The coefficient value is smaller than , but they are formally consistent with representing the same value, and consistent with the -value obtained in our isothermal calculations in §3.1. In fact, our new fit gives a value that is in agreement with the theoretical expectation for the turbulent driving, namely for purely solenoidal driving as applied here.

The exponent , which is fixed to in Equation (2), but allowed to vary in our new fit function, Equation (7), clearly shows that an almost cubic dependence on provides a better fit to the data. We find a best-fit value of in our simulations with , leading to a new form of the density variance – Mach number relation for gas,


This result presents an interesting question: is the density variance – rms Mach number relation of non-isothermal adiabatic turbulence no longer a quadratic relation, compared to the well-studied isothermal case? We will now explore whether the same/similar holds for and then address this questions in the discussion of Section 4, by comparing to a theoretical model of the relation.

3.4 Monatomic gas:

Figure 5: Evolution of simulations with and (red, small arrows), 200 (purple, normal arrows) and 400 (blue, large arrows). The arrows indicate the direction of time evolution.

In the second case we look at a monatomic EOS, Equation (5) with . As in the diatomic case we performed three separate simulations with different driving amplitude , and plot their evolutionary curves in Fig. 5 (again with arrows indicating the continuous time evolution to smaller and smaller and ). A fit with the standard relation, Equation (2), to the combined data set gives with a goodness of fit .

Given the same effect occurs as for , we fit our new power-law function, Equation (7), to data from each of the three simulations and to the combined set of data points, and compare the goodness of fit values with those from the fit to Equation (2). The results are summarised in Table LABEL:y166_fits.

Fit function
0.92 0.99
0.91 1.0
0.89 1.0
All data 0.90 0.99
Table 3: Statistical fit parameters for different functions , for . Notes. denotes the goodness-of-fit parameter. A value of indicates a perfect fit to the given data.

Once again we see a significantly better fit to the power law with exponent . We find that the driving coefficient , as for the case, indicating that the physics of the driving is indeed contained in the value of the parameter, while the fact that we deal with non-isothermal turbulence is reflected in a steeper power-law exponent , compared to the isothermal case. All three simulations fit exponents very close to cubic (), with the combined data for fitting a power law of the form


Figure 6: Exponent in the density variance – Mach number relation for different values of the adiabatic index . The data points show our simulation measurements and the dotted line is a linear fit with .

3.5 Summary of the results

In summary, both the and cases yield turbulent driving coefficients that are all consistent with the theoretical expectation for purely solenoidal driving of the turbulence, and consistent with the -values found for isothermal turbulence ().

The exponent of the relation, however, is significantly steeper with for and for , compared to the isothermal case, where provides the best fit to the data. Thus, we see that the exponent in the density variance – Mach number relation is -dependent.

In order to provide a heuristic relation that describes the behaviour of in adiabatic gases, we fit a linear function to the dependence of on in Fig. 6. It is reasonable that will continuously increase with . Thus, we chose the simplest function to approximate our data from to , i.e., a linear interpolation. The result is . The actual dependence might be somewhat different, but to determine a better fit would require us to measure for a range of values, in much smaller steps . This is beyond the scope of the paper, but we can already provide a new improved functional form of the density variance – Mach number relation that takes the adiabatic index into account. The best-fit function is given by


which is the central result of the paper. Equation (10) naturally simplifies to the well-studied isothermal case () given by Equation (2), but also approximately covers cases with , up to .

4 Discussion

In §3.3–§3.5, we found that the density variance – Mach number relation for adiabatic gases with and respectively, deviates significantly from the isothermal case, Equation (2). We quantified the discrepancy by fitting an alternative function, Equation (7), to the data, with the power-law exponent as a free fit parameter. We find that the power law provides excellent fits, with power-law exponents increasing with from for the isothermal case () to for . A heuristic function was obtained to provide a new relation that takes the dependence on into account, given by Equation (10).

We now compare this results to a recent theoretical model for the density variance – Mach number relation, in order to explain the differences of our adiabatic case to the isothermal and polytropic cases. The detailed derivation of the relation can be found in Molina et al. (2012) and Federrath & Banerjee (2015), where this relation has been explored for magnetized isothermal and polytropic gases respectively, with the latter representing a special case, where the pressure and temperature of the gas are both uniquely related to the density via


We emphasize that this is different from the adiabatic EOS, Equation (5), used here, in that the pressure depends on both density and temperature, . Thus, for any given value of density, the pressure can vary depending on the temperature, while a polytropic EOS will give only a single value of for a given input .

The basic idea of the theoretical model is to relate the density jump in a single shock to the ensemble of shocks/compressions in a turbulent medium. For that purpose, Padoan & Nordlund (2011) and Molina et al. (2012) applied the equations of mass, momentum and energy conservation across a shock,


to derive the density contrast between the pre-shock gas (denoted with index 0 and on the left-hand side of the equations) and the post-shock gas (no index; right-hand side of the shock jump equations). The equation of state, Equation (5), enters through the pressure and specific internal energy in these expressions. Combining these equations leads to the well-known Rankine-Hugoniot shock jump conditions as the solution for the density jump across the shock (Rankine, 1870; Hugoniot, 1887; Shull & Draine, 1987),


Note that we have already introduced the geometrical -parameter, because these shock jump conditions only apply to the plane-parallel component of the shock, parametrized by the parallel component of the sonic Mach number (Molina et al., 2012; Federrath & Banerjee, 2015). Following the detailed derivation in Molina et al. (2012), Equation (15) just needs to be inserted into the general expression for the ensemble of such shocks,


which leads to the density variance – Mach number relation,


We immediately see that this new -dependent density variance – Mach number relation reduces to the isothermal case, Equation (2), if we set . In the adiabatic case, however, , which leads to our theoretical prediction as a function of , given by Equation (17).

Figure 7: Combined density variance – Mach number relation plot, showing all our simulation data with (red crosses), (green boxes), and (blue diamonds). The dashed lines are our theoretical prediction given by Equation (17) with the respective values of (in the same colour as the simulation data and labelled on each curve).

Fig. 7 shows the theoretical prediction given by Equation (17) for different values of together with the isothermal solution and together with our simulation data for and . We see that the new relation qualitatively follows the trend of a slightly steeper rise with increasing for low Mach number. It also predicts that at high Mach number, the density variance saturates at lower for increasing . This is reasonable, because the density jump across shocks reduces significantly with increasing as derived in Equation (15). This regime needs to be tested in follow-up simulations that reach higher Mach numbers. However, the problem is that the adiabatic heating increases with increasing Mach, such that it quickly counteracts the effect of an increased driving amplitude.

Despite these reasonable qualitative trends produced by our new theoretical relation, Equation (17), we see that the actual simulation data with still follow a somewhat steeper curve at low Mach number, in the plane. We speculate that this discrepancy arises, because the theoretical model does not contain any information about the temporal evolution of the gas, in particular about its temperature changes along the evolutionary curve.

However, we can qualitatively argue that any shock will immediately experience the temperature and pressure increase associated with the adiabatic compression. This will reduce the density jump significantly, such that the density variance will be smaller with increasing almost immediately when these shocks are about to form, (e.g., the theoretical limit for is ). Thus, shocks in high- gas will be significantly reduced and so will be the statistical variance of density fluctuations, . The important point here is that this process almost instantaneously reduces the density variance.

At the same time, each shock dissipates energy, locally increasing the temperature and internal energy of the ideal adiabatic gas. This increases the sound speed, but at first only locally in the shocks, which leads to a decreasing pre-shock sonic Mach number over time, resulting in the time evolution (shown as arrows) in Figs. 4 and 5. We can thus qualitatively understand the time dependence and resulting relations for . The relation is steeper, because responds almost instantaneously to the local pressure and temperature increase in shocks, while the Mach number reduction is delayed, because the sound speed increases only in the post-shock gas, while our theoretical Equation (17) is based on the large-scale, volume-weighted pre-shock Mach number.

Figure 8: Slices of the normalized density (top left), temperature (top right), pressure (bottom left), and specific entropy (bottom right) at for , , and a numerical resolution of . It is evident that gas at a given density can have a wide range of temperatures and pressures, unlike a polytropic EOS where and are unique functions of the density only.

In order to substantiate our findings, we show density, temperature, pressure, and entropy slices of our highest-resolution simulation with grid cells and adiabatic in Fig. 8. We see two important points. First, the pressure and temperature are not unique functions of the density, but for a given density, the gas can have a range of temperatures and pressures, as implied by Equation (5). This is substantiated by the entropy slice shown in the bottom right panel of Figure 8, which is not uniform, demonstrating that the gas is not barotropic, but that the pressure depends on density and temperature. There is clearly viscous heating, which is primarily due to shocks in the regime, while turbulent dissipation (eddy viscosity) becomes a more important heating mechanism when . Quantifying both contributions is beyond the scope of this paper. Second, the adiabatic heating primarily occurs in the post-shock gas. The rise of the internal energy does not immediately reduce the global post-shock Mach number, but slowly diffuses to large scales, before it affects , leading to the steeper-than-isothermal relations we found for .

Figure 9: Density–temperature correlation PDFs for our adiabatic turbulence simulations with , and . The left-hand panel shows the results at , while the right-hand panel is for . The distributions are spread by more than an order of magnitude for typical gas densities, around the average correlations (shown as white lines) with . The continuous heating of the gas indicated by the overall rise in temperature at compared to primarily occurs in the post-shock gas, while the Mach number entering Equation (17) only applies to the pre-shock gas. This produces the steeper dependence of on that we find in our central result, Equation (10).

Finally, Fig. 9 shows density–temperature correlation probability density functions (PDFs). It is evident that for any given density, there is a wide range of temperatures and that the heating indeed primarily occurs in the densest gas, i.e., in the post-shock regions. Also note the continuous increase in the overall temperature of the gas between (left-hand panel) and (right-hand panel), which leads to the slowly decreasing over time.

5 Summary and conclusion

We performed hydrodynamical simulations of supersonic and subsonic turbulence, employing an ideal equation of state (EOS) with adiabatic indices (nearly isothermal), (diatomic molecular gas), and (monatomic gas). Section 3 provided a detailed analysis of the density variance – Mach number relation, , which is a key ingredient for theoretical models of the star formation rate and the initial mass function. Unlike previous studies of purely isothermal and polytropic turbulence, we find that an ideal gas EOS leads to a steeper dependence of the density variance on the rms sonic Mach number . We find a new combined approximate relation of the form given by Equation (10) for low Mach numbers, , which reduces to the well known isothermal solution for the special case , but also covers cases . We argue that the steeper-than-isothermal dependence for is a result of the local heating of the gas in post-shock regions, with the global sonic pre-shock Mach number in the relation being affected later in the evolution. This is because the turbulent driving keeps increasing the internal energy reservoir, leading to an ever increasing global sound speed in adiabatic gases. This is in stark contrast to isothermal and polytropic gases, where the sonic Mach number reaches a statistical steady state rather than continuously decreasing.

We derived a theoretical model, Equation (17), for the relation in Section 4, which is based on the Rankine-Hugoniot shock jump conditions and provides reasonable fits to all our data. It furthermore predicts a saturation of for , which is not yet in reach by numerical simulations. Such a saturation is reasonable for , given the fact that adiabatic shocks always have a finite jump in density, while isothermal shocks can theoretically have an infinitely large jump in density across the shock. Both Equation (10) for and Equation (17) for naturally simplify to the standard isothermal relation, Equation (2) for .

We conclude that changes in the adiabatic exponent can introduce important modifications in the density variance – Mach number relation and we provide an approximation of that behaviour in Equation (10). However, we emphasize that the real ISM is a mixture of atomic and molecular phases and that the effective EOS is determined by a complex balance of heating and cooling processes, which in turn depend on the chemical evolution and exposure to interstellar and local radiation fields (e.g., from massive stars). Thus, our systematic study sheds some light on the dependence of turbulent density fluctuations on the thermodynamics and composition of interstellar gas and more detailed studies including realistic heating and cooling are required to make further progress.

We hope that this work provides a more general understanding of the density variance – Mach number relation in the ISM. This is especially true for the warm, atomic part of the ISM, where the gas is clearly not isothermal and may be approximated with an adiabatic equation of state with , which was not covered by previous density variance – Mach number relations in the literature. Our new relations in this paper attempt to cover this regime and do seem to approximately do so, as tested with the set of simulations presented here. We hope that the new relations will provide useful generalisations of the previous (purely isothermal) relations, which are key ingredients for theoretical models of star formation (see e.g., the review by Padoan et al., 2014).


We thank the referee, Sam Falle, for his timely and constructive report, which improved the paper significantly. C.F. acknowledges funding provided by the Australian Research Council’s Discovery Projects (grants DP130102078 and DP150104329). We gratefully acknowledge the Jülich Supercomputing Centre (grant hhd20), the Leibniz Rechenzentrum and the Gauss Centre for Supercomputing (grant pr32lo), the Partnership for Advanced Computing in Europe (PRACE grant pr89mu), and the National Computational Infrastructure (grants n72 and ek9), supported by the Australian Government. This work was further supported by resources provided by the Pawsey Supercomputing Centre with funding from the Australian Government and the Government of Western Australia.


  • Audit & Hennebelle (2005) Audit, E., & Hennebelle, P. 2005, A&A, 433, 1
  • Audit & Hennebelle (2010) —. 2010, A&A, 511, A76+
  • Beetz et al. (2008) Beetz, C., Schwarz, C., Dreher, J., & Grauer, R. 2008, Physics Letters A, 372, 3037
  • Bland-Hawthorn et al. (2007) Bland-Hawthorn, J., Sutherland, R., Agertz, O., & Moore, B. 2007, ApJL, 670, L109
  • Burkhart & Lazarian (2012) Burkhart, B., & Lazarian, A. 2012, ApJ, 755, L19
  • Chabrier et al. (2014) Chabrier, G., Hennebelle, P., & Charlot, S. 2014, ApJ, 796, 75
  • Cooper et al. (2009) Cooper, J., Bicknell, G., Sutherland, R., & Bland-Hawthorn, J. 2009, ApJ, 703, 330
  • Elmegreen (2008) Elmegreen, B. G. 2008, ApJ, 672, 1006
  • Elmegreen & Scalo (2004) Elmegreen, B. G., & Scalo, J. 2004, ARAA, 42, 211
  • Eswaran & Pope (1988) Eswaran, V., & Pope, S. B. 1988, CF, 16, 257
  • Federrath (2013a) Federrath, C. 2013a, MNRAS, 436, 1245
  • Federrath (2013b) —. 2013b, MNRAS, 436, 3167
  • Federrath & Banerjee (2015) Federrath, C., & Banerjee, S. 2015, MNRAS, 448, 3297
  • Federrath & Klessen (2012) Federrath, C., & Klessen, R. S. 2012, ApJ, 761, 156
  • Federrath et al. (2008) Federrath, C., Klessen, R. S., & Schmidt, W. 2008, ApJL, 688, L79
  • Federrath et al. (2010) Federrath, C., Roman-Duval, J., Klessen, R. S., Schmidt, W., & Mac Low, M. 2010, A&A, 512, A81
  • Federrath et al. (2010) Federrath, C., Roman-Duval, J., Klessen, R. S., Schmidt, W., & Mac Low, M.-M. 2010, A&A, 512, A81
  • Ferrière (2001) Ferrière, K. M. 2001, Rev. Mod. Phys., 73, 1031
  • Fischera & Dopita (2005) Fischera, J., & Dopita, M. 2005, ApJ, 619, 340
  • Gazol & Kim (2013) Gazol, A., & Kim, J. 2013, ApJ, 765, 49
  • Glover et al. (2010) Glover, S. C. O., Federrath, C., Mac Low, M., & Klessen, R. S. 2010, MNRAS, 404, 2
  • Hennebelle & Audit (2007) Hennebelle, P., & Audit, E. 2007, A&A, 465, 431
  • Hennebelle & Chabrier (2008) Hennebelle, P., & Chabrier, G. 2008, ApJ, 684, 395
  • Hennebelle & Chabrier (2011) —. 2011, ApJ, 743, L29
  • Hopkins (2013a) Hopkins, P. F. 2013a, MNRAS, 430, 1653
  • Hopkins (2013b) —. 2013b, MNRAS, 430, 1880
  • Hugoniot (1887) Hugoniot, P. H. 1887, Journal de l’Ecole Polytechnique, 57, 3
  • Kitsionas et al. (2009) Kitsionas, S., Federrath, C., Klessen, R. S., et al. 2009, A&A, 508, 541
  • Konstandin et al. (2012) Konstandin, L., Girichidis, P., Federrath, C., & Klessen, R. S. 2012, ApJ, 761, 149
  • Kritsuk et al. (2007) Kritsuk, A., Norman, M., Padoan, P., & Wagner, R. 2007, ApJ, 665, 416
  • Kritsuk & Norman (2002) Kritsuk, A. G., & Norman, M. L. 2002, ApJ, 569, L127
  • Kritsuk et al. (2011) Kritsuk, A. G., Nordlund, Å., Collins, D., et al. 2011, ApJ, 737, 13
  • Krumholz (2014) Krumholz, M. R. 2014, Physics Reports, in press (arXiv:1402.0867)
  • Krumholz & McKee (2005) Krumholz, M. R., & McKee, C. F. 2005, ApJ, 630, 250
  • Mac Low & Klessen (2004) Mac Low, M.-M., & Klessen, R. S. 2004, Rev. Mod. Phys., 76, 125
  • McKee & Ostriker (2007) McKee, C. F., & Ostriker, E. C. 2007, ARAA, 45, 565
  • Molina et al. (2012) Molina, F., Glover, S., Federrath, C., & Klessen, R. S. 2012, MNRAS, 423, 2680
  • Padoan et al. (2014) Padoan, P., Federrath, C., Chabrier, G., et al. 2014, Protostars and Planets VI, 77
  • Padoan et al. (1997) Padoan, P., Jones, B., & Nordlund, A. 1997, ApJ, 474, 730
  • Padoan & Nordlund (2011) Padoan, P., & Nordlund, Å. 2011, ApJ, 730, 40
  • Passot & Vazquez-Semadeni (1998) Passot, T., & Vazquez-Semadeni, E. 1998, Phys. Rev. E, 58, 4501
  • Price & Federrath (2010) Price, D. J., & Federrath, C. 2010, MNRAS, 406, 1659
  • Price et al. (2011) Price, D. J., Federrath, C., & Brunt, C. M. 2011, ApJ, 727, L21
  • Rankine (1870) Rankine, W. J. M. 1870, Royal Society of London Philosophical Transactions Series I, 160, 277
  • Schmidt et al. (2009) Schmidt, W., Federrath, C., Hupp, M., Kern, S., & Niemeyer, J. C. 2009, A&A, 494, 127
  • Seifried et al. (2011) Seifried, D., Banerjee, R., Klessen, R. S., Duffin, D., & Pudritz, R. E. 2011, MNRAS, 417, 1054
  • Seon (2012) Seon, K.-I. 2012, ApJ, 761, L17
  • Shull & Draine (1987) Shull, J. M., & Draine, B. T. 1987, in Astrophysics and Space Science Library, Vol. 134, Interstellar Processes, ed. D. J. Hollenbach & H. A. Thronson, Jr., 283–319
  • Sutherland (2010) Sutherland, R. 2010, Ap&SS, 327, 173
  • Vázquez-Semadeni (1994) Vázquez-Semadeni, E. 1994, ApJ, 423, 681
  • Wada & Norman (2001) Wada, K., & Norman, C. A. 2001, ApJ, 547, 172
  • Wagner & Bicknell (2011) Wagner, A., & Bicknell, G. 2011, ApJ, 728, 29
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description