ABSTRACT

Direct numerical simulations (DNS) of non-rotating and rotating turbulent channel flow were conducted. The data base obtained from these DNS simulations was used to investigate the prominent coherent structures involved in the turbulence generation cycle. Predictions from three theoretical models concerning the formation and evolution of sublayer streaks, three-dimensional hairpin vortices and propagating plane waves were validated using visualizations from the present DNS data. Quadrant analysis was used to determine a phase shift between the fluctuating streamwise and wall-normal velocities as a characteristic of turbulence production in the suction region at a low rotation number.

Keywords:

Direct Numerical Simulation, Rotating Turbulent Flows, Turbulence Structures

1. Introduction

The scientific field of turbulence has posed long-standing challenges to researchers due to the inherent chaotic and irregular motions which define turbulent flows. Since turbulent flows, not laminar flows, are predominantly found in nature and with the prevalence of rotation-dependent machinery in engineering, a physical understanding of turbulent flow in these systems is a necessity for engineering and scientific analyses. The fundamental structures involved in the turbulence sustenance cycle for simple turbulent channel flow, defined as a turbulent flow field possessing only one mean flow gradient, have been well-documented in experimental [1] and DNS studies [2]. Near the channel walls, interactions between streamwise vortices and streamwise elongated sublayer streaks result in ejections of the latter structures from the near-wall region. The streaks break down into smaller instabilities which are swept back towards the wall and reform into vortices and streaks. This cyclical process is referred to as the bursting cycle and the mechanism through which turbulence is sustained [3]. Despite the lack of a constitutive theoretical model for governing the overall structure of the turbulence sustenance cycle, theoretical models have been proposed for selected components of the cycle. Using direct numerical simulation (DNS), it is instructive to assess the validity of the theoretical model predictions. The present work examines the fundamental structures of turbulence with comparison to the theoretical models proposed by [4] , [5] and [6] for the formation and evolution of such structures with the DNS data base.

When turbulent channel flow is subject to rotation in the spanwise direction, rotation-induced body forces (Coriolis, centrifugal) generate secondary flows. Regions known as the suction and pressure sides emerge, distinguished by reduced and elevated turbulence levels, respectively [7]. In the present work, a comprehensive investigation of the effects of rotation on channel turbulence including general turbulence statistics, higher-order statistics, energy budgets and coherent structures, was conducted. In the pressure region, the streaky and vortical structures associated with the turbulence sustenance cycle persist and are the primary structures which contribute towards the energetic ejection and sweep events near the channel walls [2] [4]. In the suction region, rotational forces induce flow re-laminarization and hence the turbulence sustenance cycle disappears along with the coherent structures which define it. With the absence of the turbulence generation cycle in the suction region, a sudden shift in phase relationship between streamwise and wall-normal fluctuating velocity was identified and examined using quadrant analyses. Additional aspects of the laminarized flow regime were also investigated including characteristical structures of laminar-to-turbulence transition and new observations into the mechanisms of turbulence production are provided.

The present work also provides an examination of higher-Reynolds number effects on rotational turbulence through direct comparison of lower and higher-Reynolds number simulations subject to spanwise rotation, specifically in regards to the characteristics of the coherent structures which contribute to turbulence production. Consideration of higher-Reynolds number effects is necessary due to the prevalence of higher-Reynolds number flows in engineering applications such as gas turbine blade and rotating turbomachinery design [8]. This research has been conducted by members of Prof. Sedat Biringen’s computational fluid dynamics group within the Ann and H.J. Smead Department of Aerospace Engineering Sciences at the University of Colorado-Boulder.

2. Numerical Method and Case Descriptions

The time-dependent, three-dimensional incompressible Navier-Stokes and energy equations were numerically integrated in a doubly periodic (in x and z-directions) channel flow using a fractional step method [9] [10]. With all spatial coordinates non-dimensionalized by the channel half-height $\delta$ and velocities by the laminar centerline velocity $u_c$ , the Navier-Stokes and energy equations read (in conservative form)

$\frac{\partial u_i}{\partial x_i}=0$ (1)

$\frac{\partial u_i}{\partial t}+\frac{\partial u_i{u}_j}{\partial x_j}=-\frac{\partial p}{\partial x_i}+\frac{1}{R{e}_c}\frac{{\partial }^2{u}_i}{\partial x_j\partial x_j}-R{o}_c{\epsilon }_{ijk}\frac{{\Omega }_j}{\Omega }{u}_k$ (2)

$\frac{\partial \theta }{\partial t}+\frac{\partial \left(\theta u_j\right)}{\partial x_j}=\frac{1}{R{e}_{c}Pr}\frac{{\partial }^2\theta }{\partial x_i\partial x_j}$ (3)

where $R{e}_c=u_c\delta /\nu$ , $Pr=\nu /\kappa$ , $\nu$ is the kinematic viscosity, $\kappa$ is the thermal diffusivity and the vector $u=\langle u\mathrm{,}v\mathrm{,}w\rangle$ is composed of three velocity components in the x (streamwise), y (wall-normal) and z (spanwise) directions, respectively. The rotation number (or Rossby number) is defined as

$R{o}_c=2\Omega \delta /u_c$

where $\Omega$ is the spanwise angular rotation vector. p is the non-dimensional effective pressure ( $p=p_0-\left(1/8\right)R{o}_c^2{r}^2$ ) which combines the static pressure ( $p_0$ ) and centrifugal force, r represents the nondimensional distance away from the axis of rotation. Also, t is non-dimensional time and $\theta$ is non-dimensional temperature ( $\theta =\left[T-T_L\right]/\left[T_U-T_L\right]$ ) with $T_U$ and $T_L$ representing dimensional temperatures on the upper and lower walls, respectively. Equation (2) is uncoupled from Equation (3) and buoyancy effects are neglected. The Prandtl number ( $Pr$ ) was kept constant at 0.71.

The flow geometry of the present DNS is shown in Figure 1. The velocity flow fields are assumed to be statistically homogeneous in the streamwise (x) and spanwise (z) directions allowing periodic boundary conditions in those directions. No-slip conditions are imposed on the rigid channel walls ( $y=0,2$ ). An isothermal boundary condition on the dimensionless temperature $\theta$ ( ${\theta }_U=0.5$ and ${\theta }_L=-0.5$ ) was imposed at the upper and lower walls. Temperature coupling was initiated once the turbulence field had reached a quasi steady-state solution and a linear distribution for the temperature field was assumed as an initial condition in the wall-normal direction. The time integration

scheme is a modified semi-implicit Adams-Bashforth/Crank-Nicolson method, which is second order accurate in time and fourth order accurate in space. More information on the code scheme and verification may be found in [9] and [10].

The full listing of simulations and their corresponding domain lengths and grid resolutions are found in the case descriptions (Table 1). Four simulation cases A-D were conducted for Rossby numbers $R{o}_b=0$ , 0.2, 0.5 and 0.9. An additional higher-Reynolds number simulation case E was conducted for $R{o}_b=0.2$ to examine higher-Reynolds number effects on turbulence production with rotation. For simulation cases A-D, the Reynolds number based on the laminar centerline velocity, $R{e}_c=8000$ , was kept constant which resulted in a friction Reynolds number $R{e}_{\tau }=200$ for the no-rotation case ( $R{o}_b=0$ ). The low Reynolds number was chosen such that reasonable comparisons could be made between the present DNS and other DNS studies of rotating turbulent channel flow with similar Reynolds numbers [7] [11]. The asymmetric velocity distributions due to the rotational effects decreased the value of $R{e}_{\tau }$ which was calculated as an average between the two walls. Simulations were performed at constant mass flux which resulted in a bulk Reynolds number $R{e}_b=u_b\delta /\nu$ , where $u_b$ is the mean bulk velocity defined as

$u_b=\frac{1}{2\delta }{\displaystyle {\int }_0^{2\delta }}\overline{u}\text{d}y$ (4)

where $\overline{u}$ denotes a plane-averaged quantity. The domain lengths for the simulations were selected as $L_x=4\text{π}\delta$ , $L_y=2\delta$ and $L_z=2\text{π}\delta$ such that two-point spatial autocorrelations in the streamwise and spanwise directions converged to zero at the largest separations.

Mesh independence was established by designating two high-resolution cases with grid numbers $nx\times 2ny\times nz$ for case A ( $R{o}_b=0$ ) and $2nx\times ny\times 2nz$ for case D ( $R{o}_b=0.9$ ). For both cases, the distributions of mean velocity, Reynolds stresses and turbulent kinetic energy budgets compared very favorably to those of the original simulations, demonstrating the selected meshes of the present DNS cases were mesh invariant. The grid spacing is also comparable to other DNS studies of spanwise-rotating turbulent channel flow [7] [11].

Superscript + refers to nondimensionalization by the friction velocity, $u_{\tau }=\sqrt{\nu {\partial \overline{u}/\partial y|}_{\text{wall}}}$ or friction temperature, $T_{\tau }=\left(\kappa /u_{\tau }\right){\partial \theta /\partial y|}_{\text{wall}}$ . The global friction velocity for the rotational cases is denoted as $u_{\tau }=\sqrt{u_{\tau s}^2/2+u_{\tau p}^2/2}$ where $u_{\tau s}$ and $u_{\tau p}$ are the local friction velocities at the suction and pressure walls, respectively; an equivalent calculation of the global friction temperature was used for the thermal statistics. Unless otherwise specified, all coordinate directions are non-dimensionalized by $\delta$ . For all cases, the governing equations were integrated until both the friction Reynolds number $R{e}_{\tau }$ and friction temperature converged availing a sufficiently long time window ( $t^{+}=t{u}_{\tau }^2/\nu \ge 1000$ ) to calculate statistics.

The data base used in the present work is the same one used in [12] and [13] and more information regarding the verification of simulation cases A-E via comparison to other experimental and computational studies may be found there.

In this section, visualizations from turbulent channel flow simulation case A are used to substantiate the theoretical model predictions proposed by [4] [5] [14].

The Landahl model proposes the formation of sublayers streaks, or streamwise elongated $u^{\prime }$ structures, is a consequence of algebraic instabilities commonly found in turbulent flows [15]. The streaks are generated by the continuous linear growth of these algebraic instabilities in the streamwise (x) direction, necessitated by a linear temporal growth of total streamwise momentum which continues indefinitely until viscous forces impede growth. [4] used the conditional sampling technique of variable interval time averaging (VITA) to obtain flow visualizations demonstrating the temporal structural evolution suggested by the theoretical model. The model proposed two structure classes which formed from the original algebraic instabilities: symmetrical and asymmetrical structures which correspond with oblique defomation angles in the spanwise direction of $\theta =0^{\circ }$ and $\theta =5^{\circ }$ , respectively. As the symmetrical structures did not demonstrate streamwise elongation over time, [4] related the asymmetrical structures to sublayer streaks. The asymmetrical structures demonstrated a consistent pattern dominated by a high-speed structure side-by-side in the spanwise direction with a low-speed structure. Once these structures elongated, an irregular wavy appearance was observed consistent with the oscillatory motion of streaks [1].

The present DNS results are also used to evaluate the coherent structure formation scheme proposed in the experimental study by [5] which used Particle Image Velocimetry (PIV) to examine a low-pressure-turbine blade flow regime. In contrast to the study by [4] which studied the initial development and evolution of sublayer streaks, [5] examined the formation of three-dimensional coherent structures which accompanied streak breakdown into turbulence. It was proposed that the breakdown of elongated sublayer streaks in the near-wall region induced three-dimensional vortical structures which manifested on the streak flanks. These large-scale structures were related to hairpin vortices [16] , characterized by spanwise vorticity on the top and wall-normal vorticity on the bottom legs. Smaller observed vortical structures, such as vorticity tubes, were proposed to be residuals of the hairpin vortices which contributed to the sinuous motion of the streaks. The breakdown of streaks and these vortical structures generate high velocity fluctuations and lead to the formation of other large-scale turbulence structures.

The study of [14] postulated the existence of secondary instabilities, mainly propagating plane waves, which serve as a trigger for the interactions between these primary turbulence structures. [17] was one of the first to qualitatively analyze the role of these secondary instabilities to the transition process and later studies of [18] [19] supported this mechanism by discovering these secondary instabilities travel obliquely to the streamwise direction and contributed to flow de-stabilization. The flow control studies by [20] and [21] additionally affirmed the proposal by demonstrating that the phase randomization of a small subset of propagating modes, primarily within the energy-containing scales, reduced turbulent kinetic energy and drag by significant amounts. The present work aims to substantiate the low-Reynolds number results ( $R{e}_{\tau }=120$ ) of [6] for a turbulent flow field at $R{e}_{\tau }=200$ and demonstrate the existence of these propagating waves.

2.1. Landahl Model

The Landahl model qualitatively examined the formation and evolution of sublayers streaks in the turbulence system cycle and used the variable interval time averaging (VITA) method to predict structural characteristics suggested by the theoretical model. From the present DNS results in case A, maps of fluctuating streamwise velocity $u^{\prime }$ are obtained in the x-z plane of peak turbulent kinetic energy production ( $y^{+}=15$ ). To obtain Figure 2(a), the variable interval spatial averaging (VISA) method [22] was applied with an averaging length of 200 wall units in the streamwise and spanwise directions. In accordance with the VISA method, a detection criterion was used to isolate islands of high local $u^{\prime }$ variance and the space-time position of these islands was tracked to visualize the temporal evolution of sublayer streaks. These VISA-educed structures are compared to the asymmetrical structures obtained from the modeled VITA results of [4] and numerical VISA results of [23] in Figure 2(b) and Figure 2(c), respectively.

All three patterns demonstrate a structural inclination of approximately $\theta =5^{\circ }$ in accordance with the theory of [4]. The spanwise lengths of the streak patterns are also shown to be very similar at approximately $z^{+}=100$ , corresponding to the accepted mean spacing between sublayer streaks [23]. The numerical and modeled results display a comparable streamwise streak length of approximately $x^{+}=500$ , demonstrating the three streak patterns are in a similar stage of development. The present DNS results in Figure 2(a) show more similarity with the theoretical model of [4] , demonstrating an abrupt onset of the

disturbed region and more discernable wavy structure in correspondence with the modeled results in Figure 2(b).

Similarities with the theoretical model of [4] are also observed in the temporal evolution of VISA-educed structures from the present DNS in Figure 3. The familiar spanwise array of aligned high-speed and low-speed streaks is observed and this pattern demonstrates the expected streamwise advection and elongation for increasing $t^{+}$ . Similar to the modeled asymmetrical structures, the present DNS results show a consistent oscillatory shape with an increasing amount of inflection points over time.

2.2. Lengani-Simoni Model

The [5] model proposed three-dimensional coherent structures accompanied the sublayer streaks of turbulence sustenance cycle. These three-dimensional vortical structures manifested as vorticity tubes on the streak flanks or hairpin vortices which envelop the streak. In Figure 4, a three-dimensional contour representation of the near-wall coherent structures is shown for a single low-speed sublayer streak near the bottom channel wall ( $y=0$ ) for the periodic channel. The blue streaky structure denotes high levels of negative streamwise fluctuating

velocity and the yellow hairpin vortex is composed from high levels of combined spanwise ( ${{\omega }^{\prime }}_z$ ) and wall-normal vorticity ( ${{\omega }^{\prime }}_y$ ). The smaller yellow vortical structures aligned with the streak in the spanwise direction denote high levels of streamwise vorticity ( ${{\omega }^{\prime }}_x$ ) and the vorticity field is filtered such that the coherent structures are isolated from one another.

Both views of the streak in Figure 4(a) and Figure 4(b) demonstrate good agreement with the proposed coherent structure organization by [5]. The sinuous sublayer streak is accompanied by a large hairpin vortex, characterized by spanwise and wall-normal vorticity, which envelops a portion of the streak. The smaller vortical structures are observed to coincide with the sinuous, wavy motion of the sublayer streak in Figure 4(b). These numerical results support the proposed development of three-dimensional vortical structures induced by sinuous sublayer streaks close to breakdown.

2.3. Sirovich Model

To visualize the spanwise-propagating plane waves proposed in [6] , the method of principal orthogonal decomposition was applied to a periodic channel. Principal orthogonal decomposition (POD), also known as the Karhunen-Loeve decomposition, is a procedure for extracting the coherent motions from two-point velocity correlations which contain the most energy [24]. Further detailed in [25] , this method is based on the decomposition of the velocity field

$u_i={\displaystyle \sum _{q=1}^{ny}}{a}^q{\varphi }_i^q$ (5)

where $a^q$ and ${\varphi }_i^q\left(y\mathrm{,}z\right)$ are the basis-function coefficients and basis functions (eigenfunctions), respectively. q is the quantum number which refers to two-point separations in the inhomogeneous direction y. In addition, the basis-function coefficients correspond to their respective eigenvectors via

$a^n=\frac{1}{2\delta }{\displaystyle {\int }_0^{2\delta }}{u}_i{\varphi }_i^q\text{d}y$ (6)

with satisfaction of the orthonormality condition for the eigenfunctions. The average mean energy of the velocity field is defined by

$E_{ii}=\frac{1}{2\delta }{\displaystyle {\int }_0^{2\delta }}\frac{1}{2}{u}_i^2\text{d}y$ (7)

and through substitution of $u_i$ in Equation (5), the contribution of energy from various N modes is shown through the partial sum

$E_N={\displaystyle \sum _{q=1}^N}\frac{1}{2}{a}^q{}^2$ (8)

Using the above relations, an eigenvalue problem is created using the two-point autocorrelation tensor $R_{ij}$

${\displaystyle {\int }_0^{2\delta }}{R}_{ij}{\varphi }_j^q\text{d}y=\lambda {\varphi }_i^q$ (9)

and if a structure contributes energy to the Reynolds stress tensor, it will dominate the two-point correlation statistics and manifest in the POD [25]. For the periodic channel with periodicity in the streamwise (x) and spanwise (z) directions, the three-dimensional two-point correlation tensor is [26]

$R_{ij}\left(r_x\mathrm{,}y\mathrm{,}{y}^{\prime }\mathrm{,}{r}_z\mathrm{,}t\right)=\langle u_i\left(x\mathrm{,}y\mathrm{,}z\mathrm{,}t\right)u_j\left(x+r_x\mathrm{,}{y}^{\prime }\mathrm{,}z+r_z\mathrm{,}t\right)\rangle$ (10)

where $r_x$ and $r_z$ represent the two-point separations $x-x^{\prime }$ and $z-z^{\prime }$ in the streamwise and spanwise directions, respectively. The brackets $\langle \rangle$ denote ensemble averaging in time and the homogeneous x and z directions.

For multi-dimensional POD analysis and application to three-dimensional turbulent channel flow, it is fitting to convert the two-point correlation tensor $R_{ij}$ into the spectral density correlation tensor ${\Phi }_{ij}$ [26]

${\Phi }_{ij}\left(k_x\mathrm{,}y\mathrm{,}{y}^{\prime }\mathrm{,}{k}_z\right)=\frac{1}{4{\text{π}}^2}{\displaystyle \int }{\displaystyle \int }{\text{e}}^{-i{k}_x{r}_x-i{k}_z{r}_z}{R}_{ij}\left(r_x\mathrm{,}y\mathrm{,}{y}^{\prime }\mathrm{,}{r}_z\right)\text{d}{r}_x\text{d}{r}_z$ (11)

such that the flow field may be expressed as a function of streamwise ( $k_x$ ) and spanwise wavenumbers ( $k_z$ ). For all wavenumber combinations, a $\Phi$ matrix of dimensions $3ny\times 3ny$ is assembled

$\Phi =\left[\begin{array}{ccc}{\Phi }_{11}& {\Phi }_{12}& {\Phi }_{13}\\ {\Phi }_{21}& {\Phi }_{22}& {\Phi }_{23}\\ {\Phi }_{31}& {\Phi }_{32}& {\Phi }_{33}\end{array}\right]$ (12)

Hence, a new eigenvalue problem ( $A\varphi =\lambda \varphi$ ) is generated where a, $\varphi$ and $\lambda$ represent the integrated spectral density correlation tensor, corresponding eigenfunctions and eigenvalues, respectively. To approximate the integral of $\Phi$ in the wall-normal y direction, the weighting function matrix D is calculated using the following trapezoidal numerical approximation [27]

${\displaystyle {\int }_0^{2\delta }}\Phi \text{d}y={\displaystyle \sum _{i=1}^{ny}}{\omega }_i{\Phi }_i=\frac{1}{2}{\displaystyle \sum _{i=2}^{ny}}\left(y_i-y_{i-1}\right)\left[\Phi \left(y_{i-1}\right)+\Phi \left(y_i\right)\right]$ (13)

where ${\Phi }_i$ is the value of $\Phi$ at a discrete grid point and ${\omega }_i$ is the corresponding weight function. To apply standard numerical eigenproblem solution techniques, it is required that the matrix-valued function a, or $\Phi D$ in the eigenproblem, be symmetric. This is accomplished through the following convolution

$\Phi D\varphi =\sqrt{D}\Phi \sqrt{D}\varphi =\lambda \varphi$ (14)

Once the D and $\Phi$ matrices are created, the resultant eigenproblem is solved to return a system of eigenvalues and eigenfunctions for various modal combinations: $\left(k_1,k_3,q\right)$ .

For analysis, instantaneous fluctuating velocity fields were collected from simulation case A for a large time window. For an accurate and relevant comparison to the results of [6] regarding the wavenumber ranges, these fields were interpolated onto an approximately equivalent domain: $L_x=5\delta$ , $L_y=2\delta$ and $L_z=\text{π}\delta$ . Using the previously outlined procedure for periodic POD analysis, eigenvalues and eigenfunctions were obtained for a large number of modes.

In Table 2, the top fifteen energetic modes and their corresponding energy fractions for simulation case A are compared to the results of [6] for a low Reynolds number $R{e}_{\tau }=125$ . m and n refer to indexes for the streamwise ( $k_1=2\text{π}m/L_x$ ) and spanwise wavenumbers ( $k_3=2\text{π}n/L_z$ ), respectively. It is demonstrated that nonpropagating modes, despite being a small fraction of the total number of modes, are the most energetic. Table 2 also demonstrates identical propagating modes which possess the highest energetical content such as the (1,3,1) and (1,2,1) modes, and that a small range of spanwise wavenumbers (n = 2 - 4) captures most of the highly energetic propagating modes. This finding supports the flow control design of [21] , which found that randomizing a small range of inertial scales reduced turbulent drag by significant amounts.

To visualize the presence and interactions of these propagating plane waves, [14] derived a frequency and corresponding wave speed for these structures. Using the methodology of [14] to calculate the wave speed for the most energetic propagating modes from Table 2, we plot a normal speed locus for a discrete number of points for the periodic channel in Figure 5. It is demonstrated that most waves do propagate at an oblique angle to the streamwise direction in accordance with the expectation of [6].

[6] additionally proposed a relationship between these plane wave envelopes and other energetic turbulence structures which advected with the mean flow velocity in the near-wall region, mainly the streaky structures and three-dimensional vortices discussed in the previous theoretical models, which form the bursting process. As the energy contained within the propagating modes is relatively small compared to the non-propagating modes (Table 2), these plane wave modes were proposed by [14] to be a triggering mechanism for bursting events with the non-propagating modes providing the energy cascade necessary for the bursts to occur. Through a decomposition of the Reynolds shear stress into separate contributions from the non-propagating and propagating modes, [6] demonstrated the presence of both modes were necessary for strong turbulent activity. These discoveries bode well for future flow control efforts as the present work and [14] have collectively shown that most of the turbulent kinetic energy is contained to a small range of modes.

3. Rotational Turbulence

3.1. Introduction

In this section, the DNS data base from simulation cases A ( $R{o}_b=0$ ), B ( $R{o}_b=0.2$ ), C ( $R{o}_b=0.5$ ) and D ( $R{o}_b=0.9$ ) are examined for effects of rotational forces on turbulence over a wide range of rotation rates. In spanwise-rotating turbulent channel flow, the Coriolis force acts in the wall-normal direction, resulting in asymmetry across the channel and the creation of two distinct flow regimes: the pressure and suction regions. In the pressure region of the channel, secondary flow circulation and high levels of turbulence are present and in the suction region, re-laminarization of the regime results in low levels of turbulence.

3.2. General Turbulence and Thermal Statistics

Rotational effects on the mean velocity and temperature profiles are shown for simulation cases A-D in Figure 6 and Figure 7, respectively. In Figure 6, the mean velocity profile is symmetric about the channel centerline ( $y=1$ ) for case A ( $R{o}_b=0$ ). With system rotation, the mean velocity distributions become asymmetric as the flow regime is separated into the pressure and suction regions. For rotational cases B-D, a laminar-like (parabolic) profile is observed near the suction wall ( $y=0$ ) which is characteristic of the suppressed turbulence in the suction region. As the flow progressively relaminarizes with increasing rotation number, the suction region expands. Near the pressure wall ( $y=2$ ), a constant gradient of $-2\Omega$ is shown in the mean velocity profiles at all rotation rates which is consistent with previous DNS results in [7] and [11] In Figure 7, the thickness of the thermal diffusive sublayer, characterized by large

near-wall temperature gradients, is broader near the suction wall than the pressure wall for rotating simulation cases B-D. As the rotation number increases, the size of the diffusive layer in the suction region increases and the mean temperature profile shifts towards the pressure wall [28].

In Figure 8, the Reynolds stress distributions are shown for simulation cases A-D. The suppressed amplitudes near the suction wall compared to those near the pressure wall are displayed for all figures in the rotational cases. In Figure 8(a), ${\overline{u^{\prime }{u}^{\prime }}}^{+}$ amplitudes are initially increased for case B ( $R{o}_b=0.2$ ) near the pressure wall. However, the pressure wall amplitudes gradually decrease with increasing rotation number in cases C ( $R{o}_b=0.5$ ) and D ( $R{o}_b=0.9$ ). In Figure 8(b) and Figure 8(c), ${\overline{v^{\prime }{v}^{\prime }}}^{+}$ and ${\overline{w^{\prime }{w}^{\prime }}}^{+}$ amplitudes are shown to monotonically increase with rotation number in the pressure region. In Figure 8(d), ${\overline{u^{\prime }{v}^{\prime }}}^{+}$ remain relatively consistent in the pressure region for the rotational cases.

To illustrate the effects of rotation on wall shear stress and heat transfer, the dimensionless friction Reynolds ( $R{e}_{\tau }$ ) and Nusselt ( $Nu$ ) numbers for both channel walls are provided for the present simulation cases A-D in Table 3. The introduction of rotation is shown to initially decrease $R{e}_{\tau }$ on the suction wall while increasing $R{e}_{\tau }$ on the pressure wall. At higher rotation numbers, these $R{e}_{\tau }$ trends are shown to significantly weaken or even reverse in case D ( $R{o}_b=0.9$ ) for $R{e}_{\tau }$ on the pressure wall. These results correspond well with

[7] which demonstrated that $R{e}_{\tau }$ on both walls trended towards convergence at high rotation numbers until the eventual full re-laminarization of the flow regime. For the Nusselt number, a dimensionless number generally used to represent surface heat transfer, increasing system rotation is shown to continually decrease $Nu$ on both channel walls. This trend was similarly observed in [28].

A strong correlation between streamwise velocity and temperature fluctuations is also observed. In Figure 9(a), contours of streamwise fluctuating velocity ( $u^{\prime }$ ) are shown in the crossflow plane for simulation case A ( $R{o}_b=0$ ). In the near-wall regions for both channel walls, alternating pairs of structures (sublayer streaks) are observed which contribute to the turbulence generation cycle but similar structures are also seen in the corresponding contours of fluctuating temperature ( ${\theta }^{\prime }$ ) in Figure 9(b). Structures associated with positive $u^{\prime }$ correlate with structures with negative ${\theta }^{\prime }$ and vice versa. This comparison of turbulence and thermal structures is extended to simulation case C ( $R{o}_b=0.5$ ) in Figure 10(a) and Figure 10(b), respectively. Re-laminarization of the suction region is demonstrated to suppress significant fluctuations of both quantities near the suction wall. Near the pressure wall, arrays of coherent $u^{\prime }$ and ${\theta }^{\prime }$ structures are observed. In Figure 10(a), the sublayer streaks are elongated in the wall-normal direction, corresponding to increased $v^{\prime }$ for case C compared to case A (Figure 8(b)) in the pressure region. The thermal structures exhibit similar characteristics and continue to correlate with the sublayer streaks in the presence of rotational effects.

3.3. Higher-Order Statistics

The third-moment of a fluctuating velocity component normalized by the cube of the root-mean-square (r.m.s.) velocity component is known as the skewness:

$S\left({u^{\prime }}_i\right)=\frac{\overline{{u^{\prime }}_i{u^{\prime }}_i{u^{\prime }}_i}}{{\left(\overline{{u^{\prime }}_i{u^{\prime }}_i}\right)}^{3/2}}$ (15)

The skewness quantifies the asymmetry of a variable’s probability density function (PDF) distribution about its mean and measures extreme events occurring in a velocity field. For example, positive skewness indicates large amplitude positive fluctuations have a greater likelihood for occurrence than negative

fluctuations of similar strength. The flatness (F), also known as the kurtosis, represents the fourth-order moment of a fluctuating velocity component normalized by the square of its corresponding Reynolds stress component:

$F\left({u^{\prime }}_i\right)=\frac{\overline{{u^{\prime }}_i{u^{\prime }}_i{u}^{\prime }{u^{\prime }}_i{}_i}}{{\left(\overline{{u^{\prime }}_i{u^{\prime }}_i}\right)}^2}$ (16)

The flatness is the measure of a variable’s peakedness and represents the frequency at which extreme events occur as a deviation from the Gaussian distribution ( $F=3$ ). For example, a high value of flatness ( $F>3$ ) indicates relatively large values at the edges of the PDF distribution and a higher concentration directly around the mean. In [13] , proper grid resolution for higher-order statistics is demonstrated through comparison of higher-order statistics obtained from the DNS results with the experimental measurements of [29] and [30].

In Figure 11, skewness and kurtosis distributions are shown for simulation cases A-D. Asymmetry for the rotational cases is observed in all profiles. In Figure 11(a) and Figure 11(c), rotation is shown to significantly increase values of $S\left(u^{\prime }\right)$ and $S\left(v^{\prime }\right)$ near the suction wall and decrease those values in the pressure region. In Figure 11(e), $S\left(w^{\prime }\right)$ remains zero throughout the channel

regardless of rotation number. In the kurtosis distributions, rotation is not shown to significantly alter the profile values in the pressure region even at high rotation numbers. In the suction region, there is a significant increase of kurtosis values for cases B ( $R{o}_b=0.2$ ) and C ( $R{o}_b=0.2$ ) from the non-rotational case but a decrease for case D ( $R{o}_b=0.9$ ).

3.4. Energy Budgets

It is instructive to examine the various energy budget components in order to discover which force dynamics are primarily affected by rotational forces. The Reynolds stress equation displays the intercomponent energy transfer [31]

$\frac{D\overline{{u^{\prime }}_i{u^{\prime }}_j}}{Dt}=P_{ij}+{\Pi }_{ij}-{\epsilon }_{ij}+C_{ij}+D_{ij}^T$ (17)

with the terms on the right-hand side of Equation (17) representing, respectively, the production ( $P_{ij}$ ), pressure-strain ( ${\Pi }_{ij}$ ), dissipation ( ${\epsilon }_{ij}$ ), Coriolis ( $C_{ij}$ ) and diffusion terms ( $D_{ij}^T$ ).

In the present work, the production, Coriolis and pressure-strain budgets are investigated due to their high contribution level compared to the other budget terms and correspondence with the turbulence generation cycle [1]. These quantities are expressed in tensor form in Equations (18), (19) and (20).

$P_{ij}=-\overline{{u^{\prime }}_i{u^{\prime }}_k}\frac{\partial {\overline{U}}_j}{\partial x_k}-\overline{{u^{\prime }}_j{u^{\prime }}_k}\frac{\partial {\overline{U}}_i}{\partial x_k}$ (18)

$C_{ij}=2{\Omega }_k\left(\overline{{u^{\prime }}_i{u^{\prime }}_m}{ϵ}_{mjk}-{ϵ}_{imk}\overline{{u^{\prime }}_m{u^{\prime }}_j}\right)$ (19)

${\Pi }_{ij}=-\left(\overline{{u^{\prime }}_i\frac{\partial p^{\prime }}{\partial x_j}+{u^{\prime }}_j\frac{\partial p^{\prime }}{\partial x_i}}\right)$ (20)

and the total (summation of all tensor components) distributions are shown in Figures 12(b)-(d), respectively.

In Figure 12(a), the turbulent kinetic energy (k) distributions for simulation cases A-D are shown. The expected suppression of k-amplitudes in the suction region is demonstrated for case B ( $R{o}_b=0.2$ ), but a significant peak near the suction wall continues to persist despite the elimination of the turbulence sustenance cycle; no near-wall peak is observed for cases C ( $R{o}_b=0.5$ ) and D ( $R{o}_b=0.9$ ). In Figure 12(b), rotational forces are demonstrated to suppress the amplitudes of production (P), the primary contributor to the near-wall k peak in Figure 12(a) for case A ( $R{o}_b=0$ ). A small near-wall contribution of P is shown for case B ( $R{o}_b=0.2$ ) but is significantly diminished compared to case A ( $R{o}_b=0$ ); no near-wall P contributions are observed for cases C ( $R{o}_b=0.5$ ) and D ( $R{o}_b=0.9$ ). Hence for case B ( $R{o}_b=0.2$ ), the observed peak in Figure 12(a) is demonstrated to be primarily composed of contributions from the other budget terms with a small supplement from P, demonstrating a fundamental alteration to the dominant processes which contribute towards turbulence production which manifests at low rotation numbers. At higher rotation numbers, re-laminarization mechanisms suppress turbulence production in the near-wall region of the suction side to negligible amounts.

In Figure 12(c), the amplitudes of the pressure-strain distribution are shown to be suppressed for the rotational cases with the exception of the region next to the pressure wall ( $y=2$ ). In this region, the pressure-strain budget amplitudes monotonically increase with rotation number. In Figure 12(d), the Coriolis force is shown to possess significant positive contributions in the suction region for the rotational cases, therefore supplementing the suction side amplitudes seen in Figure 12(a). However, the Coriolis force has significant negative contributions in the pressure region, counteracting the production contribution increases shown in Figure 12(b).

3.5. Coherent Turbulence Structures

In this section, visualizations of various coherent structures are extracted from DNS cases A-C to ascertain the role of these structures to turbulence production amidst contributions from rotational forces. Although the bursting cycle is maintained in the pressure region [11] , rotation-induced structures are introduced to the flow regime. Maps of fluctuating wall-normal vorticity ( ${{\omega }^{\prime }}_y$ ) for full simulation cases A-C at horizontal x-z sections near the suction wall ( $y=0.06$ ) are shown in Figure 13. Remarkable coherence of ${{\omega }^{\prime }}_y$ structures is observed. With ${{\omega }^{\prime }}_y$ dominated by $\partial u^{\prime }/\partial z$ , the elongated structures of ${{\omega }^{\prime }}_y$ shown in Figure 13(a) for case A ( $R{o}_b=0$ ) may be considered the sidewalls of sublayer

streaks as these structures are composed of significant concentrations, and therefore variations, of streamwise fluctuating velocity $u^{\prime }$ [32].

With the introduction of spanwise rotation in case B ( $R{o}_b=0.2$ ), the wall-normal vorticity field in Figure 13(b) demonstrates a significant reduction in the number of ${{\omega }^{\prime }}_y$ structures near the suction wall. The re-laminarization of the suction region also results in characteristics of reverse transition such as dense pockets of vortical structures known as turbulent spots to be observed [33] [34]. The turbulent spots appear intermittently and demonstrate an inclination to the streamwise (x) direction between 15˚ and 30˚, which corresponds with the 30˚ inclination observed in the transitional flow study by [34]. For case C ( $R{o}_b=0.5$ ) in Figure 13(c), the intermittent ${{\omega }^{\prime }}_y$ structures appear even less frequently compared to case B due to further increased re-laminarization of the suction region.

In Figure 14, maps of fluctuating spanwise vorticity ${{\omega }^{\prime }}_z$ are shown at x-y sections for simulation cases A-C. For case A ( $R{o}_b=0$ ) in Figure 14(a), spanwise vorticity is observed to be densely concentrated at both channel walls. For cases B ( $R{o}_b=0.2$ ) and C ( $R{o}_b=0.5$ ) in Figure 14(b) and Figure 14(c), respectively, the strength of ${{\omega }^{\prime }}_z$ is augmented at the pressure wall ( $y=2$ ). Near the suction wall ( $y=0$ ), large-scale vortical structures known as high shear layers are detached from the wall in the rotational cases. The wrinkling of the high shear layers is a fundamental characteristic of transition and this feature is highlighted in a smaller x-y cross-section for case B in Figure 15 [35]. The concave distortion of the high shear layer at $x=9$ is the result of a secondary instability below the shear layer.

System rotation also generates rotation-induced structures in the pressure region in the form of a spanwise array of longitudinal roll cells, known as Taylor-Gortler vortices [11]. Time-averaged spanwise and wall-normal velocity vectors from case C ( $R{o}_b=0.5$ ) are shown for a y-z crossflow section at $x=2\text{π}$ in Figure 16(a). These large-scale structures arise from the secondary flows advected by the Coriolis force and are responsible for flow circulation throughout the pressure region. The corresponding secondary flow streamlines are shown in Figure 16(b) and two full pairs of counter-rotating vortices are observed.

The three-dimensional structure of the Taylor-Gortler vortices is shown across the entire channel in Figure 17 using contours of time-averaged spanwise and wall-normal velocity from case C ( $R{o}_b=0.5$ ). The roll cells appear as streamwise-elongated cylindrical structures which persist throughout the pressure region of the channel. The study by [11] showed the number of roll cell pairs increased with increasing rotation number although the wall-normal length of the circulation region is reduced from progressive re-laminarization. The study by [11] showed the number of roll cell pairs increased with increasing rotation number although the wall-normal length of the circulation region is reduced from progressive re-laminarization.

3.6. Quadrant Analysis

Quadrant analysis of the Reynolds shear stress, which divides $\overline{u^{\prime }{v}^{\prime }}$ into four quadrants according to the signs of fluctuating streamwise ( $u^{\prime }$ ) and wall-normal ( $v^{\prime }$ ) velocity where $v^{\prime }>0$ signifies motion away from the wall, provides important information on contributions to TKE production [36]. In Equation (18), the Reynolds shear stress is shown to have a significant contribution towards the production (P) term. For non-rotating flow, the second and fourth quadrant events dominate the near-wall region of peak TKE production and at the location of peak production ( $y^{+}\approx 12$ ), the contributions from both events are approximately equal [1]. The second quadrant event, ${\left(\overline{u^{\prime }{v}^{\prime }}\right)}_2$ ( $u^{\prime }<0$ and $v^{\prime }>0$ ), contains the motion attributed with ejections of low-speed fluid away from the wall. The fourth quadrant event, ${\left(\overline{u^{\prime }{v}^{\prime }}\right)}_4$ ( $u^{\prime }>0$ and $v^{\prime }<0$ ), contains the motion attributed to an inrush of high-speed fluid into the wall region. In the present work, quadrant analysis is used to compare cases A ( $R{o}_b=0$ ) and B ( $R{o}_b=0.2$ ) in regards to the relationship between $u^{\prime }$ and $v^{\prime }$ whose outer product forms the dominant contribution towards the near-wall P peaks observed in Figure 12(b).

For simulation case A ( $R{o}_b=0$ ), Figure 18 shows maps of significant quadrant events in an x-z cross-section at $y=0.06$ , the location of peak production in Figure 12(b). The bursting event criterion recommended by [37]

$\frac{\left|u^{\prime }{v}^{\prime }\right|}{u_{rms}{v}_{rms}}\ge H$ (21)

was used to define significant quadrant events, where H is the threshold level. This H value was set to 1 in the present work and chosen to be similar to the selected threshold values in the bursting event studies by [3] and [38]. A significantly larger number of powerful second (Q2) and fourth (Q4) quadrant events are observed in comparison to the small number of powerful first (Q1) and third (Q3) quadrant events; the number of significant Q2 and Q4 events is also approximately equal. In Figure 18(b), the streamwise-elongated structures correspond with the sublayer streaks of the turbulence generation cycle.

Figure 19 shows the spatial x-z distributions of significant quadrant events for simulation case B ( $R{o}_b=0.2$ ) at $y=0.14$ , the corresponding location of peak production in Figure 12(b). The number of regions with significant Q2 and Q4 activity is significantly diminished compared to case A ( $R{o}_b=0$ ) in Figure 18. In Figure 19(b) and Figure 19(d), structures resembling “turbulent spots” [34] , A characteristical structure of transitional turbulence, are observed. These structures affirm the undergoing process of reverse transition in the suction region as rotational forces drive the flow regime from turbulence towards re-laminarization.

The spatially-averaged (in x and z) quadrant contributions are also shown in Figure 20(a) and Figure 20(b) as a function of nondimensional time (t) for cases A ( $R{o}_b=0$ ) and B ( $R{o}_b=0.2$ ), respectively. In Figure 20(a), the temporal variations are shown to be random and small for case A ( $R{o}_b=0$ ), indicative of the consistent cycle of turbulence production in the near-wall region of turbulent channel flow. Figure 20(a) also quantifies the significantly higher contributions from Q2 and Q4 to turbulence production than Q1 and Q3. In Figure 20(b), the time history for case B ( $R{o}_b=0.2$ ) displays significant changes from Figure 20(a). For all four quadrants, large, quasi-periodic temporal variations of $\overline{u^{\prime }{v}^{\prime }}$ values are observed over long time scales; the motions are also simultaneous. Although the amplitudes are heavily suppressed compared to case A ( $R{o}_b=0$ ), significant bursts of energetical activity are demonstrated at peak values. The differences between the contributions from Q2 and Q4 to those from Q1 and Q3 are also significantly diminished. At peak amplitude, Q4 events are observed to be the largest contributor to turbulence production.

To analyze the phase differences between the streamwise ( $u^{\prime }$ ) and wall-normal ( $v^{\prime }$ ) fluctuating velocities for cases A ( $R{o}_b=0$ ) and B ( $R{o}_b=0.2$ ), the temporal distributions of $\overline{u^{\prime }{v}^{\prime }}$ are decomposed using Fourier transforms. The separate $u^{\prime }$ and $v^{\prime }$ signals in Figure 20(a) and Figure 20(b) are transformed into spectral signals and a dot product is performed between the two signals. The phase component is extracted from the dot product and visualized in the polar plots in Figure 21. 0˚ and 180˚ refer to the signals being completely in-phase and out-of-phase, respectively. In Figure 21, a massive shift of the phase contours between the simulation cases A ( $R{o}_b=0$ ) and B ( $R{o}_b=0.2$ ) is observed for all quadrants. Hence a remarkable alignment between streamwise

and wall-normal fluctuating velocity, which manifests as the two quantities being almost completely in-phase or out-of-phase for all four quadrants, is demonstrated in the suction region for case B ( $R{o}_b=0.2$ ). This alignment is also corroborated by the quasi-periodic motions in Figure 20(b) which occur concurrently for all four quadrants.

4. Rotational Turbulence: Higher-Reynolds Number Effects

In engineering applications, turbulent flows often encounter higher-Reynolds number effects. Hence it is prudent to investigate how turbulence production and the corresponding turbulence structures are different in higher-Reynolds number flows compared to lower-Reynolds number flows. The higher-Reynolds number simulation case E ( $R{e}_{\tau }=406$ ) is examined and compared to the lower-Reynolds number simulation case B ( $R{e}_{\tau }=197$ ).

In Figure 22(a) and Figure 22(b), the mean velocity and temperature distributions are shown for simulation cases B ( $R{e}_{\tau }=197$ ) and E ( $R{e}_{\tau }=406$ ), respectively. In Figure 22(a), a higher-Reynolds number is shown to significantly decrease the amplitudes of the mean velocity distribution although other distribution characteristics such as the shape and slope in the pressure region are preserved. The decreased amplitudes were expected due to the significant increase of $R{e}_{\tau }$ (and subsequently $u_{\tau }$ ) in simulation case E ( $R{e}_{\tau }=406$ ), which affected the scaling of the distribution. In Figure 22(b), the mean temperature distribution for case E ( $R{e}_{\tau }=406$ ) is shown to be significantly less asymmetric than the distribution for case B ( $R{e}_{\tau }=197$ ) and resembles the distribution for the no-rotation case A shown in Figure 7.

In Figure 23, the turbulent kinetic energy and energy budget (production, pressure-strain and Coriolis) distributions for case B ( $R{e}_{\tau }=197$ ) and case E ( $R{e}_{\tau }=406$ ) are displayed. In Figure 23(a), case E demonstrates characteristics of higher-Reynolds number flows such as the shifts of both near-wall peaks towards the channel walls and the increasing amplitude of the pressure region peak in the case E ( $R{e}_{\tau }=406$ ) distribution compared to the case B ( $R{e}_{\tau }=197$ ) profile. In Figure 23(b), the near-wall peak in the pressure region for the case E ( $R{e}_{\tau }=406$ ) distribution also increases in amplitude and shifts towards the channel wall compared to the case B ( $R{e}_{\tau }=197$ ) profile.

In Figure 23(c), the amplitudes of the pressure-strain distribution for case E are shown to be very similar to that of case B with exception of the region near the pressure wall. Similar to the effects of high-rotation numbers, a higher-Reynolds number increases the amplitude of the pressure-strain budget near the pressure wall significantly. In Figure 23(d), the near-wall peaks of the Coriolis energy budget distribution are also shown to shift towards the channel walls with a higher-Reynolds number, similar to the production budget, although there is less amplitude increase compared to the other energy budget distributions.

It is also imperative to look at how turbulence structures in both the pressure and suction regions are altered by higher-Reynolds number effects. In Figure 24, x-z planar contours of fluctuating streamwise velocity are shown for simulation cases B ( $R{e}_{\tau }=197$ ) and E ( $R{e}_{\tau }=406$ ) in the region near the suction wall ( $y=0.05$ ). In Figure 24(b), the number of powerful fluctuations comprising the turbulent spots have increased significantly for simulation case E ( $R{e}_{\tau }=406$ ) and the shape of the overall structure is much more defined.

In Figure 25, contours of fluctuating streamwise velocity are shown for simulation cases B ( $R{e}_{\tau }=197$ ) and E ( $R{e}_{\tau }=406$ ) in the region near the pressure

wall ( $y=1.95$ ). Similar to Figure 24(b), Figure 25(b) demonstrate that higher-Reynolds number effects cause turbulence structures to become smaller and more numerous. Although the overall characterization of the flow field is not changed, the specific characteristics of the turbulence structures are different for higher-Reynolds number flows, specifically the diminution and number magnification of the structures.

In Figure 26, instantaneous x-y maps of fluctuating spanwise vorticity ( ${{\omega }^{\prime }}_z$ ) are shown at $z=\text{π}$ for simulation cases B ( $R{e}_{\tau }=197$ ) and E ( $R{e}_{\tau }=406$ ). In the regions near both channel walls, the overall turbulence structure is preserved: although the number of structures has increased, elongated slanted structures and clusters of high vorticity concentrations are seen near the suction and pressure walls, respectively. In the center of the channel however, vorticity has clearly increased in the case of simulation case E ( $R{e}_{\tau }=406$ ) perhaps resulting from the alterations to the flow dynamics near both channel walls.

5. Conclusions

In summary, the theoretical model predictions proposed by [4] and [5] for the coherent structures of sublayer streaks and accompanying vortical structures were validated through comparison with the present DNS results. For the [4] model, the appearance and evolution of the sublayer streaks obtained using the VISA method from the DNS data corresponded very well with the modeled VITA results. For the [5] model, snapshots of a sublayer streak and its surrounding vortical structures were obtained from both the periodic and spatial DNS simulations which matched the pictorial representation of the theoretical model. In agreement with the results of [6] for a lower-Reynolds number, non-propagating modes were found to possess the highest energetical content and a normal speed locus was generated in order to visualize the interaction and movement of propagating plane waves. With most of the turbulent kinetic energy contained to a small range of modes, this validation study bodes well for future flow control work. Since propagating waves have been established in both periodic and spatial models despite the lack of a streamwise wavenumber in the spatial channel, flow control designs in future work should aim to inhibit these particular modes to reduce turbulent kinetic energy and drag.

An examination of the coherent structures which contribute to turbulence production for various rotation rates is of great importance to understand the effects of rotational forces on turbulence. In the present work, a comprehensive investigation involving low and high-order statistics, coherent structures, budget analysis and quadrant analysis was conducted. Quadrant analysis was used in the suction region to elucidate a phase shift in the relationship between the streamwise and wall-normal fluctuating velocities, a significant finding as these components form the principal contribution to turbulence production and are involved in the ejection and sweep events which dominate the near-wall turbulence generation cycle. This investigation of the contributions of coherent structures to turbulence production was extended for a higher-Reynolds number rotational simulation in which notable characteristical changes to various turbulence structures were identified.

Acknowledgements

This work was supported in part by the Air Force Office of Scientific Research under grant number No. FA9550-15-0495.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Hsieh, A.S. and Biringen, S. (2019) Effects of Rotation on Turbulence Production. Journal of Applied Mathematics and Physics, 7, 298-330. https://doi.org/10.4236/jamp.2019.72024

References