Issue 
A&A
Volume 528, April 2011



Article Number  A145  
Number of page(s)  8  
Section  Interstellar and circumstellar matter  
DOI  https://doi.org/10.1051/00046361/201015661  
Published online  17 March 2011 
Vorticity production through rotation, shear, and baroclinicity
^{1} Nordita, AlbaNova University
Center, Roslagstullsbacken
23, SE10691 Stockholm, Sweden
email: fadiesis@gmail.com
^{2} Department of Astronomy, AlbaNova
University Center, Stockholm University, 10691
Stockholm,
Sweden
Received: 31 August 2010
Accepted: 14 February 2011
Context. In the absence of rotation and shear, and under the assumption of constant temperature or specific entropy, purely potential forcing by localized expansion waves is known to produce irrotational flows that have no vorticity.
Aims. Here we study the production of vorticity under idealized conditions when there is rotation, shear, or baroclinicity, to address the problem of vorticity generation in the interstellar medium in a systematic fashion.
Methods. We use threedimensional periodic box numerical simulations to investigate the various effects in isolation.
Results. We find that for slow rotation, vorticity production in an isothermal gas is small in the sense that the ratio of the rootmeansquare values of vorticity and velocity is small compared with the wavenumber of the energycarrying motions. For Coriolis numbers above a certain level, vorticity production saturates at a value where the aforementioned ratio becomes comparable with the wavenumber of the energycarrying motions. Shear also raises the vorticity production, but no saturation is found. When the assumption of isothermality is dropped, there is significant vorticity production by the baroclinic term once the turbulence becomes supersonic. In galaxies, shear and rotation are estimated to be insufficient to produce significant amounts of vorticity, leaving therefore only the baroclinic term as the most favorable candidate. We also demonstrate vorticity production visually as a result of colliding shock fronts.
Key words: turbulence / dynamo / magnetohydrodynamics (MHD) / ISM: bubbles / galaxies: magnetic fields / galaxies: ISM
© ESO, 2011
1. Introduction
Turbulence in the interstellar medium (ISM) is believed to be driven by supernova explosions. Such events inject sufficient amounts of energy to sustain turbulence with rms velocities of ~10 km s^{1} and correlation lengths of up to 100 pc (Beck et al. 1996). Simulations of such events can be computationally quite demanding, because the bulk motions tend to be supersonic and the flows involve strong shocks in the vicinity of individual explosion sites, as was seen early on in twodimensional simulations (Rosen & Bregman 1995). Nevertheless, such simulations are able to reproduce a number of physical phenomena such as the observed volume fractions of hot, warm, and cold gas (Rosen et al. 1996; Korpi et al. 1999a), the statistics of pressure fluctuations (Mac Low et al. 2005), the effects of the magnetic field (de Avillez & Breitschwerdt 2005), and even dynamo action (Gressel et al. 2008; Gissinger et al. 2009; Hanasz et al. 2009). These simulations tend to show the development of significant amounts of vorticity, which is at first glance surprising. Indeed, each supernova drives the gas radially outward and can roughly be described by radial expansion waves. In such a description, turbulence is forced by the gradient of a potential that consists of a timedependent spherical blob at random locations. Obviously, such a forcing is irrotational, so no vorticity is produced.
Earlier work of Mee & Brandenburg (2006) showed that under isothermal conditions only the viscous force can produce vorticity and that this becomes negligible in the limit of large Reynolds numbers or small viscosity. In principle, vorticity can also be amplified akin to the dynamo effect by the ∇ × (u × ω) term, which is analogous to the induction term in dynamo theory, where ω plays the role of the magnetic field. However, neither this nor the viscosity effect were found to operate – even at numerical resolutions of up to 512^{3} meshpoints. This disagreed with subsequent simulations by Federrath et al. (2010), who solved the isothermal inviscid Euler equations with irrotational forcing using the Flash Code. They found significant vorticity generation in proximity to shocks where some kind of effective numerical viscosity must have acted.
Given that under isothermal conditions, only viscosity can lead to vorticity production, one must ask whether numerical viscosity or effective viscosity needed to stabilize numerical codes might have contributed to the production of vorticity in some of the earlier works. Indeed, it is possible that the directional operator splitting used in the Flash Code may have been responsible for spurious vorticity generation in the work of Federrath et al. (2010; Rosner, priv. comm.). On the other hand, when cooling and heating functions are included to perform more realistic simulations of the ISM, vorticity could be produced by the baroclinic term. Furthermore, even in the isothermal case, in which the baroclinic term vanishes, vorticity could be produced if there is rotation and/or shear.
The baroclinic term results from taking the curl of the pressure gradient term and is proportional to the cross product of the gradients of pressure and density. This term can play an important role when the assumptions of isothermality or adiabaticity are relaxed. Indeed, the baroclinic term can also be written as the cross product of the gradients of entropy and temperature. This formulation highlights the need for nonideal effects, because in the absence of any other heating or cooling mechanisms, the entropy is just driven by viscosity. Again, it is not obvious that in the absence of additional heating and cooling much vorticity can be produced. On the other hand, it is clear that viscous heating must be significant even in the limit of vanishing viscosity, because the velocity gradients can be very large, especially in shocks. Of course, the assumption about additional heating and cooling is not realistic for the interstellar medium and will need to be relaxed. Finally, there are the effects of rotation and shear, that can contribute to the production of vorticity even in the absence of baroclinicity.
The goal of this paper is to study the relative importance of the individual effects that contribute to vorticity production. It is then advantageous to restrict oneself to simplifying conditions that allow one to identify the governing effects. An important simplification is the restriction to weakly supersonic conditions so that shocks and other sharp structures can still be resolved with just a uniform and constant viscosity. We also neglect the effects of stratification which can only indirectly contribute to vorticity production. In fact, a constant gravitational acceleration drops out when taking the curl. Only in the nonisothermal and nonisentropic case can gravity contribute to vorticity production by enhancing the effect of the baroclinic term. We begin with a preliminary discussion and a qualitative analysis of the important terms in the vorticity equation.
2. Preliminary considerations
We recall that in the absence of baroclinicity, rotation, and shear, the curl of the evolution equation of the velocity is given by (see, e.g., Mee & Brandenburg 2006) (1)where ν is the kinematic viscosity (assumed constant) and G_{i} = 2S_{ij}∇_{j}lnρ is a part of the viscous force that has nonvanishing curl even when the flow is purely irrotational. Here, (2)is the traceless rate of strain matrix, and commas denote partial differentiation. The G term breaks the formal analogy with the induction equation. It is convenient to express the resulting rms vorticity in terms of the typical wavenumber k_{ω} of vortical structures which we define as (3)We monitor the ratio k_{ω}/k_{f}, where k_{f} is the adopted nominal forcing wavenumber. In Mee & Brandenburg (2006), the resulting vorticity, expressed in terms of the ratio k_{ω}/k_{f}, was found to be zero within error bars. This result is compatible with the idea that the ν∇ × G term in Eq. (1) is insignificant for vorticity production. By contrast, in vortical turbulence and at moderate values of the Reynolds number, k_{ω}/k_{f} is found to be of the order of unity (Brandenburg 2001), although one should expect a mild increase proportional to the square root of the Reynolds number as this number increases.
2.1. Rotation
Rotation leads to the addition of the Coriolis force, 2Ω × u, in the evolution equation for the velocity. Taking the curl, we obtain the vorticity Eq. (1) with two additional terms, both proportional to Ω, so we have (4)where the dots denote the other terms in Eq. (1) that we discussed already. In order to estimate the production of vorticity, one could derive an evolution equation for the enstrophy density, , by multiplying the righthand side of Eqs. (1) and (4) by ω, and use a closure assumption for the resulting triple correlations. However, it is then difficult to obtain a useful prediction for ω_{rms}, because the righthand side of such an equation would necessarily be proportional to ω and would therefore vanish, unless ω_{rms} was different from zero to begin with. Instead, we estimate ω_{rms} by computing the rms value of ∂ω/∂t and replacing it by ω_{rms}/τ_{Ω}, where τ_{Ω} is a typical time scale of the problem. This leads to (5)where and ∇_{ ∥ } denote derivatives in the directions perpendicular and parallel to the rotation axis and is the velocity vector perpendicular to the rotation axis. Using Cartesian coordinates where Ω points in the z direction, we have (6)We expect τ_{Ω} to be comparable to the turnover time, τ = (u_{rms}k_{f})^{1}. We expect the rms values of the velocity derivative term in Eq. (6) to be comparable to the rms velocity and some inverse length scale. Typically, one would expect it to be proportional to u_{rms}k_{f}, although, again, there can be an additional dependence on the square root of the Reynolds number. However, for fixed Reynolds number, and not too rapid rotation, we expect ω_{rms} to increase linearly with the Coriolis number, i.e., (7)Thus, we expect k_{ω}/k_{f} = St_{Ω} Co, where we have defined an effective rotational Strouhal number, (8)We regard this as a fit parameter that will emerge as a result of the simulations. We have here introduced the quantity , where is given by the ratio of the velocity gradient terms divided by u_{rms}k_{f}. However, for larger values of Co there may be a departure from a linear dependence between k_{ω}/k_{f} and Co. (We note that, apart from a possible 4π factor, the Coriolis number is just the inverse Rossby number.) One aim of this paper is therefore to verify this dependence from simulations and to determine empirically the value of τ_{Ω}.
2.2. Shear
In the presence of linear shear with , the evolution equation for the departure from the mean shear attains additional terms, . This implies a dependence of ω_{rms} on S, analogous to the Ω dependence discussed above. In components form, this means that (9)which is quite similar to Eq. (5), except that in the penultimate term in angular brackets the indices are now interchanged, i.e. we now have u_{z,y} instead of u_{y,z}. In analogy to τ_{Ω}, we define τ_{S} as a typical time scale of the problem and we expect it to be again related to the turnover time τ. The O(xu′′) denotes the presence of additional terms that are proportional to x and to second derivatives of u. However, when adopting the shearing box approximation with shearingperiodic boundaries (Goldreich & LyndenBell 1965; Wisdom & Tremaine 1988), each point in the xy plane is statistically equivalent. We would therefore not expect there to be a systematic x dependence, which corresponds to the assumption of Galilean invariance that is sometimes used in the study of turbulent transport coefficients in linear shear flows (Sridhar & Subramanian 2009). We will postpone the possibility of additional terms until later. Since we expect τ_{S} to be comparable to τ = (u_{rms}k_{f})^{1}, the rms vorticity should be proportional to the shear parameter, (10)although for large values of  Sh  we may expect departures from a linear dependence. Determining this dependence is another aim of this paper. Again, a linear dependence is characterized by the values of τ_{S} and , where, in analogy with the previous case with rotation, the ratio is given by the derivative term in Eq. (9), normalized by u_{rms}k_{f}. A convenient nondimensional measure of the value of is what we call the shear Strouhal number, (11)which can be determined provided there is a range in Sh over which ω_{rms} increases linearly with Sh.
The study of vorticity production by rotation and shear is quite independent of thermodynamics and can in principle be studied even in the incompressible case. However, in the present paper we study this effect in the weakly compressible case of low Mach numbers and under the assumption of an isothermal equation of state, where the baroclinic term vanishes.
2.3. Baroclinicity
As mentioned in the introduction, the baroclinic term, proportional to ∇ρ × ∇p, emerges when taking the curl of the pressure gradient term, ρ^{1}∇p. This term can also be written as (12)where h and s are specific enthalpy and specific entropy, respectively, and T is the temperature. Thus, we have (13)In order to study the effect of the baroclinic term, it is useful to look at the dependence of the mean angle θ between the gradients of s and T, defined via (14)An important aspect is then to study first the dependence of the rms values of the gradients of s and T. We can do this by looking at a onedimensional model where, of course, θ = 0.
Next, we need to determine θ from threedimensional simulations. The hope is then that we can express baroclinic vorticity production in the form (15)On dimensional grounds we expect the product of (∇T)_{rms} and (∇s)_{rms} to be of the order of , and so a possible ansatz would be (16)where we have subsumed the scalings of (∇T)_{rms} and (∇s)_{rms} in that of an effective baroclinic Strouhal number .
An important issue is the fact that viscous heating leads to a continuous increases of the temperature. As a result, the sound speed changes and it becomes then impossible to study the behavior of the system in a steady state. In order to avoid this inconvenience, we add a volume cooling term that is nonvanishing when the local sound speed c_{s} is different from a given target value, c_{s0}. Thus, in the presence of finite thermal diffusivity χ, and with a cooling term governed by a cooling time τ_{cool}, our entropy equation takes the form (17)where c_{s} is the adiabatic sound speed. We assume a perfect gas so that , where γ = c_{p}/c_{v} = 5/3 for a monatomic gas, and c_{p} and c_{v} are the specific heats at constant pressure and constant volume, respectively. The value of τ_{cool} can have an influence on the results, so we need to consider different values. We express τ_{cool} in terms of c_{s0} and k_{f}, and define the nondimensional quantity St_{cool} = τ_{cool}c_{s0}k_{f}.
3. The model
In this paper we solve the continuity equation for the density ρ, (18)together with the momentum equation for the velocity u, (19)where is the advection operator with respect to the sum of turbulent flow u and laminar shear flow , p is the pressure, φ is the forcing potential, and (20)is the viscous force, where S was defined in Eq. (2). The forcing potential is given by (21)where x = (x,y,z) is the position vector, x_{f}(t) is the random forcing position that changes abruptly after a time interval Δt, R is the radius of the Gaussian, and N is a nondimensional factor proportional to Δt^{ − 1/2}. This ensures that the amplitude of the correlation function of φ is independent of Δt. Thus, we choose . Since N is nondimensional, the prefactor φ_{0} has the same dimension as φ, which is that of velocity squared. We consider two forms for the time dependence of x_{f}. First, we take x_{f} such that the forcing is δcorrelated in time. In that case, Δt is equal to the length of the time step δt. Alternatively, we choose a finite forcing time δt_{force} that defines the interval during which x_{f} remains constant, after which the forcing changes again abruptly. Thus, (22)is equal to δt in the δcorrelated case or equal to δt_{force} in the case of finite correlation time.
The work of Mee & Brandenburg (2006) showed that the peak of the energy spectrum depends on the radius R of the Gaussian. Indeed, the Fourier transform of exp(− r^{2}/R^{2}) is also a Gaussian with , where (23)In the following we use this as our definition of k_{f} and check a posteriori that this is close to the position of the peak of the energy spectrum. In the following, we characterize our simulations in terms of the ratio k_{f}/k_{1}, and consider values between 2 and 10.
We use the Pencil Code^{1}, which is a nonconservative, highorder, finitedifference code (sixth order in space and third order in time) for solving the compressible hydrodynamic and hydromagnetic equations. We adopt nondimensional variables by measuring speed in units of a reference sound speed, c_{s0}, and length in units of 1/k_{1}, where k_{1} is the smallest wavenumber in the periodic domain. This implies that the nondimensional size of the domain is (2π)^{3}.
In order to study the effects of rotation and shear, we ignore entropy effects and restrict ourselves to an isothermal equation of state with constant sound speed c_{s}. This means that ρ^{1}∇p reduces to , which has vanishing curl. Here, is the relevant enthalpy in the isothermal case. On the other hand, in order to study the effects of baroclinicity, we do need to allow the entropy to vary, so we also need to solve Eq. (17), and study the dependence of k_{ω}/k_{f} on the Mach number, (24)In order to characterize the degree of turbulence, we define the Reynolds number based on the energycarrying scale, corresponding to the typical wavenumber where the spectrum peaks, i.e. (25)For vortical turbulence, this definition is known to be a good measure of the ratio of the resulting turbulent viscosity divided by the molecular diffusivity (Yousef et al. 2003). The two numbers, Ma and Re, can be varied by changing ν and/or the strength of the forcing. In all cases we use χ = ν. Another input parameter is the forcing Strouhal number (26)which is zero for δcorrelated forcing and equal to about 0.3 in cases with finite correlation time. These are also the values used by Mee & Brandenburg (2006).
In the following we also consider kinetic energy and enstrophy spectra, E_{K}(k) and E_{ω}(k), respectively. They are normalized such that (Lesieur 1990) (27)where and are kinetic energy and enstrophy, respectively. For comparison we also consider spectra of enthalpy, E_{h}(k), which are normalized such that .
Throughout this paper we assume periodic boundary conditions, except that in the presence of shear we employ shearingperiodic boundary conditions where the x direction is periodic with respect to positions in y that shift with time, i.e. (28)where f represents any one of our four dependent variables (u,ρ). This boundary condition was first proposed by Goldreich & LyndenBell (1965) and has been routinely used in local simulations of accretion disk turbulence (Hawley et al. 1995). Note, however, that recent work of Regev & Umurhan (2008) and Bodo et al. (2008) called attention to the possibility of problems with the shearing sheet approximation when the size of the perturbations is large. In somewhat weaker form, this problem also applies to a nonshearing periodic box. Indeed, we shall kept this in mind when interpreting some of the results presented below.
Fig. 1 Timeaveraged kinetic energy and enthalpy spectra for two values of the Coriolis number for Re = 25 and St_{force} = 0.4. The two straight lines give the slopes − 2 and − 3, respectively. In both cases we have k_{f}/k_{1} = 4. 
4. Results
We begin by studying the effect of rotation. In Fig. 1 we plot timeaveraged kinetic energy and enthalpy spectra, E_{K}(k) and E_{h}(k), respectively. Note that rotation has a tendency to move the peak of E_{K}(k) to the left of the nominal value of k_{f}. However, at the Reynolds number of 25 shown here, there is no inertial range, but in all cases, the energy spectra show a clear viscous dissipation range, suggesting that these runs are sufficiently well resolved. At somewhat larger Reynolds number or smaller forcing wavenumber, earlier work of Mee & Brandenburg (2006) began to show a short k^{2} subrange. Such a slope is predicted for shock turbulence (Kadomtsev & Petviashvili 1973), and it has also been seen in the irrotational component of transonic turbulence (Porter et al. 1998).
In Fig. 2 we plot the dependence of k_{ω}/k_{f} on Co for three groups of runs: group 1 with Re = 15, k_{f}/k_{1} = 10, St_{force} = 0.3; group 2 with Re between 25 and 130, k_{f}/k_{1} = 4, St_{force} = 0.4; and group 3 with Re = 30 and 150, k_{f}/k_{1} = 2, St_{force} = 0. In all cases we use 128^{3} mesh points, average the results over at least 200 turnover times and, in some cases, even several thousand turnover times. It turns out that for St_{force} ≠ 0 a linear relationship between k_{ω}/k_{f} and Co is a good approximation for , where k_{ω}/k_{f} ≈ 0.03 Co, i.e. St_{Ω} = 0.03. Furthermore, we see from Table 1 that the normalized velocity derivative terms are all about 0.5, so the root of the sum of their squares is slightly larger than unity, corresponding to . For Co > 10 the value of k_{ω}/k_{f} seems to saturate at about unity.
Fig. 2 Dependence of k_{ω}/k_{f} on Co for three groups of runs: group 1 with Re = 15, k_{f}/k_{1} = 10, St_{force} = 0.3; group 2 with Re between 25 and 130, k_{f}/k_{1} = 4, St_{force} = 0.4; and group 3 with Re = 30 and 150, k_{f}/k_{1} = 2, St_{force} = 0. 
Fig. 3 Dependence of k_{ω}/k_{f} on Re for Co ≈ 1, k_{f}/k_{1} = 2, and St_{force} = 0 (δcorrelated forcing). 
Rootmeansquared values of components of the velocity derivative tensor, normalized by u_{rms}k_{f}, as well as the three diagonal components of the ⟨ u_{i}u_{j} ⟩ tensor for 4 values of Co.
Fig. 4 Timeaveraged enstrophy spectra, E_{ω}(k), compared with k^{2}E_{K}(k), for Re = 25, St_{force} = 0.4, and three values of the Coriolis number. The curves of k^{2}E_{K}(k) are close together and overlap for Co = 0.01 (dotted) and 0.15 (dashed), so it becomes a single dashdotted line. The k^{3} slope is shown for comparison. In all three cases we have k_{f}/k_{1} = 4. 
Fig. 5 Dependence of k_{ω}/k_{f} on Sh for Re ≈ 40 and k_{f}/k_{1} = 2 and δcorrelated forcing (St_{force} = 0). Different resolutions are shown to give similar results. At small values of  Sh  , comparisons with St_{force} = 0.3 (keeping k_{f}/k_{1} = 2) or k_{f}/k_{1} = 10 (keeping St_{force} = 0) are also shown. 
A similar result is also found for St_{force} = 0, except that there remains a spurious contamination of vorticity even for small values of Co, a limit in which we expect to observe no vorticity production. By varying the value of Re, while keeping Co ≈ 1 fixed, we see that k_{ω}/k_{f} asymptotes to zero for sufficiently small values of Re; see Fig. 3. This suggests that there can easily be spurious vorticity generation, possibly due to marginally sufficient resolution. The possibility of spurious vorticity is indeed verified by Fig. 4, where we compare enstrophy spectra, E_{ω}(k), with k^{2}E_{K}(k). Note that for large values of Co, the enstrophy spectrum decays like k^{3}. However, for smaller values of Co the level of enstrophy at the mesh scale remains approximately unchanged and is thus responsible for the spurious vorticity found above for small values of Co and not too small values of Re. Nevertheless for larger values of Co, the production of vorticity is an obvious effect of rotation in an otherwise potential velocity field, and it is most pronounced at large length scales, as can also be seen in Fig. 4.
Next, we study the dependence of the ratio k_{ω}/k_{f} on shear; see Fig. 5. We use a resolution of 64^{3} or 128^{3} mesh points, average the results over at least 200 turnover times and, in cases of lower resolution, over several thousand turnover times. It turns out that in the presence of shear, some level of helicity production can never be avoided – even in the limit of small Sh. Again, this appears spurious and demonstrates the general sensitivity of vorticity generation on resolution effects. An additional problem is of course the finite size of the shearing box (Regev & Umurhan 2008; Bodo et al. 2008), which may be responsible for spurious vorticity generation. On the other hand, there is vorticity generation even for large scaleseparation ratios, k_{f}/k_{1} = 10; see the dashdotted line in Fig. 5. This suggests the possibility of a more general problem that would not go away even in the limit of small eddies and small values of Sh . Nevertheless, there is a clear rise of k_{ω}/k_{f} when Sh > 0.1, which is in agreement with our expectations outlined in Sect. 2.2. However, the slope in this relationship is rather steep, St_{S} ≈ 6. The velocity derivative terms are only slightly larger than in the case with rotation, corresponding to ; see also Table 2. Tentatively, this suggests that for comparable values of Co and Sh, τ_{S} ≫ τ_{Ω}. On the other hand, given that even for small values of Sh there is spurious vorticity generation, we cannot be certain that the results are reliable for larger ones either. The case with shear must therefore remain somewhat inconclusive.
Finally, we consider the possibility of vorticity generation by the baroclinic term. In a preparatory step we study first the dependence of the product (∇T)_{rms}(∇s)_{rms} on both Ma and Re in a onedimensional model. In all cases we vary the strength of the forcing amplitude in the range for different values of viscosity and cooling time. As we increase the value of φ_{0}, the Reynolds number increases for a given value of the viscosity. For small values of φ_{0}, the Mach number also increases linearly, where the ratio of Ma/Re increases with increasing viscosity. However, for larger values of φ_{0} there is saturation and Ma no longer increase with φ_{0}.
Furthermore, in the range where Ma still increases linearly with φ_{0}, the rms value of the entropy gradient increases, but it also saturates when Ma saturates. The rms value of the temperature gradient, however, decreases with increasing values of φ_{0}, but this seems to be a special property of the onedimensional model that is not borne out by the threedimensional simulations where it stays approximately constant.
Remarkably, the results are fairly independent of the cooling time, except that the break point where (∇s)_{rms} saturates occurs for smaller values of φ_{0} as we increase the cooling time; see Fig. 6. This break point is also related to the point where the Mach number saturates, as can be seen from Fig. 7. However, for longer cooling times there can be longer transients, making it difficult to obtain good averages. Therefore we focus, in the rest of this paper, on the case of shorter cooling times using St_{cool} = 0.2. Another remarkable result is that the normalized value of (∇T)_{rms}(∇s)_{rms} is always of the order of about 10^{3}, independent of resolution, cooling time, or the value of the viscosity.
Fig. 6 Dependence of the rms values of temperature and entropy on φ_{0} for ν/c_{s}R = 1 and St_{cool} = 0.2 (top panel), 0.6 (middle panel), and 2 (bottom panel). 
Fig. 7 Dependence of Mach number and the rms value of entropy on φ_{0} for ν/c_{s}R = 1 and St_{cool} = 0.2, 0.6, and 2 (solid, dotted, and dashed line types, respectively). 
Fig. 8 Images of T, s, (∇T × ∇s)_{z}, and normalized vertical vorticity for a twodimensional run with δt_{force}c_{s0}/R = 0.1 at an instant shortly before the second expansion wave is launched (top row), and shortly after the second expansion wave is launched (second and third row). Note the vorticity production from the baroclinic term in the second and third row, while in the top row, (∇T × ∇s)_{z} and ω_{z} are just at the noise level of the calculation. Even under our weakly supersonic conditions shock surfaces are well localized and the zones of maximum production of vorticity appear to be those in which the fronts encounter each other. Here we have used , ν = χ = 0.1c_{s0}R, with 512^{2} mesh points. Only the inner part of the domain is shown. 
Fig. 9 Dependence of Ma, Re, the rms values of ∇s and ∇T, the angle θ between them, as well as k_{ω}/k_{f}, on for ν/c_{s}R = 1. 
Most of the basic features of the onedimensional model are also reproduced by two and threedimensional calculations. Twodimensional simulations have the advantage of being easily visualized and are therefore best suited for illustrating vorticity production by the baroclinic term. In Fig. 8 we demonstrate that vorticity production is associated with the interaction between the fronts of different expansion waves. In this example we chose δt_{force}c_{s0}/R = 0.1, so the first expansion wave is launched at t = 0 and the second one at tc_{s0}/R = 0.1. The top row of Fig. 8 shows that at tc_{s0}/R = 0.09, i.e. just before launching the second expansion wave, the baroclinic term and the vorticity are still just at the noise level. At that time the most pronounced feature is the discontinuity between the Gaussian expansion waves in the periodic domain. This leads to negligibly weak vorticity, and no baroclinic term. At tc_{s0}/R = 0.11, the effect of the second expansion wave becomes noticeable in visualizations of both (∇T × ∇s)_{z} and ω_{z}, while our visualizations of T and s barely show the second expansion wave. At tc_{s0}/R = 0.14, the first expansion wave is clearly no longer circular, which is obviously associated with the second expansion wave that is now quite pronounced in our visualizations of both T and s.
In order to have a more accurate quantitative determination of vorticity production, we now consider threedimensional models. In Fig. 9 we show the dependence of various quantities on φ_{0} for St_{cool} = 0.2 and ν/c_{s}R = 1. In all cases we use 128^{3} mesh points and average the results over between 20 and 70 turnover times. Note that here . Given that Re depends inverse proportionally on ν/Rc_{s0}, we can also write Re ≈ 0.05 φ_{0}R/c_{s0}ν. The Mach number saturates at about Ma = 3, and the rms value of the entropy gradient increases up until this point. Given that the rms value of the temperature gradient also stays approximately constant, we find a weak increase of (∇T)_{rms}(∇s)_{rms}. The value of (∇T × ∇s)_{rms} is always found to be a certain fraction below this value, resulting in a typical baroclinic angle of about 45 degrees; see the third panel of Fig. 9. Finally, the amount of vorticity production in terms of k_{ω}/k_{f} is about 0.3 for . For smaller values, on the other hand, there is an approximately linear increase with .
Fig. 10 Timeaveraged enstrophy spectra, E_{ω}(k) (thick lines), compared with k^{2}E_{K}(k) (thin lines below the corresponding thick lines), for the threedimensional baroclinic case with (dashed), 100 (dotted), and 500 (solid lines). The k^{2} slope is shown for comparison. In all three cases we have k_{f}/k_{1} = 4. 
The possibility of spurious vorticity is easily eliminated in this case by looking at enstrophy spectra; see Fig. 10, where we compare E_{ω}(k) with k^{2}E_{K}(k). All spectra fall off rapidly with increasing k. Thus, even though the initial vorticity generation occurred evidently at the smallest available scales, once the flow becomes fully developed, most of the enstrophy resides at scales equal to or larger than the driving scale. Furthermore, the spectra of E_{ω}(k) and k^{2}E_{K}(k) are close together, suggesting that the vorticity is close to its maximal value.
5. Applications
The level of vorticity that is produced in the usual case of solenoidal forcing of the turbulence is such that k_{ω}/k_{f} ≈ 1 (see, e.g., Brandenburg 2001). For turbulence whose driving force has finite correlation time (St_{force} = 0.3, for example), and small values of Re, we have k_{ω}/k_{f} = O(1) when . However, for larger values of Re, the turbulence becomes vortical already for smaller values of Co. Comparing with the galaxy, we have Ω ≈ 10^{15} s^{1}, u_{rms} = 10 km s^{1}, and an estimated correlation length of about 70 pc, so k_{f} = 3 × 10^{20} cm, so Co = 0.07, which is rather small. Thus, rotation may not be able to produce sufficient levels of vorticity. Given that in galaxies with flat rotation curves, S ≈ −Ω, shear should not be very efficient either. However, the Mach numbers are undoubtedly larger than unity in the interstellar medium, so this should lead to values of k_{ω}/k_{f} ≈ 0.3, which is the saturation value found in Fig. 9. Given that one of the reasons for studying the production of vorticity is the question of dynamo action, we should point out that such values of k_{ω}/k_{f} are large enough for the smallscale dynamo. Largescale dynamo action should be possible in galaxies as well, because of their large length scales, but it suffers from the wellknown problem of a small growth rate. It then remains difficult to explain largescale magnetic fields in very young galaxies (Beck et al. 1996).
The question of vorticity generation is also important in studies of the very early Universe, where phase transition bubbles are believed to be generated in connection with the electroweak phase transition (Kajantie & KurkiSuonio 1986; Ignatius et al. 1994). Here the equation of state is that of a relativistic fluid, p = ρc^{2}/3, where c is the speed of light. Thus, there is no baroclinic term and no obvious source of vorticity. However, the relativistic equation of state may be modified at small length scales, but it is not clear that this can facilitate significant vorticity production.
6. Conclusions
The present work has demonstrated that vorticity production is actually quite ubiquitous once there is rotation, shear, or baroclinicity. This implies that the assumption of potential flows as a model for interstellar turbulence might be of academic interest and could only be realized under special conditions of weak forcing, weak rotation, and no shear. In galaxies, however, the shear and Coriolis number are well below unity, leaving only the baroclinic term as a viable candidate for the production of vorticity. This agrees with early work of Korpi et al. (1999b), who analyzed the production terms in supersonic, supernovadriven turbulence quantitatively; see also Glasner et al. (1997), who showed that on long enough time scales significant vorticity can also be produced for subsonic flows. We have also observed how vorticity is mainly produced close to shock front encounters. This motivates a more detailed investigation of these zones as the next step in the study of vorticity generation in the interstellar medium. It should also be pointed out that the baroclinic term corresponds to the battery term in the induction equation Kulsrud et al. (1997); Brandenburg & Subramanian (2005). Thus, when studying the possibility of dynamo action, this battery term provides an intrinsic and well defined seed for the dynamo and should therefore be included as well.
Acknowledgments
We thank an anonymous referee for making a number of useful suggestions that have improved our paper. We acknowledge the allocation of computing resources provided by the Swedish National Allocations Committee at the Center for Parallel Computers at the Royal Institute of Technology in Stockholm and the National Supercomputer Centers in Linköping. This work was supported in part by the European Research Council under the AstroDyn Research Project 227952 and the Swedish Research Council grant 62120074064.
References
 Beck, R., Brandenburg, A., Moss, D., Shukurov, A., & Sokoloff, D. 1996, ARA&A, 34, 155 [NASA ADS] [CrossRef] [Google Scholar]
 Bodo, G., Mignone, A., Cattaneo, F., Rossi, P., & Ferrari, A. 2008, A&A, 487, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brandenburg, A. 2001, ApJ, 550, 824 [NASA ADS] [CrossRef] [Google Scholar]
 Brandenburg, A., & Subramanian, K. 2005, Phys. Rep., 417, 1 [NASA ADS] [CrossRef] [Google Scholar]
 de Avillez, M. A., & Breitschwerdt, D. 2005, A&A, 436, 585 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Federrath, C., RomanDuval, J., Klessen, R. S., Schmidt, W., & Mac Low, M.M. 2010, A&A, 512, A81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Gissinger, C., Fromang, S., & Dormy, E. 2009, MNRAS, 394, L84 [NASA ADS] [Google Scholar]
 Glasner, A., Livne, E., & Meerson, B. 1997, Phys. Rev. Lett., 78, 2112 [NASA ADS] [CrossRef] [Google Scholar]
 Goldreich, P., & LyndenBell, D. 1965, MNRAS, 130, 125 [NASA ADS] [CrossRef] [Google Scholar]
 Gressel, O., Elstner, D., Ziegler, U., & Rdiger, G. 2008, A&A, 486, L35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hanasz, M., OtmianowskaMazur, K., Kowal, G., & Lesch, H. 2009, A&A, 498, 335 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hawley, J. F., Gammie, C. F., & Balbus, S. A. 1995, ApJ, 440, 742 [NASA ADS] [CrossRef] [Google Scholar]
 Ignatius, J., Kajantie, K., KurkiSuonio, H., & Laine, M. 1994, Phys. Rev. D, 49, 3854 [NASA ADS] [CrossRef] [Google Scholar]
 Kadomtsev, B. B., & Petviashvili, V. I. 1973, Sov. Phys. Dokl., 18, 115 [NASA ADS] [Google Scholar]
 Kajantie, K., & KurkiSuonio, H. 1986, Phys. Rev. D, 34, 1719 [NASA ADS] [CrossRef] [Google Scholar]
 Korpi, M. J., Brandenburg, A., Shukurov, A., Tuominen, I., & Nordlund, Å. 1999a, ApJ, 514, L99 [NASA ADS] [CrossRef] [Google Scholar]
 Korpi, M. J., Brandenburg, A., Shukurov, A., & Tuominen, I. 1999b, in Interstellar Turbulence, ed. J. Franco, & A. Carramiñana (Cambridge University Press), 127 [Google Scholar]
 Kulsrud, R. M., Cen, R., Ostriker, J. P., & Ryu, D. 1997, ApJ, 480, 481 [NASA ADS] [CrossRef] [Google Scholar]
 Lesieur, M. 1990, Turbulence in Fluids (Dordrecht: Martinius Nijhoff Publishers) [Google Scholar]
 Mac Low, M.M., Balsara, D. S., Kim, J., & de Avillez, M. A. 2005, ApJ, 626, 864 [NASA ADS] [CrossRef] [Google Scholar]
 Mee, A. J., & Brandenburg, A. 2006, MNRAS, 370, 415 [NASA ADS] [Google Scholar]
 Porter, D. H., Woodward, P. R., & Pouquet, A. 1998, Phys. Fluids, 10, 237 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Regev, O., & Umurhan, O. M. 2008, A&A, 481, 21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Rosen, A., & Bregman, J. N. 1995, ApJ, 440, 634 [NASA ADS] [CrossRef] [Google Scholar]
 Rosen, A., Bregman, J. N., & Kelson, D. D. 1996, ApJ, 470, 839 [NASA ADS] [CrossRef] [Google Scholar]
 Sridhar, S., & Subramanian, K. 2009, Phys. Rev. E, 80, 066315 [NASA ADS] [CrossRef] [Google Scholar]
 Wisdom, J., & Tremaine, S. 1988, AJ, 95, 925 [NASA ADS] [CrossRef] [Google Scholar]
 Yousef, T. A., Brandenburg, A., & Rüdiger, G. 2003, A&A, 411, 321 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
All Tables
Rootmeansquared values of components of the velocity derivative tensor, normalized by u_{rms}k_{f}, as well as the three diagonal components of the ⟨ u_{i}u_{j} ⟩ tensor for 4 values of Co.
All Figures
Fig. 1 Timeaveraged kinetic energy and enthalpy spectra for two values of the Coriolis number for Re = 25 and St_{force} = 0.4. The two straight lines give the slopes − 2 and − 3, respectively. In both cases we have k_{f}/k_{1} = 4. 

In the text 
Fig. 2 Dependence of k_{ω}/k_{f} on Co for three groups of runs: group 1 with Re = 15, k_{f}/k_{1} = 10, St_{force} = 0.3; group 2 with Re between 25 and 130, k_{f}/k_{1} = 4, St_{force} = 0.4; and group 3 with Re = 30 and 150, k_{f}/k_{1} = 2, St_{force} = 0. 

In the text 
Fig. 3 Dependence of k_{ω}/k_{f} on Re for Co ≈ 1, k_{f}/k_{1} = 2, and St_{force} = 0 (δcorrelated forcing). 

In the text 
Fig. 4 Timeaveraged enstrophy spectra, E_{ω}(k), compared with k^{2}E_{K}(k), for Re = 25, St_{force} = 0.4, and three values of the Coriolis number. The curves of k^{2}E_{K}(k) are close together and overlap for Co = 0.01 (dotted) and 0.15 (dashed), so it becomes a single dashdotted line. The k^{3} slope is shown for comparison. In all three cases we have k_{f}/k_{1} = 4. 

In the text 
Fig. 5 Dependence of k_{ω}/k_{f} on Sh for Re ≈ 40 and k_{f}/k_{1} = 2 and δcorrelated forcing (St_{force} = 0). Different resolutions are shown to give similar results. At small values of  Sh  , comparisons with St_{force} = 0.3 (keeping k_{f}/k_{1} = 2) or k_{f}/k_{1} = 10 (keeping St_{force} = 0) are also shown. 

In the text 
Fig. 6 Dependence of the rms values of temperature and entropy on φ_{0} for ν/c_{s}R = 1 and St_{cool} = 0.2 (top panel), 0.6 (middle panel), and 2 (bottom panel). 

In the text 
Fig. 7 Dependence of Mach number and the rms value of entropy on φ_{0} for ν/c_{s}R = 1 and St_{cool} = 0.2, 0.6, and 2 (solid, dotted, and dashed line types, respectively). 

In the text 
Fig. 8 Images of T, s, (∇T × ∇s)_{z}, and normalized vertical vorticity for a twodimensional run with δt_{force}c_{s0}/R = 0.1 at an instant shortly before the second expansion wave is launched (top row), and shortly after the second expansion wave is launched (second and third row). Note the vorticity production from the baroclinic term in the second and third row, while in the top row, (∇T × ∇s)_{z} and ω_{z} are just at the noise level of the calculation. Even under our weakly supersonic conditions shock surfaces are well localized and the zones of maximum production of vorticity appear to be those in which the fronts encounter each other. Here we have used , ν = χ = 0.1c_{s0}R, with 512^{2} mesh points. Only the inner part of the domain is shown. 

In the text 
Fig. 9 Dependence of Ma, Re, the rms values of ∇s and ∇T, the angle θ between them, as well as k_{ω}/k_{f}, on for ν/c_{s}R = 1. 

In the text 
Fig. 10 Timeaveraged enstrophy spectra, E_{ω}(k) (thick lines), compared with k^{2}E_{K}(k) (thin lines below the corresponding thick lines), for the threedimensional baroclinic case with (dashed), 100 (dotted), and 500 (solid lines). The k^{2} slope is shown for comparison. In all three cases we have k_{f}/k_{1} = 4. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.