Issue 
A&A
Volume 562, February 2014



Article Number  A112  
Number of page(s)  20  
Section  Interstellar and circumstellar matter  
DOI  https://doi.org/10.1051/00046361/201322482  
Published online  17 February 2014 
Supersonic turbulence in 3D isothermal flow collision
^{1}
École Normale Supérieure, Lyon, CRAL, UMR CNRS 5574, Université de
Lyon,
69364
Lyon Cedex 7,
France
email:
doris.folini@enslyon.fr
^{2}
Swiss National Supercomputing Center, CSCS Lugano,
6900
Lugano,
Switzerland
Received:
11
August
2012
Accepted:
17
December
2013
Large scale supersonic bulk flows are present in a wide range of astrophysical objects, from Ostar winds to molecular clouds, galactic sheets, accretion, or γray bursts. Associated flow collisions shape observable properties and internal physics alike. Our goal is to shed light on the interplay between large scale aspects of such collision zones and the characteristics of the compressible turbulence they harbor. Our model setup is as simple as can be: 3D hydrodynamical simulations of two headon colliding, isothermal, and homogeneous flows with identical upstream (subscript u) flow parameters and Mach numbers 2 < M_{u} < 43. The turbulence in the collision zone is driven by the upstream flows, whose kinetic energy is partly dissipated and spatially modulated by the shocks confining the zone. Numerical results are in line with expectations from selfsimilarity arguments. The spatial scale of modulation grows with the collision zone. The fraction of energy dissipated at the confining shocks decreases with increasing M_{u}. The mean density is ρ_{m} ≈ 20ρ_{u}, independent of M_{u}. The root mean square Mach number is M_{rms} ≈ 0.25M_{u}. Deviations toward weaker turbulence are found as the collision zone thickens and for small M_{u}. The density probability function is not lognormal. The turbulence is inhomogeneous, weaker in the center of the zone than close to the confining shocks. It is also anisotropic: transverse to the upstream flows M_{rms} is always subsonic. We argue that uniform, headon colliding flows generally disfavor turbulence that is at the same time isothermal, supersonic, and isotropic. The anisotropy carries over to other quantities like the density variance – Mach number relation. Lineofsight effects thus exist. Structure functions differ depending on whether they are computed along a lineofsight perpendicular or parallel to the upstream flow. Turbulence characteristics generally deviate markedly from those found for uniformly driven, supersonic, isothermal turbulence in 3D periodic box simulations. We suggest that this should be kept in mind when interpreting turbulence characteristics derived from observations. Our simulations show that even a simple model setup results in a richly structured interaction zone. The robustness of our findings toward more realistic setups remains to be tested.
Key words: shock waves / turbulence / hydrodynamics / ISM: kinematics and dynamics / gammaray burst: general / binaries: close
© ESO, 2014
1. Introduction
Shockbound interaction zones are ubiquitous in astrophysics. They occur, for example, in binary star systems (Stevens et al. 1992; Myasnikov & Zhekov 1998; Parkin & Pittard 2010), the jets of highenergy objects (Panaitescu et al. 1999; Kaiser et al. 2000; Fan & Wei 2004), or galaxy formation and cosmology (Anninos & Norman 1996; Kang et al. 2005; Agertz et al. 2009; Klar & Mücket 2012). The idea that such zones may also play a crucial role in the context of molecular clouds and star formation stimulated extensive research in recent years, on their large scale properties as well as on the turbulence they harbor. The present paper stresses the combined view, the interplay between large scale aspects of the collision zone, and the characteristics of the turbulence in its interior.
Substantial progress in understanding supersonic turbulence and its characteristics has come from 3D periodic box simulations. These simulations meanwhile cover an impressive range of physics from pure isothermal hydrodynamics to the inclusion of magnetic fields, radiative cooling, selfgravity, chemistry, or relativistic flows (Stone et al. 1998; Mac Low et al. 1998; Mac Low 1999; Padoan & Nordlund 1999; Boldyrev et al. 2002b,a; Kritsuk & Norman 2004; Padoan et al. 2004; Kritsuk et al. 2006, 2007a, 2011a,b; Gazol et al. 2007; Padoan et al. 2007; Federrath et al. 2008; Schmidt et al. 2009; Federrath et al. 2010; Glover et al. 2010; Gazol & Kim 2010; Price et al. 2011; Seifried et al. 2011; Konstandin et al. 2012; Molina et al. 2012; Downes 2012; Gazol & Kim 2013; Federrath & Klessen 2013). Turbulence in such simulations is either left to decay or is forced. Typical ways of forcing are continuous energy input in each grid cell or occasional energy input at randomly selected grid cells, to mimic energy input into the interstellar medium by, for example, collimated or spherical outflows or supernova explosions (Li & Nakamura 2006; de Avillez & Breitschwerdt 2007; Wang et al. 2010; Moraghan et al. 2013). Only rather recently have analytical results been presented that demonstrate the existence of an intermediate scaling range and derive scaling relations for compressible turbulence(Aluie 2011; Galtier & Banerjee 2011; Banerjee & Galtier 2013).
For isothermal hydrodynamics, on which the present paper concentrates, the 3D periodic box studies reveal the crucial role of the form of the energy input for the characteristics of the turbulence. The wave length at which the energy is put into the system sets the largest spatial scale of the turbulent structures (Mac Low 1999; BallesterosParedes & Mac Low 2002). Purely compressible forcing leads to much sharper density contrasts and wider density PDFs than purely solenoidal forcing (Federrath et al. 2009, 2010). Likewise, stronger density contrasts result if more energy is put into the system (Mac Low 1999).
Less explored is the turbulence harbored within flow collision zones and its interplay with the overall properties of such zones. It has long been known that isothermal headon colliding flows are unstable (e.g., Vishniac 1994; Blondin & Marks 1996). The role of additional physics, like radiative cooling, selfgravity, or the inclusion of magnetic fields, has been investigated since by a number of studies (Walder & Folini 1998, 2000b; Hennebelle & Pérault 1999; Audit & Hennebelle 2005; Folini & Walder 2006; VázquezSemadeni et al. 2006, 2007; Hennebelle & Audit 2007; Hennebelle et al. 2008; Inoue & Inutsuka 2009; Heitsch et al. 2009; Folini et al. 2010; VázquezSemadeni et al. 2010; Klessen & Hennebelle 2010; Ntormousi et al. 2011; Gong & Ostriker 2011; Zrake & MacFadyen 2011; Inoue & Inutsuka 2012; Zrake & MacFadyen 2012, 2013; Inoue & Fukui 2013). Of interest in the context of the present paper is the observation that the very early evolution of the collision zone can be much more violent (Heitsch et al. 2006; VázquezSemadeni et al. 2006) than the situation later on, when a still unstable but less violent, statistically stationary situation develops (Walder & Folini 2000b; Folini & Walder 2006; Audit & Hennebelle 2010). It is this later stage of the evolution we are interested in, and there in particular in the characteristics of the supersonically turbulent interior of the collision zone and its interplay with the confining shocks, which modulate the energy input that drives the turbulence.
The paper is a follow up of Folini & Walder (2006, FW06), a 2D study of headon colliding isothermal flows. FW06 pointed out the coupling between the turbulence within the collision zone and the confining shocks of this zone. For the upstream Mach numbers considered (5 < M_{u} < 90), the fraction of upstream kinetic energy that passes the confining shocks without getting dissipated (thus is available for driving the turbulence in the collision zone) scales with the root mean square Mach number of the turbulence. In addition, FW06 derived approximate selfsimilarity relations for collision zone mean quantities and checked them against numerical results.
Here we extend this earlier study by considering isothermal headon colliding flows in 3D instead of 2D and by investigating further turbulence characteristics. In particular, the following questions shall be addressed. What level of turbulence, what root mean square Mach number M_{rms}, can be reached within the interaction zone for a given upstream Mach number M_{u}? How homogeneous and isotropic is the turbulent interior of the interaction zone? What do established characteristics of supersonic turbulence look like, in particular the density PDF, the density varianceMach number relation, and the structure functions? Finally, we contemplate on potential implications of our results for real astrophysical objects, for their observation, and on how additional physics may alter the presented results.
We stress that the setup studied here is highly idealized. In reality, there exist no strictly isothermal, uniform, or precisely headon colliding flows of equal Mach number and density. The study highlights, however, the wealth of phenomena that can result already from such a setup – and, likewise, indicates characteristics beyond reach. And it offers a somewhat complementary, thus potentially interesting, view on supersonic turbulence, as compared to the more wide spread (and often equally idealized) 3D periodic box perspective.
The structure of the paper is as follows. In Sect. 2 we present the numerical method and the physical simulation setup. Results follow in Sect. 3 and are discussed in Sect. 4. Summary and conclusions are given in Sect. 5
2. Physical model and numerical method
The physical model and numerical method used are the same as in FW06 but extended now to 3D. We refer to this paper for a detailed description and restrict ourselves here to the most important facts and information specific to the 3D simulations.
2.1. Physical model problem
We consider a 3D, planeparallel, infinitely extended (yzdirection), isothermal, shock compressed slab. Two high Machnumber flows, oriented parallel (left flow, subscript l) and antiparallel (right flow, subscript r) to the xdirection, collide headon. The resulting highdensity interaction zone we denote by CDL for “cold dense layer” to remain consistent with notation used in previous papers (Walder & Folini 1996, 1998; Folini & Walder 2006). We investigate this system within the frame of Euler equations, together with a polytropic equation of state. For the polytropic exponent, we choose γ = 1.000001. This value guarantees that jump conditions and wave speeds of a Mach90 shock are within 0.01 per cent of the isothermal values.
We only consider symmetric settings, where the left and right colliding flow have identical density and velocity (subscript u for upstream): ρ_{l} = ρ_{r} ≡ ρ_{u} and  v_{l}  =  v_{r}  ≡ v_{u}.
We express velocities in units of the isothermal sound speed a, densities in terms of the upstream density ρ_{u}, and lengths in units of Y, the yextent of the computational domain we used. This artificial choice is necessary as there is no natural timeindependent length scale to the problem (see Sect. 3). The yextent of the computational domain is the same in all our simulations and is identical with the zextent, so Y = Z.
List of performed simulations.
2.2. Numerical method
Our results are computed with the hydrodynamic code from the AMAZE codepackage (Walder & Folini 2000a; Folini et al. 2003; Melzani et al. 2013). We use the multidimensional highresolution finitevolumeintegration scheme developed by Colella (1990) on the basis of a Cartesian mesh. In all our simulations we use a version of the scheme that is (formally) second order accurate in space and time for smooth flows. We combine this integration scheme with the adaptive mesh algorithm by Berger (1985). A coarse mesh is used for the upwind flows, a finer mesh for the CDL. The meshes adapt automatically to the increasing spatial extension of the CDL.
As described in FW06, we have our CDL moving in positive xdirection at a speed of about Mach 20–40 to avoid alignment effects of strong shocks, a well known problem not particular to our scheme (Colella & Woodward 1984; Quirk 1994; Jasak & Weller 1995). We rely on the MILES approach to mimic the highwavenumber end of the inertial subrange (e.g., Boris et al. 1992; Porter & Woodward 1994), but see also Garnier et al. (1999) and Domaradzki (2010) for a discussion of limitations.
2.3. Numerical settings of simulations
The y and zextent of the computational domain (directions perpendicular to the upstream flows) are identical in all our simulations, Y = Z. The xextent is 200 times larger. Boundary conditions in xdirection are “supersonic inflow”. Periodic boundary conditions are used in y and zdirection.
The runs performed differ in their upstream Machnumber M_{u}, with 2 ≲ M_{u} ≲ 45, and in their discretization, with 128 or 256 cells in y and zdirection on the finest level of refinement, which only covers the CDL. Runs are labled with M_R, where M is the upwind Machnumber and R indicates the refinement: 128 cells (R = 1) or 256 cells (R = 2). The runs are listed in Table 1.
At t = 0, no CDL is present in our simulations. Random density perturbations of up to 2% are added initially to the flow left of the interface separating both flows, up to a distance Y/2 from the interface. This initialization results in a fast development of turbulence in the CDL without imposing any artificial mode.
Fig. 1 Different geometry of confining shocks. In 2D (periodic in ydirection, infinite in zdirection), the confining shocks resemble corrugated sheets (left panel). In 3D (periodic in y and zdirection) they evoke the idea of cardboard egg wrappings (right panel). Both simulations have M_{u} = 21.7. Color coded is density (log10(ρ); particles/cm^{3} in 2D; g/cm^{3} in 3D) and, right panel only, the driving efficiency f_{eff} at the shock surface. 
2.4. Data analysis
The relevant quantity for the time evolution of CDL quantities is the average xextension of the CDL, ℓ_{cdl} ≡ V/(YZ) with V the volume of the CDL (see Sect. 3 or the 2D case in FW06). The average CDL extension does, however, not need to be a strictly monotonically increasing function of time, as the CDL compression fluctuates slightly with time.
A proxy to ℓ_{cdl} that does not show such fluctuations can be defined through the monotonically increasing CDL mass. Expressing ℓ_{cdl} in terms of the CDL mass m_{cdl} and mean density ρ_{m} as ℓ_{cdl} = m_{cdl}/(ρ_{m}Y^{2}), where we used Y = Z, we can replace ρ_{m} by ρ_{u} and normalize with Y to obtain a monotonically increasing, dimensionless quantity, ℓ ≡ m_{cdl}/(ρ_{u}Y^{3}).
All simulations were advanced till at least ℓ = 12, some much further. Except for the lowest Mach number simulations, ℓ ≈ 12 corresponds to within 10% to a spatial extent of Y/2, i.e., the xextension of the CDL is about half the domain size in yzdirection. Data analysis is done offline, on previously dumped data sets. Temporal averages, unless otherwise stated, are taken over ℓ ∈ [11,12]. This interval comprises at least three data sets.
Part of the data analysis was done with VisIt^{1}, combined with python scripting. We developed a reader for VisIt that can cope with our block structure adaptive mesh refinement data, stored in hdf5 format. VisIt proved particularly useful for the slice wise analysis of the data (see Sect. 3.2.3) and for the numerical determination of the driving efficiency. The latter involves local decomposition of the upstream flow into components perpendicular and parallel to the confining shocks and integration over these shocks (see Appendix A). We use VisIt to create an isosurface representative of the confining shock, defined by a density threshold slightly greater than ρ_{u}, and then exploit VisIts capabilities of providing surface normals to isosurfaces. An example of the reconstructed shock surface, colored in driving efficiency f_{eff} (see Eq. (9)), is shown in Fig. 1, right panel.
3. Results
The section is structured along two main perspectives: the view on the CDL as an entity and the view on its interior structure. In Sect. 3.1, we concentrate on approximate scaling relations for some CDL mean quantities in 1D to 3D. In Sect. 3.2 we focus on the anisotropy and inhomogeneity of the CDL interior. 2D slices are used to analyze how mean quantities change with distance from the CDL center. Section 3.3 is dedicated to density PDFs, Sect. 3.4 to structure functions. The physical results are complemented by considerations on numerical resolution in Sect. 3.5.
3.1. CDL Mean quantities: approximate scaling relations
Within the frame of Euler equations, i.e., in the absence of viscous dissipation, and making some further assumptions as detailed below, approximate selfsimilar scaling relations for CDL mean quantities can be derived and expressed in terms of upstream flow parameters. The relations are qualitatively different in 1D as compared to 2D and 3D (see FW06). In particular, in 1D the mean density increases with the upstream Mach number squared. In 2D/3D the CDL mean density is expected to be limited, the root mean square Mach number should increase linearly with the upstream Mach number.
3.1.1. Analytical scaling relations for CDL mean quantities
Looking first at the 1D case and denoting the density and velocity of the CDL by ρ_{1d} and v_{1d}, those of the left and right upwind flows by ρ_{i} and v_{i} (i = l,r), and the isothermal sound speed by a, it can be shown that (e.g., Courant & Friedrichs 1976) in the rest frame of the CDL The approximation holds for large Machnumbers. A relation between characteristic length and time scales of the solution, the selfsimilarity variable κ_{1d}, can be obtained as the ratio between the spatial extension ℓ_{1d} of the CDL and the time τ needed to accumulate the corresponding column density N_{1d}: and using (5)For strong shocks and symmetric settings (l = r) the above relations reduce to and κ_{1d} = 2a/M_{u}.
For the 2D case, FW06 derived scaling relations for five CDL mean quantities: mean density ρ_{m}; root mean square Mach number M_{rms}; driving energy , i.e., the kinetic energy density entering the CDL per time and unit length in the ydirection; , the energy density dissipated per time within an average column of length ℓ_{cdl}; the change per time of the average kinetic energy energy density contained within such an average column, .
Extension to 3D is straightforward. The only derivation we repeat here (see Appendix A) is that of the driving energy , which involves an integral over the confining shocks. Part of the total (left plus right flow) upwind kinetic energy flux density, , is thermalized at these shocks. The remaining part, , drives the turbulence in the CDL. In analogy with the 2D case, we write , where the driving efficiency f_{eff} is a function of the upwind Machnumber alone. The formal dependencies, involving upstream quantities ρ_{u} and M_{u}, as well as parameters β_{i} and η_{i} (to be determined), are as in the 2D case: Here, β_{1} = 0 and β_{2} = 1 from analytical considerations, whereas β_{3} as well as η_{i}, i = 1,2,3, can be determined only from numerical simulations and turn out to be different in 2D and 3D. The 3D numerical results, to be presented in Sect. 3.1.2, confirm β_{1} = 0 and β_{2} = 1, and yield β_{3} ≈ −1.17, η_{1} ≈ 25, η_{2} ≈ 0.2, and η_{3} ≈ 5 as best fit values.
The above relations imply that the mean density of the CDL does not depend on the upstream Mach number M_{u}, in strong contrast to the 1D case where . The CDL root mean square Mach number should increase linearly with M_{u}. Also: the larger M_{u}, the larger f_{eff}, i.e., the more energy passes the confining shocks of the CDL unthermalized
We stress the four basic assumptions made in deriving the above relations: a) a simple Machnumber dependency of ρ_{m}, v_{rms}, and f_{eff}; b) CDL density and velocity are uncorrelated; c) we have high Machnumbers in the sense that or ; d) v_{diss} ∝ v_{rms}. For a discussion of the validity of these assumptions we refer to FW06.
3.1.2. 3D scaling relations: numerical results
The analytical relations in Eqs. (6) to (12) basically state that CDL mean quantities can be expressed in terms of upstream quantities and that the corresponding relations do not change with time/CDL size. Comparing the analytical predictions with our numerical results thus means checking for two things: do the relations hold across simulations, i.e., across widely varying upstream Mach numbers, and do the relations hold independent of time/CDL size.
We first focus on the comparison of analytical and numerical results with regard to the first question, the validity of the analytical relations across different upstream Mach numbers. In doing so, we tacitly assume that βi and η_{i} do not change with time/CDL size.
Looking first at ρ_{m}, Fig. 2, top panel, and ignoring for the moment the obvious drift with increasing CDL thickness, and concentrating instead on CDL thicknesses 6 ≤ ℓ ≤ 12, the density compression ratios can be seen to lie mostly in a range of 20 ≤ η_{1} ≤ 25. This narrow range is remarkable given the range of upstream Mach numbers (2.7 ≤ M_{u} ≤ 43.4) and the associated 1D compression ratios ( 7 < ρ_{1d}/ρ_{u} < 1900). The only clear exceptions here are simulation R2_2 and R4_2, whose upstream Mach numbers allow at most for 1D compression factors of about 7 and 16, respectively, roughly the values attained also in the 3D case (dark blue and dark green curves in Fig. 2). The analytically predicted value of β_{1} = 0 is confirmed by the numerical results in that the scatter among the different curves in Fig. 2, top panel, cannot be improved by any other choice for β_{1}.
For the root mean square Mach number (Fig. 2, middle panel), we find η_{2} = M_{rms}/M_{u} = 0.25 to within about 10% for all simulations except those with very low upstream Mach number (M_{u} ≤ 7). As in the case of density, changing the predicted dependence of M_{rms} on M_{u} does not reduce the scatter of the curves, unless low Mach number simulations are included, in which case β_{2} = 1.1 or even 1.2 gives a smaller spread among curves than the predicted value of β_{2} = 1. Using Eq. (7), the above range for η_{2} implies for a value in the range from 13 to 20, somewhat smaller but roughly consistent with the numerically obtained density compression ratios discussed above. Better consistency is achieved if η_{2} is determined using Eq. (11) instead of Eq. (7). From Fig. 2, bottom panel, and neglecting again low Mach number simulations, we find η_{2} ≈ 0.2, which corresponds to η_{1} ≈ 25. We ascribe the ambiguity in the numerical value of η_{2} to slight violations of the underlying assumptions in deriving Eqs. (6) to (12), in the case of low Mach number simulations especially the assumption that M_{rms} ≫ 1.
Fig. 2 Mean density ρ_{m} (top panel), root mean square Mach number M_{rms} (middle panel), and (bottom panel) as function of the monotonically increasing CDL size ℓ for simulations R*_2. Individual curves denote M_{u} = 2 (dark blue), 4 (dark green), 5 (yellow), 7 (purple), 8 (light blue), 11 (black), 16 (green), 22 (red), 27 (pink), 33 (cyan), and 43 (magenta). 
Fig. 3 Driving efficiency of 3D slabs, same simulations and color coding as in Fig. 2. Top panel: comparison of numerical results with expected scaling relation β_{3}log _{10}(M_{u}) ∝ log _{10}(1 − f_{eff}). Each star denotes the maximum driving efficiency of one simulation. Bottom panel: evolution of f_{eff} as function of CDL size. 
Numerical results regarding the remaining scaling relations, those involving the driving efficiency f_{eff} and associated parameters η_{3} and β_{3}, are shown in Fig. 3. The top panel relates the maximum driving efficiency achieved in each simulation with its upstream Mach number. For M_{u} ≥ 5, an excellent linear relation ship between log(M_{u}) and log(1 − f_{eff}) is found with slope β_{3} = −1.17 ± 0.02 and η_{3} = 5.
Changing perspective and looking at the evolution with CDL thickness (or time) in Fig. 2, all simulations show a persistent and clear decrease of the root mean square Mach number with time and an associated increase of the density compression ratio, in clear disagreement with the expected scaling relations. In units M_{rms}/M_{u}, the decrease is about the same for all simulations, around 0.02 to 0.03 (between about 10% and 20%) for a doubling of the CDL thickness ℓ from 6 to 12. Relative changes are more pronounced for simulations with smaller M_{u}.
The same effect was observed in the 2D case and there attributed to the increasing role of numerical dissipation with increasing CDL size (FW06). We favor the same explanation given there also for the 3D case. The volume (or mass) fraction of the CDL that is subsonic remains about constant with time (Fig. 4). In absolute terms, the subsonic parts of the CDL grow as the CDL grows. Energy loss in these parts, dominated by numerical dissipation/viscous dissipation in terms of the MILES approach, grows in proportion – while . The effect leads to an overall decay of CDL turbulence. It is more severe for smaller M_{u}, as a larger fraction of the CDL is subsonic (Fig. 4), thus in violation of the assumptions made in deriving the analytical scaling relations in Sect. 3.1.1.
Similar drifts with time exist for the driving efficiency, Eq. (9) and Fig. 3 (bottom panel). However, considering again the period 6 ≤ ℓ ≤ 12, the drifts with time of f_{eff} can be to either larger or smaller values. Apparently, the larger M_{u}, the larger the CDL has to grow before the driving efficiency reaches its peak value.
In summary, we conclude from the above findings that the numerical results essentially confirm the validity of the selfsimilar analytical relations Eqs. (6) to (12) across different upstream Mach numbers, although the numerical results display a systematic drift with increasing CDL size, which we attribute to the effect of numerical dissipation/viscous dissipation.
3.2. Inhomogeneity and anisotropy of CDL
The driving of the turbulence within the CDL is highly anisotropic and inhomogeneous: energy is injected only at the two confining shocks of the CDL and only in the form of flows directed parallel or antiparallel to the xaxis. The amount of kinetic energy entering the CDL, as well as its decomposition into compressible and solenoidal components, depends on the orientation of the confining shock with respect to the upstream flow and thus changes with position along the confining shock. It is thus to be expected that the turbulence within the CDL is neither homogeneous, nor isotropic. Viewing angle effects in observational data are a logical consequence.
3.2.1. Mach number anisotropy
From Table 1 it can be taken that the CDL root mean square Mach number is highly anisotropic. Root mean square Mach numbers in direction transverse to the upstream flow (yzplane, M_{rms, ⊥}) are subsonic in all our experiments, with concrete Mach numbers ranging from about 0.2 to 0.8. In the direction parallel to the upstream flow, supersonic rootmeansquare Mach numbers M_{rms, ∥} are found from R7_2 onward. The ratio M_{rms, ∥}/M_{rms, ⊥} ranges from around 2 to over 11.
Indeed, analytical considerations show the difficulty associated with substantially deviating isothermal headon colliding flows, a necessary condition to reach isotropic conditions, while retaining supersonic flow conditions. To illustrate the point, look at a uniform flow hitting an oblique shock. Let α denote the angle between the flow direction and the surface normal of the shock. Decomposing the upstream and downstream flow velocities, v_{u} and v_{d}, into components parallel and perpendicular to the shock, assuming M_{u} ≫ 1, and applying isothermal shock jump conditions , one can write These relations can be recast in the form of the Mach number reduction factor F ≡ v_{d}/v_{u} and the deviation angle β between the upstream and the downstream flow direction, where the approximations hold for large M_{u} and angles α not too close to either 0 or π/2. From Eqs. (15) and (16) it is apparent that large deviation angles β come at the price of substantial reduction factors F, at least for single shock passages. Multiple “grazing” shock passages (large α) allow, in principle, for large total deviation angles β_{tot} while retaining an overall reduction factor F_{tot} close to one. For n subsequent shock passages at the same angle α one has β_{tot} = n(π/2 − α) and F_{tot} = (sinα)^{n}. The situation is illustrated in Fig. 5 for β_{tot} ≥ π/2. For example, two subsequent shock passages at α = π/4 allow for a total deviation of π/2 while retaining about half of the initial upstream Mach number (F_{tot} ≈ 0.5), while seven passages at an angle of close to 80° are needed to retain 80% of M_{u}. From Fig. 5 it can also be taken that already one shock passage at an angle α < 30° reduces M_{u} by at least 75%.
Fig. 5 Deviating a highly supersonic, isothermal flow by (at least) β_{tot} = π/2 through n (right yaxis) subsequent shock passages at angle α (xaxis; in degrees; angle between flow direction and surface normal of shock) results in a total reduction of the Mach number by a factor F_{tot} (left yaxis). 
It is this later situation, even a few shock passages at small angles α, which we believe to be responsible for the low transverse Mach numbers in our 3D slabs, despite the generally rather “grazing” angles at the confining shocks.
3.2.2. Density variance – Mach number relation
This anisotropy of M_{rms} carries over to the density variance – Mach number relation, (17)From 3D periodic box simulations of driven, isothermal, supersonic, hydrodynamical turbulence a value of b ∈ [1/3,1] is found, depending on the nature of the forcing, purely compressible (b = 1) or purely solenoidal (b = 1/3) (Padoan et al. 1997; Passot & VázquezSemadeni 1998; Federrath et al. 2010).
Using Eq. (17) for our CDL with M_{rms} the total rootmeansquare Mach number, we find b ∈ [0,1/3] (see Table 1), a range hardly overlapping with 3D periodic box results. Higher bvalues are obtained for the CDL if only the transverse component M_{rms, ⊥} is considered: around 0.6 to 0.7 for small M_{u} and up to 0.85 for large M_{u}. In the context of turbulence in flow collision zones, the density variance – Mach number relation thus seems rather indicative of the viewing angle (lineofsight along M_{rms, ∥} or M_{rms, ⊥}) than of the driving of the turbulence.
Fig. 6 2D slice (yzplane) averages as function of (scaled) xcoordinate for R22_2 for different times (ℓ = 12,17,23,29, and 34, for the blue, cyan, green, magenta, and red curves, respectively). Top panel: fractional 2D CDL area F(x). Bottom panel: root mean square Mach number in units of M_{u}. For details see text, Sect. 3.2.3. 
3.2.3. Role of distance to confining shocks: 2D slices
CDL quantities are not only anisotropic, but they also differ between regions close to the confining shocks, the only place where energy is injected into the CDL, and areas close to the central plane of the CDL. What “close to the confining shocks” really means may be debated. One interpretation could be to study flow properties as a function of distance from the confining shock. This we will not attempt here, as our spatial resolution is likely too coarse anyway to capture the details of such “boundary layer flows”.
Instead, we study how averages over 2D slices (yzplane) of the CDL vary with the position of the slice from the CDL center. Technically speaking, we take 2D slices (yzplane) at different positions along the xaxis, identify those cells within a given slice that belong to the CDL (and not the upstream flow) by checking for M_{rms} < M_{u}, and compute average quantities for each such 2D CDL area A(x). Introducing the fractional 2D CDL area F(x) = A(x)/(YZ), with YZ the yzcrosssection of the computational domain, slices close to the CDL center are characterized by F(x) = 1, whereas slices close to outermost tips of the confining shocks have F(x) ≈ 0. To plot the resulting 2D averages, we use scaled xunits, which we define by mapping the core part of the CDL with F(x) > 0.9 onto the interval [–1, 1] in scaled xunits. We expect this mapping to make sense because of (approximate) selfsimilarity. This is confirmed by Fig. 6, top panel: the curves for the fractional 2D CDL area F(x) as a function of scaled xunits (simulation R22_2) fall essentially on top of each other, although the CDL size increases from ℓ ≈ 12 to ℓ ≈ 34. Similar figures are obtained for other simulations.
Denoting parts where F(x) < 1 as “CDL boundary layer” and parts with F(x) = 1 as “CDL core region”, Fig. 6, top panel, further shows that the relative widths of the “CDL boundary layer” as compared to the “CDL core region” remain essentially constant throughout the simulation. The boundary layer does not increase with respect to the core region or vice versa, the spatial extension of each part obeys approximate selfsimilarity again. Within the “CDL boundary layer”, M_{rms} generally decreases with time for any given fixed relative xposition (Fig. 6, bottom panel). A possible explanation for the decrease could be that as the CDL gets wider, individual wiggles of the confining shocks extend farther and comprise larger volumes (“get fatter”), the surface to volume ratio decreases. In a 2D cut half way between F(x) = 0.9 and F(x) ≈ 0, the CDL patches get larger for larger CDL thickness ℓ. The distance from the patch center to its boundary – the confining shock where energy is injected – gets larger. The patch area (where turbulence must be driven) divided by the patch circumference (where energy is injected) gets larger.
Fig. 7 2D slices, average quantities for different runs (color coding as in Fig. 2) at CDL size ℓ ≈ 12. From top to bottom: fractional CDL area F(x), root mean square Mach number M_{rms}, scaled root mean square Mach number M_{rms}/M_{u}, mean density ρ_{m}, and b = σ(ρ)/(ρ_{m}M_{rms}). 
Passing from one simulation at different times to comparison of different simulations for the same CDL size, as shown in Fig. 7 for different quantities, reveals again the prominent role of the upstream Mach number. Larger upstream Mach numbers result in larger “CDL boundary layers”, relative to the “CDL core region” (Fig. 7, first panel). The root mean square Mach number in the “CDL core region” is slightly lower than the CDL mean (around 15% to 20% of M_{u}, second and third panel of Fig. 7). The density distribution is more peaked around the CDL center than the M_{rms} distribution (Fig. 7, fourth panel). Consistent with CDL mean values of Sect. 3.1, the largest densities are reached for intermediate upstream Mach numbers, not for the highest. A tentative interpretation here is that low upstream Mach numbers, on the one hand, cannot generate the necessary compression, whereas high upstream Mach numbers generate too much turbulence, leading again to somewhat lower mean densities. Numerical limitations (resolution) to maximum compression ratios may also play a role (see Sect. 3.5). Finally, the density varianceMach number relation shows less and less dependence on the slice position with increasing M_{u} (Fig. 7, fifth panel)
In summary, the analysis in terms of 2D slices illustrates the inhomogeneity of the CDL turbulence, the CDL core region is less turbulent than the CDL boundary layer. The relative weight of core and boundary layer is independent of the CDL extension.
3.3. Density probability distribution function
The probability distribution function (PDF) of turbulent density fluctuations is of key interest in connection with star formation and as such extensively covered in the literature (e.g., Padoan & Nordlund 2002; Krumholz & McKee 2005; Hennebelle & Chabrier 2011). A number of studies suggest that the density PDF of isothermal, isotropic, and homogeneous supersonic turbulence is lognormal, i.e., s = ln(ρ/ρ_{m}) follows a Gaussian distribution (e.g., VazquezSemadeni 1994; Passot & VázquezSemadeni 1998). The PDF may deviate from a lognormal distribution if additional physics is included, e.g., selfgravity, or if subsonic root mean square Mach numbers result from solenoidal forcing (Federrath et al. 2009; Kitsionas et al. 2009).
The PDFs of ln(ρ/ρ_{m}) we find for our different runs also clearly deviate from a lognormal distribution (Fig. 8, top panel). Each PDF is the average over three snapshots of the CDL with seizes ℓ ∈ [11,12].
To the right of their maxima, the PDFs can be approximately captured by a lognormal distribution, but typically have a fat tail (Fig. 8, bottom panel). The fit (dotted line) was obtained by first mirroring the right half of the PDF at its peak, then determining the FWHM w of the resulting distribution, from which σ can be obtain using w = 2(2ln(2))^{1/2}σ and assuming that the distribution is normal.
To the left of their maxima, the PDFs show a more intricate structure. For runs with M_{rms} < 1, the peak value decreases sharply before a plateau region is formed. This plateau region is likely associated with grid cells close to or even within the confining shocks of the CDL, which on numerical grounds are about 3 cells wide. Except for this plateau region, the shape of the PDFs at these low Mach numbers resembles density PDFs found by Kitsionas et al. (2009) for 3D periodic box simulations with purely solenoidal forcing at subsonic root mean square Mach numbers. As we go to higher Mach numbers, the decline of the PDF to the left of its peak value becomes shallower, reaching about a constant slope for runs with M_{u} > 11. The plateau region of the low Mach number runs turns into a steep slope for the high Mach number runs.
Fig. 8 Top panel: volumeweighted PDFs of s = ln(ρ/ρ_{m}) for the different runs, color coding as in Fig. 2. Each curve is an average over three data sets between ℓ ≈ 11 and ℓ ≈ 12. Bottom panel: centrally (at peak value of PDF) mirrored PDF with Gaussian fitting (dotted line). 
In summary, despite having purely isothermal conditions and no additional physics, like for example selfgravity, the density PDFs of our CDL always remain far from a lognormal distribution, even if the turbulence is highly supersonic (large M_{rms}).
3.4. Structure functions
Fig. 9 A 5% random perturbation of Z_{p} (30 000 realizations, at the example of Z_{p} from Boldyrev 2002, p = 1 to p = 5) causes smearing of codimensions C derived from Eq. (20) (left panel). A double peak distribution results if C is computed from Eq. (21) (middle panel). The analytical value C = 1 is shown as a red line, green lines in the left panel indicate where the 2/3 largest values in the histogram are located. In the right panel, the 2D histogram (log10 contours, spacing 0.5, spanning three orders of magnitude) of codimensions from perturbed Z_{p} (500 000 realizations) is overplotted on the codimension C = Δ/(1 − β) (color coded) in the Δβ plane (Eq. (21)). Two sharp peaks can be distinguished, at Δ ≈ 0.4 and β ≈ 0.1 (C ≈ 0.44), as well as at Δ ≈ 1.0 and β ≈ 0.5 (C ≈ 1.9). The red dot is the analytical value by Boldyrev (2002), Δ = 2/3, β = 1/3, C = 1. 
Since the seminal work of Kolmogorov (1941, K41), velocity structure functions have become an established tool for the characterization of isotropic, homogeneous, incompressible, and stationary turbulence. By contrast, the understanding of anisotropic, inhomogeneous, compressible, nonstationary, or bound turbulence is only in its infancy (Kurien & Sreenivasan 2000; Biferale et al. 2008; Byrne et al. 2011; Hansen et al. 2011; Biferale et al. 2012; Grauer et al. 2012). Characterizing our CDL data in terms of structure functions thus is rather exploratory. We nevertheless consider presentation of associated results worthwhile, to make a few first steps in this direction, to compare our setup with results from 3D periodic box simulation, and to illustrate the potential role of lineofsight effects, the latter also in view of published, observation based structure functions for molecular clouds whose physical interpretation continues to be debated (e.g., Padoan et al. 2003; Gustafsson et al. 2006; HilyBlant et al. 2008).
3.4.1. Working with low resolution CDL data
Anticipating quantitative values given in Sect. 3.4.2, here we first want to carve out two implications of the low resolution of our CDL data: the compelling use of extended selfsimilarity and why we determine the codimension of the most dissipative structure by using the one parameter expression of Eq. (20) instead of the two parameter expression of Eq. (21). We also introduce some notation used later on.
The structure function of order p is a scalar function of distance r, defined as (18)where ⟨ ... ⟩ denotes the average over all positions x within the sample and over all directed distances r from each position. Structure functions are termed longitudinal () if only orientations u ∥ r enter the average, and perpendicular () if only u ⊥ r enters. In the inertial range, both can be well approximated by power laws with scaling exponents ζ_{p}(19)The numerical determination of ζ_{p} is difficult if the inertial range is small, as is the case for our data (see Sect. 3.4.2). The extended selfsimilarity (ESS) hypothesis (Benzi et al. 1993) can remedy the situation to some extent. It basically states that the selfsimilarity of the velocity field, if it exists, extends beyond the usual inertial range, into the dissipation range. This extended range can be explored by considering not absolute scaling exponents ζ_{p} but ratios Z_{p} = ζ_{p}/ζ_{3}.
ESS scaling exponents Z_{p} can be related to the codimension C = 3 − D of the most dissipative structures (dimension D), (20)For 1D vortex filaments C = 2 (She & Leveque 1994), for sheetlike structures C = 1 (Boldyrev 2002). Intermediate values of C are found in numerical simulations (Padoan et al. 2004). From multiphase 3D periodic box simulations Kritsuk & Norman (2004) find D = 2.3, which is, as pointed out by the authors, surprisingly close to observationally determined fractal dimensions of molecular clouds (about 2.3, e.g., Elmegreen & Falgarone 1996; RomanDuval et al. 2010). Roughly speaking, velocity structure functions link the dimension of the most dissipative turbulent structures to (observable) velocity differences.
Behind the oneparameter (C) view of Eq. (20) is the twoparameter (Δ and β) expression for Z_{p} by Dubrulle (1994), (21)with C = Δ/(1 − β). For Δ = 0, the K41 law is recovered. For incompressible turbulence, She & Leveque (1994) postulated Δ = 2/3 and with β = 2/3 got C = 2. For supersonic turbulence, Boldyrev (2002) kept Δ = 2/3 but chose β = 1/3, giving C = 1. There is, however, no reason why Δ = 2/3 should also apply in the supersonic case, and why one may collapse in this way the twoparameter model of Eq. (21) into the oneparameter model of Eq. (20)^{2}.
Using Eq. (21) to determine C thus seems preferential. However, in an application like ours, another factor comes into play: Eq. (21) is much more sensitive to uncertainties in the Z_{p} than Eq. (20). Figure 9 illustrates the point. We start with theoretical values Z_{p} from Boldyrev (2002) for p = 1 to p = 5. We randomly perturb each Z_{p} by at most 5% (our accuracy requirement for Z_{p}, see Appendix B.2), thus obtaining 30 000 perturbed data sets. For each data set, we determine C as best fit to either Eqs. (20) or (21) (best fit = smallest root mean square error between analytical Z_{p}, p = 1..5, and Z_{p} of perturbed data set). If Eq. (20) is used, the 5% perturbation (uncertainty) of Z_{p} results in a smearing of C, with C mostly in a range 0.9 to 1.1 (Fig. 9, left panel). Using Eq. (21) results in a distribution with two peaks at C ≈ 0.44 and C ≈ 1.9 (Fig. 9, middle panel). A similar dichotomy was reported by Kritsuk et al. (2007b) and by Schmidt et al. (2009). Seen in the Δβplane (Fig. 9, right panel) the 5% uncertainty in Z_{p} translates into a whole band of (Δ, β) pairs and associated codimensions (Δ ≈ 0.4 and β ≈ 0.1 for C ≈ 0.44; Δ ≈ 1.0 and β ≈ 0.5 for C ≈ 1.9).
Fig. 10 “Full in” case for simulations R*_2 at ℓ ≈ 12, color coding as in Fig. 2. Top panel: structure functions (upper curves) and (lower curves), curves vertically shifted to start from one common point. Stars connected by lines indicate data used to determine ESS scaling exponents. Also indicated are slopes and as obtained by Kritsuk et al. (2007a) for 3D periodic box simulations. Bottom panel: fit (solid lines) of ESS scaling exponent as slope of versus . Filled symbols denote data used for the fit, empty symbols denote all data available. 
In view of these results, whose deeper discussion is beyond the scope of this paper, we decided to use the more robust one parameter expression Eq. (20) to compute the codimension C from the ESS exponents Z_{p}.
3.4.2. Computing structure functions, ESS scaling exponents, and codimensions from CDL data
To compute S_{p}(r) we use data at ℓ ≈ 12. Only points x within the CDL may contribute, the number of points contributing to S_{p}(r) therefore decreases with r. When a directed distance r leaves the CDL it is forbidden to reenter. The computation, detailed in Appendix B, is repeated in four different ways. Spherical averaging in Eq. (18) (“unified” case); spherical averaging but retaining only pairs (x, r) for which all distances r lie still within the CDL (“full in” case); along a direction (lineofsight, LOS) either perpendicular or parallel to the upstream flow. We stress that no density variations or radiative transfer effects are taken into account.
Given the small inertial range of our structure functions, apparent in Fig. 10, top panel, we make use of the ESS hypothesis and compute linear fits to Z_{p}, i.e., to log 10(S_{p}) versus log 10(S_{3}). For the “full in” case, the quality of the fits is illustrated in Fig. 10, bottom panel, for longitudinal structure functions and p = 2. In Appendix B, the quality of the fits is illustrated for all cases considered (“full in”, “unified”, different LOSs) for Z_{5} (Fig. B.1) and computational details are given. From Z_{p}, p = 1 to 5, a best estimate for the codimension C is derived for each run (using Eq. (20)). Our demand that the uncertainty of the linear fit to Z_{p} is smaller than 5% translates roughly into an uncertainty of ± 0.1 for C (see Sect. 3.4.1 and Fig. 9).
3.4.3. CDL structure functions: case “full in”
The “full in” case is the most likely to bare similarities, if any, with 3D periodic box simulations. This because our demand that all 20 rays lie inside the CDL (see Sect. 3.4.2) rather disfavors contributions from the vicinity of the confining shocks with their (narrow) wiggles, where turbulence is most inhomogeneous and anisotropic (see Sect. 3.2.3).
Looking at the third order structure functions in Fig. 10, top panel, five observations may be made. First, our numerical resolution (256 cells in yzdirection) hardly allows for the formation of a proper inertial range. The approximately linear dependence on r covers less than one order of magnitude in r. Scaling exponents ζ_{p} (see Eq. (19)) thus cannot be determined directly with sufficient quality. Second, this rapid leveling off seems not related to changing statistics with distance r. It takes place at distances r where the sample size has not yet decreased substantially (see Fig. A.1, top panel). Third, larger values of M_{u} result in shallower slopes of log 10(S_{3}) versus log 10(r) that level off earlier. Fourth, the scatter of third order structure functions due to changing M_{u} is larger for than for , even at small radii. typically also levels off earlier than . Fifth, slopes and from 3D periodic box simulations at M_{rms} = 6 (Kritsuk et al. 2007a), black lines in Fig. 10, appear to form an upper limit to the slab turbulence studied here, which covers a range of 0.33 < M_{rms} < 11.5 (see Table 1).
ESS scaling exponents are given in Table 2. Noteworthy are the following points. Longitudinal structure functions are consistent with codimension C ≈ 1 or D ≈ 2, i.e., indicating shocks as the most dissipative structures. They may show a tendency toward larger codimension C for larger M_{u}. The tendency stems to a good part from Z_{5}, all other Z_{p} show no clear dependence on M_{u}. The fits to Z_{5} look reasonable (see Fig. B.1, top left), but firm conclusions have to wait for better resolved simulation data. Transverse structure functions display a systematic dependency on M_{u}: the codimension decreases with increasing M_{u} from C ≈ 1.2 for M_{u} = 2 to C ≈ 0.5 for M_{u} = 33. The robustness of this decrease is supported by the fact that the related increase in Z_{1} is clearly larger than the 5% accuracy limit imposed when determining individual ESS exponents (see Appendix B).
It is interesting to note that the notion of universal scaling (see e.g., Schmidt et al. 2009) would imply especially Z_{2} to increase with M_{rms} or, equivalently, with M_{u}. Our transverse structure functions are compatible with this expectation. For the longitudinal structure functions, our data are inconclusive. One may speculate that the increase of Z_{5} with M_{u} rather points toward a decrease of Z_{2} – a behavior incompatible with universal scaling.
The above findings are robust against increasing CDL size. For R22_2, the simulation which we integrated longest, we computed ESS scaling exponents for CDL sizes ℓ = 12,17,23,29, and 34. The resulting maximum spread of Z_{p} is less than 7% for p = 1 and less than 5% for p = 2 to p = 5. Similar results were obtained for two other simulations, R11_2 and R33_2.
Longitudinal (top) and transverse (middle) ESS scaling exponents Z_{p} and best estimates for codimension C, case “full in” for the different runs at ℓ ≈ 12.
3.4.4. CDL structure functions: beyond the “full in” case
The CDL being not strictly periodic and its turbulence asymmetric and inhomogeneous, structure functions are expected to depend on how the averages in Eq. (18) are taken. Apart from the “full in” case presented above, we looked at three more cases: “unified” and LOSs parallel and perpendicular to the upstream flows (see Sect. 3.4.2). For the perpendicular LOS we tested that results do not depend on whether the y or zdirection is analyzed, thus only the ydirection is shown. Third order structure functions are shown in Fig. 11, ESS scaling exponents are given in Table B.1.
The “unified” case by and large resembles the “full in” case. The structure functions look similar (Figs. 10 and 11, top panel), the quality of ESS fits is comparable (Fig. B.1, columns one and two). Codimensions tend to be smaller in the “unified” case (0.9 ≤ C ≤ 0.4), indicating a dimension D > 2 for the most dissipating structures. This seems plausible as the “unified” case covers the CDL more completely than the “full in” case, especially with regard to the highly compressible regions close to the confining shocks. As in the “full in” case, there is a tendency of the codimension of longitudinal (transverse) structure functions to increase (decrease) with increasing M_{u}.
Structure functions taken along a LOS parallel or perpendicular to the upstream flow look widely different among themselves (Fig. 11, middle and bottom panel) and with respect to both the “unified” and “full in” case. First, the spread induced by M_{u} is larger for the perpendicular LOS than for the parallel LOS. We speculate that this may be related to our demand that a ray must not reenter the CDL (see Sect. 3.4.2), which gains relevance as M_{u} increases, and implies that regions close to the confining shocks tend to be neglected as r increases. Second, longitudinal structure functions from a parallel LOS do not level off but keep increasing over the distance considered. This seems plausible as with increasing distance the velocity difference in Eq. (18) approaches the difference between the two post shock flow velocities at each of the confining shocks.
Best ESS fits are of reasonable quality for transverse structure functions and low to intermediate M_{u} (Fig. B.1). The concrete values of ESS scaling exponents given in Table B.1 (middle and right columns) reveal the crucial role of the orientation of the LOS with respect to the upstream flow: transverse ESS scaling exponents indicate codimensions C > 1 for a LOS parallel to the upstream flow and C < 1 for a perpendicular LOS. For a parallel LOS there is again a tendency toward smaller codimension with increasing M_{u}.
Best ESS fits for the longitudinal structure functions are typically of mediocre quality and difficult to interpret. While we give the corresponding data in the appendix (Table B.1), we renounce at their discussion. We only mention two features whose further analysis may be of interest, but is beyond the scope of the present paper. First, for the parallel LOS best fits are typically obtained at intermediate values of (Fig. B.1, top row, third panel), which translates to intermediate values of r. Second, for perpendicular LOSs there is a steeper (shallower) linear slope for larger (smaller) values of (Fig. B.1, top row, fourth panel). As M_{u} increases, the shallower slope gradually vanishes and only the steeper slope prevails.
In summary, our results indicate that the LOS plays a crucial role for the CDL structure functions and derived quantities, notably the codimension C.
Fig. 11 Structure functions (solid) and (dashed) for the case “unified” (top panel), as well as computed along the xdirection (middle panel) and ydirection (bottom panel), i.e., parallel and perpendicular to the upstream flow, for simulations R*_2 at ℓ = 12. Color coding as in Fig. 2. Also indicated are slopes and as obtained by Kritsuk et al. (2007a) from 3D periodic box simulations. Curves are vertically shifted, to start from one common point. 
3.5. Numerical resolution
It is clear from the literature (e.g., Kritsuk et al. 2007a; Federrath et al. 2010) that the present simulations are at the brink of the necessary resolution to get physically meaningful results in the context of isothermal supersonic turbulence. Computational costs limit concrete numerical resolution studies.
An obvious resolution dependence that can be taken from Table 1, and that was already found for the 2D case (FW06), is that lower resolution results in higher M_{rms} for the same M_{u} and for the same CDL size (ℓ ≈ 12). This although there is a tendency toward smaller driving efficiencies (thus less energy input into the CDL) at lower resolutions. We can only speculate about potential reasons for this observed behavior, as to why we find higher resolution to result in lower M_{rms} (and not higher M_{rms}, as in the subsonic case), despite increased energy input.
One potential explanation we can think of, and which we put forward already in the 2D case: in supersonic turbulence dissipation is dominated by shocks and with increasing resolution more shocks are resolved, thus dissipation increases, thus M_{rms} decreases. This mechanism we would expect to act primarily on the parallel flow component, as only M_{rms, ∥} ≫ 1.
Another possible explanation that comes to mind could be that finer resolution tends to promote coupling of parallel and transverse directions. The latter thus may be exploited more efficiently for energy dissipation, predominantly viscous dissipation as M_{rms, ⊥} < 1. Total energy dissipation (parallel plus transverse directions) thus may increase with increasing resolution. Further potential explanations may exist.
The first idea, enhanced dissipation in shocks, is attractive as it acts predominantly on M_{rms, ∥}, which dominates M_{rms} by far. The second idea, enhanced viscous dissipation in transverse directions, fits well with our data in that M_{rms, ∥}/M_{rms, ⊥} decreases with increasing resolution, indicating better coupling of transverse and parallel directions. Within the frame of our model set up it is not easy to proof that any of these two ideas are indeed at work, admittedly leaving us for the time being with an unsatisfying situation.
It would be interesting to check whether a similar dependence of M_{rms} on resolution exists for driven turbulence in 3D periodic box simulations. This could be done by monitoring the amount of energy that has to be injected to maintain a prescribed level of M_{rms}. If the first of the above potential explanations were true, we would expect that for higher resolution more energy has to be injected per time to maintain a prescribed level of turbulence. By contrast, the second of the above potential explanations may hardly manifest itself in the 3D periodic box case, as driving there is typically isotropic and thus coupling among directions is less of an issue than in the present study. The authors are not aware of any dedicated study in this direction.
Besides M_{rms}, other quantities show some resolution dependence as well. Clearly increasing with increasing resolution is σ(ρ)/ρ_{m} (last column of Table 1). As ρ_{m} shows no clear tendency, this increase is an increase in density variance with resolution. The increase of peak densities with increasing resolution is well known (Hennebelle & Audit 2007; Kitsionas et al. 2009; Federrath et al. 2010) and expected on numerical grounds.
In summary, we expect our results to essentially hold qualitatively but to change somewhat quantitatively if we were to go to higher resolution. In particular, we would expect somewhat lower M_{rms} for the same M_{u} and a larger density variance.
4. Discussion
4.1. 3D slabs versus 2D slabs
The relevance of dimensionality for headon supersonic collision of isothermal flows was already pointed out in Sect. 3.1: while compression of the collision zone scales as in 1D, mean densities of the CDL are largely independent of M_{u} in 2D and 3D. Comparing the 3D numerical results presented here with corresponding 2D results (FW06) further shows that for the same upstream conditions the interaction zone is more turbulent in 3D than in 2D: M_{rms} is about 25% of M_{u} in 3D, but only about 20% in 2D. For the same upstream conditions, driving efficiency is considerably higher in 3D than in 2D (see Fig. 12) and the power law dependence of the driving efficiency on the upstream Mach number is steeper: β_{3} = −0.7 in 2D but β_{3} = −1.17 in 3D.
We attribute the larger driving efficiency in 3D (a larger fraction of the upstream kinetic energy density traverses the confining shocks of the CDL unthermalized) to the different shock geometry. In 2D, with the CDL extending to infinity in zdirection, the confining shocks resemble a corrugated sheet. In 3D, they look more like a cardboard egg wrapping. Regions where the confining shocks are perpendicular to the upstream flow, thus maximum thermalization occurs, are linelike in the 2D case but pointlike in the 3D case (see Fig. 1).
Numerical resolution we exclude as an explanation, as 2D simulations for M_{u} = 22 and three different resolutions yield a range of f_{eff} ∈ [0.55,0.62] (FW06), whereas the corresponding 3D value is f_{eff} = 0.84. This for about the same xextension of the CDL and about the same coverage (up to a factor of two) of this extension in terms of grid cells.
Fig. 12 Role of dimensionality for driving efficiency f_{eff} as function of upstream Mach number M_{u}. Shown are simulations R*_2 (red stars) and corresponding values for 2D slabs (dark red diamonds, simulations R*_0.2.4 in FW06). 
4.2. 3D slabs versus 3D periodic boxes
The main difference of the simulations presented here to 3D periodic box simulations lies in the driving of the turbulence – and resulting consequences. In 3D periodic box models, energy is typically injected in a controlled, rather homogeneous and isotropic way, although concrete realizations differ. The 3D slab model may be regarded as the other extreme: energy input into the turbulent interaction zone only occurs at the two confining shocks and this in a comparatively uncontrolled way concerning both the amount of energy and the form of the driving (compressible/solenoidal and driving wave length). The kinetic energy flow entering the interaction zone is predominantly directed parallel to the upstream flows. Its absolute amount is modulated by the wiggling of the confining shocks, whose spatial scale in turn correlates with the spatial extension of the interaction zone.
The different forcing is reflected in the character of the turbulence: isotropic, homogeneous, and stationary in the case of 3D periodic box models, anisotropic, inhomogeneous, and with some time evolution for 3D flow collision zones. We even argue (Sect. 3.2.1) that for headon colliding, isothermal flows it is far from trivial, if not impossible, to reach isotropy while retaining M_{rms} ≫ 1. As a corollary, we see some danger in neglecting the concrete large scale driving of the turbulence and replace it by an isotropic random driving instead. For uniform isothermal, headon colliding flows as studied here, the anisotropy and inhomogeneity of the turbulence are a direct consequence of the large scale driving. It seems plausible that a similar imprint of the large scales on the small scales is also present in more complex flows, e.g., accretion onto compact objects (Walder et al. 2010). This is not to say that isothermal isotropic supersonic turbulence cannot exist. But we caution that its realization may not be straight forward, at least not in the context of large scale colliding flows.
Density PDFs deviate from lognormal already in 3D periodic box simulations of driven isothermal turbulence if the driving has a substantial compressible component (Schmidt et al. 2009; Federrath et al. 2010). It thus seems plausible that our density PDFs deviate from lognormal. However, the deviations we find are more pronounced than those seen in compressibly forced 3D periodic box simulations. A potential reason could be the inhomogeneity of the turbulence in the collision zone. The more violent turbulence close to the confining shocks may be responsible not only for high compressions but also for pronounced low density voids, as visible in the 2D case, Fig. 1, left panel. Density PDFs for individual 2D slices may shed more light on this issue but were beyond the scope of the present paper. For molecular clouds, if they indeed result from colliding flows, the highdensity tail we find suggests that current star formation models based on lognormal density distributions may severely underestimate high mass star formation, even more than already postulated in Schmidt et al. (2009) on the basis of their simulations.
Regarding ESS scaling exponents, Table 3 shows that our values lie well within the range of published values. They basically range from the theoretical model by She & Leveque (1994) (C = 2) to numerical 3D periodic box results by Schmidt et al. (2009). The very low values for Z_{4} and Z_{5} found in the later study are, however, never reached in our data. Nevertheless, our data are much closer to these values than any other of the listed data.
The range of Z_{p} values we observe results, on the one hand, from the dependence of Z_{p} on M_{u}. The dependence is particularly clear for transverse structure functions. Here, Z_{2} increases with M_{u}, in line with the notion of universal scaling. Longitudinal structure functions show rather inconclusive results in this respect, possibly contradicting expectations from universal scaling. Whether isothermal turbulence, in flow collision zones or 3D periodic box simulations, indeed displays nonuniversal scaling properties, possibly due to an additional degree of freedom related to the large scale forcing as suggested by Schmidt et al. (2009), remains to be seen.
The other reason for the large range of Z_{p} are LOS effects, a point we take up again in Sect. 4.3.2. Here we only note that parallel and perpendicular LOSs (rows 3 and 4 in Table 3) display nearly complementary ranges in terms of Z_{p}, and that it is mainly the parallel LOS which is responsible for Z_{p} values close to theoretical values by She & Leveque (1994). This kind of LOS effects is likely absent in 3D periodic box simulations.
Comparing ESS scaling exponents Z_{p}, p = 1 to p = 5, following Table 4 in Schmidt et al. (2009).
4.3. Toward real objects?
The physical model examined here is clearly too simple for direct comparison with real world observations. Nevertheless, we consider it worthwhile to contemplate on potential realworld implications of two results, putting them at the same time into perspective by speculating on effects of some neglected physics. Finally, we reverse the perspective and ask for which classes of real objects 3D slab studies may provide useful physical insight.
4.3.1. Driving the turbulence: high M_{rms}, high M_{u}?
The first point we want to discuss more thoroughly is the rather poor conversion from M_{u} to M_{rms}, the latter being only about 10% to 20% of the former in the core region of the CDL. For molecular clouds, where observations indicate M_{rsm} ≈ 5–50 (Zuckerman & Evans 1974; McKee & Ostriker 2007; Klessen 2011; Polychroni et al. 2012), and assuming that molecular clouds result from large scale flow collisions, this would imply M_{u} ≈ 40–500 within the frame of our model.
Several questions may be asked. How could such fast flows be produced? Could inclusion of additional physics within the frame of colliding flows improve the conversion from M_{u} to M_{rsm}? Are driving mechanisms other than colliding flows compelling? We will mostly dwell on the second question.
The inclusion of radiative cooling may or may not improve conversion from M_{u} to M_{rsm}. Simulations where the cooling limit, and thus the lowest temperature within the CDL, is equal to the upstream flow temperature show less turbulence (smaller M_{rms}) than isothermal simulations with the same upstream Mach number (Walder & Folini 2000b; Folini et al. 2010). If significant cooling layers form at the confining shocks, the associated large thermal postshock pressure will tend to straighten the confining shocks. A positive feedback results: the more straight the confining shocks, the more upstream kinetic energy is thermalized at these shocks, the higher the postshock temperature, the stronger the corresponding thermal pressure and the longer the cooling time for the shockheated gas. And the less kinetic energy is available for driving the turbulence.
The situation may be different during the very early collision phase, when the CDL is still very thin, and especially if the post shock temperature falls within a temperature range where cooling instabilities occur, i.e., where cooling becomes more and more efficient as the temperature decreases. The early (thin) CDL as a whole becomes very unstable under such conditions (Walder & Folini 2000b; Heitsch et al. 2005, 2011; Folini et al. 2010). The situation may also change if the cooling limit is decidedly below the temperature of the upstream flow (Pittard et al. 2005). The CDL then can cool to temperatures well below the upstream flow temperature, the sound speed within the CDL decreases and, consequently, the Mach number increases. Except for this last possibility, it seems rather unlikely that inclusion of radiative cooling improves the conversion from M_{u} to M_{rsm}. Isothermal conditions rather set an upper limit.
Other physics of potential relevance for M_{rms} (and beyond) but not covered by our model include magnetic fields, which may increase the CDL turbulence by additional energy input through reconnection or decrease the turbulence as the field guides the flow direction. More specifically for molecular clouds, inclusion of selfgravity may augment the turbulence within the cloud. MHD waves generated by moving accreting cores may stir the turbulence (Folini et al. 2004). Turbulence may also be driven through jets, winds, and radiation from young or forming stars.
All possibilities just listed to augment M_{rms} have in common that they drive the turbulence from within the CDL. For molecular clouds this may seem in contradiction with observation based results favoring driving on the scale of the cloud (Ossenkopf & Mac Low 2002; Heyer et al. 2006; Brunt et al. 2009; RomanDuval et al. 2011). These results are also somewhat challenging with regard to our findings. For our CDL, the wiggling of the confining shocks modulates the incoming flows. It seems plausible that the spatial scale of this modulation affects, or even sets, the driving scale of the turbulence. If so, this would imply an energy injection scale that is smaller than the CDL size, as the spatial scale of this wiggling is smaller than the spatial extension of the CDL. The scale of the wiggling increases, however, in proportion to the CDL size as the CDL grows (see FW06). Also, our simulations suggest that spatial modulation of the incoming flow comprises not one scale but some continuum of scales.
A speculative solution to the above, somewhat conflicting, results and demands for molecular clouds could be that various driving mechanisms and associated driving wave lengths coexist: external large scale and internal small scale driving. The latter would add to M_{rms}, thus would help to overcome the poor conversion from M_{u} to M_{rms} and the associated demand for very high Mach number flows. It also would help to make the turbulence more isotropic, something we have shown is difficult to achieve within the frame of isothermal colliding flows. Large scale forcing would equally add to M_{rms} and leave a characteristic largescaledriving imprint in the turbulence.
4.3.2. Viewing angle effects?
The second finding that in our opinion deserves discussion is the strong anisotropy of CDL velocities. It means that lineofsight velocities of an observer, and any derived quantities, will depend on the viewing angle of the CDL. The same CDL can display supersonic or subsonic lineofsight velocities, depending on whether the lineofsight of the observer is oriented parallel or perpendicular to the incoming flows. The density Mach number relation suffers from the same problem. If taken as indicative of the nature of the driving (b = 1/3 for solenoidal driving and b = 1 for compressible driving, from isothermal 3D periodic box simulations, see Federrath et al. 2008), the viewing angle will codecide whether one concludes for the same CDL that its driving is rather compressible or solenoidal.
Equally affected are ESS scaling exponents Z_{p}. For transverse structure functions and derived codimensions we have shown that results depend crucially on the adopted viewing angle: codimensions C > 1 (C < 1) are obtained for a lineofsight parallel (perpendicular) to the upstream flow. As stressed before, these results are preliminary in that they do not take into account radiative transfer, density effects, or projection effects. That such factors matter is well established in the literature (Stutzki et al. 1998; Brunt & Mac Low 2004; Sánchez et al. 2005; Ossenkopf et al. 2006; Esquivel et al. 2007; Brunt et al. 2009; Federrath et al. 2010). Further analysis of our data in this direction is envisaged but beyond the scope of the present paper.
Despite this cautionary remark, it is tempting to dwell a moment longer on the subject and add a concrete illustration of our point that lineofsight effects should be considered as one more factor when interpreting observational data. As an example, consider simulation R5_2, whose transverse structure functions are generally well behaved. In terms of Z_{p} we find Z_{1} = 0.38 and Z_{5} = 1.49 when looking at the CDL parallel to the upstream flow, but Z_{1} = 0.45 and Z_{5} = 1.37 when the lineofsight is perpendicular to the upstream flow. In the first case, with values roughly similar to observation based results for the Polaris Flare (HilyBlant et al. 2008), we deduce a codimension of C = 1.5 and may conclude that the forcing is predominantly solenoidal and, possibly, that the turbulence is not too compressible. If we look at the same data but from a different viewing angle, now perpendicular to the upstream flow, we deduce a codimension of C = 0.9 and may conclude that the turbulence is highly compressible and the forcing rather compressible than solenoidal.
4.3.3. 3D slabs: a useful concept for real objects?
In a number of astrophysical objects, large scale colliding flows play a key role, i.e., flows with a dominant bulk velocity component on top of any additional flow structure. In fact, the present study was motivated by numerical simulations of entire objects by the authors (Folini & Walder 2000; Dumm et al. 2000; Folini & Walder 2002; Harper et al. 2005; Walder et al. 2005, 2008; Georgy et al. 2013). We advocate that insight into the physics of such collision zones can be gained from 3D slab studies. The present study, for example, shows that a richly structured interaction zone exists even for rather simple physics and although the colliding flows are essentially reduced to their bulk properties. Of course, real colliding flows likely harbor more structure and physics. The relative importance of such additional complexity may be gradually addressed within the context of 3D slabs.
Taking the perspective of astrophysical objects, 3D slab studies may be useful for the wind collision zone in massive binaries. Concrete questions range from the Xray and nonthermal emission of such zones (e.g., Pittard & Parkin 2010; Zhekov 2012; De Becker & Rauw 2013) to dust formation near periastron passage, presumably in the collision zone and despite the strong stellar radiation fields (Tuthill et al. 1999; Marchenko et al. 1999; Cherchneff et al. 2000; Folini & Walder 2002; Dougherty et al. 2005; Williams et al. 2012). Such studies clearly require more physics than covered in this paper. Also, massive stellar winds are likely clumped (e.g., Owocki et al. 1988; Moffat et al. 1988), which can affect the turbulence in the collision zone (e.g., Walder & Folini 2002; Pittard 2007; Pittard et al. 2009, 2010). Another large class of objects where 3D slabs might provide some insight are internal collisions in jets, from young stars (e.g., Flower et al. 2003; Cunningham et al. 2006, 2011) to highenergy objects, like gamma ray burst (Mimica et al. 2007; Zitouni et al. 2008; Bošnjak et al. 2009; Mimica & Aloy 2010; Granot 2012). In the later case, however, relativistic effects are likely relevant. A third context where 3D slab studies might be useful is star formation on galactic scales. Turbulence in the wake of flow collision is one suggested explanation (e.g., Gabor & Bournaud 2014) for the observed delay and potential suppression of star formation for z > 2 (e.g., Bouché et al. 2010; Weinmann et al. 2012). The case of 3D slabs and molecular clouds we briefly touched in Sect. 4.3.1. A multiphase medium is likely relevant for such studies (e.g., Hennebelle & Pérault 1999; Hennebelle et al. 2007; Gray & Scannapieco 2013). Also, the flows probably already carry (turbulent) structure, whose effect on the collision zone remains to be clarified. In fact, we consider it a most interesting question whether gradual inclusion of such and other effects allows to recover the 3D periodic box view on molecular cloud turbulence (see Sect. 4.2) in the context of 3D slabs. We are not aware of any such studies.
5. Summary and conclusions
We have performed 3D simulations of headon colliding isothermal flows with upstream Mach numbers in the range 2 ≤ M_{u} ≤ 43 and presented a first analysis of the turbulence in the resulting interaction zone, the CDL (cold dense layer). We find that the characteristics of the CDL turbulence deviate markedly from what is obtain in corresponding 3D periodic box simulations.
As in 2D, approximate selfsimilar scaling relations also hold in 3D for the mean density, root mean square Mach number, driving efficiency, and energy dissipation in the CDL, albeit with different numerical constants and a different (steeper) scaling exponent for the driving efficiency. For the same M_{u}, the driving efficiency is larger in 3D than in 2D because of the different shock geometry (“egg wrapping” in 3D, “corrugated sheet” in 2D).
Density PDFs are not lognormal but show fat tails at high densities and a “twopowerlaw” composite for low densities.
The CDL turbulence is inhomogeneous and anisotropic. Root mean square Mach numbers are around 15% (CDL center) to 25% (CDL average) of M_{u}. They are always larger parallel to the upstream flow than in perpendicular direction, where they never reach supersonic values. The density variance – Mach number relation inherits this anisotropy. Isotropization of the flow hardly takes place and we argue, in fact, that isotropization is hard to achieve while retaining substantially supersonic root mean square Mach numbers. Even in the (small scale) central region of the CDL the turbulence remains anisotropic, even here the imprint of the large scale driving is still felt.
The viewing angle of the CDL is shown to play a prominent role for different turbulence characteristics like ESS scaling exponents or the density variance – Mach number dependence. ESS scaling exponents of transverse structure functions imply for the most dissipative structures a dimension D < 2 if computed for a lineofsight parallel to the upstream flow, but D > 2 if taken along a lineofsight perpendicular to the upstream flow. We suggest to keep this in mind when interpreting observational data.
Our results show that even our very simple model setup results in a very richly structured, turbulent interaction zone. The physical system studied thus opens a perspective on supersonic turbulence that is rewarding and complementary to 3D periodic box simulations. Under what conditions both perspectives may be united we consider an interesting question for the future. The finite size of the interaction zone is a challenge with regard to data analysis and interpretation of results. But this finite size makes the data also potentially very interesting, as real objects are of finite size as well, and therefore likely incorporate boundary effects. Moreover, bulk flows in real astrophysical objects likely carry some structure (are inhomogeneous) and require a more elaborate physical description than the isothermal approach taken here. The assumption of headon collision equally needs to be relaxed and the role of numerical resolution for supersonic turbulence needs investigation. The robustness of our findings against these factors remains to be tested. Nevertheless, we consider the present study a useful contribution for the physical understanding of complex real objects.
For a good discussion, see Schmidt et al. (2009), who follow Frisch (1995) and Pan et al. (2008).
Acknowledgments
The authors would like to thank an anonymous referee for his/her constructive criticism. RW and DF acknowledge support from the French Stellar Astrophysics Program PNPS. We acknowledge support from the Pôle Scientifique de Modélisation Numérique (PSMN) and from the Grand Équipement National de Calcul Intensif (GENCI), project number x2012046960. The research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework (FP7/20072013)/ERC grant agreement no. 320478.
References
 Agertz, O., Teyssier, R., & Moore, B. 2009, MNRAS, 397, L64 [Google Scholar]
 Aluie, H. 2011, Phys. Rev. Lett., 106, 174502 [NASA ADS] [CrossRef] [Google Scholar]
 Anninos, P., & Norman, M. L. 1996, ApJ, 460, 556 [NASA ADS] [CrossRef] [Google Scholar]
 Audit, E., & Hennebelle, P. 2005, A&A, 433, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Audit, E., & Hennebelle, P. 2010, A&A, 511, A76 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 BallesterosParedes, J., & Mac Low, M. 2002, ApJ, 570, 734 [NASA ADS] [CrossRef] [Google Scholar]
 Banerjee, S., & Galtier, S. 2013, Phys. Rev. E, 87, 013019 [NASA ADS] [CrossRef] [Google Scholar]
 Benzi, R., Ciliberto, S., Tripiccione, R., et al. 1993, Phys. Rev. E, 48, 29 [Google Scholar]
 Berger, M. J. 1985, Lect. Appl. Math., 22, 31 [Google Scholar]
 Biferale, L., Lanotte, A. S., & Toschi, F. 2008, Physica D, 237, 1969 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Biferale, L., Musacchio, S., & Toschi, F. 2012, Phys. Rev. Lett., 108, 164501 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Blondin, J. M., & Marks, B. S. 1996, New Astron., 1, 235 [Google Scholar]
 Boldyrev, S. 2002, ApJ, 569, 841 [NASA ADS] [CrossRef] [Google Scholar]
 Boldyrev, S., Nordlund, Å., & Padoan, P. 2002a, ApJ, 573, 678 [NASA ADS] [CrossRef] [Google Scholar]
 Boldyrev, S., Nordlund, Å., & Padoan, P. 2002b, Phys. Rev. Lett., 89, 1102 [Google Scholar]
 Boris, J., Grinstein, F., Oran, E., & Kolbe, R. 1992, Fluid. Dynam. Res, 10, 199 [Google Scholar]
 Bouché, N., Dekel, A., Genzel, R., et al. 2010, ApJ, 718, 1001 [NASA ADS] [CrossRef] [Google Scholar]
 Bošnjak, Ž., Daigne, F., & Dubus, G. 2009, A&A, 498, 677 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Brunt, C. M., & Mac Low, M.M. 2004, ApJ, 604, 196 [NASA ADS] [CrossRef] [Google Scholar]
 Brunt, C. M., Heyer, M. H., & Mac Low, M.M. 2009, A&A, 504, 883 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Byrne, D., Xia, H., & Shats, M. 2011, Phys. Fluids, 23, 095109 [NASA ADS] [CrossRef] [Google Scholar]
 Cherchneff, I., Le Teuff, Y. H., Williams, P. M., & Tielens, A. G. G. M. 2000, A&A, 357, 572 [NASA ADS] [Google Scholar]
 Colella, P. 1990, J. Comput. Phys., 87, 171 [NASA ADS] [CrossRef] [Google Scholar]
 Colella, P., & Woodward, P. R. 1984, J. Comput. Phys., 54, 174 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Courant, R., & Friedrichs, K. O. 1976, Supersonic flow and shock waves (New York: SpringerVerlag) [Google Scholar]
 Cunningham, A. J., Frank, A., & Blackman, E. G. 2006, ApJ, 646, 1059 [NASA ADS] [CrossRef] [Google Scholar]
 Cunningham, A. J., Klein, R. I., Krumholz, M. R., & McKee, C. F. 2011, ApJ, 740, 107 [NASA ADS] [CrossRef] [Google Scholar]
 de Avillez, M. A., & Breitschwerdt, D. 2007, ApJ, 665, L35 [NASA ADS] [CrossRef] [Google Scholar]
 De Becker, M., & Rauw, F. 2013, A&A, 558, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Domaradzki, J. A. 2010, Int. J. Comput. Fluid Dyn., 24, 435 [CrossRef] [Google Scholar]
 Dougherty, S. M., Beasley, A. J., Claussen, M. J., Zauderer, B. A., & Bolingbroke, N. J. 2005, ApJ, 623, 447 [NASA ADS] [CrossRef] [Google Scholar]
 Downes, T. P. 2012, MNRAS, 425, 2277 [NASA ADS] [CrossRef] [Google Scholar]
 Dubrulle, B. 1994, Phys. Rev. Lett., 73, 959 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Dumm, T., Folini, D., Nussbaumer, H., et al. 2000, A&A, 354, 1014 [NASA ADS] [Google Scholar]
 Elmegreen, B. G., & Falgarone, E. 1996, ApJ, 471, 816 [NASA ADS] [CrossRef] [Google Scholar]
 Esquivel, A., Lazarian, A., Horibe, S., et al. 2007, MNRAS, 381, 1733 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Fan, Y. Z., & Wei, D. M. 2004, ApJ, 615, L69 [NASA ADS] [CrossRef] [Google Scholar]
 Federrath, C., & Klessen, R. S. 2013, ApJ, 763, 51 [NASA ADS] [CrossRef] [Google Scholar]
 Federrath, C., Klessen, R. S., & Schmidt, W. 2008, ApJ, 688, L79 [NASA ADS] [CrossRef] [Google Scholar]
 Federrath, C., Klessen, R. S., & Schmidt, W. 2009, ApJ, 692, 364 [NASA ADS] [CrossRef] [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]
 Flower, D. R., Le Bourlot, J., Pineau des Forêts, G., & Cabrit, S. 2003, MNRAS, 341, 70 [NASA ADS] [CrossRef] [Google Scholar]
 Folini, D., & Walder, R. 2000, Ap&SS, 274, 189 [NASA ADS] [CrossRef] [Google Scholar]
 Folini, D., & Walder, R. 2002, in Interacting Winds from Massive Stars, eds. A. F. J. Moffat, & N. StLouis, ASP Conf. Ser., 260, 605 [Google Scholar]
 Folini, D., & Walder, R. 2006, A&A, 459, 1 (fW06) [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Folini, D., Walder, R., Psarros, M., & Desboeufs, A. 2003, in Stellar Atmosphere Modeling, eds. I. Hubeny, D. Mihalas, & K. Werner, ASP Conf. Ser., 288, 433 [Google Scholar]
 Folini, D., Heyvaerts, J., & Walder, R. 2004, A&A, 414, 559 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Folini, D., Walder, R., & Favre, J. M. 2010, in Numerical Modeling of Space Plasma Flows, Astronum2009, eds. N. V. Pogorelov, E. Audit, & G. P. Zank, ASP Conf. Ser., 429, 9 [Google Scholar]
 Frisch, U. 1995, Turbulence (Cambridge University Press) [Google Scholar]
 Gabor, J. M., & Bournaud, F. 2014, MNRAS, 437, L56 [NASA ADS] [CrossRef] [Google Scholar]
 Galtier, S., & Banerjee, S. 2011, Phys. Rev. Lett., 107, 134501 [NASA ADS] [CrossRef] [Google Scholar]
 Garnier, E., Mossi, M., Sagaut, P., Comte, P., & Deville, M. 1999, J. Comput. Phys., 153, 273 [NASA ADS] [CrossRef] [Google Scholar]
 Gazol, A., & Kim, J. 2010, ApJ, 723, 482 [NASA ADS] [CrossRef] [Google Scholar]
 Gazol, A., & Kim, J. 2013, ApJ, 765, 49 [NASA ADS] [CrossRef] [Google Scholar]
 Gazol, A., Kim, J., VázquezSemadeni, E., & Luis, L. 2007, in SINS – Small Ionized and Neutral Structures in the Diffuse Interstellar Medium, eds. M. Haverkorn, & W. M. Goss, ASP Conf. Ser., 365, 154 [Google Scholar]
 Georgy, C., Walder, R., Folini, D., Bykov, A., Marcowith, A., & Favre, J. M. 2013, A&A, 559, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Glover, S. C. O., Federrath, C., Mac Low, M.M., & Klessen, R. S. 2010, MNRAS, 404, 2 [NASA ADS] [Google Scholar]
 Gong, H., & Ostriker, E. C. 2011, ApJ, 729, 120 [NASA ADS] [CrossRef] [Google Scholar]
 Granot, J. 2012, MNRAS, 421, 2467 [NASA ADS] [CrossRef] [Google Scholar]
 Grauer, R., Homann, H., & Pinton, J.F. 2012, New J. Phys., 14, 063016 [NASA ADS] [CrossRef] [Google Scholar]
 Gray, W. J., & Scannapieco, E. 2013, ApJ, 768, 174 [NASA ADS] [CrossRef] [Google Scholar]
 Gustafsson, M., Brandenburg, A., Lemaire, J. L., & Field, D. 2006, A&A, 454, 815 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hansen, C. E., McKee, C. F., & Klein, R. I. 2011, ApJ, 738, 88 [NASA ADS] [CrossRef] [Google Scholar]
 Harper, G. M., Brown, A., Bennett, P. D., et al. 2005, AJ, 129, 1018 [NASA ADS] [CrossRef] [Google Scholar]
 Heitsch, F., Burkert, A., Hartmann, L. W., Slyz, A. D., & Devriendt, J. E. G. 2005, ApJ, 633, L113 [NASA ADS] [CrossRef] [Google Scholar]
 Heitsch, F., Slyz, A. D., Devriendt, J. E. G., Hartmann, L. W., & Burkert, A. 2006, ApJ, 648, 1052 [NASA ADS] [CrossRef] [Google Scholar]
 Heitsch, F., Stone, J. M., & Hartmann, L. W. 2009, ApJ, 695, 248 [NASA ADS] [CrossRef] [Google Scholar]
 Heitsch, F., Naab, T., & Walch, S. 2011, MNRAS, 415, 271 [NASA ADS] [CrossRef] [Google Scholar]
 Hennebelle, P., & Audit, E. 2007, A&A, 465, 431 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hennebelle, P., & Chabrier, G. 2011, ApJ, 743, L29 [NASA ADS] [CrossRef] [Google Scholar]
 Hennebelle, P., & Pérault, M. 1999, A&A, 351, 309 [NASA ADS] [Google Scholar]
 Hennebelle, P., Audit, E., & MivilleDeschênes, M.A. 2007, A&A, 465, 445 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hennebelle, P., Banerjee, R., VázquezSemadeni, E., Klessen, R. S., & Audit, E. 2008, A&A, 486, L43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Heyer, M. H., Williams, J. P., & Brunt, C. M. 2006, ApJ, 643, 956 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 HilyBlant, P., Falgarone, E., & Pety, J. 2008, A&A, 481, 367 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Inoue, T., & Fukui, Y. 2013, ApJ, 774, L31 [NASA ADS] [CrossRef] [Google Scholar]
 Inoue, T., & Inutsuka, S.I. 2009, ApJ, 704, 161 [NASA ADS] [CrossRef] [Google Scholar]
 Inoue, T., & Inutsuka, S.I. 2012, ApJ, 759, 35 [NASA ADS] [CrossRef] [Google Scholar]
 Jasak, H., & Weller, H. 1995, internal Report, CFD research group, Imperial College, London [Google Scholar]
 Kaiser, C. R., Sunyaev, R., & Spruit, H. C. 2000, A&A, 356, 975 [NASA ADS] [Google Scholar]
 Kang, H., Ryu, D., Cen, R., & Song, D. 2005, ApJ, 620, 21 [NASA ADS] [CrossRef] [Google Scholar]
 Kitsionas, S., Federrath, C., Klessen, R. S., et al. 2009, A&A, 508, 541 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Klar, J. S., & Mücket, J. P. 2012, MNRAS, 423, 304 [NASA ADS] [CrossRef] [Google Scholar]
 Klessen, R. S. 2011, in EAS Pub. Ser. 51, eds. C. Charbonnel, & T. Montmerle, 133 [Google Scholar]
 Klessen, R. S., & Hennebelle, P. 2010, A&A, 520, A17 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kolmogorov, A. 1941, Akademiia Nauk SSSR Doklady, 30, 301 (K41) [Google Scholar]
 Konstandin, L., Girichidis, P., Federrath, C., & Klessen, R. S. 2012, ApJ, 761, 149 [NASA ADS] [CrossRef] [Google Scholar]
 Kritsuk, A. G., & Norman, M. L. 2004, ApJ, 601, L55 [NASA ADS] [CrossRef] [Google Scholar]
 Kritsuk, A. G., Norman, M. L., & Padoan, P. 2006, ApJ, 638, L25 [NASA ADS] [CrossRef] [Google Scholar]
 Kritsuk, A. G., Norman, M. L., Padoan, P., & Wagner, R. 2007a, ApJ, 665, 416 [NASA ADS] [CrossRef] [Google Scholar]
 Kritsuk, A. G., Padoan, P., Wagner, R., & Norman, M. L. 2007b, in Turbulence and Nonlinear Processes in Astrophysical Plasmas, eds. D. Shaikh, & G. P. Zank, AIP Conf. Ser., 932, 393 [Google Scholar]
 Kritsuk, A. G., Nordlund, Å., Collins, D., et al. 2011a, ApJ, 737, 13 [NASA ADS] [CrossRef] [Google Scholar]
 Kritsuk, A. G., Norman, M. L., & Wagner, R. 2011b, ApJ, 727, L20 [NASA ADS] [CrossRef] [Google Scholar]
 Krumholz, M. R., & McKee, C. F. 2005, ApJ, 630, 250 [NASA ADS] [CrossRef] [Google Scholar]
 Kurien, S., & Sreenivasan, K. R. 2000, Phys. Rev. E, 62, 2206 [NASA ADS] [CrossRef] [Google Scholar]
 Li, Z.Y., & Nakamura, F. 2006, ApJ, 640, L187 [NASA ADS] [CrossRef] [Google Scholar]
 Mac Low, M.M. 1999, ApJ, 524, 169 [NASA ADS] [CrossRef] [Google Scholar]
 Mac Low, M.M., Klessen, R., & Burkert, A. 1998, Phys. Rev. Lett., 80, 2754 [NASA ADS] [CrossRef] [Google Scholar]
 Marchenko, S. V., Moffat, A. F. J., & Grosdidier, Y. 1999, ApJ, 522, 433 [NASA ADS] [CrossRef] [Google Scholar]
 McKee, C. F., & Ostriker, E. C. 2007, ARA&A, 45, 565 [NASA ADS] [CrossRef] [Google Scholar]
 Melzani, M., Winisdoerffer, C., Walder, R., et al. 2013, A&A, 558, A133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mimica, P., & Aloy, M. A. 2010, MNRAS, 401, 525 [NASA ADS] [CrossRef] [Google Scholar]
 Mimica, P., Aloy, M. A., & Müller, E. 2007, A&A, 466, 93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Moffat, A. F. J., Drissen, L., Lamontagne, R., & Robert, C. 1988, ApJ, 334, 1038 [NASA ADS] [CrossRef] [Google Scholar]
 Molina, F. Z., Glover, S. C. O., Federrath, C., & Klessen, R. S. 2012, MNRAS, 3020 [Google Scholar]
 Moraghan, A., Kim, J., & Yoon, S.J. 2013, MNRAS, 432, L80 [NASA ADS] [CrossRef] [Google Scholar]
 Myasnikov, A. V., & Zhekov, S. A. 1998, MNRAS, 300, 686 [NASA ADS] [CrossRef] [Google Scholar]
 Ntormousi, E., Burkert, A., Fierlinger, K., & Heitsch, F. 2011, ApJ, 731, 13 [NASA ADS] [CrossRef] [Google Scholar]
 Ossenkopf, V., & Mac Low, M.M. 2002, A&A, 390, 307 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ossenkopf, V., Esquivel, A., Lazarian, A., & Stutzki, J. 2006, A&A, 452, 223 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Owocki, S. P., Castor, J. I., & Rybicki, G. B. 1988, ApJ, 335, 914 [NASA ADS] [CrossRef] [Google Scholar]
 Padoan, P., & Nordlund, Å. 1999, ApJ, 526, 279 [NASA ADS] [CrossRef] [Google Scholar]
 Padoan, P., & Nordlund, Å. 2002, ApJ, 576, 870 [Google Scholar]
 Padoan, P., Nordlund, A., & Jones, B. J. T. 1997, MNRAS, 288, 145 [NASA ADS] [CrossRef] [Google Scholar]
 Padoan, P., Boldyrev, S., Langer, W., & Nordlund, Å. 2003, ApJ, 583, 308 [NASA ADS] [CrossRef] [Google Scholar]
 Padoan, P., Jimenez, R., Nordlund, Å., & Boldyrev, S. 2004, Phys. Rev. Lett., 92, 191102 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Padoan, P., Nordlund, Å., Kritsuk, A. G., Norman, M. L., & Li, P. S. 2007, ApJ, 661, 972 [NASA ADS] [CrossRef] [Google Scholar]
 Pan, L., Wheeler, J. C., & Scalo, J. 2008, ApJ, 681, 470 [NASA ADS] [CrossRef] [Google Scholar]
 Panaitescu, A., Spada, M., & Mészáros, P. 1999, ApJ, 522, L105 [NASA ADS] [CrossRef] [Google Scholar]
 Parkin, E. R., & Pittard, J. M. 2010, MNRAS, 406, 2373 [NASA ADS] [CrossRef] [Google Scholar]
 Passot, T., & VázquezSemadeni, E. 1998, Phys. Rev. Lett., 58, 4501 [Google Scholar]
 Pittard, J. M. 2007, ApJ, 660, L141 [NASA ADS] [CrossRef] [Google Scholar]
 Pittard, J. M., & Parkin, E. R. 2010, MNRAS, 403, 1657 [NASA ADS] [CrossRef] [Google Scholar]
 Pittard, J. M., Dobson, M. S., Durisen, R. H., et al. 2005, A&A, 438, 11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pittard, J. M., Falle, S. A. E. G., Hartquist, T. W., & Dyson, J. E. 2009, MNRAS, 394, 1351 [NASA ADS] [CrossRef] [Google Scholar]
 Pittard, J. M., Hartquist, T. W., & Falle, S. A. E. G. 2010, MNRAS, 405, 821 [NASA ADS] [Google Scholar]
 Polychroni, D., Moore, T. J. T., & Allsopp, J. 2012, MNRAS, 422, 2992 [NASA ADS] [CrossRef] [Google Scholar]
 Porter, D. H., & Woodward, P. R. 1994, ApJS, 93, 309 [NASA ADS] [CrossRef] [Google Scholar]
 Price, D. J., Federrath, C., & Brunt, C. M. 2011, ApJ, 727, L21 [NASA ADS] [CrossRef] [Google Scholar]
 Quirk, J. J. 1994, Int. J. Numer. Methods Fluids, 18, 555 [NASA ADS] [CrossRef] [Google Scholar]
 RomanDuval, J., Jackson, J. M., Heyer, M., Rathborne, J., & Simon, R. 2010, ApJ, 723, 492 [NASA ADS] [CrossRef] [Google Scholar]
 RomanDuval, J., Federrath, C., Brunt, C., et al. 2011, ApJ, 740, 120 [NASA ADS] [CrossRef] [Google Scholar]
 Sánchez, N., Alfaro, E. J., & Pérez, E. 2005, ApJ, 625, 849 [NASA ADS] [CrossRef] [Google Scholar]
 Schmidt, W., Federrath, C., Hupp, M., Kern, S., & Niemeyer, J. C. 2009, A&A, 494, 127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Seifried, D., Schmidt, W., & Niemeyer, J. C. 2011, A&A, 526, A14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 She, Z.S., & Leveque, E. 1994, Phys. Rev. Lett., 72, 336 [Google Scholar]
 Stevens, I. R., Blondin, J. M., & Pollock, A. M. T. 1992, ApJ, 386, 265 [NASA ADS] [CrossRef] [Google Scholar]
 Stone, J. M., Ostriker, E. C., & Gammie, C. F. 1998, ApJ, 508, L99 [NASA ADS] [CrossRef] [Google Scholar]
 Stutzki, J., Bensch, F., Heithausen, A., Ossenkopf, V., & Zielinsky, M. 1998, A&A, 336, 697 [NASA ADS] [Google Scholar]
 Tuthill, P. G., Monnier, J. D., & Danchi, W. C. 1999, Nature, 398, 487 [NASA ADS] [CrossRef] [Google Scholar]
 VazquezSemadeni, E. 1994, ApJ, 423, 681 [NASA ADS] [CrossRef] [Google Scholar]
 VázquezSemadeni, E., Ryu, D., Passot, T., González, R. F., & Gazol, A. 2006, ApJ, 643, 245 [NASA ADS] [CrossRef] [Google Scholar]
 VázquezSemadeni, E., Gómez, G. C., Jappsen, A. K., et al. 2007, ApJ, 657, 870 [NASA ADS] [CrossRef] [Google Scholar]
 VázquezSemadeni, E., Colín, P., Gómez, G. C., BallesterosParedes, J., & Watson, A. W. 2010, ApJ, 715, 1302 [NASA ADS] [CrossRef] [Google Scholar]
 Vishniac, E. T. 1994, ApJ, 428, 186 [NASA ADS] [CrossRef] [Google Scholar]
 Walder, R., & Folini, D. 1996, A&A, 315, 265 [NASA ADS] [Google Scholar]
 Walder, R., & Folini, D. 1998, A&A, 330, L21 [NASA ADS] [Google Scholar]
 Walder, R., & Folini, D. 2000a, in Thermal and Ionization Aspects of Flows from Hot Stars: Observations and Theory, eds. H. J. G. L. M. Lamers, & A. Sapar, ASP Conf. Ser., 281 [Google Scholar]
 Walder, R., & Folini, D. 2000b, Ap&SS, 274, 343 [NASA ADS] [CrossRef] [Google Scholar]
 Walder, R., & Folini, D. 2002, in Interacting Winds from Massive Stars, eds. A. F. J. Moffat, & N. St.Louis, ASP Conf. Ser., 260, 595 [Google Scholar]
 Walder, R., Burrows, A., Ott, C. D., Livne, E., Lichtenstadt, I., & Jarrah, M. 2005, ApJ, 626, 317 [NASA ADS] [CrossRef] [Google Scholar]
 Walder, R., Folini, D., & Shore, S. N. 2008, A&A, 484, L9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Walder, R., Folini, D., Favre, J. M., & Shore, S. N. 2010, in Numerical Modeling of Space Plasma Flows, Astronum2009, eds. N. V. Pogorelov, E. Audit, & G. P. Zank, ASP Conf. Ser., 429, 173 [Google Scholar]
 Wang, P., Li, Z.Y., Abel, T., & Nakamura, F. 2010, ApJ, 709, 27 [NASA ADS] [CrossRef] [Google Scholar]
 Weinmann, S. M., Pasquali, A., Oppenheimer, B. D., et al. 2012, MNRAS, 426, 2797 [NASA ADS] [CrossRef] [Google Scholar]
 Williams, P. M., van der Hucht, K. A., van Wyk, F., et al. 2012, MNRAS, 420, 2526 [NASA ADS] [CrossRef] [Google Scholar]
 Zhekov, S. A. 2012, MNRAS, 422, 1332 [NASA ADS] [CrossRef] [Google Scholar]
 Zitouni, H., Daigne, F., Mochkovich, R., & Zerguini, T. H. 2008, MNRAS, 386, 1597 [NASA ADS] [CrossRef] [Google Scholar]
 Zrake, J., & MacFadyen, A. 2011, in ASP Conf. Ser. 1358, eds. J. E. McEnery, J. L. Racusin, & N. Gehrels, 102 [Google Scholar]
 Zrake, J., & MacFadyen, A. I. 2012, ApJ, 744, 32 [NASA ADS] [CrossRef] [Google Scholar]
 Zrake, J., & MacFadyen, A. 2013, ApJ, 763, L12 [NASA ADS] [CrossRef] [Google Scholar]
 Zuckerman, B., & Evans, II, N. J. 1974, ApJ, 192, L149 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Derivation of analytical expression for the driving efficiency
Fig. A.1 Sample size versus radial distance, normalized to the sample size at the smallest radius for cases “full in” (first panel), “unified” (second panel), parallel LOS (third panel), and perpendicular LOS (fourth panel) for the different runs, all at ℓ ≈ 12. Color coding as in Fig. 2. 
Part of the total (left plus right) upwind kinetic energy flux density, , is thermalized at the shocks confining the CDL. The remaining part, , drives the turbulence in the CDL. In analogy with the 2D case, we assume that and ℱ_{ekin,u} are related by a function of the upwind Machnumber alone, (A.1)We call the function f_{eff} the driving efficiency. An expression for f_{eff} can be derived by using the jump conditions for strong, oblique shocks, (A.2)The subscript d denotes downstream quantities, right after shock passage, the subscripts ⊥ and ∥ denote flow components perpendicular and parallel to the shock, and α is the absolute value of the angle between the xaxis and the normal to the shock. Using Eq. (A.2) we obtain where the integrals over S_{l,r} and YZ_{l,r} run over the left and right shock and where it was used that dscos(α) = dydz. We further used the shock jump conditions given in Eq. (A.2), that , and that . The term is omitted as the shocks we observe in our simulations fulfill for the most part. Analogous to the 2D case we thus obtain (A.5)where we used the midpoint rule.
Appendix B: Structure functions from CDL data: computational details and additional results
Appendix B.1: Computing structure functions from CDL data
To compute S_{p}(r) from our slab data, we consider data of only one time step, at ℓ ≈ 12. To average over all positions x, we step through all grid points in the CDL. To average over all directed distances r, we cast rays at each position x into 20 directions, given by the surface normals of an icosahedron centered at x: (± 1, ± 1, ± 1), (0, ± 1/ϕ, ± ϕ), (± 1/ϕ, ± ϕ,0), and (± ϕ,0, ± 1/ϕ) with . We tested the robustness of our results against rotation of this icosahedron. We step along each ray and compute the longitudinal and transverse contribution to S_{p}(r). Each ray is covered by 23 logarithmically spaced steps from 3 Δx to 100 Δx, Δx being the grid cell size. Much larger distances (beyond 128) make no sense because of the periodic domain. Only points within the CDL may contribute to S_{p}(r). A ray that once left the CDL is forbidden to reenter the CDL at a larger r.
As for a given position x and distance r some rays may already reach outside the CDL while others still lie within, we distinguish two cases for each pair (x, r): all 20 rays are still
within the CDL (“full in” case) or at least one ray has left the CDL (“some out” case). The “full in” case ascertains that only symmetric contributions (all 20 directions) enter the mean in Eq. (18). The price to pay is a rather rapidly decreasing and nonuniform sample size: only points x close to the center of the CDL contribute to S_{p}(r) for large r. Considering “full in” and “some out” contributions in the mean in Eq. (18) (“unified” case) samples the CDL as completely as possible, but introduces some asymmetry in the computation of the mean in Eq. (18): some directions are missing for some x and r, especially for large r. The different sample sizes as function of distance r are illustrated in Fig. A.1. Third order structure functions are shown in Figs. 10 and 11.
While structure functions based on spherical averaging are of theoretical interest, observation based structure functions are necessarily based on lineofsight data. This motivated us to compute structure functions also along individual one dimensional directions, parallel or perpendicular to the upstream flow. The direction of computation may be interpreted as a lineofsight (LOS), although we stress that no density variations or radiative transfer effects are taken into account. Also, our demand that a ray having left the CDL may not reenter it, implies for perpendicular LOSs that regions close to the (wiggled) confining shocks contribute to S_{p}(r) only for small r, whereas large distances are dominated again by contributions from the central part of the CDL. For a parallel LOS, the velocity difference in Eq. (18) approaches, as r increases, the difference between the two post shock flow velocities at each of the confining shocks.
Appendix B.2: Computation of ESS scaling exponents
Given the small inertial range of our structure functions, apparent in Fig. 10, top panel, we make use of the ESS hypothesis and compute linear fits to Z_{p}, i.e., to log 10(S_{p}) versus log 10(S_{3}). We restrict the analysis to distances where the third order structure functions still increase, i.e., we demand that log 10(S_{3}(r_{i})) > 0.99 log 10(S_{3}(r_{i − 1})) for i > 4. We further require that at least eight data points enter the linear fit (corresponding to roughly half a magnitude in r). Within these limits, we search for the best linear fit (of Z_{p}) in terms of relative standard error. Finally, we reject even the best fit if its relative standard error exceeds 5%. While the choice of 5% is arbitrary, we found the results to remain essentially unchanged if we repeated the analysis with a value of 10%. For the “full in” case, the quality of the fits is illustrated in Fig. 10, bottom panel, for longitudinal structure functions and p = 2. In Appendix B, the quality of the fits is illustrated for all cases considered (“full in”, “unified”, different LOSs) for Z_{5} (Fig. B.1). From Z_{p}, p = 1 to 5, and Eq. (20) a best estimate for the codimension C is derived for each run (minimal RMS error). Note that typically more than the required 8 points enter each ESS fit and that the best fit is found for the smallest distances r. This situation is typical for nearly all the situations considered. The only exceptions are ESS for longitudinal structure functions computed along one dimensional lineofsights (see Sect. 3.4.4).
Fig. B.1 Fits (solid lines) of ESS scaling exponents Z_{5} as the slope of log 10(S_{5}) versus log 10(S_{3}) for longitudinal (top row) and transverse (bottom row) structure functions for the cases “full in” (first column), “unified” (second column), parallel LOS (third column), and perpendicular LOS (fourth column). Filled symbols denote data used for the fit, empty symbols denote all data available. Color coding as in Fig. 2. 
Longitudinal (top) and transverse (bottom) ESS exponents Z_{p} and best estimates for codimension C for the case “unified” (left columns) and lineofsights parallel (middle columns) and perpendicular (right columns) to the upstream flow, for the different runs at ℓ ≈ 12.
All Tables
Longitudinal (top) and transverse (middle) ESS scaling exponents Z_{p} and best estimates for codimension C, case “full in” for the different runs at ℓ ≈ 12.
Comparing ESS scaling exponents Z_{p}, p = 1 to p = 5, following Table 4 in Schmidt et al. (2009).
Longitudinal (top) and transverse (bottom) ESS exponents Z_{p} and best estimates for codimension C for the case “unified” (left columns) and lineofsights parallel (middle columns) and perpendicular (right columns) to the upstream flow, for the different runs at ℓ ≈ 12.
All Figures
Fig. 1 Different geometry of confining shocks. In 2D (periodic in ydirection, infinite in zdirection), the confining shocks resemble corrugated sheets (left panel). In 3D (periodic in y and zdirection) they evoke the idea of cardboard egg wrappings (right panel). Both simulations have M_{u} = 21.7. Color coded is density (log10(ρ); particles/cm^{3} in 2D; g/cm^{3} in 3D) and, right panel only, the driving efficiency f_{eff} at the shock surface. 

In the text 
Fig. 2 Mean density ρ_{m} (top panel), root mean square Mach number M_{rms} (middle panel), and (bottom panel) as function of the monotonically increasing CDL size ℓ for simulations R*_2. Individual curves denote M_{u} = 2 (dark blue), 4 (dark green), 5 (yellow), 7 (purple), 8 (light blue), 11 (black), 16 (green), 22 (red), 27 (pink), 33 (cyan), and 43 (magenta). 

In the text 
Fig. 3 Driving efficiency of 3D slabs, same simulations and color coding as in Fig. 2. Top panel: comparison of numerical results with expected scaling relation β_{3}log _{10}(M_{u}) ∝ log _{10}(1 − f_{eff}). Each star denotes the maximum driving efficiency of one simulation. Bottom panel: evolution of f_{eff} as function of CDL size. 

In the text 
Fig. 4 Supersonic mass fraction (solid) and volume fraction (dotted), color coding as in Fig. 2. 

In the text 
Fig. 5 Deviating a highly supersonic, isothermal flow by (at least) β_{tot} = π/2 through n (right yaxis) subsequent shock passages at angle α (xaxis; in degrees; angle between flow direction and surface normal of shock) results in a total reduction of the Mach number by a factor F_{tot} (left yaxis). 

In the text 
Fig. 6 2D slice (yzplane) averages as function of (scaled) xcoordinate for R22_2 for different times (ℓ = 12,17,23,29, and 34, for the blue, cyan, green, magenta, and red curves, respectively). Top panel: fractional 2D CDL area F(x). Bottom panel: root mean square Mach number in units of M_{u}. For details see text, Sect. 3.2.3. 

In the text 
Fig. 7 2D slices, average quantities for different runs (color coding as in Fig. 2) at CDL size ℓ ≈ 12. From top to bottom: fractional CDL area F(x), root mean square Mach number M_{rms}, scaled root mean square Mach number M_{rms}/M_{u}, mean density ρ_{m}, and b = σ(ρ)/(ρ_{m}M_{rms}). 

In the text 
Fig. 8 Top panel: volumeweighted PDFs of s = ln(ρ/ρ_{m}) for the different runs, color coding as in Fig. 2. Each curve is an average over three data sets between ℓ ≈ 11 and ℓ ≈ 12. Bottom panel: centrally (at peak value of PDF) mirrored PDF with Gaussian fitting (dotted line). 

In the text 
Fig. 9 A 5% random perturbation of Z_{p} (30 000 realizations, at the example of Z_{p} from Boldyrev 2002, p = 1 to p = 5) causes smearing of codimensions C derived from Eq. (20) (left panel). A double peak distribution results if C is computed from Eq. (21) (middle panel). The analytical value C = 1 is shown as a red line, green lines in the left panel indicate where the 2/3 largest values in the histogram are located. In the right panel, the 2D histogram (log10 contours, spacing 0.5, spanning three orders of magnitude) of codimensions from perturbed Z_{p} (500 000 realizations) is overplotted on the codimension C = Δ/(1 − β) (color coded) in the Δβ plane (Eq. (21)). Two sharp peaks can be distinguished, at Δ ≈ 0.4 and β ≈ 0.1 (C ≈ 0.44), as well as at Δ ≈ 1.0 and β ≈ 0.5 (C ≈ 1.9). The red dot is the analytical value by Boldyrev (2002), Δ = 2/3, β = 1/3, C = 1. 

In the text 
Fig. 10 “Full in” case for simulations R*_2 at ℓ ≈ 12, color coding as in Fig. 2. Top panel: structure functions (upper curves) and (lower curves), curves vertically shifted to start from one common point. Stars connected by lines indicate data used to determine ESS scaling exponents. Also indicated are slopes and as obtained by Kritsuk et al. (2007a) for 3D periodic box simulations. Bottom panel: fit (solid lines) of ESS scaling exponent as slope of versus . Filled symbols denote data used for the fit, empty symbols denote all data available. 

In the text 
Fig. 11 Structure functions (solid) and (dashed) for the case “unified” (top panel), as well as computed along the xdirection (middle panel) and ydirection (bottom panel), i.e., parallel and perpendicular to the upstream flow, for simulations R*_2 at ℓ = 12. Color coding as in Fig. 2. Also indicated are slopes and as obtained by Kritsuk et al. (2007a) from 3D periodic box simulations. Curves are vertically shifted, to start from one common point. 

In the text 
Fig. 12 Role of dimensionality for driving efficiency f_{eff} as function of upstream Mach number M_{u}. Shown are simulations R*_2 (red stars) and corresponding values for 2D slabs (dark red diamonds, simulations R*_0.2.4 in FW06). 

In the text 
Fig. A.1 Sample size versus radial distance, normalized to the sample size at the smallest radius for cases “full in” (first panel), “unified” (second panel), parallel LOS (third panel), and perpendicular LOS (fourth panel) for the different runs, all at ℓ ≈ 12. Color coding as in Fig. 2. 

In the text 
Fig. B.1 Fits (solid lines) of ESS scaling exponents Z_{5} as the slope of log 10(S_{5}) versus log 10(S_{3}) for longitudinal (top row) and transverse (bottom row) structure functions for the cases “full in” (first column), “unified” (second column), parallel LOS (third column), and perpendicular LOS (fourth column). Filled symbols denote data used for the fit, empty symbols denote all data available. Color coding as in Fig. 2. 

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.