3.1. Velocity derivatives

In axisymmetric homogeneous stratified turbulence there are five statistically independent local velocity derivatives. The probability density functions (p.d.f.s) of four are plotted in the first two rows of figure 1, with the notation $u_{i,j}$ meaning that $i$ is the velocity component and $j$ is the derivative direction, with index 3 denoting the direction parallel to gravitational acceleration. Each point represents the mean value for a bin with 16 777 216 samples. The r.m.s. of the fluctuations of the values in each bin about the bin mean provides an estimate of the statistical uncertainty that can be plotted as a horizontal error bar. In figure 1, the error bars are comparable to the linewidth and are not visible (cf. figure 3). The moments of the velocity derivatives are given in table 2.

Figure 1.

The p.d.f.s of the velocity derivatives, each normalised by its own mean and standard deviation.

Table 2.

The variance ${\it\sigma}^{2}$ , skewness ${\it\mu}_{3}$ and kurtosis ${\it\mu}_{4}$ of the velocity derivatives.

Considering first the reference case R4, the longitudinal and transverse velocity derivatives are visibly different in figure 1, as expected. The variances are almost exactly isotropic, the skewness of the longitudinal derivatives is approximately −0.6, and the kurtoses indicate tails of the p.d.f.s heavier than Gaussian, all in agreement with Ishihara et al. (Reference Ishihara, Kaneda, Yokokawa, Itakura and Uno2007) and the references therein.

Using case R4 as a reference, we now consider the stratified cases, starting with F1, which has the highest buoyancy Reynolds number. As with R4, the behaviours of the longitudinal and transverse derivatives differ but there is little difference between the two longitudinal or the two transverse quantities. When $\mathit{Re}_{b}$ is reduced to 48 in case F2, several changes in the distributions become apparent. First, for all but $u_{1,3}$ the probability of values near zero increases so that the distributions, particularly $u_{3,1}$ , have sharp peaks. Second, the two longitudinal distributions become markedly different. This is even more obvious in case F3, for which the distribution of $u_{1,3}$ is similar to that in the isotropic case, whereas $u_{3,1}$ exhibits a sharp peak near zero. The probability of very-small-magnitude transverse derivatives of $u_{3}$ is almost an order of magnitude greater than that for $u_{1}$ . These characteristics of the p.d.f.s are due to stratified flows forming quasi-horizontal vortices, sometimes referred to as ‘pancake eddies’ (e.g. Lin & Pao Reference Lin and Pao1979), in which $u_{1,3}$ is dominated by shearing between these eddies. As the stratification increases, the pancakes get thinner (Billant & Chomaz Reference Billant and Chomaz2001) and the magnitudes of $u_{1,3}$ increase (see table 2 and Hebert & de Bruyn Kops Reference Hebert and de Bruyn Kops2006b , figure 5), but the shape of $P(u_{1,3})$ does not change very much. On the other hand, the shape of $P(u_{3,1})$ is affected by the strength of the stratification because stratified turbulence tends to form in patches, which can be seen for a lower-resolution realisation of F3 in Almalkie & de Bruyn Kops (Reference Almalkie and de Bruyn Kops2012b , figure 1). The regions of relative quiescence contribute to the peak of $P(u_{3,1})$ while the regions of more intense turbulence contribute to the tails, and this effect is more pronounced for stronger stratification.

Further analyses are made using the data in table 2. It is observed that the moments of the velocity derivatives depend on $\mathit{Re}_{b}$ and approach their isotropic values as $\mathit{Re}_{b}$ increases. The importance of the departure from isotropy is subjective, but we can make some general observations to put it into perspective. To aid in doing so, some of the moments of the derivatives averaged over spherical volumes of radius $r$ are shown in figure 2. These locally filtered quantities are computed in terms of Fourier series according to the unpublished derivation of Sykora, which is reproduced in the appendix of Almalkie & de Bruyn Kops (Reference Almalkie and de Bruyn Kops2012a ). Because the velocity fields in our simulations are computed in terms of finite Fourier series, the filtered quantities are exact for any averaging dimension $r$ up to one-half of the size of the domain.

Figure 2.

Moments of locally averaged velocity derivatives for selected cases scaled by $15{\it\nu}/\langle {\it\epsilon}\rangle$ . To aid in curve identification, the symbols on the variance plots go from top to bottom in the legend in the same order as they go from top to bottom in the plot. In the online version of the paper, the symbol and line colours are the same and are unique for each derivative.

From the tabulated and graphical data for the variances, it is apparent that at low $\mathit{Re}_{b}$ there are three distinct groups of derivatives, whereas at higher $\mathit{Re}_{b}$ the data suggest two groups as in the isotropic case. The skewnesses of the longitudinal velocity derivatives are negative because the flows are statistically stationary: viscous effects continually destroy vorticity, and so there must be generation by vortex stretching, which implies negative skewness of the longitudinal derivatives (Taylor Reference Taylor1938). The skewnesses of $u_{1,2}$ and $u_{2,1}$ averaged over large $r$ are far from zero because of how the large scales of the flow are forced. Even as $r\rightarrow 0$ , though, these skewnesses include contributions from the large scales, which is consistent with them apparently asymptoting to non-zero values.

3.2. Scalar derivatives

The scalar field ${\it\phi}$ is the temperature in the stratified simulations and a passive scalar in simulation R4. The three spatial derivatives are ${\it\phi}_{,1}$ , ${\it\phi}_{,2}$ and ${\it\phi}_{,3}$ , where the first two are in the horizontal directions and statistically interchangeable. The p.d.f.s of each derivative are plotted as figure 3, and their moments are given in table 3.

Figure 3.

The p.d.f.s of the scalar derivatives, each normalised by its own mean and standard deviation. The error bars are defined in § 3.1.

Table 3.

The variance ${\it\sigma}^{2}$ , skewness ${\it\mu}_{3}$ and kurtosis ${\it\mu}_{4}$ of the scalar derivatives; $D$ is the molecular diffusivity.

First, let us check the reference case R4 against the data of Yeung et al. (Reference Yeung, Xu and Sreenivasan2002) and Donzis & Yeung (Reference Donzis and Yeung2010) for unity Schmidt number and various Reynolds numbers, and more generally against the data referenced in Warhaft (Reference Warhaft2000). From the cited data, we expect the variance of ${\it\phi}_{,3}$ to be slightly higher than that in the directions perpendicular to the mean scalar gradient, the skewness ${\it\phi}_{,3}$ to be positive whereas the other skewnesses are nearly zero, and the kurtoses to indicate tails heavier than Gaussian and not strongly dependent on the direction of the derivative. It is concluded that case R4 is behaving as expected and so is a useful reference case against which to compare the stratified cases.

Next, we consider case F1. We recall that this case has Froude number of order unity and is highly non-isotropic at large scales. Nevertheless, the moments of the p.d.f.s of the local velocity derivatives are only slightly affected by buoyancy. As the buoyancy Reynolds number is reduced (cases F2 and F3), the main effect of stratification is to increase the kurtoses of all of the scalar derivatives. The same is observed for the velocity derivatives and will be considered further in the context of the local dissipation rates in § 5. Another effect of low buoyancy Reynolds number is that the vertical derivatives of the scalar are larger than the horizontal derivatives. This is expected since, in the conceptual configuration of horizontal pancake-like structures shearing against each other, the vertical derivatives will be dominant.

The bottom rows of figures 1 and 3 show how the p.d.f.s of the velocity and scalar respectively change with buoyancy Reynolds number. In each pair of plots, the horizontal derivatives are on the left and the vertical derivatives on the right. For velocity, changing $\mathit{Re}_{b}$ has a dramatic effect on the shape of the p.d.f. of the horizontal derivative and a minor effect on that of the vertical derivative, as discussed in § 3.1. For the scalar, the effect is minor for all derivatives. Our use of the term ‘minor’ presumes a focus on the portion of the p.d.f.s where the probability density is high, and it is noted that the logarithmic vertical axis emphasises the tails of the p.d.f.s, which are important in some applications (Sreenivasan & Antonia Reference Sreenivasan and Antonia1997). It is also noted that the horizontal axes for the scalar derivatives extend to 15 standard deviations from the mean, again because the tails are of interest to some readers. Therefore we conclude that, even with unity Prandtl number, stratification has a stronger effect on the derivatives of velocity than on those of the scalar.

The non-zero skewness of the scalar derivative parallel to the mean gradient, first observed more than 40 years ago (Gibson, Stegen & McConnell Reference Gibson, Stegen and McConnell1970), indicates that the scalar is not isotropic at small scales. In his review article, Warhaft explains the mechanism that causes the non-zero skewness and concludes that the isotropy assumption underlying KOC scaling cannot be satisfied for a scalar having a mean gradient (Warhaft Reference Warhaft2000). Nevertheless, consistency with KOC scaling for flows with a passive scalar and a mean gradient has proven useful (Yeung et al. Reference Yeung, Donzis and Sreenivasan2005). We conclude that, despite evidence of anisotropy in the data, it is worthwhile to examine KOC scaling for the stratified cases, which we do in the next section.

4.1. Behaviour of the velocity fields

The most commonly known classical scaling hypotheses are those presented in Kolmogorov (Reference Kolmogorov1941) and Oboukhov (Reference Oboukhov1941a ,Reference Oboukhov b ). With the assumption that $\langle {\it\epsilon}\rangle$ and the kinematic viscosity ${\it\nu}$ determine the form of the velocity statistics in the equilibrium range, and the assumption of local isotropy, the longitudinal and transverse velocity spectra in the inertial range take the forms

and

and the second-order longitudinal velocity structure function takes the form

The textbook values for the constants are (Pope Reference Pope2000)

The spectra and structure function are plotted in figures 4 and 5. Full advantage of homogeneity and the statistical interchangeability of the two horizontal directions was taken in computing these quantities. Considering first the longitudinal spectra, KOC scaling collapses all the cases to a single curve in the far dissipation range. This is not obvious with the linear scaling of the ordinate in the plot of the spectra, but is evident if logarithmic scaling is used, and it is clear from the plot of the second-order structure functions. Furthermore, the third-order structure functions collapse to the relation $S/(15)^{3/2}(r/{\it\eta})^{2}$ for small $r/{\it\eta}$ , where $S$ is the velocity derivative skewness. The behaviours of the spectra and structure functions in the dissipation range suggest that ${\it\eta}$ is the appropriate length scale for studying KOC scaling rather than, say, a length scale based on the sum of $\langle {\it\epsilon}\rangle$ and $\langle {\it\chi}\rangle$ . We consider mean dissipation rates further in § 4.4.

Figure 4.

Longitudinal and transverse velocity spectra and their ratios.

Figure 5.

Second- and third-order longitudinal velocity structure functions. The line segments for $r{\it\eta}<0.6$ are not data but rather the theoretical relationships annotated on the plots.

Next, we consider the longitudinal spectra at wavenumbers smaller than the dissipation range. All the cases have a region with a plateau, so one might conclude that the scaling holds but with $C_{1}$ and $C_{1}^{\prime }$ dependent on the buoyancy Reynolds number. This is consistent with conclusions that might be drawn from figure 9 in Lindborg (Reference Lindborg2006). A closer look, though, reveals that the plateaus do not occur in the same ranges of ${\it\kappa}_{1}{\it\eta}$ for all cases. In cases F1 and R4, the plateau is in the vicinity of ${\it\kappa}_{1}{\it\eta}\approx 0.01$ , which is consistent with the laboratory data, simulations and models for unstratified turbulence (Saddoughi & Veeravalli Reference Saddoughi and Veeravalli1994; Gotoh, Fukayama & Nakano Reference Gotoh, Fukayama and Nakano2002; Meyers & Meneveau Reference Meyers and Meneveau2008). In cases F2 and F3, the plateaus occur at higher values of ${\it\kappa}_{1}{\it\eta}$ . The transverse spectra lead us to the same conclusions for cases F1 and R4, but for this quantity plateaus are not evident for cases F2 and F3.

To investigate further, let us consider the bottleneck effect, which appears as a bump in the near dissipation range of compensated spectra. In case R4, the bump is clearly visible at ${\it\kappa}_{1}{\it\eta}\approx 0.1$ , as expected (Yeung & Zhou Reference Yeung and Zhou1997; Meyers & Meneveau Reference Meyers and Meneveau2008; Donzis & Sreenivasan Reference Donzis and Sreenivasan2010). Using the same data sets as in this paper (but at lower resolution), Almalkie & de Bruyn Kops (Reference Almalkie and de Bruyn Kops2012b ) conclude that stratification impedes downscale transfer of kinetic energy, which is consistent with the smaller bump for case F1 relative to case R4. In cases F2 and F3, the bottleneck effect is not apparent in the longitudinal spectra, but is visible in the transverse spectra. The plateau regions in the longitudinal spectra for cases F2 and F3, to the extent that they might exist, are in the range of ${\it\kappa}_{1}{\it\eta}$ where the bottleneck effect is expected.

A consequence of Kolmogorov’s hypotheses is $E_{22}/E_{11}=4/3$ in the power-law range. This ratio is plotted in figure 4(c). The ratio equality is satisfied for cases R4 and F1 in the same range of wavenumbers for which plateaus appear in the spectra. For case F3, the ratio equality is not satisfied over any range of wavenumbers. For case F2, the results are not definitive; the ratio might suggest a short scaling range in the vicinity of ${\it\kappa}_{1}{\it\eta}\approx 0.04$ .

Putting together observations about the plateaus, the bottleneck effect and the ratio of the spectra, we arrive at the same conclusion as Bartello & Tobias (Reference Bartello and Tobias2013), namely that KOC scaling of spectra, if it occurs, does so only if there is sufficient scale separation between the Ozmidov and Kolmogorov length scales, i.e. if the buoyancy Reynolds number is high enough. Whereas it is concluded in Bartello & Tobias (Reference Bartello and Tobias2013) that $\mathit{Re}_{b}>O(10)$ is required for KOC scaling, the current data show that with $\mathit{Re}_{b}=48$ the scaling region is obscured by the bottleneck effect. Cases F2 and F3 suggest that when $\mathit{Re}_{b}$ is not high, a region with ${\it\kappa}_{1}^{-5/3}$ scaling may appear in the near dissipation range due to the competing effects of the bottleneck tending to cause a positive slope and dissipation effects tending to cause a steeply negative slope.

Lastly in our discussion of velocity spectra, we examine the numerical values of $C_{1}$ and $C_{1}^{\prime }$ . Since case F1 is close to locally isotropic and exhibits plateaus in $E_{11}$ and $E_{22}$ in the range of ${\it\kappa}{\it\eta}$ typical for flows that are consistent KOC scaling, we might expect it to exhibit the universal values for $C_{1}$ and $C_{1}^{\prime }$ . Indeed, the textbook constants of 0.49 and 0.65 fit the data extremely well.

Turning now to velocity structure functions, let us consider figure 5. Consistency with KOC scaling will be exhibited by a plateau in the second-order structure function at 2 and in the third-order structure function at $4/5$ . The isotropic case satisfies both relations almost exactly. This is an important criterion because it indicates whether a simulation has sufficient dynamic range to not just show ${\it\kappa}_{1}^{-5/3}$ scaling in the spectra but also to satisfy other relationships that follow from that scaling, which simulations with lower dynamic range do not do (Yeung et al. Reference Yeung, Donzis and Sreenivasan2005). We recall that these simulations are forced and note that Antonia & Burattini (Reference Antonia and Burattini2006) show that the $4/5$ law is exhibited at lower dynamic range in forced than in decaying turbulence. The two most strongly stratified cases do not appear to satisfy either of the scaling relationships for structure functions. For the second-order structure function, case F1 almost plateaus at 2.0, but it appears that slightly more scale separation is needed for the plateau to be clear. The third-order structure functions, though, for the stratified cases have pronounced positive slopes in the range of $r/{\it\eta}$ where case R4 shows KOC scaling.

4.2. Behaviour of the scalar fields

The one-dimensional horizontal scalar spectra and corresponding structure functions are shown in figure 6. The horizontal reference lines at 0.4 and 1.6 on the spectra and structure function plots respectively are based on experimental data (Sreenivasan & Kailasnath Reference Sreenivasan and Kailasnath1996). The horizontal reference line at 2.0 on the structure functions marks the plateau in simulation data of Yeung et al. (Reference Yeung, Xu and Sreenivasan2002) for comparable Reynolds and Schmidt numbers. Also shown is the asymptotic behaviour at small $r/{\it\eta}$ (Yeung et al. Reference Yeung, Xu and Sreenivasan2002),

Figure 6.

One-dimensional scalar spectra (a) and the second-order scalar structure functions (b).

The scalar spectra and structure functions confirm the conclusions drawn from the velocity data. For the unstratified case, there is a plateau in the spectrum around ${\it\kappa}{\it\eta}\approx 0.01$ at the ordinate value reported in the literature, and the bottleneck effect is pronounced. The data for case F1 are broadly consistent with the unstratified case, although the bottleneck effect is weaker. For cases F2 and F3 there is no bottleneck effect; there is a plateau in case F2 where the bottleneck effect is expected, and no plateau in case F3. In the plot of the structure functions, the theoretical asymptotic behaviour at low $r/{\it\eta}$ is evident for the cases R4 and F1; for case F3, the curve does not quite reach the asymptotic slope. This behaviour is consistent with lack of scale separation, that is, low buoyancy Reynolds number. It is also consistent with insufficient small-scale resolution to reveal the asymptotic behaviour, but this explanation is discounted since the values of $r/{\it\eta}$ shown on the plot are in the far dissipation range and the expected asymptotic behaviours are observed for the second- and third-order velocity structure functions. It is concluded, therefore, that effects of low buoyancy Reynolds number are observed far into the dissipation range of the scalar when $\mathit{Re}_{b}=13$ and not when $\mathit{Re}_{b}=220$ .

4.3. Yaglom relationship

Yeung et al. (Reference Yeung, Donzis and Sreenivasan2005) comment that the most important relationship in physical space that pertains to KOC scaling is that derived by Yaglom (Reference Yaglom1949) for the mixed velocity–scalar structure function:

The appropriately scaled structure function is plotted as figure 7. As with the scalar structure functions, the expected asymptotic behaviour at low $r/{\it\eta}$ is observed for the unstratified case and at high $\mathit{Re}_{b}$ . For case F3, effects of stratification are evident in this statistic even at length scales in the far dissipation range. For all of the cases the mixed structure function peaks at approximately $2/3$ , although even case R4 does not show a distinct plateau; figure 3 in Yeung et al. (Reference Yeung, Donzis and Sreenivasan2005) indicates that an even higher Reynolds number is required for a clear plateau to be apparent.

Figure 7.

The velocity–scalar structure function. Line types are the same as in figure 6.

4.4. Scaled dissipation rates

A common premise about turbulence is that $\langle {\it\epsilon}\rangle$ is dependent only on the size and energy of the large-scale motions in a flow (Tennekes & Lumley Reference Tennekes and Lumley1972; Sreenivasan Reference Sreenivasan1984). This is not part of the KOC hypotheses, but it is in the same vein, and so we perform the analysis on our four data sets. As noted in § 4.1, stronger stratification somewhat reduces the absolute dissipation rate of kinetic energy in simulations when the large horizontal length scales of the horizontal velocities are forced to be the same. Now let us consider the scaled dissipation rate

Here, $L$ is an outer scale, evaluated as the average of the longitudinal integral scales of the two horizontal velocities, and $u_{rms}$ is the r.m.s. velocity.

For a variety of flow configurations, Sreenivasan (Reference Sreenivasan1998) report $\mathscr{D}$ in the range of approximately 0.5–1.0 at sufficiently high Reynolds number, but the definitions of $L$ and $u_{rms}$ vary in the literature when the flow is not isotropic. Nevertheless, it is useful to consider the values of $\mathscr{D}$ for the current simulations, as shown in table 4. The table also lists the values of

where ${\it\epsilon}_{p}$ is the rate at which available potential energy is converted to background potential energy by molecular diffusion. The value of $\mathscr{D}$ for R4 is very close to the theoretical estimate of Rotta (Reference Rotta1972), but is higher than reported for a number of isotropic homogeneous simulations (Kaneda et al. Reference Kaneda, Ishihara, Yokokawa, Itakura and Uno2003). The value for forced simulations, though, tends to be higher than that for decaying simulations (Sreenivasan Reference Sreenivasan1998). For the stratified cases, $\mathscr{D}$ is near the low end of the range reported by Sreenivasan (Reference Sreenivasan1998), and there is significant spread between the stratified and unstratified cases. However, $\mathscr{D}_{tot}$ is approximately the same for all four of our data sets. This is consistent with the rate of kinetic energy transfer out of the large scales being determined only by the large scales and independent of whether it eventually dissipates as kinetic or potential energy.

Table 4.

Scaled dissipation rates. The ranges indicate the maximum and minimum values obtained using the two horizontal integral length scales available from assuming statistical interchangeability of the horizontal velocity components.

6.1. Intermittency exponents

While Kolmogorov and Oboukhov assumed the lognormal distribution expressed in (6.1) and (6.2), the theory of intermittency exponents was soon generalised so as not to require that assumption. In particular, Novikov & Stewart (Reference Novikov and Stewart1964) postulated that the second-order moment of ${\it\epsilon}_{r}$ has a universal form in the inertial range:

Additional history on the topic is summarised by Sreenivasan & Kailasnath (Reference Sreenivasan and Kailasnath1993), who also consider several approaches for computing ${\it\mu}$ . Given the foregoing, ${\it\mu}$ and ${\it\mu}_{{\it\phi}}$ are of interest in stratified turbulence even though the lognormal assumption is not valid for the range of $\mathit{Re}_{b}$ we are considering.

The textbook value for ${\it\mu}$ with $r$ in the inertial range is ${\it\mu}\approx 0.25$ (Sreenivasan & Antonia Reference Sreenivasan and Antonia1997; Pope Reference Pope2000). Mydlarski & Warhaft (Reference Mydlarski and Warhaft1998) report ${\it\mu}_{{\it\phi}}\approx 0.2$ for grid turbulence at $\mathit{Re}_{{\it\lambda}}=100$ with a slight upward trend with increasing $\mathit{Re}_{{\it\lambda}}$ . Muschinski (Reference Muschinski2004) reports ${\it\mu}_{{\it\phi}}\approx 0.3$ in the stable nocturnal boundary layer. Sreenivasan & Kailasnath (Reference Sreenivasan and Kailasnath1993) suggest ${\it\mu}_{{\it\phi}}\approx 0.35$ based on Sreenivasan, Antonia & Danh (Reference Sreenivasan, Antonia and Danh1977) and Prasad, Meneveau & Sreenivasan (Reference Prasad, Meneveau and Sreenivasan1988). It should be noted that Almalkie & de Bruyn Kops (Reference Almalkie and de Bruyn Kops2012a ) report that ${\it\mu}$ computed from ${\it\epsilon}_{r}$ estimated from a single velocity derivative will be higher than ${\it\mu}$ computed from the true ${\it\epsilon}_{r}$ .

Figure 12.

(a) The scale dependence of the second-order moments of ${\it\epsilon}_{r}$ . (b) The local scaling exponents. The Taylor microscale marked at $r/{\it\eta}=42$ is within $\pm 2$ of the true value for all cases (see table 1).

For the current data sets, the exponents are evaluated from the scale dependence of the variances of the filtered dissipation rates (cf. Sreenivasan & Kailasnath Reference Sreenivasan and Kailasnath1993), which are plotted as figures 12 and 13. Figures 12(a) and 13(a) show the variances of the locally averaged dissipation rates. Figures 12(b) and 13(b) are the local scaling exponents computed by taking the derivative of the data in figures 12(a) and 13(a). The local scaling exponent for $\langle {\it\epsilon}_{r}^{2}\rangle$ with $r$ in the inertial range is ${\it\mu}$ , and similarly for ${\it\mu}_{{\it\phi}}$ . The plots reveal no broad range of $r$ with power-law scaling, and the local scaling exponents decrease monotonically with increasing $r$ for $r$ in the vicinity of  ${\it\lambda}$ . Similar behaviour is observed in high-Reynolds-number experimental measurements and atmospheric observations of Praskovsky & Oncley (Reference Praskovsky and Oncley1997) and so is presumably not an artefact of the limited dynamic range of the simulations. Since there is no broad range of $r$ where the local scaling exponent is constant, all we can conclude is that ${\it\mu}=0.25\pm 0.05$ and ${\it\mu}_{{\it\phi}}=0.35\pm 0.1$ in the vicinity of ${\it\lambda}$ , which is consistent with the published values cited at the beginning of this section.

Figure 13.

(a) The scale dependence of the second-order moments of  ${\it\chi}_{r}$ . (b) The local scaling exponents. The Taylor microscale marked at $r/{\it\eta}=42$ is within $\pm 2$ of the true value for all cases (see table 1).