1 Introduction

The closure modelling techniques are a fundamental theoretical tool for understanding and investigating the nonlinear behaviors of turbulent flows. The quasi-normal assumption, first introduced by Chou ^[1], is the most famous and successful technique for obtaining closures in the spectral space. Subsequent developments of the quasi-normal assumption led to the eddy-damping quasi-normal Markovian (EDQNM) closure ^[2-3], which is usually regarded as one of the most accurate analytical models of turbulence. There are also some investigations on relaxing the quasi-normal assumption to consider the fourth-order accumulations ^[4-5], and investigations on the possibility of using the restricted Euler (RE) assumption to replace the quasi-normal assumption ^[6].

However, the quasi-normal assumption only provides a description in the spectral space which introduces links between second-order and fourth-order spectral correlations. This restricts its implementation in complex turbulent flows which are difficult to be described in the spectral space.

In the present study, we attempt to apply a closure to the fourth-order structure functions in the physical space, and discuss the possibility to reproduce the correct scaling behaviors. In particular, we discuss the pressure effects in the next-order structure function equations. The present approach is expected to inspire future studies on closure models in the physical space for describing complex turbulent flows.

2 The model

We consider a statistically stationary homogeneous isotropic turbulence field such that for any statistical quantity, there is ∂_t〈·〉 = 0, where t is time, and 〈〉 is ensemble average. The classical Kolmogorov equation can be simply written as ^[7-9]

(1)

where υ is the kinematic viscosity, r is the two-point distance, and ϵ is the dissipation rate. D₁₁(r)=〈 (u₁(re₁) -u₁(0))²〉 and D₁₁₁(r)=〈 (u₁(re₁) -u₁(0))³〉 are the second-order and third-order longitudinal velocity structure functions, respectively, with e₁ being the unit vector in the x₁-axis. Clearly, this equation shows the relation between the second-order and third-order structure functions.

The relation between the third-order and fourth-order structure functions can be found in Ref. [10]. Neglecting the non-stationary term ∂_t D₁₁₁ and assuming the first-order Taylor approximation that Z₁₁₁ ≈ 6C in Ref. [10], we obtain

(2)

where D₁₁₁(r)=〈 (u₁(re₁) -u₁(0))⁴〉 and D₁₁₂₂(r)=〈 (u₁(re₁) -u₁(0))² (u₂(re₁) -u₂(0))²〉 are the fourth-order longitudinal and cross velocity structure functions, respectively, and T₁₁₁(r) is a velocity-pressure correlation function.

In order to close the fourth-order structure functions in Eq. (2), we employ the extended scale similarity (ESS) theory ^[11-12] which describes the relation between the fourth-order and second-order longitudinal structure functions, and the quasi-normal assumption which considers the relation between the fourth-order longitudinal and cross structure functions ^[6].

The ESS theory implies that in a wide range, the fourth-order and second-order longitudinal structure functions can be expressed by a constant scale-similarity fractal scaling. From Table 2 of Ref. [12], this scaling can be estimated as p=1.28/0.70 ≈ 1.83. McComb et al. ^[13] showed that the ESS theory is valid at almost all scales in homogeneous isotropic turbulence (see Fig. 1 of Ref. [13]), which directly supports the present closure idea. Therefore, the following relation can be obtained:

(3)

where c_p is a constant (see Subsection 4.1 for discussion).

Moreover, there are quite few studies on different fourth-order structure functions. Therefore, the quasi-normal assumption might be the only choice in the present model. From Eq. (3) in Ref. [6], by using the isotropy condition and the scaling in the inertial range (Subsection 4.3 for details), we can obtain

(4)

Although we have not found any direct support for this relation, there is numerical evidence that implies the rationality of using the quasi-normal assumption to consider the relation between the fourth-order longitudinal and cross structure functions. An example was given in Fig. 4(a) of Ref. [6]. The quasi-normal (i.e., the "Gaussian" state) approximation works better than the RE approximation in many situations in equilibrium turbulence and low Mach turbulence. However, in order to consider the non-equilibrium procedures, one needs to consider the Gaussian-RE line ^[6] instead of the quasi-normal assumption in Eq. (4). The Gaussian-RE line involves an additional parameter on the non-equilibrium property ^[14]. In the present paper, we only focus on the statistically stationary turbulence where the flow cannot be non-equilibrium. Thus, the quasi-normal assumption is acceptable. Actually, as will be discussed in Subsection 4.3, different formulations for the cross structure function do not qualitatively affect the results.

From Eqs. (1) -(4), we can finally obtain the following equation:

(5)

This equation keeps only the second-order longitudinal velocity structure function rather than the higher-order moments. However, there is still a velocity-pressure correlation term T₁₁₁ on the right-hand side which cannot be simply neglected. There is currently no convincing theory for closing this terms, while existing models often require complicated phenomenological assumptions ^[15-17]. In order to simplify this model, we note the observations from Figs. 1 and 2 of Ref. [10] that T₁₁₁ is approximately proportional to both and at all scales. Therefore, by using the quasi-normal assumption in Eq. (4) and the isotropic condition, we assume either a proportional relation

(6)

with α being a model constant, or a proportional relation

(7)

with β being another model constant.

We will prove that from a scaling law, the values of α and β can be analytically obtained. Suppose the classical 2/3 scaling in the inertial range, i.e., . From Eqs. (5) and (7), we can finally obtain

(8)

which yields α=(p-1) /(p+3) for large r (corresponding to the inertial scales). A value of α=0.1718 can be derived using the value of p=1.83 obtained previously by the ESS theory. Similarly, from Eqs. (5) and (7), we can finally obtain

(9)

which yields β=(p-1) /4 for large r. The model constant of β=0.2075 is calculated with the same value of p=1.83 by the ESS theory.

Combining Eqs. (6) and (5) leads to the α formulation,

(10)

while substituting Eq. (7) into Eq. (5) leads to the β formulation,

(11)

Equation (10) or (11) then defines a closure model on velocity structure functions in isotropic turbulence.

Similarly, if we consider the anomalous scaling laws ^[18] instead of the 2/3 scaling in the inertial range, the values of α and β can also be calculated. For example, supposing D₁₁∝ r^0.7 will yield the constants α=(0.7p-2/3) /(0.7p+2) and β=3(0.7p-2/3) /8. When p=1.8, they are, respectively, α=0.187 2 and β=0.230 4.

3 Results

Both Eqs. (10) and (11) are nonlinear. Thus, it is difficult to find analytical general solutions to them. Instead, we present the corresponding numerical solutions in this section to validate these models. Some typical parameters are defined as initial conditions such as

(12)

Here, the value of D₁₁(0.01) is arbitrarily defined, while the values of D₁₁^'(0.01) and D₁₁^"(0.01) are defined to guarantee that D₁₁(r)∝ r² at very small r. The α formulation and β formulation are tested independently. We also perform two groups of simulations to compare the results obtained by different scaling laws for D₁₁. The parameters c_p=0.01 and υ=0.1 are artificially defined. Results are shown in Fig. 1, where the scaling exponent is calculated as n(r)= r D₁₁^'/D₁₁ ^[9]. Clearly, in each case, the scaling exponent is 2 at small scales and approaches the expected scaling (i.e., 2/3 and 0.7, respectively) at large scales. This behavior is qualitatively in agreement with the literature ^{[9, 19]}, but in the present study, all curves show oscillations, which contradict the monotonous values in the literature. These oscillations are perhaps caused by the inappropriate assumptions on the pressure term, and will be discussed in Subsection 4.2. Also, in each group, the α and β formulations yield similar results, illustrating the consistency between the two formulations.

4 Discussion 4.1 Model parameters

In the previous section, we have shown that the α formulation and β formulation produce similar results, while the values of α and β are directly related to the asymptotic scaling exponent in the inertial range. In this section, the sensitivities to other model parameters, such as the ESS coefficient c_p and the viscosity υ, will be discussed. For brevity, we will use the α formulation with the classical Kolmogorov scaling r^2/3 in the following discussion.

From the ESS theory, the self similarity between structure functions is fractal (i.e., p=1.83 instead of p=2). We emphasize that this does not mean that the Galilean invariance is broken, since all phenomena in the ESS theory must be observed by using the root-mean-square velocity u^' for nondimensionalization. This implies that although formally c_p is fractal, the turbulence dissipative scale is always considered to determine the coefficient c_p. The sensitivity of the model to different values of c_p is shown in Fig. 2. Clearly, when we replace the original length scale r with the fractal length scale , the scaling exponents are in excellent coincidence with each other. Because in Eqs. (10) and (11), the first two terms (which are dominant for large r) have the spatial dimension c_p r^2p+1, replacing r with leads to the vanishing of c_p. In brief, the ESS coefficient c_p represents a spatial scaling with consideration of the turbulence dissipative scale.

The role of υ in the present model is similar to that of c_p, which corresponds to a spatial scaling. From Eqs. (10) and (11), it is clear that c_p~ υ^-2. Hence, in Fig. 3, we can use the length scale υ^-1/pr to achieve coincidences among different cases.

4.2 Corrections to the model of T₁₁₁

As discussed in the previous section, the oscillating phenomenon is unphysical. This may be caused by the inappropriate assumption on the pressure term T₁₁₁. In Section 2, we assumed that the coefficients α and β are constant, which are determined by an asymptotic analysis for large r. However, for small r, there should be another asymptotic limit and thus a transition against r for the coefficients α and β. When r is small, from the Taylor expansion, we have D₁₁(r)∝ r². Then, from Eqs. (5) and (6), we can obtain

(13)

which differs from the inertial-scale results α_l=(p-1) /(p+3) and β_l=(p-1) /4. In order to introduce a transition between these values, we introduce an exponential function, and the transitional α(r) writes

(14)

with c₁ being constant. The results with c₁=0.1 are plotted in Fig. 4 in comparison with the constant α formulation. A better asymptotic behavior is observed with a reduced oscillating amplitude. This suggests the possibility of using a transition model in the closure model, and that the present closure model has potential to be improved by using an appropriate approximation for the pressure term T₁₁₁.

4.3 On the closure of the fourth-order cross structure function

In the present closure model, we use the quasi-normal assumption to close the fourth-order cross structure function. Here, we present more details about this derivation and discuss the role of the Gaussian-RE line.

In the quasi-normal assumption, there are relations,

(15)

Equation (4) can be derived in the inertial range by applying the classical scaling r^2/3 under the isotropy condition . However, this equation is not true in the dissipative range where D₁₁∝ r² and D₂₂=2D₁₁ ^[20], and there will be ^[6]

(16)

Moreover, we can also use the RE assumption to close the fourth-order cross structure function. The RE assumption yields ^[6]

(17)

In the inertial range, if we employ the scaling law D₁₁₁∝ r^q with q ≈1.28 ^[18], there is

(18)

In the dissipative range, the Taylor expansion D₁₁₁∝ r⁴ leads to

(19)

Therefore, there are four different relations between D₁₁₁₁ and D₁₁₂₂, i.e., Eqs. (4), (16), (18), and (19), derived under different assumptions and scales, respectively. All these formulations with various related transitional models are tested, but the unphysical oscillations are not eliminated (for brevity, these results are not plotted). Therefore, we remark that the relation between D₁₁₁₁ and D₁₁₂₂ is not a dominant factor in the present closure.

5 Conclusions

The closure of turbulence field is a longstanding fundamental problem, while most closure models are introduced in the spectral space. In this paper, we present a preliminary attempt to close the fourth-order structure functions via the second-order longitudinal structure function. In this model, the Kolmogorov equation and the next-order equation by Hill and Boratav are used as the basic equations, while both the ESS theory and the quasi-normal assumption for fourth-order moments are employed. In addition, a linear model for the pressure-velocity correlation term is introduced.

The present closure model successfully reproduces the asymptotic scalings, which can be classical r^2/3 or anomalous r^0.7, for both small and large scales in statistically stationary homogeneous isotropic turbulence. Both the α and β formulations lead to similar results. Due to the fractal formulations, the parameters c_p and υ are normalized by the root-mean square velocity u^', implying the relative scales to the dissipative scale. The distance r can be re-scaled by .

The oscillations in the results are unphysical and may be caused by the inappropriate assumption of the pressure-velocity correlation. A transition α formulation with artificial interpolation is shown to improve the results, suggesting that the pressure models have potential to be improved by an appropriate consideration of the pressure-velocity correlations.

Being different from the existing closures in the spectral space, the proposed model has the potential to be employed in complex flows in the physical space. Further improvements of this model are expected to take into account the mean velocity or non-stationarity with appropriate relaxation of the homogeneity condition. Also, in addition to the longitudinal components considered in the present paper, the summation of D_ⅱ which describes anisotropic complex flows could be another possibility.