*Article* **Numerical Investigation of the Turbulent Wake-Boundary Interaction in a Translational Cascade of Airfoils and Flat Plate**

#### **Xiaodong Ruan 1,\*, Xu Zhang 1, Pengfei Wang 2, Jiaming Wang <sup>1</sup> and Zhongbin Xu <sup>3</sup>**


Received: 24 July 2020; Accepted: 27 August 2020; Published: 31 August 2020

**Abstract:** Rotor stator interaction (RSI) is an important phenomenon influencing performances in the pump, turbine, and compressor. In this paper, the correlation-based transition model is used to study the RSI phenomenon between a translational cascade of airfoils and a flat plat. A comparison was made between computational results and experimental results. The computational boundary layer velocity is in reasonable agreement with the experimental velocity. The thickness of boundary layer decreases as the RSI frequency increases and it increases as the fluid flows downstream. The spectral plots of velocity fluctuations at leading edge *x*/*c* = 2 under RSI partial flow condition *f* = 20 Hz and *f* = 30 Hz are dominated by a narrowband component. RSI frequency mainly affects the turbulence intensity in the freestream region. However, it has little influence on the turbulence intensity of boundary layer near the wall. A secondary vortex is induced by the wake–boundary layer interaction and it leads to the formation of a thickened laminar boundary layer. The negative-vorticity wake also facilitates the formation of a thickened boundary layer while the positive-vorticity wake has a similar effect, like a calmed region which makes the boundary layer thinner.

**Keywords:** rotor stator interaction; boundary layer; secondary vortex

#### **1. Introduction**

Rotor stator interaction (RSI) is a common phenomenon in pumps and turbomachines, especially in the case of multistage compressors. Many works have been carried out to investigate the influence of RSI on energy calculation [1], pressure fluctuation [2], internal fluid flow [3], noise [4], and machine vibration for turbomachines [5]. The unsteady interaction has an important influence on the hydrodynamic performance of turbomachines in many aspects such as blade loading, stage efficiency, heat transfer and noise generation [6]. On the one hand, the wake generated by the trailing edge of moving blades is periodically cut by the stationary blades, and pressure fluctuation is created, pressure fluctuation on the blade may lead to fatigue failure due to periodic load variation. On the other hand, the unsteady interaction changes the local pressure distribution, the inverse pressure gradient, and viscous effects may lead to the separation of the boundary layer where a back-flow region is generated, which makes the energy loss increase and the efficiency decrease. Therefore, in order to improve the efficiency and control the blade load, it makes sense to carry out an investigation of the mechanism by which RSI influences the boundary layer and how the boundary layer develops with the interaction.

Rotor stator interaction might be divided into two different interaction mechanisms, namely potential flow interaction and wake interaction. In general, wake interaction ranges farther downstream, but if the impellers and guide vanes are closely spaced, both mechanisms will occur simultaneously [7]. However, if the impellers and guide vanes are closely spaced, the fluid-induced force will be a very high level which will quickly lead to blade failure [8], so a close space design between rotor and stator is usually avoided. As a result, the primary mechanism of rotor stator interaction is wake interaction. Wakes generated by upstream impeller trailing edge convey low-speed fluids, and wakes are transported to the downstream where they interact with the boundary layer of the downstream guide vane in form of travelling wave [9]. This interaction will directly affect the development of the boundary layer of downstream blades. Although many investigations on boundary layer profiles have been made both theoretically and experimentally [10–12], the nature of the mechanism of wake–boundary interaction and its impact on the unsteadiness of boundary layer, including velocity profiles, turbulence intensity distribution, shear stress in different phases, and boundary layer structure, is not well understood.

In recent decades, there have been many theories on rotor stator interaction [13–15], but two procedures on rotor stator interaction are widely accepted. One was proposed by Kubota [16] and Franke [17] based on the diametrical mode theory. In this procedure, the periodic flow induced by the wake effects could be represented by the Fourier series, and rotor stator interaction can be calculated by the Fourier series superposition, whereby the diametrical mode number k could indicate the number of high or low-pressure fluid regions of the specific frequency component in the circumferential direction. This procedure can indicate the relative amplitude of each harmonic frequency components with different rotor/stator blade number [18]. The other procedure was proposed by Rodriguez [19] based on the sequence of interactions. This approach enables us to regard the vibration as the result of a modulation in the amplitudes of the interactions and rather than a consequence of an exciting diametrical mode number as the first theory proposed. Furthermore, this procedure could help to determine the origin of the harmonics, thus providing guidance for obtaining lower amplitude at specific harmonics or to avoid the generation of harmonics [20]. However, due to the complexity of the internal flow in turbomachines, the interaction is not only a consequence of potential-wake flow interaction but also a consequence of wake–boundary interaction. These two procedures neglect the unsteady wake–boundary interaction and naturally are both lack of engineering cases to corroborate its approach.

The combination of experiment and numerical calculation is an important method to study unsteady flow field [21–23]. To precisely compute the boundary layer flow, the choice of numerical methodology is very important. Sanders [24] made a comparison between computational fluid dynamics (CFD) predictions using Reynolds-averaged Navier–Stokes equations (RANS) and experimental results of unsteady wake effects on a highly loaded low-pressure turbine blade. And the results showed that the RANS method, especially the shear stress transport (SST) model, had a good prediction on the unsteady L1A blade flow field. Due to the shortcomings of the RANS method, the LES method has also recently been used to study the interaction of wake boundary layers if more accurate results are desired. Sarkar [25] used large eddy simulation (LES) to investigate the influence of wake structures on the boundary layer evolution at the suction side of a high-lift low-pressure turbine blade, and his research illustrated that in addition to the kinematics of wakes, the high- pressure oscillation and the curl of the separated shear layer along nearly half of the suction surface were determined by the length of the convective wake. Direct numerical simulation (DNS) is acknowledged to be the most perfect numerical method in accordance with experiments, because it is a direct solution to the N-S equation. However, due to its high spatial and temporal resolution, it is computationally intensive, time-consuming, and highly dependent on computer memory. Currently using DNS to simulation is still in the exploratory stage. Wissink [26] deliberated the effect of wake turbulence intensity on transition in a compressor cascade by using DNS, and the results indicated that the separation was intermittently inhibited due to the intense wake reaction as the periodically passing wakes tried to induce turbulent spots upstream of separation point. Particle Image Technology (PIV) has high measurement accuracy due to non-contact and small interference, so it has become a hot topic in recent fluid mechanics measurement research. Jia [27] conducted a PIV measurements of rotor boundary layer development with the effects of wakes generated by inlet guide vane (IGV), the research found that for the stronger wake, the interaction between the unsteady wake and the separated boundary layer on the suction surface of the airfoil was dominated by the high turbulence of wake and negative jet behavior of the wake. As the wake strength increased, the separation was almost totally inhibited. The different boundary layer states also affected the development of disturbances conversely.

From the above previous studies, we know that it is essential to study the wake generation and its effects on the boundary layer development. Before attempting the extremely complex three-dimensional and rotating flow on a turbomachine, the present study based on a simpler setup can be used as a preliminary orientation for the turbomachinery blade. This setup is abstracted from the blade row interaction problem of turbomachinery, the advantages are the geometrical complexity and pressure gradients are removed and it is convenient for calculations. The present study used the experimental data in the literature [28], which has convenient measurements and accurate data in the boundary layer. Comparisons were made between the experimental data and the CFD results to ensure the numerical simulation is reliable. The main objective of this work is to investigate the mechanism of the wake effects on boundary layer development, especially the RSI frequency effects on the boundary layer and its time–space variation.

#### **2. Model and Methods**

#### *2.1. Computational Domain and Mesh Generation*

The experimental model of turbomachine stage in the literature [28] was represented by a low speed wind tunnel with an unsteady wake generating rig. Although this model is a simple test rig, the model is quite representative of general realistic. When the turbomachinery impellers are expanded along the circumferential direction, the impeller movements can be considered as vertical movements, as shown in Figure 1. The stator blade was replaced by a flat plate in the test section, and the rotor blade was represented by a translational cascade of airfoils (NACA0024). The moving airfoils generate periodic wake disturbances and spread to downstream, so the wake reacts with the flat plate when it passes through the test section.

**Figure 1.** A global view of the computational domain, (**left**) 3D schematic of geometry (spanwise× 10), (**right**) plane view with boundary conditions.

The entire fluid domain consists of the moving cascade of airfoils and the stationary flat plate, as shown in Figure 1. 3D calculation was carried out in this study. In order to find a balance between the quality of simulation and computational cost, a geometric optimized 3D mesh with spanwise-reduced domain was set up to perform these numerical simulations. As a guide [29], the spanwise depth of the calculated domain was set as *LS* = 0.16*c*. Block-structured with O-grid around the airfoil and flat plate was adopted for the mesh designs. The γ-*Re*θ*<sup>t</sup>* transition SST model was used in the transient numerical simulation. Therefore, the grid cells near the wall were encrypted with a cell stretching

factor equal to 1.2 from the walls to ensure the value y+ less than 1 at the first node from the solid walls, as is required by the SST turbulence model [30,31]. The maximum number of layers is 25, and the wall y+ distribution is shown in Figure 2. It can be observed that the average value of y+ is 0.2 for airfoil and 0.4 for plate. A snapshot of mesh encryption is shown in Figure 3, it shows that the mesh of leading edge and trailing edge of airfoil and the boundary layer are encrypted. The quality of mesh is larger than 0.7 and 0.8 for airfoil domain and plate domain correspondingly, thus a good mesh quality is guaranteed. CFD simulations were performed with commercial software ANSYS Fluent 16.0 [32] which was run on an HP Z840 high-performance workstation using a total of 32 cores with 64GB memory for the unsteady calculations. These simulations required approximately 500 core-hours for all unsteady simulations on the Workstation with 3.1GHz Intel Xeon(R) CPU E5-2687w v3 (2 CPU).

**Figure 2.** Wall *y*+ distribution of airfoil blade (**left**) and flat plate (**right**).

**Figure 3.** Snapshot of boundary layer mesh encryption (**left**: airfoil blade, **right**: flat plate).

#### *2.2. Theoretical Methods*

The viscous sublayer flow can be directly solved by the SST k-ω turbulence model, The SST k-ω model combined with turbulent shear stress transmission and is quite accurate in the prediction of flow separation near the wall under an inverse pressure gradient [33], however, this model has some defects in the prediction of the transition phenomenon [34]. Based on an empirical formula, the γ-*Re*θ*<sup>t</sup>* transition SST model proposed by Langtry et al. [35] is more accurate in the prediction of the location of the transition point. This model takes advantage of the local parameters provided by a RANS solver rather than employing an additional boundary layer solver. The γ-*Re*θ*<sup>t</sup>* transition SST model [36] is coupled with another two transport equations based on the SST k-ω model, and it does not influence the flow filed of fully turbulent region. So, theoretically, the transition SST model has higher prediction accuracy for the velocity distribution in the near-wall region. Many studies [37–39] have confirmed the accuracy of this model.

The transport equations of turbulent kinetic energy and specific dissipation rate are:

$$\frac{\partial(\rho k)}{\partial t} + \frac{\partial(\rho u\_l k)}{\partial \mathbf{x}\_l} = P\_k + D\_k + \frac{\partial}{\partial \mathbf{x}\_l} \Big[ (\mu + \sigma\_k \mu\_l) \frac{\partial \mathbf{k}}{\partial \mathbf{x}\_l} \Big] \tag{1}$$

$$\frac{\partial(\rho\omega)}{\partial t} + \frac{\partial(\rho u\_l \omega)}{\partial \mathbf{x}\_l} = P\_{\omega} + D\_{\omega} + \frac{\partial}{\partial \mathbf{x}\_l} \Big[ (\mu + \sigma\_{\omega} \mu\_l) \frac{\partial \omega}{\partial \mathbf{x}\_l} \Big] \tag{2}$$

where *Pk* and *P*<sup>ω</sup> are the production terms, *Dk* and *D*<sup>ω</sup> are the destruction terms, *ui* is the flow velocity, ρ is the density of fluid, μ is the dynamic viscosity of the fluid, σ*<sup>k</sup>* and σω are the model coefficients, *k* is the turbulent kinetic energy, and ω is the specific dissipation rate.

The γ-*Re*θ*<sup>t</sup>* transition model can be written as:

$$\frac{\partial(\rho \gamma)}{\partial t} + \frac{\partial(\rho u\_l \gamma)}{\partial \mathbf{x}\_l} = P\_{\gamma} + E\_{\gamma} + \frac{\partial}{\partial \mathbf{x}\_l} \left[ (\mu + \frac{\mu\_l}{\sigma\_f}) \frac{\partial \gamma}{\partial \mathbf{x}\_l} \right] \tag{3}$$

$$\frac{\partial(\rho \vec{\text{Re}}\_{\partial t})}{\partial t} + \frac{\partial(\rho u\_i \vec{\text{Re}}\_{\partial t})}{\partial \mathbf{x}\_i} = P\_{\partial t} + \frac{\partial}{\partial \mathbf{x}\_i} \left[ \sigma\_{\partial t} (\mu + \mu\_t) \frac{\partial \mathbf{\hat{R}} \mathbf{c}\_{\partial t}}{\partial \mathbf{x}\_i} \right] \tag{4}$$

where transport Equation (3) is about intermittency γ and Equation (4) is about the local transition onset momentum thickness Reynolds number *Re*˜ <sup>θ</sup>*t*. *P*<sup>γ</sup> and *E*<sup>γ</sup> are the production term and the destruction term of the transport equation for intermittency γ, respectively. *P*θ*<sup>t</sup>* is the source item of the transport equation for the local transition onset momentum thickness Reynolds Number *Re*˜ <sup>θ</sup>*t*, μ is the dynamic viscosity coefficient, μ*<sup>t</sup>* is the eddy viscosity, σ*<sup>f</sup>* and σθ*<sup>t</sup>* are the constant coefficients.

#### *2.3. Numerical Scheme*

In order to get the accurate boundary layer velocity distribution, a reasonable boundary condition setting is very important. The computational domain is based on the parameters in the experiment. As the previous section states, the model is a representative of general realistic turbomachinery. The flow conditions were chosen based on the typical demand conditions. The "velocity inlet" was used in the tunnel inlet with velocity = 3 m/s according to the experiments. The turbulence parameters were specified according to the turbulence intensity and hydraulic diameter of the inlet. The "pressure outlet" was set in the outlet of the wind tunnel with pressure = 0 Pa, and the turbulence parameters were specified according to turbulence intensity and viscosity ratio. Besides, a smooth no-slip wall condition was specified for the solid wall boundaries of the airfoil surface, the flat plate wall, and the wind tunnel walls. The upper and lower boundaries were specified to be "periodic" boundary conditions. The front face and back face in the spanwise direction were set as symmetrical. The rationality of the periodic boundary is explained as follows: the rotor part itself is periodic, and the upper and lower areas of the stator side are fully expanded to match the rotor area. The extension of the stator area does not affect the stator flow field. The periodic boundary is added to the virtual extension to match the rotor area for calculation. The rotor-stator interface was set to periodic repeats; the detailed information of boundary condition settings is shown in Figure 1.

The RSI frequency *f* is defined as the ratio of the vertical velocity of moving airfoil *Ur* to the distance S between two adjacent airfoils. The undisturbed free-stream velocity, *U*0, is set to 3 m/s for all simulation, and the vertical velocity of the rotor, *Ur*, ranged between 2 m/s and 4 m/s, so the RSI frequency *f* ranged from 20~40 Hz. The Reynolds number is 2.37 <sup>×</sup> 104 based on the half height of wind tunnel. Different rotor velocity corresponds to different incidence angles. The incidence angle is the largest at *f* = 20 Hz, so it is under partial flow condition, while the incidence angle at *f* = 40 Hz is the smallest, so in this case it is under full flow condition. The turbulence intensity in the wind tunnel is 0.7% in experiments, so a turbulence intensity of 0.7% is also applied at the velocity inlet during the calculation. The steady calculation is performed before the transient simulation was carried out, and the results of steady calculation are taken as the initial condition of the transient simulation. The unsteady flow is calculated by a second-order upwind spatial discretization with a second-order implicit velocity formulation. And a pressure-based solver is applied, the SIMPLE pressure-velocity coupling method is used for the numerical simulations. The convergence criteria of 1 <sup>×</sup> 10−<sup>5</sup> is used for the residual of continuity equations, momentum equations, and turbulence variables. The time step is set to 1 <sup>×</sup> 10−<sup>4</sup> <sup>s</sup> so that there are 250, 333, and 500 time steps in a blade-passing cycle for *Ur* = 4 m/s, 3 m/s and 2 m/s correspondingly. To avoid the influence of unstable effects at the beginning, time-averaging process and analysis are made after the computation is conducted for at least 5 cycles. The velocity is recorded at points that are at a distance of *x*/*c* = 2, 6, 10, 14 from the leading edge of the flat plate.

#### **3. Validation of the Numerical Model**

Verification of mesh dependency was carried out to ensure the mesh has enough resolution. 3 different mesh resolutions were tested including coarse mesh case 3 (2.4 <sup>×</sup> 10<sup>5</sup> elements), standard mesh case 2 (3.6 <sup>×</sup> <sup>10</sup><sup>5</sup> elements) and fine mesh case 1 (4.8 <sup>×</sup> <sup>10</sup><sup>5</sup> elements), and the results are compared with experimental results as shown in Figure 4. Ultimately, the fine mesh case 1 was used for current computations. This proved that the numerical methodology used in this study is feasible.

**Figure 4.** Boundary layer velocity profiles calculated by three different mesh sizes (steady condition).

In order to validate the accuracy of the numerical method, different turbulence models including Re- normalization group (RNG) k-ε turbulence model which has a quite high accuracy in predicting vortex flow, SST k-ω turbulence model which has the ability to predict the near wall flow, and transition SST model were chosen to predict the boundary layer velocity profiles at rotor-stator interaction frequency *f* = 30 Hz. The experimental velocity in the literature [28] was measured by a hot wire probe. *U*<sup>0</sup> is the velocity in the freestream, the thickness of boundary layer which is donated as δ<sup>99</sup> is the distance from the wall of the plate to the region where *U* = 0.99*U*0. For transient flow field calculations, the arithmetic mean of the fluid velocity at different time steps can be expressed as

$$\mathcal{U}\_{mcm} = \frac{\sum\_{i=1}^{N} u\_i}{N} \tag{5}$$

where *Umean* is time-averaged mean velocity, *ui* is the transient velocity, and *N* is the calculated time steps.

The time-averaged predicting results in Figure 5a shows the mean velocity distribution in the boundary layer. Comparing the calculated results, the RNG k-ε turbulence model overestimates the velocity in the free stream flow by approximately 10%, the flow continuity was checked at each streamwise position, and the error was less than 1%. The reason why the RNG k-ε turbulence model overestimates the velocity in the free-flow region by 10% is that the RNG k-ε turbulence model underestimated the velocity in other flow regions. The results showed that the RNG k-ε turbulence model severely underestimated the velocity in the large velocity gradient region. Results calculated by the SST k-ω turbulence model and the transition SST model showed a more precise velocity prediction which has no more overestimation in the freestream region and less error within the large velocity gradient region. Figure 5b shows the boundary layer thickness distribution at different plate streamwise locations. The numerical results of the transition SST model coincide with the experimental results better, the boundary layer thickness ranges from 0.68 *y*/*Ls* to 1.31 *y*/*Ls* calculated by transition SST turbulence model at RSI frequency *f* = 30 Hz. As the fluid is transported downstream, the thickness of boundary layer increases gradually. Figure 6 is the time history contrast of velocity calculated by different turbulence model and from experimental results at *x*/*c* = 2, *f* = 20 Hz. All of the computational results can capture the pulsation in velocity. However, it is evident that the RNG k-ε turbulence model excessively overestimated the velocity. The SST k-ω turbulence model and the transition SST model show a similar prediction in the transient velocity. Ultimately, comparing the three predicting results with the experimental results, it can be concluded that the transition SST turbulence model is more accurate, and thus it is reasonable to use this turbulence model in this paper.

(**a**) Time-averaged mean velocity distribution of boundary layer at different streamwise positions.

(**b**) Boundary layer thickness ߜଽଽ at different streamwise positions.

**Figure 5.** Comparison of boundary layer distribution between results calculated by different turbulence model and experimental results.

**Figure 6.** Time history of instantaneous velocity calculated by different turbulence model and experimental results at *x*/*c* = 2, *y*/*Ls* = 0.25, *f* = 20 Hz.

#### **4. Results and Discussions**

#### *4.1. Velocity Profiles in the Boundary Layers*

The statistical time-averaging process is made with results obtained from multiple periodic cycles. Figure 7a is the comparison of time-averaged mean velocity distributions calculated by different rotor-stator interaction frequency. It can be clearly observed that the velocity profiles of boundary layer are similar in the streamwise direction on the plate. Firstly, the velocity increases rapidly near the wall, and then the velocity growth rate decreases. When the wall distance *y* is greater than the thickness of the boundary, the velocity maintains the maximum and remains unchanged. For the three simulations with different RSI frequency, the velocity gradient in the boundary layer decreases from the leading edge of the plate at *x*/*c* = 2 to the trailing edge of plate at *x*/*c* = 14. This is mainly because the boundary layer of the flat plate has not been fully developed, the thickness of the boundary layer is still increasing as shown in Figure 7b, so its velocity gradient decreases with the increase of the streamwise distance. Moreover, the RSI frequency range investigated in this study shows that the significant fluctuation of time-average mean velocity distribution is a function of RSI frequency in unsteady flow. It can be observed that the velocity gradient in the boundary layer at *f* = 40 Hz is more significant than the velocity gradient at *f* = 20 Hz. This is mainly due to the decrease of incidence angle at *f* = 40 Hz, which results in a smaller scale wake vortex at the trailing edge of airfoil, and the thickness increasement of boundary layer caused by wake–boundary interaction decreases as shown in Figure 7b. The detailed boundary layer thickening mechanism will be discussed in later section. Figure 7b also shows the thickness of the laminar boundary layer and turbulent boundary layer of the flat plate according to Equations (6) and (7) [40]. It can be observed that the thickness of boundary layer is mostly between the thickness of laminar boundary layers and turbulent boundary layers due to the influence of the rotor wakes.

$$\delta\_1(\mathbf{x}) = 5 \sqrt{\frac{\nu \mathbf{x}}{\mathcal{U}\_{\infty}}} \tag{6}$$

$$\delta\_2(\mathbf{x}) = 0.37 \left(\frac{\nu}{l L\_{\infty}}\right)^{1/5} \mathbf{x}^{4/5} \tag{7}$$

(**a**) Time-averaged mean velocity distribution of boundary layer for different RSI frequency at different streamwise position.

(**b**) Boundary layer thickness ߜଽଽ varies with the streamwise position for different RSI frequency.

**Figure 7.** Comparisons of boundary layer calculated by different rotor-stator interaction frequency.

Figure 8 is the transient velocity fluctuation in the streamwise direction and wall-normal direction. The spanwise fluctuations are not shown as the amplitude of spanwise fluctuations are far smaller than for the streamwise and wall-normal fluctuations. The reason why the spanwise fluctuation is so small is that the movement of airfoil blade causes the fluid to pulsate in the streamwise direction and wall-normal direction, and there is no source of disturbance in the spanwise direction except for turbulent random pulsation, which is far less than the periodic fluctuations. The transient results show that as the wake impinges on the boundary layer of flat plate, the velocity in the boundary layer shows obvious periodic fluctuations. The fluctuation amplitudes of *ux* and *uy* at the leading edge *x*/*c* = 2 are far greater than that at trailing edge *x*/*c* = 14. The amplitude of fluctuation velocity decreases with the increase of streamwise direction, this is mainly due to the fluctuation energy is viscously dissipated when the fluctuation spreads downstream. And it is evident that the fluctuation amplitude of *ux* is

much higher than that of *uy* . For the flat plate, the velocity fluctuation amplitude in the streamwise direction and the wall-normal direction depends on two aspects. On the one hand, the installation angle of the front airfoil blade has an effect on the velocity fluctuation amplitude. The smaller the angle between the airfoil chord and the plate, the higher the fluctuation amplitude in wall-normal direction generated by the vertical movement of the airfoil blade. On the other hand, the velocity fluctuation depends on the inflow attack angle of the plate. The smaller the inflow attack angle, the lower the absolute amplitude of the velocity fluctuation in the wall-normal direction. Obviously, in the case of present study, it is mainly due to the small inflow attack angle of the plate, which causes the amplitude of the velocity fluctuation in wall-normal direction to be much lower than that of the streamwise direction. For the three RSI frequencies, the amplitude of fluctuation velocity at *f* = 20 Hz is the highest while the amplitude of fluctuation velocity at *f* = 40 Hz is the lowest. This is because the inflow attack angle of airfoil blade at *f* = 20 Hz is larger than that of *f* = 40 Hz, when the inflow attack angle of the airfoil blade is relatively large, the flow separation is liable to occur at the boundary layers, and large scale of vortex shedding is apt to form, which will lead to a high amplitude of velocity fluctuation.

(**a**) Instantaneous velocity fluctuation at *f* = 20 Hz, left: ݑ௫ <sup>ᇱ</sup> , right: ݑ௬ ᇱ .

(**b**) Instantaneous velocity fluctuation at *f* = 30 Hz, left: ݑ௫ <sup>ᇱ</sup> , right: ݑ௬ ᇱ .

#### **Figure 8.** *Cont.*

(**c**) Instantaneous velocity fluctuation at *f* = 40 Hz, left: ݑ௫ <sup>ᇱ</sup> , right: ݑ௬ ᇱ .

**Figure 8.** Instantaneous velocity fluctuation at *y*/*Ls* = 0.25 over 8 passing cycle: (**a**) *f* = 20 Hz (**b**) *f* = 30 Hz (**c**) *f* = 40 Hz.

Due to the isotropic properties of turbulent flow in free flow, the fluctuation amplitude of *ux* , *uy* , *uz* should be approximately equal. However, because the specific motion of the rotor airfoil blades imposes different degrees of perturbation on the fluid in different directions, the wake conveys this perturbation into the boundary layer of the flat plate, which makes the fluctuation amplitude in the boundary layer of the flat plate vary greatly. Further perceptions into the growth of the fluctuations in the flow streamwise direction is revealed by spectral analysis. Figure 9 shows the spectral plots of velocity fluctuations at the leading edge *x*/*c* = 2 and the trailing edge *x*/*c* = 14. It can be observed that the fluctuation amplitudes of *ux* , *uy* at *x*/*c* = 2 with RSI frequency *f* = 20 Hz are dominated by a narrowband component centered around 20 Hz. However, the fluctuation amplitude of *uz* does not exist a peak at a certain frequency. The structure of the fluctuation is dominated by an increasingly broadband spectrum. This is because the movement of rotor airfoil induces a relatively large disturbance of the fluid in the streamwise and wall-normal directions, but the disturbance in the spanwise direction is minimal. At the trailing edge *x*/*c* = 14, the frequency component with high peak amplitude disappears in the three directions, indicating that the high peak frequency component is dissipated. The structure of fluctuations becomes increasingly chaotic, resulting in the disappearance of narrowband component. The velocity fluctuation is dominated by an increasingly broadband spectrum. At *f* = 30 Hz, the spectral plots of velocity fluctuation are similar to the case of *f* = 20 Hz. When RSI frequency increases to *f* = 40 Hz, the velocity fluctuations at the leading edge *x*/*c* = 2 and trailing edge *x*/*c* = 14 are both dominated by an increasingly broadband spectrum. This is mainly induced by the decrease of the frequency component at *f* = 40 Hz, so there is no high-peak component in the entire fluctuation spectral plots. As the frequency of RSI increases, the difference of velocity fluctuation amplitude with *ux* , *uy* , and *uz* at the leading edge of the plate *x*/*c* = 2 also becomes larger. However, the amplitude of fluctuation at the trailing edge *x*/*c* = 14 has little to do with the variation of RSI frequency.

#### *4.2. Turbulence Intensity in the Boundary Layer*

Many experimental results have confirmed that the velocity of wake is quite low compared with the surrounding fluid, the wake is a negative jet with high turbulence [41]. When the convective wake with increased turbulence intensity permeates the laminar boundary layer on the plate surface, it disturbs the velocity in the boundary layer as the wake conveys high turbulence kinetic energy. The turbulence intensity in the plate boundary layer is investigated and the results are shown in Figure 10. Turbulence intensity profiles are similar for the three different frequencies. The turbulence intensity in the boundary layer increases first and then decreases. A peak value exists in the boundary layer. In the viscous sublayer, the viscosity stress is dominated for fluid motion. The viscous damping has a strong inhibitory effect on the tangential velocity fluctuation and normal velocity fluctuation. As the distance from the wall increased, the turbulent Reynolds stress increases rapidly, so the turbulence intensity increases rapidly in the buffer layer. This indicates that the buffer layer is the primary turbulence energy producer; it provides the needed energy for the maintenance of wall turbulence. The outer region is fully developed turbulent flow. The turbulence intensity dissipated as the wall distance increase, and finally, it maintains a relatively low level in the outer region. The turbulence intensity profile is quite different from others for *f* = 20 Hz at *x*/*c* = 2. The profile is much fuller than others. At *f* = 20 Hz, the incidence angle of fluid around the airfoil is the largest, so flow separation occurs in the forward suction surface of airfoil, thus larger scale wake vortexes which contains high turbulence level are generated, wake vortex shedding impacts on the leading edge of the plate. The large scale of wake vortex makes the turbulence intensity profile much fuller. For the three RSI frequencies, the value of turbulence intensity is the highest at *f* = 20 Hz, nearly double the value at *f* = 40 Hz in the outside region. However, their peak value has a small difference. This indicates that the RSI frequency has little influence on the turbulence intensity in the near-wall boundary layer; the turbulence intensity near the wall might be affected by the wall shape. RSI frequency mainly affects the turbulence intensity in the freestream flow region.

(**a**) Spectral plots of velocity fluctuation at *f* = 20 Hz, left: ȁݑ௫ <sup>ᇱ</sup> |, right: ȁݑ௬ <sup>ᇱ</sup> ȁ.

(**b**) Spectral plots of velocity fluctuation at *f* = 30 Hz, left: ȁݑ௫ <sup>ᇱ</sup> |, right: ȁݑ௬ <sup>ᇱ</sup> ȁ.

(**c**) Spectral plots of velocity fluctuation at *f* = 40 Hz, left: ȁݑ௫ <sup>ᇱ</sup> |, right: ȁݑ௬ <sup>ᇱ</sup> ȁ.

**Figure 9.** Spectral plots of velocity fluctuation at *y*/*Ls* = 0.25 with different RSI frequency: (**a**) *f* = 20 Hz (**b**) *f* = 30 Hz (**c**) *f* = 40 Hz.

**Figure 10.** Turbulence intensity distribution in the boundary layer at different streamwise position.

Figure 11a–c is the time-space diagram of turbulence intensity calculated by different RSI frequency. The horizontal axis is the plate length; the ordinates are the time normalized by the passing period. Figure 11d is the schematic diagram of turbulence intensity at *f* = 40 Hz which outlines the boundary layer states. At this moment, it is necessary to introduce the boundary layer states described by Halstead et al. [42]. Due to the relative movement of blade rows, the wake from upstream airfoils periodically impacts the surface of the downstream plate, thus make the turbulence intensity in the plate periodically increase and decrease. The turbulence spots are initiated at some locations downstream from the leading edge by the wake disturbance. Region C in low turbulence intensity is a laminar region between wakes. Region A was used to represent the turbulence initiated by the wake. It represents a "wake-induced transitional strip". The boundary becomes turbulent when transition is completed within the region A, the strip continues to propagate to the trailing edge of the plate as a "wake-induced turbulent strip" labeled as E, in which turbulence intensity is decreasing along the streamwise direction. The turbulence spots are followed by a calmed region B. The high shear stress in calmed region B is considered to be one of the characteristics that can be effectively inhibit flow separation and delay transition. The region D followed the calmed region B represents the transition between wakes, and it finally becomes a turbulent region F. The calmed region extends to *x*/*c* = 5, 2.2, 2 from the leading edge for *f* = 20 Hz, *f* = 30 Hz, and *f* = 40 Hz correspondingly. Although the calmed region at *f* = 20 Hz is extended, it also could be found that in two adjacent periods, the area of the calmed region in one period is reduced a lot. Moreover, this result can be observed at *f* = 30 Hz as well. The calmed region in two adjacent periods is not equal until the RSI frequency increases to *f* = 40 Hz. This is mainly because the size of wake at *f* = 20 Hz and *f* = 30 Hz is quite large, so the trailing edge of the last wake on the boundary interacts with the leading edge of the next wake. As a result, it makes the laminar flow between wakes turns into turbulent flow in advance.

#### *4.3. Convection of Wake and Wake-Induced Structure*

By means of analyzing the distributions of the phase-averaged results, the whole process of turbulent wake–boundary layer interaction in a wake cycle can be precisely revealed. Figures 12–14 show the phase-averaged results of selected moments during one passing cycle with different RSI frequency. Three representative phases in a passing cycle are selected to show the wake movement of wakes and the interaction with the boundary layer. The first column is the velocity distribution in three selected phases, the second column is the turbulence kinetic energy distribution, and the third column is the z vorticity distribution.

**Figure 11.** Time-Space diagram of turbulence intensity at *y*/*Ls* = 0.25 for different frequency (**a**) *f* = 20 Hz (**b**) *f* = 30 Hz (**c**) *f* = 40 Hz (**d**) schematic of boundary layer states at *f* = 40 Hz.

(**a**) Contours at phase 1, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

(**b**) Contours at phase 2, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

**Figure 12.** *Cont.*

(**c**) Contours at phase 3, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

**Figure 12.** Distributions of phase-averaged results at three selected phases at *f* = 20 Hz.

(**a**) Contours at phase 1, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

(**b**) Contours at phase 2, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

(**c**) Contours at phase 3, left: velocity, middle: turbulence kinetic energy, right: z vorticity. **Figure 13.** Distributions of phase-averaged results at three selected phases at *f* = 30 Hz.

(**a**) Contours at phase 1, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

(**b**) Contours at phase 2, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

(**c**) Contours at phase 3, left: velocity, middle: turbulence kinetic energy, right: z vorticity.

**Figure 14.** Contours of phase-averaged results at three selected phases at *f* = 40 Hz.

As many experimental results have shown [43], the velocity distributions show that the wake is a negative jet with high turbulence kinetic energy. Moreover, it can be observed from the vorticity distribution results in Figures 12–14 that the wake is opposite around the centerline, the vorticity from the pressure surface is positive while the vorticity from the suction surface is negative. The vorticity distributions on both sides of the airfoil are asymmetrical due to the different flow on the pressure surface and suction surface at the trailing edge of airfoil. Comparing the phase-averaged results in different RSI frequency, it is found that the rotor wake is twisted and produces a pocket vortex at the end of wake vortex and a thickened boundary layer at the front of wake vortex. The thickened boundary layer has more concentrated vorticity and higher turbulence kinetic energy than the phase-averaged results of the wake. From the velocity distribution in Figures 12–14, it is observed that there are a high-speed region and a low-speed region in the surrounding area of the thickened boundary layer, and those flow structures move downstream along the plate boundary together.

When the wake vortex is close enough to the wall, the pressure gradient is strong enough to make the local fluid separate, thus a bubble containing vorticity of inverse sense to main vortex is produced (Figure 12a). Propagating downstream, this bubble grows rapidly to the separation point where it departs from the wall and evolves into a secondary vortex (red dotted line), which is supplied by a vortex sheet from the separation point (Figure 12b). However, the secondary vortex at *f* = 30 Hz and *f* = 40 Hz has been always attached to the wall and has no separation. This is mainly because the wake vortex strength at *f* = 30 Hz and *f* = 40 Hz is weakened compared with that of *f* = 20 Hz. According to the theory of Doligalski [44], the fractional convection rate α increases as the absolute strength of the vortex is decreased. For α < 0.75, in the region behind the vortex, a rapid lifting of the boundary-layer streamlines and strong growth of boundary-layer occurs. The fractional convection rate α < 0.75 at *f* = 20 Hz while α > 0.75 at *f* = 30 Hz and *f* = 40 Hz. So, the secondary vortex induced by wake vortex keep attached to the plate surface. What's more, the phase-averaged results could also explain the length difference of the calmed region shown in Figure 11. The calmed region of *f* = 20 Hz is obviously longer than that of *f* = 30 Hz and *f* = 40 Hz. This is because the tail of positive vorticity wake (wake between red dotted line and blue dotted line) interacts with the boundary layer like a calmed region, so the calmed region is extended. It decreases the thickness of the boundary layer. It can be observed that the tail of positive vorticity wake interacts with the boundary layer for a long distance at *f* = 20 Hz. In contrast, the positive vorticity wake interacts with the boundary layer for a pretty short distance at *f* = 30 Hz and *f* = 40 Hz. So, the calmed region at *f* = 30 Hz and *f* = 40 Hz is shorter compared with that of *f* = 40 Hz.

Figure 15 shows the boundary layer thickness δ<sup>99</sup> at the leading edge *x*/*c* = 2 with different RSI frequency. The black dashed-dotted line represents the center of the wake, labeled as C. While the black dashed lines around the center line represent the leading and trailing edges of the wake, labeled as LE and TE. It is apparent that the thickness of boundary layer rapidly increases because of the presence of secondary vortex, the maximum thickness of boundary layer lies on the center line of secondary vortex. Because of the presence of calmed region, the thickness of boundary layer becomes thinner. However, due to the interaction of the negative vorticity vortex (Figure 13 blue dotted line), a low speed region is created on the plate surface, thus a platform appears on the decreasing curve. But there is no platform on the decreasing curve at *f* = 40 Hz, this is mainly because the negative vorticity vortex in the wake decreases in size and does not contact the plate surface. Finally, the thickness of boundary layer continues to decrease to a laminar boundary layer in the absence of wake interaction. It can be concluded that the wake itself cannot suppress separation and transition. However, the calmed region induced by the interaction between the airfoil wake and the boundary layer has the effects of suppressing separation and making the boundary layer thinner.

**Figure 15.** Boundary layer thickness δ<sup>99</sup> variation in two periods at different RSI frequency.

#### **5. Conclusions**

The wake–boundary layer interactions in the turbomachine field are key fundamental issues from a physics perspective, i.e., the boundary layer behavior and properties with the effects of rotor wakes, as well as the engineering design, such as flow separation near the boundary layer and boundary layer heat transfer. The research in this article mainly focus on the fundamental physical mechanism of wake–boundary layer interaction. In order to extensively study the wake–boundary layer interaction, the γ-*Re*θ*<sup>t</sup>* transition SST model is adopted in the transient numerical simulation of the 3D airfoil cascade and flat plate. Furthermore, the effects of different RSI frequency are investigated. The time-averaged results and transient results are compared with the experimental results, and the comparisons show that the numerical results have a reasonable agreement. This research deepens the understanding of the complex flow of wake–boundary layer interaction in turbomachinery. The results obtained in this article could provide some useful references for the design of turbomachinery. The main conclusions are as follows:

The boundary layer profiles are much fuller at RSI frequency *f* = 40 Hz. The boundary layer thickness of the flat plate decreases as the RSI frequency increases. The boundary layer thickness increases as the fluid flows downstream. As the wake impinges on the boundary layer of plate, the streamwise velocity fluctuations in the boundary layer are far greater than for the wall-normal fluctuation and spanwise fluctuation. For the three RSI frequency, the fluctuation amplitude at partial flow condition *f* = 20 Hz is the highest, while the fluctuation amplitude at full flow condition *f* = 40 Hz is the lowest. The spectral plots of *ux* , *uy* at leading edge *x*/*c* = 2 under partial flow condition *f* = 20 Hz and *f* = 30 Hz are dominated by a narrowband component. However, the amplitude of fluctuation at the trailing edge *x*/*c* = 14 has little to do with the variation of RSI frequency and it is dominated by an increasingly broadband spectrum.

The turbulence intensity in the boundary layer increases first and then decreases; a peak value exists in the buffer layer. For the three RSI frequencies, the value of turbulence intensity is the highest at partial flow condition *f* = 20 Hz, nearly double the value at design condition *f* = 40 Hz in the outside freestream region. However, their peak value has a small difference. This indicates that the RSI frequency has little influence on the turbulence intensity in the near-wall boundary layer. The turbulence intensity might be affected by the wall shape and freestream velocity. RSI frequency mainly affects the turbulence intensity in the freestream region.

A secondary vortex with high turbulence kinetic energy is produced as the rotor wake interacts with the boundary layer of plate. And a calmed region is formed behind the secondary vortex. The calmed region has an effective suppression of the boundary layer transition. The calmed region at partial flow condition *f* = 20 Hz is far longer than that for design condition *f* = 40 Hz. This is mainly because the positive vorticity wake has a similar effect as a calmed region. The tail of the positive vorticity wake interacts with the boundary layer for a long distance at partial flow condition *f* = 20 Hz while it only continues a short distance at design condition *f* = 40 Hz. In the leading edge, the thickness of boundary layer rapidly increases because of the presence of the secondary vortex. The maximum thickness of boundary layer lies on the center of secondary vortex. Then the thickness of boundary layer decreases because of the calmed region. Next, the thickness of boundary layer goes into a platform as the negative-vorticity wake interacts with the boundary layer. Finally, the boundary layer becomes a laminar boundary layer in the absence of wake interaction.

**Author Contributions:** Conceptualization, X.R.; methodology, X.Z.; software, X.Z.; validation, J.W., Z.X.; formal analysis, X.Z.; investigation, X.R.; resources, X.R.; data curation, Z.X.; writing—original draft preparation, X.Z.; writing—review and editing, P.W.; visualization, X.Z.; supervision, X.R.; project administration, X.R.; funding acquisition, X.R. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the Science Fund for Creative Research Groups of National Natural Science Foundation of China [Grant No. 51821093] and the National Basic Research Program of China [Grant No. 2015CB057301].

**Acknowledgments:** We also highly appreciate reviewers for useful comments and editors for improving the manuscript.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
