**1. Introduction**

For the commonly used three-stage high-lift device, the flow characteristic is that the leading slat and trailing-edge flaps work together to improve the lift–drag characteristics of the entire airfoil. Van Dam [1] reviewed the recent developments in aerodynamic design and analysis methods for multi-element high-lift systems in transport airplanes, including necessity, the major design objectives, and constraints. For the commonly used three-stage high-lift device, the flow characteristic is that the leading slat and trailing-edge flaps work together to improve the lift–drag characteristics of the entire airfoil. In 1970, A. Smith [2] analyzed the aerodynamic principle of the multi-element configuration, and believed that it was affected by the slat effect, circulation effect, off-the-surface pressure recovery, and boundary layer effect of each wing segment. It is well recognized that the usual function of the slot is that of a boundary layer control device, permitting highly adverse upper-surface pressure gradients to be sustained without incurring severe separation. There appear to be five primary effects of gaps, and here we speak of properly designed aerodynamic slots. According to the current classical active and passive flow control theories, the mechanism of the three-segment high-lift device can be classified into the passive control method. The passive flow control technology does not add energy to the flow, but changes the flow boundary or pressure gradient through some methods to control the flow, such as optimizing the geometric parameters of the airfoil, adding or combining flow control devices, and through fixed-point displacement or deformation under specific working conditions. Compared with active flow control technology [3], it is widely valued for its advantages of simple structure, lack of additional energy consumption or lower energy

**Citation:** Yu, J.; Mi, B. A New Flow Control Method of Slat-Grid Channel-Coupled Configuration on High-Lift Device. *Appl. Sci.* **2023**, *13*, 3488. https://doi.org/10.3390/ app13063488

Academic Editors: Francesca Scargiali, Gabriel Bugeda Castelltort and Josep Maria Bergadà

Received: 30 December 2022 Revised: 22 February 2023 Accepted: 6 March 2023 Published: 9 March 2023

**Copyright:** © 2023 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/).

consumption, and low cost. For the three-stage high-lift device, the leading-edge slat can accelerate the airflow through the slots formed with the main wing, forming a similar effect to the jet flow, which can favorably disturb the airflow of the main wing and improve its wind resistance. However, the practical application of this three-element high-lift device has many limitations. For example, the aircraft is significantly affected by the ground effect during takeoff and landing. For the high-lift device, the flaps are deflected closer to the ground [4], meaning that the impact is greater. For amphibious aircraft [5,6], flap failure may even occur due to body splash. In order to further improve the aerodynamic characteristics and practical value of the high-lift device, slats are the preferred component for optimization and improvement.

The basic principle of slat lift enhancement is that there is a gap between the slat and the leading-edge surface of the main wing. Specifically, the velocity of the airflow on the lower surface of the airfoil increases in the gap and flows to two places, respectively. Part of the airflow accelerates from the gap between the tail of the slat and the leading edge of the main wing, flows backwards with the airflow on the upper surface of the airfoil, and accelerates the flow to the convergence of the upper surface of the main wing. Moreover, the velocity of airflow in the boundary layer on the upper surface of the wing is increased to eliminate a large number of vortices [7]; the other part flows to the cavity of the concave pressure surface of the slat, forming a recirculation zone. The velocities of these two streams differ greatly; a strong free shear layer is, therefore, formed at the interface. This shear layer originates from the sharp lower trailing edge and develops downstream with the separation vortex. When the separation vortex hits the trailing edge wall of the slat, a part of the vortex will flow back into the concave cavity of the slat, interfering with the new shear layer near the lower trailing edge [8] and forming return feedback [9]. In order to improve the lift–drag characteristics, grid flow channels can be added to the slats to improve the flow field in the slats. Grid fins [10], as multi-surface control systems consisting of an external frame and an internal grid partition, originate from the lift device composed of multiple elements used in the early design of aircraft, and have the advantages of high control efficiency, a large stall angle of attack, small hinge torque, and easy folding. In particular, their lift characteristics are good, and it is not easy to stall the devices with a high angle of attack, meaning they maintain good linearity. However, due to the defects in the structure, weight, resistance and technology of multi-element lift devices under the technical conditions of the time, single-element lift devices have been dominant in aircraft design for a long time. With the development of technology, the performance characteristics of grid fins derived from multi-element lift devices are increasingly understood, and their advantages as aerodynamic stabilizing and controlling surfaces for various types of aircraft continue to be highlighted. Research from Russia, the United States, Germany, and China has strengthened the understanding of grid fins; as a result, theoretical and applied achievements have been realized [11–15].

Wang Long et al. [16] carried out a numerical simulation on four different slotting positions of wind turbine blades by means of a computational fluid dynamic method, and obtained the flow field and lift coefficient of the blades at different angles of attack. They found that, under high angles of attack, the jet formed by slots can effectively inhibit flow separation, and the aerodynamic performance was better when the slotting position was near the separation point. Hui Zenghong et al. [17] studied the aerodynamic characteristics of multi-element leading slat lift devices using the pressure measurement method and obtained a better leading slat scheme by comparing the influence of different configuration slats on the airfoil lift effect and aerodynamic efficiency. Zhang Lijun et al. [18] proposed a scheme of symmetrical slitting in the wind separation zone for the vertical axis wind wing NACA0012 studied, which delayed the separation of the airfoil boundary layer and improved the aerodynamic performance of the airfoil. Liu Zhongyuan et al. [19], taking the L1T2 three-element airfoil as the research object, proposed carrying out slot control on the leading-edge slat to slow down the boundary layer separation and reduce its impact on the wake. They focused on the slot location, the flow direction at the slot outlet, and the flow control technology in order to improve the stall characteristics of the lift enhancement device. Although there is a lot of research on the flow control of slats, most of this research focuses on the slot position and slot angle. There is relatively little optimization research on the slot layout scheme, slot width, and slot and chord angles of slotted airfoils of multi-element configuration, which limits the potential of the flow control effect of slats.

In this paper, a grid-slat flow control configuration is proposed to increase lift, reduce drag and suppress flow separation. The quadratic curve parametric modeling method is used for the slot, and the design of the experiment aims to optimize its design. The optimal design parameters of the grid-slat configuration under the given design space and flow conditions are determined through numerical simulation. Moreover, the flow control mechanism of the optimal configuration and the influence and correlation of each design parameter on the lift–drag coefficient are analyzed. It provides a new idea and method for the design of three-dimensional grid multi-element airfoil lift devices in the future.

## **2. Flow Control Optimization Strategy for Leading-Edge Slot of the Airfoil 30P30N Based on Grid Fin**

The grid-slat configuration can improve the flow field level of the slot by increasing grid flow channels. Its basic principle is that, in addition to increasing the flow path of the leading-edge slats, the slats themselves generate additional lift. More importantly, a new additional slot is formed on the basis of the original single channel, which can enhance the flow characteristics of the original slot, improve its disturbance characteristics to the flow of the main wing, and thus enhance the flow control effect. However, since there are many design parameters involved due to the addition of slots in the slats, this paper constructs an aerodynamic optimization design method for grid-slat flow control based on the DOE algorithm.

### *2.1. Design of Grid Fin–Slat Composite Configuration*

The grid-slat configuration is designed on the multi-section airfoil 30P30N of the high-lift device, and the configuration is parameterized to facilitate the transformation of the original optimization of the geometric model of the configuration into the optimization of the configuration design parameters.

Analytic geometric curves can be expressed using explicit functions with few design variables. For example, in the typical conic curve method, section control can be realized by the beginning and end coordinates of curves, as well as by the tangent slope and shape parameters. It has fewer control parameters, which is conducive to the optimal design of the aerodynamic shape. The general equation of a conic curve is: *y* <sup>2</sup> + *ax*<sup>2</sup> + *bxy* + *cx* + *dy* + *e* = 0. The conic curve can be regarded as a conic obtained by performing plane oblique cutting. By changing the angle of the tangent plane, conic curves such as circle, ellipse, parabola and hyperbola can be obtained [20].

The slot at the leading edge is separately described by two curve equations, as shown in Figure 1, curves 1-2-3 and 4-5-6, respectively. Assuming that point 2 is the origin of the curve coordinate system, where the slope of the curve is infinite, the 1-2-3 elliptic equation can be simplified as *y* <sup>2</sup> + *ax*<sup>2</sup> + *bxy* + *cx* = 0. Then, the equation can be obtained according to the coordinates and slope of the curve at point 1 and the tangent equation of the curve at point 3.

Among them, the designer gives the x coordinates of point 1, and the x and y coordinates of point 2. The method to determine the tangent equation at point 3 is shown in the figure above. First, a straight line should be made through the two points denoted as AB. The angle between AB and the horizontal direction is *k*3; then, the distance of *y*<sup>3</sup> is offset upward, where *k*<sup>3</sup> and *y*<sup>3</sup> are also given by the designer. The curve 4-5-6 is determined in the same way with the parameters *x*4, *x*5, *y*<sup>5</sup> and *k*6. In conclusion, the chord length direction is assumed to be the x axis and the direction perpendicular to the chord length is assumed to be the y axis. Considering the calculation time and cost as well as the

configuration changes, the control variables of the slot at the leading edge are shown in Table 1:


**Table 1.** Design variable and value range.

Among them, the designer gives the x coordinates of point 1, and the x and y coordinates of point 2. The method to determine the tangent equation at point 3 is shown in the figure above. First, a straight line should be made through the two points denoted as AB. The angle between AB and the horizontal direction is ଷ; then, the distance of ଷ is offset upward, where ଷ and ଷ are also given by the designer. The curve 4‐5‐6 is determined in the same The position and size of the slot inlet can be adjusted by adjusting *x*<sup>5</sup> and *x*<sup>4</sup> The degree of shrinkage can be adjusted by adjusting *x*<sup>2</sup> and *x*5. Adjusting *y*<sup>2</sup> and *y*<sup>5</sup> can change the longitudinal position of the transition between the contraction and expansion of the joint. The position of the slot outlet position can be adjusted by adjusting *x*<sup>3</sup> and *x*6. The shape of the expansion segment can be controlled by adjusting *k*3, *y*3, *k*<sup>6</sup> and *y*6.

### way with the parameters ସ, ହ, ହ and . In conclusion, the chord length direction is as‐ sumed to be the x axis and the direction perpendicular to the chord length is assumed to be *2.2. Aerodynamic Optimization Model of Grid-Slat Composite Configuration*

the y axis. Considering the calculation time and cost as well as the configuration changes, the control variables of the slot at the leading edge are shown in Table 1: **Table 1.** Design variable and value range. **Design Variable () () () () (**°**)** value range 11 േ 30% 9 േ 15% 11 േ 15% 0–6 40 േ 15% **Design Variable () () () () (**°**)** value range 16 േ 30% 9 േ 15% 11 േ 15% 0–6 35 േ 15% The position and size of the slot inlet can be adjusted by adjusting ହ and ସ The degree of shrinkage can be adjusted by adjusting ଶ and ହ. Adjusting ଶ and ହ can The DOE module in the Isight-integrated platform is adopted, and the numerical optimization of control design parameters is used to replace the optimization of grid-slat configuration geometry. It is a structured research method used to systematically study the relationship between the independent variable and the dependent variable. This method studies the effect of one or more input parameters (factors) on multiple output results (responses) by arranging experiments scientifically and reasonably and replacing comprehensive tests with partial tests. The seven steps of the DOE experiment design are target determination, process analysis, factor screening, fast approach, factorial experiment, regression experiment, and robust design. The optimized design variables are parameters in the parameterized expression of grid-slat curve configuration (*x*1, *x*2, *y*2, *y*3, *k*3, *x*4, *x*5, *y*5, *y*6, *k*6). The objective of optimization is to obtain the maximum lift–drag ratio K. The calculation formula is as follows:

$$\mathbf{C}\_{l} = \frac{L}{0.5\rho L l\_{\infty}^{2} c} \tag{1}$$

$$\mathbf{C}\_d = \frac{D}{0.5\rho l l\_\infty^2 c} \tag{2}$$

$$\mathbf{K} = \mathbf{C}\_l / \mathbf{C}\_d \tag{3}$$

configuration geometry. It is a structured research method used to systematically study

results (responses) by arranging experiments scientifically and reasonably and replacing comprehensive tests with partial tests. The seven steps of the DOE experiment design are target determination, process analysis, factor screening, fast approach, factorial

As shown in Figure 2, the optimization process of grid-slat configuration is as follows: As shown in Figure 2, the optimization process of grid‐slat configuration is as follows:

(1)

<sup>ଶ</sup> (2)

ൌ/ௗ (3)

(1) The optimization module assigns values to ten design variables for the first time according to the constraints, and Catia establishes a two-dimensional parametric model of grid-slat configuration according to the assignment; (1) The optimization module assigns values to ten design variables for the first time ac‐ cording to the constraints, and Catia establishes a two‐dimensional parametric model of grid‐slat configuration according to the assignment;

experiment, regression experiment, and robust design. The optimized design variables are parameters in the parameterized expression of grid‐slat curve configuration (ଵ, ଶ, ଶ, ଷ, ଷ, ସ, ହ, ହ, , ). The objective of optimization is to obtain the maximum

> ൌ .ହఘಮ <sup>మ</sup>

ௗ ൌ 0.5ஶ


**Figure 2.** Optimization flow chart of grid‐slat configuration. **Figure 2.** Optimization flow chart of grid-slat configuration.

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 5 of 16

lift–drag ratio K. The calculation formula is as follows:

### **3. Verification 3. Verification**

#### *3.1. Calculation Method and Verification of Aerodynamic Characteristics of the Multi‐Section Airfoil 30P30N 3.1. Calculation Method and Verification of Aerodynamic Characteristics of the Multi-Section Airfoil 30P30N*

In the present simulation, the flow field was assumed to be described by the 2D Reyn‐ olds‐averaged Navier–Stokes (RANS) equations. The turbulence model adopted was the െ SST model. This turbulence model uses the damping function in the vortex viscosity model. It is more accurate at simulating the aerodynamic characteristics and pressure dis‐ tribution of multi‐section airfoil [21,22]. This has wide applicability and good prediction performance for separation flows with high‐lift configurations. In the present simulation, the flow field was assumed to be described by the 2D Reynolds-averaged Navier–Stokes (RANS) equations. The turbulence model adopted was the *k* − *ω* SST model. This turbulence model uses the damping function in the vortex viscosity model. It is more accurate at simulating the aerodynamic characteristics and pressure distribution of multi-section airfoil [21,22]. This has wide applicability and good prediction performance for separation flows with high-lift configurations.

The െ model formula is given as: The *k* − *ω SST* model formula is given as:

$$\begin{array}{l} \mu\_{l} = \frac{\rho a\_{l}k}{\max(\mu\_{l}\omega\mu\Omega F\_{2})}\\ \frac{\partial(\rho k)}{\partial t} + \mu\_{l}\frac{\partial(\rho k)}{\partial x\_{i}} = P\_{k} - \beta\_{k}\rho\mathbf{k}\omega + \frac{\partial}{\partial x\_{i}}\left[\left(\mu\_{l} + \frac{\mu\_{l}}{\sigma\_{k}}\right)\frac{\partial \mathbf{k}}{\partial x\_{i}}\right] \\ \frac{\partial(\rho \omega)}{\partial t} + \mu\_{l}\frac{\partial(\rho \omega)}{\partial x\_{i}} = C\_{\omega i}P\_{\omega} - \beta\_{\omega}\rho\omega^{2} + \frac{\partial}{\partial x\_{i}}\left[\left(\mu\_{l} + \frac{\mu\_{l}}{\sigma\_{\omega}}\right)\frac{\partial \omega}{\partial x\_{i}}\right] + 2\rho(1 - F\_{1})\frac{1}{\tau\_{\omega 2}}\frac{1}{\omega}\frac{\partial \mathbf{k}}{\partial \overline{x\_{i}}}\frac{\partial \omega}{\partial \overline{x\_{i}}} \end{array} \tag{4}$$

The modified *k* − *ω* SST turbulence model constructs the model coefficient in the vortex viscosity coefficient of the equation as the following equation:

$$\alpha^\* = \alpha^\*\_{\infty} \left[ \frac{\alpha^\*\_0 + \frac{\text{Re}\_l}{\mathcal{R}\_k}}{1 + \frac{\text{Re}\_l}{\mathcal{R}\_k}} \right] \tag{5}$$

where Re*<sup>t</sup>* = *ρκ*/*µω* is the turbulent Reynolds number.

A numerical simulation of the 30P30N airfoil is taken as an example to verify the numerical simulation method of the two-dimensional lifting device adopted in this paper. The multi-section airfoil 30P30N of the high-lift device is widely used by CFD workers [23,24]. The deflection angle of the leading-edge slat and trailing-edge flap of this airfoil is 30◦ , the

slot width of the leading-edge slat is 2.95%, and the overhang is −2.5%. The trailing-edge flap slot width is 1.27%, and the extension is 0.25% [23], as shown in Figure 3. is 1.27%, and the extension is 0.25% [23], as shown in Figure 3.

1

1 1 2 1

2

 

(5)

, the slot width

(4)

 

 

 

*t*

1 1 2

*a k a F*

max( ; )

 

> 

where Re / *<sup>t</sup>*

**Figure 3.** Geometry profile of 30P30N.

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 6 of 16

 

 

*t*

*<sup>k</sup> u CP F t x x x xx*

vortex viscosity coefficient of the equation as the following equation:

 

  *t*

The modified െ SST turbulence model constructs the model coefficient in the

0

A numerical simulation of the 30P30N airfoil is taken as an example to verify the numer‐ ical simulation method of the two‐dimensional lifting device adopted in this paper. The multi‐ section airfoil 30P30N of the high‐lift device is widely used by CFD workers [23,24]. The de‐

of the leading‐edge slat is 2.95%, and the overhang is −2.5%. The trailing‐edge flap slot width

Re <sup>1</sup>

*R*

 

Re

*R*

*t k t k*

*i i i ii*

 

is the turbulent Reynolds number.

flection angle of the leading‐edge slat and trailing‐edge flap of this airfoil is 30°

2

*i i ki*

*i kk l*

 

*i l*

 

 

*k k <sup>k</sup> u Pk tx x x*

**Figure 3.** Geometry profile of 30P30N. The comparative test data are the wind tunnel test results of McDonnell Douglas Aerospace and NASA Langley [25,26]. According to actual wind tunnel test conditions, the attack angles of the incoming flow are െ0.17<sup>∘</sup> ∼ 21<sup>∘</sup> . The Reynolds number is 9 ൈ 10. The Mach number is 0.2. All the calculations provided here were obtained using the cou‐ pled solver of Fluent code. A second‐order upwind difference scheme was chosen for mo‐ mentum and turbulent kinetic energy, respectively. The incoming flow pressure is the standard atmospheric pressure. The mesh is generated as shown in Figure 4. For the com‐ ponent grid, we used a hybrid mesh with 40,792 cells. Wall spacing of the first layer is 1.71 ൈ 10ିଷ, which satisfies 1. The comparative test data are the wind tunnel test results of McDonnell Douglas Aerospace and NASA Langley [25,26]. According to actual wind tunnel test conditions, the attack angles of the incoming flow are −0.17◦ ∼ 21◦ . The Reynolds number is <sup>9</sup> <sup>×</sup> <sup>10</sup><sup>6</sup> . The Mach number is 0.2. All the calculations provided here were obtained using the coupled solver of Fluent code. A second-order upwind difference scheme was chosen for momentum and turbulent kinetic energy, respectively. The incoming flow pressure is the standard atmospheric pressure. The mesh is generated as shown in Figure 4. For the component grid, we used a hybrid mesh with 40,792 cells. Wall spacing of the first layer is 1.71 <sup>×</sup> <sup>10</sup>−<sup>3</sup> , which satisfies *y*+ ≤ 1. *Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 7 of 16

**Figure 4.** Computational mesh. **Figure 4.** Computational mesh.

In the process of calculation, we monitored the convergence of residuals on the one hand, and found that the energy item decreases to 10−<sup>9</sup> and other residuals except the energy item decrease to 10−<sup>4</sup> . On the other hand, the calculation results cl or cd basically do not change with the increase in iteration steps. In Figure 5a, it seems that good agreement is achieved by comparing the calculated lift coefficient of the airfoil with the wind tunnel test results and that the deviation is large only near the stall angle of attack. As shown in Figure 5c,d, the surface pressure coefficients distribution of the main wing and flap coincide very well, and the pressure coefficients of the slats are basically consistent with the experimental data.

### *3.2. Grid Independence Verification*

*3.2. Grid Independence Verification*

(**a**) (**b**) In order to ensure the accuracy of numerical calculation and reduce the calculation cost, the multi-section airfoil 30P30N is taken as an example to verify the independence of the computational grid. Three-grid generation schemes are used for simulation calculation. Table 2 shows the grid convergence verification results at an attack of 8 ◦ . The lift coefficient *Cl* converges with the increase in the number of grid cells. The drag coefficients, shown in Figure 6a, are basically the same using 2- and 3-grid generation schemes at different angles of attack. As shown in Figure 6c,d, the surface pressure coefficient distribution of the main wing and flap coincide very well. Using 2- and 3-grid generation schemes, the calculation results are similar, which meets the requirements of grid independence. Therefore, the grid division scheme of Case2 is adopted for subsequent calculations of cost.

(**c**) (**d**)

In order to ensure the accuracy of numerical calculation and reduce the calculation cost, the multi‐section airfoil 30P30N is taken as an example to verify the independence of the computational grid. Three‐grid generation schemes are used for simulation calcula‐

coefficient converges with the increase in the number of grid cells. The drag coeffi‐ cients, shown in Figure 6a, are basically the same using 2‐ and 3‐grid generation schemes at different angles of attack. As shown in Figure 6c,d, the surface pressure coefficient dis‐ tribution of the main wing and flap coincide very well. Using 2‐ and 3‐grid generation

tion. Table 2 shows the grid convergence verification results at an attack of 8°

. The lift

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 7 of 16

**Figure 4.** Computational mesh.

**Figure 5.** Verification results of 30P30N airfoil. (**a**) Lift Coefficient; (**b**) Slat Pressure Distribution [18,19]; (**c**) Main Pressure Distribution; (**d**) Flap Pressure Distribution. **Figure 5.** Verification results of 30P30N airfoil. (**a**) Lift Coefficient; (**b**) Slat Pressure Distribution [18,19]; (**c**) Main Pressure Distribution; (**d**) Flap Pressure Distribution. Case3 61,923 3.144985

**Figure 6.** Grid independence verification. (**a**) Drag Coefficient; (**b**) Slat Pressure Distribution; (**c**) Main Pressure Distribution; (**d**) Flap Pressure Distribution. **Figure 6.** Grid independence verification. (**a**) Drag Coefficient; (**b**) Slat Pressure Distribution; (**c**) Main Pressure Distribution; (**d**) Flap Pressure Distribution.

with the original airfoil at an angle of attack of 16.3°

Through optimization design, the values of design variables are shown in Table 3, and the configuration at this time is shown in Figure 7a. The grid‐slat configuration plays

.

comparison of the lift coefficient, drag coefficient and lift–drag ratio between the original airfoil and the new configuration at different angles of attack. The simulation and optimi‐ zation results show that the lift coefficient increases by 3.3%, the drag coefficient decreases by 12.7%, and the lift–drag ratio increases by 18.4% of the optimal configuration compared

**4. Optimization Example and Analysis**


**Table 2.** Grid independence verification.

### **4. Optimization Example and Analysis**

Through optimization design, the values of design variables are shown in Table 3, and the configuration at this time is shown in Figure 7a. The grid-slat configuration plays a role in increasing lift and reducing drag through calculations. Figure 7b–d show the comparison of the lift coefficient, drag coefficient and lift–drag ratio between the original airfoil and the new configuration at different angles of attack. The simulation and optimization results show that the lift coefficient increases by 3.3%, the drag coefficient decreases by 12.7%, and the lift–drag ratio increases by 18.4% of the optimal configuration compared with the original airfoil at an angle of attack of 16.3◦ . *Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 9 of 16

> **Table 3.** Value of design variable. **Table 3.** Value of design variable.


**Figure 7.** Numerical simulation comparison between original airfoil and optimal configuration. (**a**) Lift coefficient of slat; (**b**) Total lift coefficient; (**c**) Drag Coefficient; (**d**) Lift‐drag ratio. **Figure 7.** Numerical simulation comparison between original airfoil and optimal configuration. (**a**) Lift coefficient of slat; (**b**) Total lift coefficient; (**c**) Drag Coefficient; (**d**) Lift-drag ratio.

In order to further clarify the mechanism of increasing lift and reducing drag of grid‐ slat configuration, Figure 8 describes the total lift coefficients of the original airfoil and the optimized airfoil at different angles of attack. The lift coefficients, extracted separately

also improves the lift of the main wing.

In order to further clarify the mechanism of increasing lift and reducing drag of gridslat configuration, Figure 8 describes the total lift coefficients of the original airfoil and the optimized airfoil at different angles of attack. The lift coefficients, extracted separately from the slats, show that the configuration not only improves the lift of the slat itself but also improves the lift of the main wing. *Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 10 of 16

**Figure 8.** Comparisons of lift coefficients between original airfoil and optimal configuration. (**a**) Lift coefficient of slat; (**b**) Total lift coefficient. **Figure 8.** Comparisons of lift coefficients between original airfoil and optimal configuration. (**a**) Lift coefficient of slat; (**b**) Total lift coefficient.

The simulation results of the original airfoil and the optimized configuration at four typ‐ ical angles of attack ( ൌ 8° , 16.3° , 19° , 22.36° ) are compared and analyzed. The development process of wake flow is shown in Figure 9. For the original airfoil at an angle of attack of 22.36° , the wake flow area generated by the separation bubble is strongly mixed with the slot jet and the boundary layer of the main wing section, so that the rear section of the main airfoil and trailing‐edge flap are in the mixed wake flow area. After optimization, the separation of grid‐slat configuration is weakened, and the wake area is greatly reduced. The simulation results of the original airfoil and the optimized configuration at four typical angles of attack (*α* = 8 ◦ , 16.3◦ , 19◦ , 22.36◦ ) are compared and analyzed. The development process of wake flow is shown in Figure 9. For the original airfoil at an angle of attack of 22.36◦ , the wake flow area generated by the separation bubble is strongly mixed with the slot jet and the boundary layer of the main wing section, so that the rear section of the main airfoil and trailing-edge flap are in the mixed wake flow area. After optimization, the separation of grid-slat configuration is weakened, and the wake area is greatly reduced.

Superimposed contour plots of velocity magnitude and instantaneous streamlines at different angles are shown in Figure 10. Although the jet velocity has not yet reached the peak value, the velocities near the upper trailing edge are improved to varying degrees. The downstream flow-separated area near the upper trailing edge becomes a flow-attached area.

Figure 11 shows the surface pressure coefficient Cp of the original airfoil and the optimized layout under two operating conditions. For the convenience of differentiation, blue represents the original airfoil and red represents the optimal grid-slat configuration. The negative pressure value in the figure represents the suction on the upper surface of the wing. The positive pressure value represents the pressure on the lower surface of the wing. The pressure surface curve is located below the suction surface curve, indicating that the pressure on the lower surface is greater than that on the upper surface, which means the lift of the airfoil is positive [27,28]. The area enclosed by the pressure coefficient curve can roughly reflect the lift coefficient.

 (**a**) ൌ 8° As is shown in Figure 11a, at an angle of attack of 19◦ , each suction peak of optimized configuration for each wing is slightly higher than that of the original wing, while the pressure surface curve is basically the same, so the optimal grid-slat configuration has a higher lift. As can be seen in Figure 11b, the suction peak of the original airfoil decreases sharply, and the area enclosed by the pressure coefficient curve decreases. However, the pressure coefficient curve of the optimized configuration is significantly different from that of the original airfoil, indicating that stall status had not been reached. Therefore, it can be concluded that the optimized configuration improves the aerodynamic characteristics at a high angle of attack and additionally increases the stall angle of attack.

(**b**) ൌ 16.3°

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 10 of 16

Lift coefficient of slat; (**b**) Total lift coefficient.

, 16.3°

, 19°

grid‐slat configuration is weakened, and the wake area is greatly reduced.

, 22.36°

ical angles of attack ( ൌ 8°

22.36°

) are compared and analyzed. The development

(**a**) (**b**) **Figure 8.** Comparisons of lift coefficients between original airfoil and optimal configuration. (**a**)

The simulation results of the original airfoil and the optimized configuration at four typ‐

, the wake flow area generated by the separation bubble is strongly mixed with the slot jet and the boundary layer of the main wing section, so that the rear section of the main airfoil and trailing‐edge flap are in the mixed wake flow area. After optimization, the separation of

process of wake flow is shown in Figure 9. For the original airfoil at an angle of attack of

**Figure 9.** Comparison of wake flow development of the original airfoil and optimal grid‐slat con‐ figuration. **Figure 9.** Comparison of wake flow development of the original airfoil and optimal grid-slat configuration.

Superimposed contour plots of velocity magnitude and instantaneous streamlines at different angles are shown in Figure 10. Although the jet velocity has not yet reached the peak value, the velocities near the upper trailing edge are improved to varying degrees.

(**a**) ൌ 8°

tached area.

figuration.

tached area.

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 11 of 16

(**c**) ൌ 19°

(**d**) ൌ 22.36° **Figure 9.** Comparison of wake flow development of the original airfoil and optimal grid‐slat con‐

Superimposed contour plots of velocity magnitude and instantaneous streamlines at different angles are shown in Figure 10. Although the jet velocity has not yet reached the peak value, the velocities near the upper trailing edge are improved to varying degrees. The downstream flow‐separated area near the upper trailing edge becomes a flow‐at‐

(**d**) ൌ 22.36°

**Figure 10.** Contour comparison of the original airfoil and optimal grid‐slat configuration at different attack angles. **Figure 10.** Contour comparison of the original airfoil and optimal grid-slat configuration at different attack angles.

Figure 11 shows the surface pressure coefficient Cp of the original airfoil and the optimized layout under two operating conditions. For the convenience of differentiation, blue represents the original airfoil and red represents the optimal grid‐slat configuration.

the wing. The positive pressure value represents the pressure on the lower surface of the wing. The pressure surface curve is located below the suction surface curve, indicating that the pressure on the lower surface is greater than that on the upper surface, which

*Appl. Sci.* **2023**, *13*, x FOR PEER REVIEW 13 of 16

As is shown in Figure 11a, at an angle of attack of 19°

curve can roughly reflect the lift coefficient.

**Figure 11.** Comparison of pressure coefficient between the original airfoil and optimal grid‐slat con‐ figuration at different attack angle. **Figure 11.** Comparison of pressure coefficient between the original airfoil and optimal grid-slat configuration at different attack angle.

means the lift of the airfoil is positive [27,28]. The area enclosed by the pressure coefficient

configuration for each wing is slightly higher than that of the original wing, while the pressure surface curve is basically the same, so the optimal grid‐slat configuration has a higher lift. As can be seen in Figure 11b, the suction peak of the original airfoil decreases sharply, and the area enclosed by the pressure coefficient curve decreases. However, the pressure coefficient curve of the optimized configuration is significantly different from that of the original airfoil, indicating that stall status had not been reached. Therefore, it can be concluded that the optimized configuration improves the aerodynamic character‐

istics at a high angle of attack and additionally increases the stall angle of attack.

, each suction peak of optimized

#### **5. Influence Analysis of Design Parameters 5. Influence Analysis of Design Parameters**

In order to further clarify the design method of grid‐slat configuration, the influence of each parameter on lift–drag coefficient and lift–drag ratio is analyzed according to the calculation results under the different design parameters calculated in the optimization process. First, the Pearson correlation coefficient [29] is selected to calculate and analyze the correlation between design parameters, lift–drag coefficient and lift–drag ratio, as shown in Table 4. The results show that the five design parameters are significantly re‐ lated to the lift–drag ratio , where ଵ, ଷ, are positively correlated with , ଷ, are negatively correlated with . Pearson's formula [30] is: In order to further clarify the design method of grid-slat configuration, the influence of each parameter on lift–drag coefficient and lift–drag ratio is analyzed according to the calculation results under the different design parameters calculated in the optimization process. First, the Pearson correlation coefficient [29] is selected to calculate and analyze the correlation between design parameters, lift–drag coefficient and lift–drag ratio, as shown in Table 4. The results show that the five design parameters are significantly related to the lift–drag ratio *K*, where *x*1, *k*3, *k*<sup>6</sup> are positively correlated with *K*, *y*3, *y*<sup>6</sup> are negatively correlated with *K*.


, ൌ covሺ, ሻ ൌ ሺሻ െ ሺሻሺሻ **Table 4.** The correlation analysis results between design parameters and lift–drag coefficient.

\*\*. At 0.01 level (2-tailed), the correlation is significant. \*. At 0.05 level (2-tailed), the correlation is significant.

Pearson's formula [30] is:

$$\rho\_{X,Y} = \frac{\text{cov}(X,Y)}{\sigma\_X \sigma\_Y} = \frac{E(XY) - E(X)E(Y)}{\sqrt{E(X^2) - (E(X))^2} - \sqrt{E(Y^2) - (E(Y))^2}} \tag{6}$$

where cov(*X*,*Y*) is the covariance of *X* and *Y*, and *σX, σ<sup>Y</sup>* is the standard deviation of *X* and *Y*, respectively.

This shows that within the given range of variables, the shape of the slot changes with the change in five design parameters, which will greatly change the flow field. The change of *x*<sup>1</sup> results in the change of seam entry position, while the change of the other four design parameters *y*3, *y*6, *k*<sup>3</sup> and *k*<sup>6</sup> results in the change in seam shape. The design parameter *x*<sup>1</sup> is positively correlated, indicating that the lift force is better when the inlet position of the slot is closer to the right. The design parameters *k*<sup>3</sup> and *k*<sup>6</sup> are positively correlated, indicating that as the expansion degree of the slot increases, the acceleration degree of the slot to the flow is greater, and the lift force is better.

### **6. Conclusions**

In this article, a flow control method of grid-slat configuration was proposed, and its optimal design was configured using the DOE algorithm. The optimal configuration was obtained and verified through numerical simulation, and we analyzed the flow control mechanism and the influence of design parameters. Two main conclusions were drawn from this study.

(1) By slitting the leading-edge slats, the lift–drag characteristics of the multi-element configuration can be effectively improved, which shows the feasibility and effectiveness of conducting flow control through the leading-edge slats. The simulation and optimization results show that the lift coefficient increases by 3.3%, the drag coefficient decreases by 12.7%, and the lift–drag ratio increases by 18.4% of the optimal configuration compared with the original airfoil at an angle of attack of 16.3◦ . The grid-slat can effectively improve the flow field level, and the flow always remains attached.

(2) The analysis of the correlation between the design parameters and the configuration of aerodynamic performance shows that the slot design parameters are highly sensitive to the final flow control effect, which also indicates that the flow control potential of the slats is very good.

In this paper, a flow control model and design optimization method based on slat slitting have only been constructed for a single grid slot. In future work, we expect to build a more effective slot scheme with the possibility of three-dimensional applications.

**Author Contributions:** Conceptualization, B.M.; Methodology, J.Y.; Writing—original draft, J.Y.; Writing—review & editing, B.M. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by [National Natural Science Foundation of China] grant number [12202363] and [Natural Science Basic Research Program of Shaanxi] grant number [2021JQ-084].

**Institutional Review Board Statement:** Ethical review and approval were waived for this study because this article does not involve human or animal research.

**Informed Consent Statement:** Informed consent was obtained from all subjects involved in the study.

**Data Availability Statement:** The data is unavailable due to privacy.

**Acknowledgments:** The authors would like to acknowledge the support of National Natural Science Foundation of China (Grant No. 12202363) and the support of Natural Science Basic Research Program of Shaanxi (Program No. 2021JQ-084).

**Conflicts of Interest:** The authors declared no potential conflict of interest with respect to the research, authorship, and/or publication of this article.

### **References**


**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
