*Article* **Aerodynamic Efficiency Improvement on a NACA-8412 Airfoil via Active Flow Control Implementation**

**Nil Couto <sup>1</sup> and Josep M. Bergada 2,\***


**\*** Correspondence: josep.m.bergada@upc.edu; Tel.: +34-937398771

**Abstract:** The present paper introduces a parametric optimization of several Active Flow Control (AFC) parameters applied to a NACA-8412 airfoil at a single post-stall Angle of Attack (AoA) of 15◦ and Reynolds number *Re* <sup>=</sup> 68.5 <sup>×</sup> <sup>10</sup><sup>3</sup> . The aim is to enhance the airfoil efficiency and to maximize its lift. The boundary layer separation point was modified using Synthetic Jet Actuators (SJA), and the airfoil optimization was carried on by systematically changing the pulsating frequency, momentum coefficient and jet inclination angle. Each case has been evaluated using Computational Fluid Dynamic (CFD) simulations, being the Reynolds Averaged Navier–Stokes equations (RANS) turbulence model employed the Spalart Allmaras (SA) one. The results clarify which are the optimum AFC parameters to maximize the airfoil efficiency. It also clarifies which improvement in efficiency is to be expected under the operating working conditions. An energy balance is presented at the end of the paper, showing that for the optimum conditions studied the energy saved is higher than the one needed for the actuation. The paper clarifies how a parametric analysis has to be performed and which AFC parameters can be initially set as constant providing sufficient previous knowledge of the flow field is already known. A maximum efficiency increase versus the baseline case of around 275% is obtained from the present simulations.

**Keywords:** Computational Fluid Dynamics (CFD); Active Flow Control (AFC); synthetic jet; boundary layer; aerodynamic efficiency

### **1. Introduction**

The Active Flow Control (AFC) technology applied to airfoils is beginning to be of common use, then a considerable number of recent papers tackle this matter [1–3]. Parametric analysis is still the most common methodology employed to evaluate the effect of the AFC on airfoils, yet the use of optimizers is becoming popular to obtain the most appropriate AFC parameters in any given application [4,5]. Some initial studies analyzed the performance when using continuous and pulsating actuations, and it was soon discovered that pulsating flow was coupling with the boundary layer natural instabilities and so being more energetically efficient [6–8].

Airfoil aerodynamic efficiency decisively effect fuel consumption, and AFC can improve it. Briefly, the AFC technology adds or subtracts momentum to/from the main flow with the aim of interacting with the boundary layer and delaying or advancing its separation. One of the advantages over passive flow control is that it is not creating drag increase when operating outside the design conditions. Cattafesta and Sheplak [9] divided the AFC techniques into three categories. (1) Moving body actuators, they do not add or substract mass and induce local fluid motion [10]. (2) Plasma actuators, they generate high frequency ionized jets of fluid [11–14]. (3) Fluidic actuators (FA), create pulsating flow and in some configurations no moving parts are required. In some particular FA designs, the origin of the self-sustained oscillations was recently unveiled in [15–17].

**Citation:** Couto, N.; Bergada, J.M. Aerodynamic Efficiency Improvement on a NACA-8412 Airfoil via Active Flow Control Implementation. *Appl. Sci.* **2022**, *12*, 4269. https://doi.org/ 10.3390/app12094269

Academic Editor: Wei Huang

Received: 13 March 2022 Accepted: 20 April 2022 Published: 23 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 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 (https:// creativecommons.org/licenses/by/ 4.0/).

From the different sorts of FA, Zero Net Mass Flow Actuators (ZNMFA) also called Synthetic Jet Actuators (SJA), have gained particular attention in the scientific community due to their high capabilities in controlling boundary layer separation [18–21]. Another advantage of SJAs is that do not require external fluid supply.

When considering the implementation of SJAs, two parameters are particularly relevant: the non-dimensional frequency *F* <sup>+</sup> = *f C*/*U*∞, where *f* is the dimensional frequency, *C* the airfoil chord and *U*<sup>∞</sup> the free-stream velocity, and the momentum coefficient (*Cµ*). *C<sup>µ</sup>* is defined as *C<sup>µ</sup>* = (*hρjet*(*U*<sup>2</sup> *j* ) sin *θ*)/(*Cρ*∞(*U*<sup>2</sup> <sup>∞</sup>)), where *h* is the jet width, *ρjet* is the jet density and *ρ*<sup>∞</sup> is the far field one, *U<sup>j</sup>* characterizes the jet maximum velocity, the parameter *θ* defines the inclination angle of the jet versus the adjacent surface.

Amitay et al. [22] and Amitay and Glezer [23] studied the SJA frequency, position and *C<sup>µ</sup>* when applied to a symmetric airfoil and observed that when locating the actuator near to the boundary layer separation point, the momentum coefficient needed to reattach the flow was smaller. Regarding the actuation frequency, they found that when it was an order of magnitude higher than the vortex shedding one (*F* <sup>+</sup>O(10)), the flow was completely reattached.

SJAs on a NACA-0015 airfoil at *Re* <sup>=</sup> 8.96 <sup>×</sup> <sup>10</sup><sup>5</sup> were experimentally studied by Gilarranz et al. [24]. They increased the stall angle from 12◦ to 18◦ and observed that for AoA > 10◦ , the effectiveness of the actuation considerably increased. At AoA > 25◦ , considerably high actuation frequencies were needed to observe large flow variations. The same airfoil and Reynolds number was studied by You and Moin [25] using Large Eddy Simulation (LES). A lift increase of 70% was obtained when using *C<sup>µ</sup>* = 0.0123, *F* <sup>+</sup> = 1.284 and *<sup>θ</sup>* = 30.2º. The same airfoil at *Re* <sup>=</sup> 3.9 <sup>×</sup> <sup>10</sup><sup>4</sup> was experimentally studied by Tuck and Soria [26]. The optimum momentum coefficient and actuation frequencies were respectively, *C<sup>µ</sup>* = 0.0123O and *F* <sup>+</sup> = 0.7 and 1.3, being the highest frequency the most effective one. When using these forcing conditions, an improvement of the stall AoA from 10◦ to 18◦ was observed. Kitsios et al. [27] numerically studied the same airfoil and Reynolds number via LES simulations. They observed that the optimal frequencies coincided with the vortex shedding one obtained from the baseline case (*fwake*) and its first harmonic (2 *fwake*). Buchmann et al. [28] reached the same conclusion through experimental tests. Flow control using SJAs on a NACA-23012 at *Re* <sup>=</sup> 2.19 <sup>×</sup> <sup>10</sup><sup>6</sup> was numerically investigated by Kim and Kim [29], Monir et al. [30]. Maximum lift was obtained when the jet was placed nearby the boundary layer separation point and *F* <sup>+</sup> = 1. When applying SJA tangentially, clear advantages were observed than when injecting/sucking fluid perpendicular to the surface.

The separation control generated by SJA on a NACA-0025 airfoil at *Re* = 10<sup>5</sup> and AoA = 5◦ were experimentally studied by Goodfellow et al. [31]. Momentum coefficient was found to be the primary control parameter, a drag decrease of nearly 50% was obtained when *C<sup>µ</sup>* was higher than a certain threshold value. The same Reynolds number and wing profile, but at AoA = 10◦ was evaluated by Feero et al. [32]. They observed that the momentum coefficient needed *C<sup>µ</sup>* to keep the flow reattached, was an order of magnitude lower when the excitation frequencies happened to be around the vortex shedding one. The same Reynolds number and wing profile but at AoA = 12◦ was analyzed in [33], where the jet position effectiveness in controlling the flow separation was evaluated. The optimum results were obtained when locating the jet groove nearby (whether downstream or upstream) the boundary layer separation point, a small improvement was observed in the upstream location. The NACA-0018 profile at AoA = 10◦ and at *Re* = 1000, was evaluated using direct numerical simulation (3D-DNS) by Zhang and Samtaney [34]. Three non-dimensional frequencies (*F* <sup>+</sup> = 0.5, 1 and 4) were considered. The optimum frequency was *F* <sup>+</sup> = 1 although improvement was observed in all of them. In Rodriguez et al. [35], the SD7003 airfoil at *Re* <sup>=</sup> <sup>6</sup> <sup>×</sup> <sup>10</sup><sup>4</sup> and under three AoA = <sup>4</sup> ◦ , 11◦ and 14◦ was numerically studied. At AoA = 14◦ they obtained an aerodynamic efficiency increase of 124%. The same Reynolds number and wing profile but at AoA = 14◦ , and AoA = 13◦ & 16◦ was studied by [4,5], respectively. Due to the fact that in these studies optimization processes were employed, the maximum efficiency increase was obtained to be 280% and 591%, respectively.

The airfoil chosen for the present study is the NACA-8412, a Reynolds number of *Re* = 68.5<sup>×</sup> <sup>10</sup><sup>3</sup> and AoA = <sup>15</sup>◦ were considered. According to the authors knowledge no previous studies for this airfoil and Re configuration have been found in the literature. The usage of the NACA-8412 airfoil was motivated by its large curvature, which in combination with the high AoA used will produce a separation of the flow close to the leading edge. This early separation will remark the effects of introducing SJA, and enable a clearer comparison of the results with the baseline case. Furthermore, the methodology employed in the present manuscript, which is based on a previous knowledge of the flow field and consist in initially fixing some AFC parameters, is proven to be very useful in obtaining good results while saving computational resources.

The rest of the manuscript is structured as follows. The formulation of the problem, the numerical methods employed and the mesh independence study are presented in Section 2. The definition of the AFC parameters and their implementation are introduced in Section 3. The Results section is introduced in Section 4 and the summary of the work is presented in Section 5.

### **2. Numerical Method**

### *2.1. Governing Equations and Turbulence Model*

Due to the unsteadiness of the problem, three possible turbulence models arose: Direct Numerical Simulation (DNS), Large Eddy Simulations (LES) and Unsteady Reynolds Averaged Navier–Stokes (URANS). Due to the large computational power required to perform both DNS and LES simulations, URANS was selected as the turbulence model used. More specifically, the Spalart–Allmaras model was chosen to perform all the simulations involved in the present study because of its suitability for low Reynolds cases and its easy convergence when the Courant–Friedrichs–Levy number *CFL* < 1. The Navier–Stokes equations for incompressible flow, take the form:

$$\frac{\partial u\_i}{\partial x\_i} = 0 \tag{1}$$

$$\frac{\partial u\_i}{\partial t} + \frac{\partial u\_i u\_j}{\partial x\_j} = -\frac{1}{\rho} \frac{\partial p}{\partial x\_i} + \nu \frac{\partial^2 u\_i}{\partial x\_j \partial x\_j} \tag{2}$$

Substituting each variable *φ* by the average *φ*¯ and fluctuation *φ* 0 terms.

$$
\phi = \Phi + \phi'\tag{3}
$$

Being *φ*(*x*, *t*) a generic flow variable, it can be written as the sum of its mean *φ*¯(*x*, *t*) and fluctuating components *φ* 0 (*x*, *t*), being *φ*(*x*, *t*) = *φ*¯(*x*, *t*) + *φ* 0 (*x*, *t*). By substituting it in the Navier–Stokes equations and taking the time average, it results in the following Reynolds-averaged Navier–Stokes equations,

$$\frac{\partial \overline{u}\_i}{\partial x\_i} = 0 \tag{4}$$

$$\frac{\partial \overline{u}\_i}{\partial t} + \frac{\partial \overline{u}\_i \overline{u}\_j}{\partial x\_j} = -\frac{1}{\rho} \frac{\partial \overline{p}}{\partial x\_i} + \nu \frac{\partial^2 \overline{u}\_i}{\partial x\_j \partial x\_j} - \frac{\partial u\_i' u\_j'}{\partial x\_j} \tag{5}$$

where the term *u* 0 *i u* 0 *j* is the Reynolds stress tensor denoted by *Rij* and was approximated using Boussinesq hypothesis. The deviatoric part of the tensor is given as follows,

$$R\_{ij} - \frac{1}{3} R\_{kk} \delta\_{ij} = -2\nu\_t \tilde{S}\_{ij} \tag{6}$$

where *ν<sup>t</sup>* is the kinematic eddy viscosity that needs to be modeled using one of the RANS models.

As it has been previously mentioned, the turbulence model selected was the Spalart– Allmaras (S-A). The model proposed by [36] solves a single transport equation for the modified form of the turbulent kinetic energy called *ν*˜. This new parameter is identical to *ν<sup>t</sup>* except in the viscous-affected region nearby the wall.

$$\nu\_t = \tilde{\nu}\left(\frac{\chi^3}{\chi^3 + \mathcal{C}\_{v1}^3}\right) \tag{7}$$

The transport equation of *ν*˜ is given as:

$$\frac{\partial \boldsymbol{\vartheta}}{\partial t} + \boldsymbol{u}\_{j} \frac{\partial \boldsymbol{\vartheta}}{\partial \mathbf{x}\_{j}} = \frac{1}{\sigma} \left[ \frac{\partial}{\partial \mathbf{x}\_{j}} \left( (\boldsymbol{\nu} + \boldsymbol{\vartheta}) \frac{\partial \boldsymbol{\vartheta}}{\partial \mathbf{x}\_{j}} \right) + \mathbf{C}\_{\delta 2} \frac{\partial \boldsymbol{\vartheta}}{\partial \mathbf{x}\_{i}} \frac{\partial \boldsymbol{\vartheta}}{\partial \mathbf{x}\_{i}} \right] + \mathbf{C}\_{\delta 1} (1 - f\_{j2}) \tilde{\mathbf{S}} \boldsymbol{\vartheta} - \left[ \mathbf{C}\_{w1} f\_{w} - \frac{\mathbf{C}\_{\delta 1}}{\kappa^{2}} f\_{l2} \right] \left( \frac{\boldsymbol{\vartheta}}{d} \right)^{2} \tag{8}$$

where:

$$\begin{aligned} f\_{\omega} &= g \left[ \frac{1 + \mathcal{C}\_{t\omega}}{\mathcal{S}^{6} + \mathcal{C}\_{\omega 3}} \right]^{1/6}, \quad g = r + \mathcal{C}\_{\omega 2} (r^{6} - r), \quad r = \min \left[ \frac{\wp}{\mathcal{S} \kappa^{2} d^{2}}, 10 \right] \\ \mathcal{S} &= \Omega + \frac{\mathcal{V}}{\kappa^{2} d^{2}} f\_{v2}, \quad f\_{t2} = \mathcal{C}\_{t3} \exp(-\mathcal{C}\_{t4} \chi^{2}), \quad f\_{v2} = 1 - \frac{\chi}{1 + \chi f\_{v1}} \\ \mathcal{f}\_{v1} &= \frac{\chi^{3}}{\chi^{3} + \mathcal{C}\_{v1}^{3}}, \quad \chi = \frac{\mathcal{V}}{\nu} \end{aligned} \tag{9}$$

Being *d* the distance from a given point to the nearest wall, and Ω the magnitude of the vorticity. The model constants have the following default values:

$$\begin{aligned} \text{C}\_{b1} &= 0.1355, \quad \text{C}\_{b2} = 0.622, \quad \sigma = \frac{2}{3}, \quad \text{C}\_{\omega1} = \frac{\text{C}\_{b1}}{\text{\textkappa}^2} + \frac{(1 + \text{C}\_{b2})}{\sigma} \\ \text{C}\_{\omega2} &= 0.3, \quad \text{C}\_{\omega3} = 2.0, \quad \text{C}\_{v1} = 7.1, \quad \text{\kappa} = 0.4187, \quad \text{C}\_{l3} = 1.2 \\ \text{C}\_{t4} &= 0.5 \end{aligned} \tag{10}$$

### *2.2. Numerical Domain and Boundary Conditions*

The NACA-8412 airfoil's leading edge (LE) was located at the origin of the coordinate system, and the airfoil was fixed at zero degree angle with respect to the *x*-axis. The poststall angle of attack of 15◦ was introduced by decomposing the freestream velocity in the corresponding *x* and *y* components. Note that all distances have been expressed as function of the airfoil chord *C* (see Figure 1). The distance between the inlet and the LE following the *x*-axis was set at 10*C*, while the distance between the trailing edge (TE) and the outlet was set at 18*C*. The height of the domain is 16*C*, symmetrically distributed with respect to the airfoil. An almost identical computational domain was recently employed in [4].

**Figure 1.** Sketch of the computational domain. The inlet and outlet boundaries have their limits in the points defined by the letters A and B.

The inlet was defined by the line connecting points A and B (see Figure 1) in a counter-clockwise direction, while the outlet was defined by the line connecting A and B in clockwise direction. At the inlet, Dirichlet boundary conditions for velocity were set by imposing the decomposition of the freestream velocity as it has been previously mentioned. Regarding pressure, Neumann boundary conditions were used. At the outlet, Neumann boundary conditions for velocity and Dirichlet boundary conditions for pressure were imposed. Over the airfoil's surface Dirichlet boundary conditions for velocity were imposed in order to guarantee the non-slip condition, Neumann boundary conditions for pressure were considered.

The point *A* represented in Figure 1, which corresponds to one of the connecting points between the inlet and the outlet, had to guarantee that the angle formed by its tangent was larger than the angle of attack studied (AoA = 15◦ ). For the present configuration, the tangent line to the inlet curve at the point *A* and the horizontal form an angle of 36.87º, complying with the mentioned restriction. Figure 2 introduces some AFC parameters as the jet position, jet width and jet inclination angle. The International Standard Atmosphere (ISA) model at sea level was used to select the air density, the pulsating flow fluid density, as well as the fluid dynamic viscosity. A relatively low Reynolds number was selected in order to have proportionally low mesh computational requirements and computational time. The values of the physical parameters involved in the present study are summarized in Table 1.

**Figure 2.** Sketch of the groove location and injection angle over the airfoil.

**Table 1.** Physical conditions.


*2.3. Non-Dimensional Parameters*

In this section the non dimensional groups used along the paper are introduced. **Reynolds number**

$$Re = \frac{\mathcal{U}\_{\infty}\mathcal{C}}{\nu} \tag{11}$$

where *U*<sup>∞</sup> is the free-stream velocity, *C* the chord length and *ν* the fluid kinematic viscosity. **Lift coefficient**

$$\mathcal{C}\_{l} = \frac{2F\_{l}}{\rho \mathcal{C} U\_{\infty}^{2}} \tag{12}$$

where *ρ* is the free-stream fluid density, and *F<sup>l</sup>* the dimensional lift force. **Drag coefficient**

$$\mathcal{C}\_d = \frac{2F\_d}{\rho \mathcal{C} U\_\infty^2} \tag{13}$$

where *F<sup>d</sup>* is the dimensional drag force. Note that the lift and drag coefficients are a function of the angle of attack (AoA) and Reynolds number.

**Aerodynamic efficiency**

$$E = \frac{F\_l}{F\_d} = \frac{\mathcal{C}\_l}{\mathcal{C}\_d} \tag{14}$$

**Pressure coefficient**

$$\mathcal{C}\_P = \frac{P - P\_{\infty}}{\frac{1}{2}\rho U\_{\infty}^2} \tag{15}$$

The parameter *P* is the generic static pressure and *P*<sup>∞</sup> is the free-stream one. **Friction coefficient**

$$C\_F = \frac{\tau\_w}{\frac{1}{2}\rho U\_\infty^2} \tag{16}$$

where *τ<sup>w</sup>* are the wall shear stresses.

**Dimensionless wall distance**

$$y^{+} = \frac{u\_{\tau}y}{\nu} \tag{17}$$

With:

$$u\_{\tau} = \sqrt{\frac{\tau\_w}{\rho}}\tag{18}$$

*u<sup>τ</sup>* characterizes the friction velocity and *y* is the dimensional distance with respect to the wall. **Jet momentum coefficient**

$$\mathbf{C}\_{\mu} = \frac{h(\rho \mathcal{U}\_{\text{max}}^2) \sin \theta\_{\circ}}{\mathcal{C}(\rho \mathcal{U}\_{\text{\infty}}^2)} \tag{19}$$

where *h* characterizes the jet width, *Umax* is the maximum jet velocity and *θ<sup>j</sup>* stands for the jet inclination angle measured versus the airfoil surface (see Figure 2).

**Non-dimensional forcing frequency, (Strouhal number)**

$$Strouhal = F^{+} = \frac{f\mathbb{C}}{\mathcal{U}\_{\infty}}\tag{20}$$

The non-dimensional forcing frequency is often given as:

$$F\_\*^+ = \frac{f}{f\_0} \tag{21}$$

The variable *f* is a generic frequency and *f*<sup>0</sup> is the vortex shedding frequency obtained from the baseline case. As in the preset study *C* = 1 and *U*<sup>∞</sup> = 1, the relation between *F* + and *F* + <sup>∗</sup> becomes, *<sup>F</sup>* <sup>+</sup> = *F* + ∗ *f*0.

**Courant–Friedrichs–Levy number**

$$\text{CFL} = \frac{u\Delta t}{\Delta x} \tag{22}$$

where *u* is the generic fluid velocity, ∆*t* is the time step and ∆*x* the generic mesh cell length.

### *2.4. Mesh Assessment*

Due to the fact that no previous experimental or numerical studies are available in the literature, it is essential to validate the numerical model. In order to study the baseline case and the later AFC configurations, a mesh independence test was performed to guarantee independence between the mesh resolution and the simulated results. A hybrid mesh configuration, composed by a rectangular structured and a triangular unstructured mesh region, was implemented in order to capture the boundary layer while saving computational power required on the outer regions of the domain. The surrounding region of the airfoil consists on two body-fitted structured sub-meshes with a progressive growth of 1% for the one closest to the airfoil and 3% for the one surrounding the first sub-mesh. The remaining portion of the domain was filled with an unstructured mesh, the resolution was particularly high on the wake region. An overall view of the computational domain is presented in Figure 3, a close up view of the regions surrounding the airfoil for the baseline case and for the corresponding AFC mesh implementation is shown in Figure 4 and Figure 5, respectively.

**Figure 3.** Overall view of the computational domain.

**Figure 4.** Detailed view of the near-wall region of the mesh. (**a**) Close up view of the structured sub-mesh. (**b**) Close up view of the LE region structured sub-mesh.

**Figure 5.** Mesh modification for the AFC jet implementation. (**a**) Close up view of the structured sub-mesh. (**b**) Close up view of the jet area.

Four mesh candidates were designed according to four decreasing *y* <sup>+</sup> values, and consequently an increasing mesh resolution, as it can be observed in Table 2 from the number of mesh cells *Ncell* employed in the CFD models. Due to the lack of experimental data for the NACA-8412 airfoil at *Re* <sup>=</sup> 68.5 <sup>×</sup> <sup>10</sup><sup>3</sup> , the evaluation of the proper mesh could not be based on the relative error between experimental and computational data. Alternatively, the strategy followed was to study the temporal averaged *C<sup>l</sup>* and *C<sup>d</sup>* , as well as *C<sup>P</sup>* and *C<sup>F</sup>* distributions along the chord for each mesh candidate, and to compare them to observe for which *y* <sup>+</sup> the results did not significantly change. This approach was also applied to evaluate whether the relative error between the results compensated in terms of computational cost. A baseline case simulation was performed for each mesh candidate for a computational time of 30 s, which was determined based on the evaluation of the velocity and pressure residuals. The time step was of <sup>2</sup> <sup>×</sup> <sup>10</sup>−<sup>4</sup> s. In order to avoid the transient phase, the resulting time averaged for the lift *C<sup>l</sup>* and drag *C<sup>d</sup>* coefficients were computed considering the last 10 s of each simulation. The resulting averaged aerodynamic coefficients together with the relative error *ClError*, *CdError* given in percentage are presented in Table 2. Calculated *y* <sup>+</sup> values were assessed on the upper side of the airfoil due to its importance in the later AFC implementation, and its maximum value obtained for each case is also included in Table 2. The resulting pressure and friction coefficients are introduced in Figures 6 and 7, respectively.

Table 2 shows that the time averaged aerodynamic coefficients stabilize as *y* <sup>+</sup> decreases. For instance, the relative error between *y* <sup>+</sup> = 0.3 and *y* <sup>+</sup> = 1 is 14.93% for *C<sup>l</sup>* and 9.14% for *C<sup>d</sup>* . On the other hand, the relative error between *y* <sup>+</sup> = 0.3 and *y* <sup>+</sup> = 0.5 is 3.35% for *C<sup>l</sup>* and 1.97% for *C<sup>d</sup>* . As the difference among the aerodynamic coefficients for *y* <sup>+</sup> = 0.3 and *y* <sup>+</sup> = 0.5 can be considered relatively small (∆*C<sup>l</sup>* <sup>=</sup> <sup>−</sup>0.0425 and <sup>∆</sup>*C<sup>d</sup>* <sup>=</sup> 4.8 <sup>×</sup> <sup>10</sup>−<sup>3</sup> , with reference values corresponding to the *y* <sup>+</sup> = 0.3 results), and the computational time for the *y* <sup>+</sup> = 0.3 resulted in nearly the double than the one required for *y* <sup>+</sup> = 0.5, the results corresponding to the mesh with an estimated *y* <sup>+</sup> = 0.5 were considered acceptable for the present simulations. It is relevant to highlight that the calculated maximum *y* <sup>+</sup> values on the airfoil upper surface were slightly lower than the estimated ones used to determine

the first cell height, which certifies that the mesh was properly designed and ensured the condition of *y* <sup>+</sup> < 1 necessary to evaluate the boundary layer.

**Table 2.** Mesh candidates for the Mesh Independence Test.


**Figure 6.** Baseline case pressure coefficient comparison.

In order to further check the mesh to be employed, the pressure coefficient *C<sup>P</sup>* evaluated at each point of the airfoil surface and for all *y* <sup>+</sup> studied is presented in Figure 6. The *C<sup>P</sup>* distribution along airfoil lower surface showed a very similar behaviour for every mesh candidate, despite that *y* <sup>+</sup> = 1 and *y* <sup>+</sup> = 0.7 showed a slightly lower *C<sup>P</sup>* value for approximately the 80% of the chord. A more notable difference can be observed in the airfoil upper surface, where *C<sup>P</sup>* progressively grows as *y* <sup>+</sup> decreases. It can be observed how *C<sup>P</sup>* in the airfoil upper surface for *y* <sup>+</sup> = 0.3 and *y* <sup>+</sup> = 0.5 are very similar along the whole chord, only showing a minor difference between the leading edge of the airfoil and *x*/*C* = 0.3. When considering the friction coefficient, see Figure 7, no notable differences can be noticed between the *C<sup>F</sup>* distribution for the different meshes studied. The separation point is well predicted by all meshes evaluated, while small differences are observed in the reattachment point. After analyzing the *CP*, *CF*, *C<sup>l</sup>* and *C<sup>d</sup>* coefficients, we can conclude that the results obtained for the mesh with an estimated *y* <sup>+</sup> = 0.5 are representative and almost mesh independent. This mesh was selected for the baseline case evaluation as well as for the later AFC implementation simulations. In addition, selecting this mesh exploited the fact that the computational time required to complete the simulations was drastically lower than for the most dense mesh, without compromising the solution in a great extent. To conclude with the baseline case analysis, the streamlines of the averaged flow field are presented in Figure 8. A huge separation bubble is generated very close to the airfoil LE, as indicated in Figure 7, and it extends up to the vicinity of the TE, where a counter rotating vorticular structure is generated. Notice that the boundary layer separation point is observed to be at about 10% of the chord, *x*/*C* = 0.1045. The objective of the AFC implementation is to postpone the boundary layer separation, breaking this large vorticular

structure onto smaller ones, or even completely suppressing it. This will reduce the drag generated by the large vortex and increase the lift generated by the airfoil.

**Figure 7.** Baseline case friction coefficient comparison.

**Figure 8.** Averaged streamlines of the baseline case.

### **3. AFC Implementation**

### *3.1. Jet Location and Mesh Modification*

In order to evaluate the influence of AFC on aerodynamic efficiency, three design variables have been studied: the jet inclination angle *θ*, the jet momentum coefficient *C<sup>µ</sup>* and the non-dimensional forcing frequency *F* + ∗ . In the present study, the groove width *h* and location were kept constant. The design of the AFC synthetic jet begins by studying its location over the airfoil. According to [33] among others, the optimum location is just upstream of the boundary layer separation point, which can clearly be seen in Figures 7 and 8. The separation point is located at *x*/*C* = 0.1045, considerably close to the LE. This is a consequence of the high AoA used (15º) and of the geometry of the airfoil selected, as it shows a pronounced curvature near the leading edge. At the separation point *C<sup>F</sup>* begins to adopt negative values due to the adverse pressure gradient, leading to flow separation. Based on what it was defined in reference [4], the synthetic jet has been placed at a location of 0.01*C* upstream of the flow separation point, with a slot width of 0.01*C* (1% of the chord), from *x*/*C* = 0.08 to *x*/*C* = 0.09. The next step is to design the AFC velocity signal, which must be time periodic as the jet corresponds to the ZNMF typology. Thus, in the following subsections, the frequency, amplitude and jet inclination angle will need to be evaluated. The mesh was further refined in the slot region due to the quick variation of the flow properties close to the jet. For this reason, the number of partitions at the vicinity of the jet have been increased as it can be seen in Figure 5. The total number of mesh cells increased to 162,789 for the cases where AFC was considered. Regarding the instantaneous velocity along the groove it was decided to use a top hat velocity distribution.

### *3.2. AFC Velocity Input*

Synthetic jets, also called ZNMF actuators, are characterized by transporting the same quantity of fluid in the blowing and suction phases. For this reason, a sinusoidal velocity input as in Equation (23) must be designed with the objective of extracting and injecting moment from the boundary layer. An initial approach must be carried out to estimate a reference value for the velocity input frequency *f* and amplitude *A*.

$$\mathcal{U}\_{jet}(t) = A \sin(2\pi ft) \tag{23}$$

The introduction of periodic forcing is aimed to modify the airfoil vortex shedding frequency and amplitude. The vortex shedding natural frequency for the baseline case must be assessed by analyzing its frequency spectrum. Following the methodology presented in [37], a Fast Fourier Transform for the lift coefficient of the baseline case was performed. The resulting vortex shedding frequency was *f*<sup>0</sup> = 0.6965 Hz. The amplitude *A* of the velocity signal corresponds to the maximum velocity of the jet *Umax* when sucked/injected from/to the boundary layer. The parametric study presented in this paper is designed as follows: the groove position and width were kept constant, and the natural vortex shedding frequency was selected as input frequency. Then, for a given jet inclination angle, CFD simulations were performed for a set of momentum coefficients *C<sup>µ</sup>* ranging between 0.0001 ≤ *C<sup>µ</sup>* ≤ 0.1. The process was performed for three different jet inclination angles of *θ* = 20º, 30º and 40º. Via following this procedure the optimum momentum coefficient and jet injection angle (from the range chosen) were obtained. As a final step, a set of jet pulsating frequencies ranging between 0.5 ≤ *F* + <sup>∗</sup> ≤ <sup>7</sup> were evaluated for the optimum *<sup>C</sup><sup>µ</sup>* and *θ* configuration previously obtained.

### **4. Results**

Three jet inclination angles *θ* have been studied (20º, 30º and 40º). For each *θ*, twelve momentum coefficients *C<sup>µ</sup>* have been evaluated to obtain the optimal *θ* and *C<sup>µ</sup>* configuration among those studied. The rest of the parameters were kept constant. The range of *C<sup>µ</sup>* values have been selected based on previous validated studies [38].
