Anisotropic Characteristics of Turbulence Dissipation in Swirling Flow : A Direct Numerical Simulation Study

1 Institute of Nuclear and New Energy Technology and Key Laboratory of Advanced Reactor Engineering and Safety, Ministry of Education, Beijing 100084, China 2College of Mechanical and Transportation Engineering, China University of Petroleum, Beijing 102249, China 3School of Mechanical Engineering, Northwestern Polytechnical University, Xi’an 710072, China 4China Academy of Space Technology, Beijing 100094, China 5School of Aerospace, Mechanical & Manufacturing Engineering, RMIT University, Melbourne, VIC 3083, Australia


Introduction
Vortex breakdown is an intriguing and important phenomenon that occurs in a variety of natural and technological swirling flows, for example, swirling combustor, cyclone, bathtub vortex, hurricanes and tornadoes, and spiral galaxies.In general, the appearance of swirl is caused by the impartment of rotating motion upon the jet which makes the flow more complicated and has been widely studied in the literature [1][2][3][4][5][6][7][8][9][10].The scope of this study does not include summarizing existing vast literature of swirling flow investigation; thus, only a few investigations on this topic were covered.
For example, Chen and Sun [1] addressed the nonlinear 3D instability of a specific type of viscous swirling flow, the Ossen vortex, by using direct numerical simulation at Re = 5000.They considered the global optimal perturbation as the initial perturbation and characterized different flow regimes in axisymmetric cases.Wang and Chen [2] studied vortex breakdown by solving 3D unsteady Navier-Stokes equations for swirling pipe flows, including the flow structures in the bubble domain and the tails behind the breakdown vortex.Moreover, Shtern et al. [4] studied symmetry breaking in a meridional steady motion of viscous incompressible fluids using the laminar axisymmetric "vortex dynamo." They demonstrated the feasibility of a supercritical pitchfork bifurcation from an initially nonswirling flow to a steady swirling regime, as Re exceeds a critical value.In addition, with regard to model study, di Pierro and Abid [5] investigated modeled inviscid swirling flows, addressing the weak variation of axial and azimuthal velocities.They checked the asymptotic results using numerically computed growth rates of linearized Euler equations for a family of variable-density Batchelor-like vortices as base flows and so forth.
Turbulent flows are well known to have scientific and practical importance.However, the full spectrums of the length and time scales of turbulent flows are impossible to solve via computer simulation.Thus, some types of modeling for Reynolds stresses are needed to simulate high Reynolds turbulence.For example, the turbulent kinetic energy and dissipation are involved in the standard form of the twoequation Reynolds stress turbulence model based on the Boussinesq-type approximation [11]: where  is the turbulent kinetic energy,   is the timeaveraged velocity, and ]  is the isotropic turbulent viscosity.However, this turbulence model is fairly ineffective in simulating anisotropic turbulence due to the theoretical deficiency for anisotropic nature of turbulent motions, although it is widely applied in engineering.Thus, numerous studies have been carried out to further understand the physical aspects of anisotropic turbulent flows and vortex dynamics [12][13][14].
In conclusion, swirling flow is present in anisotropic turbulent flows and is a striking and intriguing case of the generation and breakdown of strong vortices.However, the characteristics of the anisotropic nature of turbulent motion and energy dissipation are not clear, especially in strong swirling turbulent jet flows.Thus, the present study carried out a numerical study on the characteristics of energy dissipation tensor in a rotating jet flow.Three swirling numbers are used, corresponding to weak, intermediate, and strong swirling levels.The relations of the components of energy dissipation corresponding to normal and shearing turbulent fluctuations are explored, including their probability density function.The specific locations in the swirling flows, corresponding to the regions of extremely anisotropic and nearly isotropic turbulent dissipation, are also presented.The locations indicate the correlation of anisotropic turbulence dissipation to the large-scale structure of vortex breakdown and the correlation of isotropic turbulence dissipation to small-scale vortices.

Numerical Description
2.1.Governing Equations.For incompressible Newtonian fluids, the governing equations can be expressed in dimensionless form as follows: where u and  are the velocity and pressure, respectively.Re =  0 ⋅ /] is the Reynolds number, in which  0 is the inflow velocity,  is the jet diameter at the inlet, and ] is the kinematic viscosity.Equations ( 2) and (3) deal with a three-dimensional timedependent flow problem without body force.To solve them, the finite difference method is applied, where the convection term is discretized by upwind compact schemes [15], and the space derivatives and pressure-gradient terms are discretized by fourth-order compact difference schemes [16], respectively.The third-order explicit schemes are used to deal with the boundary points and to maintain the global fourthorder spatial accuracy.The time stepping process is integrated by using the fourth-order Runge-Kutta schemes [17].The pressure-Poisson equation is solved to obtain pressure via the fourth-order finite difference method [18].The validation of the codes has been done in a recently published work on swirling flows [19].

Boundary Conditions.
As shown in Figure 1(a), the flow configuration contains a rectangular flow domain of 20 × 10 × 10, where  is the diameter of the jet inlet and the jets are injected from the inlet with a mean velocity  0 .The entire flow domain is discretized by 512 × 256 × 256 Cartesian mesh grids, and the mesh scale is not larger than Δ = 15.625 m.The Kolmogorov length scale is estimated in the same order of the finest mesh scale, namely,  ∼ (Δ) [19,20].According to the suggestion by Moin and Mahesh [21], the mesh scale in the same order of the Kolmogorov scale is fine enough to capture the smallest scale of turbulence.Thus, the grid testing procedure is omitted here.The time step is Δ = 0.001, and 20000 time steps are simulated for each case.The nonreflecting boundary condition is utilized for the outlet condition [22], and the sidewalls are set as nonslipping wall boundaries.For axial inflow velocity, a hyperbolic tangent profile is used within the range of | −   | < /2: where   is the center position of the jet.The reference system is centered at the inlet with   = 0.  is the inflow momentum thickness (Table 1) and the ratio of the jet width to the inflow momentum thickness is / = 20.For azimuthal inflow velocity, a polynomial expression is used; that is, where the coefficients   are listed in Table 2.The combined profiles of axial and azimuthal velocities are shown in   In addition, swirl number  is defined as the ratio of the axial flux of angular momentum to the axial momentum; that is, Based on ( 4) and ( 5), three swirl numbers are used (Table 2).
From the engineering viewpoint,  = 0.28 is considered as low swirl jets, whereas  > 0.59 is considered as strong swirling flows and  = 0.45 is intermediate.

Reynolds Stress Transport Equation.
The Reynolds-averaged Navier-Stokes equation [23] is usually used to solve the turbulent flows in industrial scales: where ⟨⋅⟩ is the time-averaged operator.⟨      ⟩ is the Reynolds stress tensor, which can be solved through the Reynolds stress transport equation and is given by where   = 2(1/ Re)⟨(   /  )(   /  )⟩ is the turbulent energy dissipation tensor.The Reynolds stress tensor is relevant to the dissipation of turbulent kinetic energy and is a symmetric positive definite tensor.Thus, it corresponds to a diagonalizable matrix, and the anisotropic nature of the turbulent kinetic energy dissipation is indicated from the eigenvalues of the diagonalizable matrix.The anisotropic tensor has three unequal eigenvalues (i.e., at least two of them are unequal).Thus, the discrepancy between the eigenvalues, or equivalently, their relation and distribution, can indicate the anisotropic nature of turbulent kinetic energy dissipation in swirling flows.

Coherent Vortex Structures.
The typical vortex structures are shown in Figure 1(a) when  1 = 0.28 and  3 = 0.59.For  1 = 0.28, the Kelvin-Helmholtz instability dominates the vortex evolution due to the existence of a shear layer between the jets and ambient fluids.In contrast, the strong rotating effect dominates the evolution of the vortex structures for  3 = 0.59.A cone-type vortex breakdown is established, and the K-H instability causes the breakdown of the "cone" feature from its terminal into the turbulent vortex streets.
Figure 2 shows that the typical cross-sectional visualizations of the vortex structures are visualized for  3 = 0.59 (referring to the jet diameter  for the region width).The early evolution of the vortex for  = 5 (Figure 2(a)) has a full ring-type structure, with a positive stream-wise velocity (  > 0) in the inner side of the vortex ring (flood by red color in Figure 2(a) for 0.2 <   / 0 < 0.8) and a negative stream-wise velocity (  < 0) in the outer side of the ring (flood by blue color in Figure 2 Based on the observations of Figure 2, the intrinsic characteristics of the vortex can be reflected by turbulent energy dissipations.However, these characteristics need to be explored to show the anisotropic/isotropic turbulent energy dissipation tensor and its correlation to vortex structures, among others.

Relation of Components of 𝜀
is a symmetric positive definite tensor, which can be diagonalized.In general, three eigenvalues and three eigenvectors can be obtained in the diagonalization process.The eigenvectors are orthogonal and denote the three characteristic directions.The corresponding eigenvalues indicate the magnitudes of dissipation in each direction.Thus, an anisotropic tensor should have three (at least two) different eigenvalues, which correspond to different levels of dissipation in the three characteristic directions.In this way, the characteristics of anisotropic turbulent motion ⟨      ⟩ can be reflected by the characteristics of the eigenvalues and eigenvectors of dissipation tensor   .  = ⟨          ⟩ was used because Re = 5000 is a constant.
First, Figure 3 illustrates the relation of the components of the tensor   .Each point in Figure 3

Distribution of the Eigenvalues.
Second, tensor   = diag(ε 1 , ε2 , ε3 ) was diagonalized, where ε1 , ε2 , and ε3 are the three eigenvalues of the tensor, which were sorted as ε1 ≥ ε2 ≥ ε3 ≥ 0 for simplicity.The distribution and relation of ε are illustrated in Figure 4 using  1 = 0.28 and  3 = 0.59 for the case study.Figure 4(a) ( 1 = 0.28)  shows that the distribution of (ε 1 , ε3 ) has a widespread area below the line of ε3 = ε1 , indicating that it has ε1 > ε3 in most regions.Similar results were obtained in Figure 4(c) for  1 = 0.59, and extreme discrepancies between ε3 and ε1 were observed, even as large as more than two or three orders.This observation indicates the extreme anisotropic characteristics of turbulent kinetic energy dissipation when the swirling level is sufficiently strong.Moreover, for (ε 2 , ε3 ) (Figures 4(b) and 4(d)), although extreme events seem to be prohibited due to the possible existence of superior limit line, the distribution area is open and enlarged as ε2 increases.In conclusion, Figure 4 shows the validation of the anisotropic characteristics of turbulent energy dissipation tensor, especially for strong swirling flows where extremely anisotropic events may occur.
Figure 5 shows the PDFs for different swirling flows.For all the cases, that is,  1 = 0.28 (Figure 5 for the same eigenvalue of ε0 and ε1 |  0 > ε2 |  0 > ε3 |  0 for the same value of PDF  0 , which means ε1 > ε2 > ε3 is true in most ranges; that is, the anisotropic turbulent dissipation occurs in most ranges. (2) A subrange of ε ∈ (  ,  ℎ ), which has (ε 1 ) = (ε 2 ) = (ε 3 ), always exists, except that (ε 1 ) for most strong swirling flows  3 = 0.59.This phenomenon indicates the existence of subrange isotropic turbulent dissipation when the swirling level is not strong.On the other hand, turbulent energy dissipation in the major characteristic direction corresponding to the largest eigenvalues ε1 is always larger than the other characteristic directions for strong swirling flows ( 3 ≥ 0.59), which indicates highly anisotropic turbulent energy dissipation characteristics in the main flow region.Moreover, Figures 5(a) to 5(c) show that the isotropic subrange can be characterized by a regression line, designated as  1 ,  2 , and  3 for  1 = 0.28,  2 = 0.45, and  3 = 0.59, respectively.The data in these subranges are illustrated in Figure 5(d) for clarity.Figure 5(d) shows a perfect linearity between log (ε) and log(ε) for the isotropic subrange, which obtained the following equation using linear regression: log ( (ε)) =  ⋅ log (ε) + , namely  (ε) =  ⋅ ε .(11) The regression coefficients are listed in Table 3, which shows that the power-law exits in the isotropic subrange of turbulent energy dissipation and that the power exponent is approximately −1 to −0.8.Moreover, the power-law occurs mainly for low turbulent energy dissipations, and the width of strong swirling flows is reduced.In other words, weak swirling flow has wider subranges of isotropic turbulent energy dissipation than the strong swirling flow and vice versa.This phenomenon indicates the significant role of swirling motion in augmenting the anisotropic characteristics of turbulence.To speak specifically, it is possible to assume that turbulent kinetic energy dissipation is more intensive and more nonlinear in strongly swirling flows than in weakly swirling flows.Therefore, the subranges of linear relations between probability distribution functions and eigenvalues of dissipation rates are reduced in strongly swirling flows compared to those in weakly swirling flows.

Correlation to Vortex Structures.
Previously mentioned results are based on statistical analysis; thus visualizing the locations of extremely anisotropic and nearly isotropic turbulence dissipation should be helpful to further understand anisotropic turbulent swirling flows.For this reason, the data of |ε 1 /ε 3 | > 100 (see Figures 4(a) and 4(c)) and −0.8 < log ε1 , log ε2 , log ε3 < 0.4 (see Figure 5(d) or equivalently 0.158 < ε1 , ε2 , ε3 < 2.5) was extracted and the locations for these data in the background of vorticities of  = 0.28 and 0.59 were shown.
Figure 6 shows that the locations of extremely anisotropic turbulence dissipation (e.g., |ε 1 /ε 3 | > 100, designated by colored spheres) are concentrated in the central region of the jet (Figure 6(a)), corresponding to strong twisted vortices when the swirl level is low.On the other hand, the locations of nearly isotropic turbulence dissipation (0.158 < ε1 , ε2 , ε3 < 2.5) are dispersed in the peripheral region of the jet (Figure 6(b)), corresponding to the region with small-scale and weak turbulence.Moreover, the locations of extremely anisotropic turbulence dissipation are concentrated in the central region of vortex breakdown (Figure 6(c)), corresponding to strong swirling large-scale vortex structure when the swirl level is large.The locations of nearly isotropic turbulence dissipation are dispersed in the peripheral region of the vortex breakdown, such as strong small-scale vortices (Figure 6(d)).
Figure 6 conclusively shows that strong anisotropic turbulence dissipation occurs concentratively in the vortex breakdown region or is closely related to the large-scale vortex structure.On the other hand, nearly isotropic turbulence dissipation occurs dispersively in the peripheral region of the strong swirling flows, that is, closely related to small-scale vortices.

Conclusion
This work was carried out to investigate the physical aspects of anisotropic turbulent motions and dissipations in swirling flows.Based on the observation and analysis of the DNS results, the following results were found. (

Figure 1 :
Figure 1: Sketch of simulation setup: typical vortex streets (a) and inflow velocity profiles (b).

Figure 1 (
Figure 1(b).In addition, no initial turbulence is introduced to show the intrinsic full evolution of coherent vortex structures and interactions.In addition, swirl number  is defined as the ratio of the axial flux of angular momentum to the axial momentum; that is,
(a) for −0.8 <   / 0 < −0.2).Central recirculation then occurs in the center of the jet, where the velocities are weakly negative.The ring structure of vortex for  = 15 (Figure2(b)) evolves into more complicated structures.The diameter of the vortex ring tube becomes larger and starts to break down, and braid vortices are then formed.The complex vortex structures indicate the existence of turbulent vortex motions and dissipation of turbulent energy.

Table 1 :
Typical real values of variables used in the numerical simulation.

Table 2 :
Coefficients of the expressions for inflow velocity   .
corresponds to the energy dissipation characteristics at each point location in the flow, with   part in the -axis and  , ̸ = part in the axis, because the dissipation tensor   has different values for different locations in the flow domain.The magnitude of  , ̸ = fluctuation is almost the same as   , despite the swirling levels.Thus the energy dissipation of ⟨

Table 3 :
Coefficients of the regression lines of  1 ,  2 , and  3 .