The aerodynamic optimization of a Gurney flap attached to the NACA 0012 airfoil in rarefied gas flow is studied in this study. In order to investigate the effects of Gurney flap geometry (height, mounting angle, and mounting location) and angle of attack on the aerodynamic characteristics of the airfoil in the rarefied gas flow, the flow fields around NACA 0012 airfoils with Gurney flap attached are simulated by using the direct simulation Monte Carlo method. Furthermore, the multi-objective aerodynamic optimization of the airfoil with Gurney flap is performed at the single design point, and the Pareto front is obtained. Based on the benchmark configuration (obtained from the Pareto front), aerodynamic optimization by controlling the mounting angle of the Gurney flap and the angle of attack is investigated over a wide range of Mach numbers (ranging from 2 to 4) and Knudsen numbers (ranging from 0.024 to 0.24) by constructing an artificial neural network model. The results show that through optimization, the lift-to-drag ratio of the airfoil is increased by up to 29.25% under the combination of different Mach numbers and Knudsen numbers, which significantly improves the aerodynamic characteristics.

## I. INTRODUCTION

The micro air vehicle (MAV) has attracted substantial attention in contemporary times due to its potential to perform various missions and low production cost compared with full-size aircraft. The aerodynamic attributes of MAVs have constituted a focal point for numerous researchers.^{1} Notably, prevalent MAVs are scaled-down renditions of traditional aircraft with chords spanning mere centimeters, once the characteristic scale is sufficiently small, the non-equilibrium, rarefied flow phenomenon becomes dominant.^{2} During the preceding decades, researchers have studied the aerodynamic characteristics of small-size wings by wind tunnel test and CFD method^{3–5} and obtained many valuable conclusions. Recently, researchers began to investigate the flow field simulation surrounding airfoils operating within rarefied gas flow conditions by using the Direct Simulation Monte Carlo (DSMC) method.^{2,6–8} They illuminated the direct correspondence between an augmented Knudsen number and an augmented drag coefficient, ultimately resulting in the deterioration of aerodynamic performance. Hence, the necessity to enhance the aerodynamic performance of airfoils in rarefied gas flow has become palpable.

As a simple and low-cost method, the Gurney flap is often used to enhance the aerodynamic attributes of airfoils.^{9,10} Amini *et al.*^{11} investigated the impacts of AOA and Gurney flap height on the aerodynamic performance of the NACA 0012 airfoil within the domain of rarefied conditions. They considered three distinct Kn numbers (0.026, 0.1, 0.26), three variable Gurney flap heights (2%, 5%, 7%), and a wide range of AOA. The findings underscored the potency of 5% and 7% Gurney flaps in enhancing the airfoil’s aerodynamic behavior. Liebeck^{12} evidenced the consequential enhancement in lift coefficient due to the Gurney flap but also revealed a proportional rise in drag when flap heights surpassed 2% of chord length. Li *et al.*^{13–15} embarked upon a detailed examination of the influence of Gurney flaps on the aerodynamic characteristics of the NACA 0012 airfoil at Ma = 0.7. They found that the Gurney flap significantly improved the aerodynamic performance of the airfoil. The results notably disclosed a remarkable 41% increase in maximum lift-to-drag ratio attained with a Gurney flap of 2% chord length. Troolin *et al.*^{16} delved into the dynamic evolution of the wake structure of the NACA 0012 airfoil, both with and without Gurney flap appendages. Employing time-resolved Particle Image Velocimetry (PIV), they analyzed the underpinning mechanisms behind the lift augmentation attributed to Gurney flaps. Naderi *et al.*^{17} studied the effect of Gurney flap height, mounting angle, and mounting location on the aerodynamic performance of NACA 0012 pitching airfoil under low Reynolds number regimes. Their outcomes affirmed the heightened aerodynamic efficacy engendered by Gurney flaps, and the flap height is the main factor affecting the flow field. In addition, Gurney flaps are also extensively used in diverse applications, such as dynamic stall,^{18} energy extraction,^{19,20} and rotor blade load control.^{21,22}

In addition to exploring the effect of various factors on aerodynamic characteristics, aerodynamic design optimization has also attracted significant scholarly attention.^{23,24} Amini *et al.*^{25} adopted the adjoint shape optimization method to minimize the drag coefficient of airfoil attached the Gurney flap in different Mach numbers, and the Gurney flap heights were 1%c and 2%c. The numerical simulation results show that the optimization process markedly diminishes the drag coefficient, while the corresponding influence on the lift coefficient remains marginal. In addition to the airfoil shape optimization, the improvement of aerodynamic performance can be achieved through optimization of Gurney flap configurations. Joo *et al.*^{26} studied the design optimization of a combined nose droop and Gurney flap, aiming to improve dynamic stall and post-stall aerodynamic behavior within rotor airfoils. The optimized combination configuration was obtained that can substantially improve aerodynamic performance. In addition, Yee *et al.*^{27} found the existence of optimal heights for Gurney flaps that concurrently minimize drag and pitching moment, while maximizing both the lift coefficient and the lift-to-drag ratio. Previous scholarly investigations suggest that aerodynamic optimization is feasible by optimizing the configuration of Gurney flaps.

The main aim of this study is to improve the aerodynamic performance of NACA 0012 airfoil with an attached Gurney flap over a wide range of Knudsen and Mach numbers by controlling the gurney flap mounting angle and angle of attack. To accomplish this aim, we first simulate the rarefied gas flow around the NACA 0012 airfoil with Gurney flaps by using the DSMC method and investigate the effects of geometrical and freestream parameters on aerodynamic characteristics under rarefied gas flow conditions. Then, to determine the benchmark airfoil for control, we use the ANN surrogate model and NSGA-II algorithm to optimize Gurney flaps. Finally, we construct an ANN surrogate model for the fast prediction of aerodynamic characteristics and the optimization through the precise control of both the Gurney flap mounting angle and angle of attack. The remainder of this paper is organized as follows: Sec. II introduces the optimization methods in detail. Then, in Sec. III, the physical problem is shown. In Sec. IV, the effect law of each parameter on aerodynamic characteristics is studied. In Sec. V, the aerodynamic optimization of Gurney flaps at single-point design and a wide range of flight conditions is introduced. Finally, concluding remarks are presented in Sec. VI.

## II. METHOD INTRODUCTION

In this section, we will introduce our comprehensive framework for aerodynamic optimization, which is shown in Fig. 1. First, the sampling method in the design space is described. Then, the DSMC method is introduced. Finally, we expound upon the surrogate model and numerical optimization algorithm.

### A. Design of experiments (DoE)

^{28}to obtain the initial sample points within the design space. A good experimental design is supposed to be “space-filling” so that the sample points are distributed as evenly as possible in the design space to obtain an accurate surrogate model. A Latin hypercube design (LHD) is a modern popular DoE method, which has been widely used in various fields. An LHD is expressed as an

*n*×

*m*matrix, in which each column is a random permutation of {1, 2, …,

*n*}. The optimal Latin hypercube design is to add an optimization algorithm based on the Latin hypercube to amplify its space filling performance. The optimal criterion we used is the maximin distance design proposed by Johnson

*et al.*,

^{29}aimed at the maximization of the minimal inter-site distance,

*x*

_{i}and

*x*

_{j},

*t*is 2 in this study.

### B. DSMC method

Knudsen number (Kn) is a dimensionless parameter to define the rarefaction degree of flow. The rarefaction regimes can be generally classified into three categories: continuum flow (Kn $<0.1$), transition flow (0.1 < Kn $<10$), and free molecular flow (Kn $>10$).^{30} There are diverse methodologies to simulate different regimes. In the case of continuum flow, the N–S equation is suitably applicable, while the Boltzmann equation can be employed with the Kn number in the full range, particularly in rarefied flow. DSMC is an efficient method to solve the Boltzmann equation, which was invented by Bird,^{31,32} and has evolved into the favored numerical approach for simulating rarefied fluids after decades of development.

This method solves the Boltzmann equation by using thousands of simulation particles, which represent plenty of real gas molecules and have physical parameters, such as position, velocity, and molecular internal energy. The average molecular diameter in rarefied flow is much smaller than the average molecular space, so molecular motion is allowed to be decoupled from molecular collision.^{33} The DSMC method mainly includes four calculation processes: moving the particles, indexing particles, simulating collisions, and sampling aimed at extracting macroscopic properties of the flow field. In each time step, the particles move in the grid, collide with other particles or walls, and redistribute energy through kinetic theory; thus, their position and velocity have changed. Finally, after reaching a certain number of sampling steps, the average macroscopic parameters of the flowfield are calculated.

The range of the Kn number investigated within this research is set as 0.024–0.24, so the DSMC stands as the fitting choice for conducting aerodynamic computations under this condition. The DSMC algorithm used is our in-house parallel DSMC solver,^{34–36} and its accuracy has been validated.^{37} Particle collisions are simulated utilizing the variable hard sphere (VHS) model in conjunction with the timeless collision (NTC) scheme. In the context of the translation-rotation (T-R) energy exchange process, the generalized Larsen–Borgnakke model^{38} is harnessed to facilitate the redistribution of rotational energy, with a designated rotational collision number of *Zr* = 10.

### C. Artificial neural networks

^{39}is employed for the establishment of a multilayer perceptron (MLP) neural network, which serves as a feedforward ANN. The MLP configuration is characterized by an initial layer with

*R*inputs, followed by

*M*− 1 hidden layers, with each layer

*m*being comprised of a certain count of neurons denoted as

*S*

^{m}and an output layer equipped with

*S*

^{M}outputs. Notably, in the context of an MLP, every layer is fully interconnected with the subsequent layer, which means that the output of each neuron within layer

*m*is input to every neuron situated in layer

*m*+ 1 (and in the initial layer inputs possess connections to all neurons). The computation of the output $xjm+1$ for each neuron

*j*in the (

*m*+ 1) layer can be facilitated through the following mathematical formulation:

*w*

_{ij}and then add the bias weight

*b*

_{j}, and the cumulative outcome is the output $njm+1$. The function

*f*(

*n*), designated as the transfer function or activation function, in this research, equals tanh(

*n*) for all neurons residing in the hidden layers and equals

*n*in the output layer.

Upon the successful construction of the neural network, we can begin to train the matched input–output data. This training process entails the iterative adjustment of the network’s learnable parameters, orchestrated to progressively minimize the error *E* existing between the model’s predictions and the actual provided outputs. Within the Pyrenn framework, the calculation of this error adopts the mean square error criterion, and the optimization is facilitated through the employment of the Levenberg–Marquardt algorithm,^{40} a second-order quasi-Newton optimization technique. The training procedure proceeds until either the predefined threshold of error or the maximum stipulated number of iterations is attained, culminating in the model can serve as an original surrogate model representing the system.

### D. Numerical optimization algorithm

The Differential Evolution (DE) algorithm is a heuristic search algorithm based on group differences. The implementation of the DE algorithm mainly includes four steps: population initialization, mutation, crossover, and selection. Price and Storn^{41} proposed various differential strategies of the DE algorithm and described as DE/x/y/z, where x and y refer to the way the basis vectors are chosen and the number of differential vectors in the mutation operation, respectively, and z refers to the way of crossover operation. The convergence speed and solution accuracy of the algorithm are influenced by the differential strategy; the differential strategy selected in this paper is DE/best/1/bin.

*f*

_{i}(

*X*) is the objective function,

*X*is the decision vector, and

*g*

_{j}(

*X*) ≥ 0 and

*h*

_{k}(

*X*) = 0 represent inequality constraints and equality constraints, respectively. Given the multidimensional nature of the objective function, attaining an unequivocal optimal solution to the multi-objective optimization problem is a challenge of considerable magnitude. As a result, the focal pursuit within the domain of multi-objective optimization primarily centers upon the identification of the Pareto optimal solution set, which is recognized more commonly as the Pareto front.

The NSGA-II applied in this paper represents a prominent member of the multi-objective genetic algorithm repertoire, which is an improved NSGA algorithm proposed by Deb and Pratap for the shortcomings of NSGA.^{42} The core operational framework of this algorithm is founded upon the application of non-dominated sorting to preserve individuals demonstrating advantageous attributes. Notably, the inclusion of rapid non-dominated sorting serves the dual purpose of streamlining algorithmic complexity. In addition, the algorithm also uses the elite strategy, which means both the parent population and the offspring population participate in the competition to produce the next generation population, thereby accelerating the population optimization speed.

## III. PROBLEM DEFINITION

In this study, first, the precision of the DSMC method developed to calculate the airfoil aerodynamic attributes is verified through comparison with the results of references. Next, univariate analysis is conducted to explore the effects and underlying principles governing the influence exerted by the Gurney flap on the aerodynamic characteristics of the NACA 0012 airfoil in the rarefied gas flow. The outcomes of univariate analysis contribute to the comprehension and analysis of subsequent optimal results and control strategies. Then, a multi-objective optimization method is used at the design point to procure an optimal combination of the lift-to-drag ratio and pitch moment coefficients to determine the benchmark Gurney flap configuration. Finally, the control optimization of the mounting angle and angle of attack for the Gurney flap over a wide range of rarefied freestream conditions is investigated based on the benchmark configuration.

In this study, the chord length (c) of the airfoil was selected to be 4 cm^{6}. The Gurney flap manifests as a planar appendage characterized by a height ranging from 1% to 5% of the chord length that is affixed proximate to or upon the trailing edge of the wing, located either perpendicularly or at an angle with respect to the airfoil’s string line. The specifics of Gurney flap implementation upon the NACA 0012 airfoil are referenceḍ,^{15} as shown in Fig. 2(a). The pivotal design parameters include three geometric parameters: Gurney flap height (h), mounting position (s), and mounting angle (Φ), and three inflow parameters: angle of attack (*α*), Ma, and Kn. One thing to note is that the Gurney flap height and mounting position are maintained constant, but the mounting angle and angle of attack can be manipulated. These parameters are selected due to experiments by Li *et al.*^{13,15} affirming their pronounced influence on the airfoil’s aerodynamic attributes, including but not limited to the lift-drag ratio and pitching moment coefficient; the impact of the Knudsen number on aerodynamic characteristics is underscored in Ref. 11.

As shown in Fig. 2(b), the computational domain consists of a semi-ellipse with long axes of 6c and short axes of 5c upstream of the airfoil as well as a downstream rectangular expanse spanning 2c. The inlet boundary is established as a freestream characterized by a velocity of 594 m/s (Mach = 2) and a temperature of 219.6 K. Concurrently, the temperature of the airfoil surface is set at 300 K. In the DSMC simulation, the air is treated as a mixture composed of 71% nitrogen and 21% oxygen. By the independence study, in the benchmark case, the number of nodes on the airfoil surface is 360, the grid number is 138 000, the time step is selected as 5 × 10^{−7} s, and the sampling times are 10 000. Similar to Ref. 11, the particle number is selected as 12. With this combination of parameters, the errors of lift, drag, and pitch moment coefficients obtained are all less than 2.5% as compared with strict benchmark cases, which satisfies the requirement of DSMC computational accuracy. The details of the time step and mesh independence study are listed in the Appendix. Simulations were performed on the center for high-performance computing of SJTU using 40 cores for ∼6 h,^{43} and the parameters of the benchmark case are set as given in Table I.

Parameter type . | Parameter value . |
---|---|

Freestream Knudsen number | 0.0245 |

Freestream pressure (pa) | 5.2209 |

Freestream Mach number | 2 |

Freestream velocity(m/s) | 594 |

Wall temperature (K) | 300 |

Time step (s) | 5 × 10^{−7} |

Sampling steps | 10 000 |

Particle per cells | 12 |

Parameter type . | Parameter value . |
---|---|

Freestream Knudsen number | 0.0245 |

Freestream pressure (pa) | 5.2209 |

Freestream Mach number | 2 |

Freestream velocity(m/s) | 594 |

Wall temperature (K) | 300 |

Time step (s) | 5 × 10^{−7} |

Sampling steps | 10 000 |

Particle per cells | 12 |

In Table II, the drag coefficients computed by DSMC for NACA 0012 airfoil at Ma = 2, Kn = 0.026, and three distinct AOA conditions are compared with the reference data.^{11} The maximum relative error between the drag coefficients computed by the DSMC we used and the results given^{11} is less than 2%, and the error is probably due to differences in the computational domain and the number of grids. In Fig. 3, to further verify the accuracy of the DSMC method we employed to simulate the flow field around the NACA 0012 airfoil, the density field (*ρ*/*ρ*_{∞}) simulated by DSMC was compared with the results obtained experimentally by Allegre and Lengrand^{6} under the conditions characterized by Ma = 2, Kn = 0.026, and AOA = 0. A concordant concurrence between the DSMC simulation result and the experimental result is conspicuously evident. Upon example verification, it is shown that the numerical computational method we used is reliable for the problem studied in this paper.

. | Present result . | Reference 11 . | . |
---|---|---|---|

AOA . | C_{D}
. | C_{D}
. | Relative error(%) . |

0° | 0.4008 | 0.3998 | 0.25 |

10° | 0.4526 | 0.4544 | 0.40 |

20° | 0.6058 | 0.6158 | 1.62 |

. | Present result . | Reference 11 . | . |
---|---|---|---|

AOA . | C_{D}
. | C_{D}
. | Relative error(%) . |

0° | 0.4008 | 0.3998 | 0.25 |

10° | 0.4526 | 0.4544 | 0.40 |

20° | 0.6058 | 0.6158 | 1.62 |

## IV. UNIVARIATE ANALYSIS

This section investigates the effect laws of design variables on aerodynamic characteristics of NACA 0012 airfoil under Kn = 0.0245 and Ma = 2. The range of design variables is referred to in Refs. 9 and 11 and defined in Table III. Figure 4(a) graphically imparts the contour of pressure for the benchmark case. Notably, the main effect of the Gurney flap on the flow field is at the airfoil’s trailing edge. The pressure contour demonstrates that a high-pressure region was created at the lower airfoil surface and the Gurney flap upstream, causing a concurrent amplification in both lift and drag coefficients.

Design variables . | GF height (%c) . | GF location (%c) . | GF mounting angle (°) . | AOA (°) . |
---|---|---|---|---|

x_{min} | 1 | 0 | 90 | 0 |

x_{max} | 5 | 8 | 150 | 50 |

Design variables . | GF height (%c) . | GF location (%c) . | GF mounting angle (°) . | AOA (°) . |
---|---|---|---|---|

x_{min} | 1 | 0 | 90 | 0 |

x_{max} | 5 | 8 | 150 | 50 |

*L*denotes the lift force,

*ρ*represents the density of the freestream,

*v*

_{∞}represents the velocity of the freestream,

*c*represents the chord length of the airfoil, and

*D*denotes the drag force. Meanwhile, the pitch moment coefficient of the airfoil holds significance in the pitch balance during flight; hence, it usually needs to be constrained within a certain range,

### A. Effect of the gurney flap height

The effect of varying Gurney flap heights on aerodynamic characteristics of the airfoil at Ma = 2, Kn = 0.0245, and AOA = 20° is shown in Fig. 5. As shown in Fig. 5(a), under benchmark conditions, the introduction of Gurney flaps effectually augments the lift coefficient in comparison to the clean airfoil, and the lift coefficient augments with the augmentation of gurney flap height. When the Gurney flap height reaches 5%c, a notable 11.3% enhancement in the lift coefficient relative to the clean airfoil is observed. As the lift coefficient increases with the height of the Gurney flap, so does the drag coefficient, which means that higher gurney flaps will result in a larger drag increment. As a specific instance, at a Gurney flap height equivalent to 5%c, the drag coefficient increases by 7% relative to the clean airfoil.

The pressure coefficient distribution on the airfoil surface is depicted in Fig. 4(b). As can be seen, for both attached and clean Gurney flap airfoils, no evident alteration is observed in the pressure distribution along the upper surface. This phenomenon is attributed to the phenomenon of flow compression transpiring in the upstream vicinity of the Gurney flap, thus causing the increase in the lift coefficient.

The effect of the Gurney flap height on the lift-to-drag ratio and pitching moment coefficient is shown in Fig. 5(b). Under benchmark conditions, Gurney flaps with a height of below 3%c manifest negligible effect on the lift-to-drag ratio, and Gurney flaps with a height of 5%c significantly enhance the lift-to-drag ratio, which is similar to the conclusion obtained in Ref. 11. Simultaneously, with an increase in the Gurney flap height, a commensurate rise is discernible in the nose-down pitching moment coefficient. These outcomes suggest that the Gurney flap’s impact is to augment the effective camber of the airfoil.

### B. Effect of the Gurney flap mounting angle

Figure 6 shows the effect of changing the Gurney flap mounting angle on airfoil aerodynamic characteristics under benchmark conditions. As can be observed in Fig. 6(a), when the angle of attack is 20°, the lift coefficient experiences an incremental augmentation followed by subsequent attenuation as the Gurney flap mounting angle increases, and the maximum value appears near 120°. As the Gurney flap mounting angle increases, the drag coefficient decreases monotonically. This trend underscores the potential of employing an inclined Gurney flap configuration to heighten the airfoil’s lift-to-drag ratio by reducing drag penalties.

Figure 4(c) depicts the distribution of the pressure coefficient on the airfoil surface when the Gurney flap mounting angle is varying. Notably, the pressure distribution along the upper surface of the airfoil is not obviously affected by the Gurney flap. The drop in the mounting angle induces a decline in the pressure magnitude along the lower surface upstream of the Gurney flap. However, downstream of the Gurney flap, with the mounting angle tilt, the location of the pressure drops sharply farther back. The lift effect needs to take into account both the upstream and downstream effects of the Gurney flap.

The impact of the Gurney flap mounting angle on the lift-to-drag ratio and pitching moment is visually shown in Fig. 6(b). Notably, the lift-to-drag ratio follows a monotonic trajectory of augmentation with the Gurney flap mounting angles increasing. The augmentation in the nose-down pitching moment coefficient on the quarter chord of the airfoil exhibits a direct proportionality to the augmentation in the lift coefficient, which is consistent with the results obtained in Ref. 15.

### C. Effect of the Gurney flap location

Figure 7 shows the effect of Gurney flaps mounting at different locations from the trailing edge of the airfoil on aerodynamic characteristics under the benchmark condition. In Fig. 7(a), the lift coefficient is maximum when the Gurney flap is mounted at the airfoil’s trailing edge (*s* = 0). As the gurney flap moves forward away from the trailing edge, its lift-enhancing efficacy is reduced. The drag coefficient of the airfoil increases sharply when the flap is mounted at the range of 4%c–8%c measured from the trailing edge. The Gurney flap moving forward a large distance causes the apparent thickness of the airfoil to increase, leading to an increase in the drag coefficient.

The effect of the Gurney flap mounting location on the lift-to-drag ratio and pitching moment coefficient is shown in Fig. 7(b). Optimal aerodynamic performance is achieved when the Gurney flap is installed at the trailing edge of the airfoil (*s* = 0), while a reduction in the lift-to-drag ratio when the flap moves forward. The pitching moment coefficient on the quarter chord consistently surpasses that of a clean airfoil. However, as the Gurney flap is mounted further forward, a notable mitigation is discernible in the nose-down pitching moment coefficient.

As can be observed in Fig. 8, the Gurney flap moving forward from the trailing edge will result in negative pressure at the rear of the flap. This indicates that there is a large recirculation flow zone, which prevents a section of the airfoil located behind the Gurney flap from generating lift, which also results in a decrease in total lift. This phenomenon is also can be seen in Fig. 4(d), which presents the distribution of pressure coefficient across the airfoil surface. There is a sharp drop in pressure downstream of the flaps. However, when the Gurney flap moves forward to a certain location, the rear area is larger than the recirculation flow zone and the ability to generate lift is restored, so the lift coefficient increases.

### D. Effect of angle of attack

Figure 9 delineates the effect of altering the angle of attack on the aerodynamic characteristics of the Gurney flap-equipped airfoil based on the benchmark case. As shown in Fig. 9(a), it is evident that the lift coefficient experiences augmentation with the increase in the angle of attack, up until the angle remains below 40°, thereafter undergoing a decrement as the angle attains ∼40°. Within the angle of attack range from 0° to 50°, the drag coefficient shows a monotonous incremental trajectory.

Figure 9(b) shows the impact of the angle of attack on the lift-to-drag ratio and pitching moment coefficient. For angle of attack values below 20°, a progressive augmentation in the lift-to-drag ratio accompanies the increase of the angle. However, once the angle of attack surpasses 20°, a decrease in the lift-to-drag ratio ensues, coinciding with a decline in overall aerodynamic performance. This observed trend can be attributed to the differential growth rates between the lift and drag coefficients: when the angle of attack remains below 20°, the lift coefficient experiences a more pronounced growth rate compared to the drag coefficient. Nevertheless, once the angle of attack surpasses 20°, the increment in the drag coefficient outpaces that of the lift coefficient, resulting in the overall reduction in the lift-to-drag ratio. It is easy to understand that the larger the angle of attack is the more inclined the airfoil is and the larger the pitching moment is. Therefore, a positive correlation manifests between the coefficient of the nose-down pitching moment and the angle of attack increment.

## V. AERODYNAMIC OPTIMIZATION

In the pursuit of enhancing the aerodynamic performance of the NACA 0012 airfoil equipped with Gurney flaps in the rarefied gas flow, this section begins with the multi-objective optimization of the Gurney flap at the design point to get the benchmark configuration. Then, we investigate the optimization by controlling the Gurney flap mounting angle and angle of attack over a wide range of Mach numbers and Knudsen numbers. It needs to be noted that the Gurney flap height and mounting location are fixed at the second optimization stage and equal to the flap height and mounting position of the benchmark configuration. In addition, the optimization results and control laws are understood and analyzed based on the univariate phenomena and principles in Sec. IV.

### A. Multi-objective optimization at design point

*δ*) of the predicted value and the true value of the testing data is defined as follows:

*f*

_{test}and

*f*

_{predict}are the true value of the objective function calculated by the DSMC method and the predicted value of ANN, respectively.

Figure 10 illustrates the comparison of 15 distinct test cases, which were obtained by the ANN and by the DSMC method. Figure 10(a) shows the comparison of the values of the lift-to-drag ratio. It can be observed that the absolute error between the predicted value of the ANN surrogate model and the simulated value of DSMC is less than 0.01. Due to the minute magnitude of the DSMC-simulated value within case 1, its *δ* still reaches 300% although it is less than 0.01 compared with the predicted value of ANN. *δ* of other cases is less than 1%. This observation can be rationalized, particularly in instances where the lift-to-drag ratio converges closely to zero. Given that such regions are not of primary interest, the surrogate model can address our concerns.

Figure 10(b) shows the comparison of the value of the pitching moment coefficients. It can be seen that the simulated value of DSMC is almost consistent with the predicted value of ANN, characterized by an average *δ* of 1.1% and a maximum *δ* of less than 3.5%. The model meets the predicted performance requirements and can be used in the subsequent optimization process.

Based on the ANN surrogate model constructed before, we use the NSGA-II multi-objective optimization algorithm for optimization, and the Pareto front is obtained as shown in Fig. 11(a). Given that the lift-to-drag ratio and pitching moment coefficient are mutually restricted, it is impossible to reach the optimum at the same time; thus, it is impossible to obtain the optimal individual with all objectives. To harmonize the performance of the lift-to-drag ratio and pitching moment coefficient, a suitable configuration is singled out from the Pareto front. The selected point is marked in Fig. 11(a), where *h* is 1.00%c, *s* is 0.54%c, *ϕ* is 150°, and AOA is 13.4°. Since the selected point is the predicted value of ANN, the results simulated by DSMC are used to verify that the relative error of the two points is less than 2%, and the accuracy meets the requirements.

Figure 11(b) shows the velocity contours of the flow field, with the uppermost representation depicting the initial configuration and the lowermost portraying the optimized configuration. It can be observed that the wake velocity loss of the optimized airfoil is smaller than that of the initial airfoil, and the total momentum loss is reduced and so is the drag coefficient. According to the univariate analysis, the lower Gurney flap height can reduce the nose-down pitching moment coefficient, but at the cost of the lift coefficient. The mounting angle of the Gurney flap is more inclined, which is conducive to reducing the pitching moment coefficient and increasing the lift-to-drag ratio. In addition, the Gurney flap is located closer to the trailing edge, which reduces lift loss from the recirculation flow zone and increases the lift-to-drag ratio. Summarily, the benchmark configuration selected from the Pareto front can strike a harmonious equilibrium between lift-to-drag ratio and pitching moment considerations so that the airfoil’s aerodynamic performance is more reasonable.

### B. Lift-to-drag ratio optimization at a wide range of Mach number and Knudsen number

^{44,45}Our study focuses on the optimization of the airfoil with the Gurney flap over a range of flight conditions, which is based on the benchmark configuration obtained by single-point optimization. The optimization is expressed in the following formulas:

Figure 12 shows the lift-to-drag ratio variation of the benchmark airfoil at a wide range of Mach numbers and Knudsen numbers. It indicates a reduction in the lift-to-drag ratio as the Mach number increases, a trend paralleling the similar response observed in relation to an augmented Knudsen number. After studying the effect of Ma and Kn on the lift-to-drag ratio of the benchmark airfoil, the next step is to consider the establishment of a wide range of aerodynamic prediction models. Central to this endeavor is the adoption of an artificial neural network (ANN) framework.

We eventually have 75 DSMC samples for training, where 50 samples were obtained by the optimal Latin hypercube sampling and 25 samples from the single factor analysis, and 20% of them are randomly selected as testing data. In terms of neural network structure, four neurons (Ma, Kn, AOA, *ϕ*) are set in the input layer, two neurons ($CLCD$ and *C*_{m}) are in the output layer, and 20 hidden neurons are set. Thirumalainambi and Bardina^{46} indicated that the number of hidden neurons is five times that of input neurons, which can predict the aerodynamic coefficient accurately.

Figure 13(a) shows the prediction accuracy of the model for the lift-to-drag ratio. Notably, the average relative error of the lift-to-drag ratio is 0.66%, while the maximum relative error remains below 4%. Figure 13(b) delineates the prediction accuracy of the model for the pitching moment coefficient. The average relative error of the pitching moment coefficient is 0.59%, with the maximum relative error again residing below the 4% threshold. The demonstrated predictive precision of the artificial neural network model meets the requirement, and it can be used as the model of aerodynamic performance prediction. The neural network model provides rapid aerodynamic analysis.

With the constraint of the pitch moment coefficient, the lift-to-drag ratio can be maximum at different flight conditions by controlling the mounting angle of the Gurney flap and angle of attack. Figure 12(b) shows the maximum lift-to-drag ratio across diverse flight conditions subsequent to the optimization of the Gurney flap’s mounting angle and the angle of attack. Offering a comparative perspective, Figs. 12(a) and 12(b) show the lift-to-drag ratio distributions before and after the control optimization, respectively. A discernible observation is that post-control optimization, the lift-to-drag ratio still exhibits a decreasing trend with increasing Knudsen number (Kn); however, with the increase in the Ma number, the amplitude of the variation of the lift-to-drag ratio decreases greatly.

With the aerodynamic coefficient value of the benchmark airfoil over a wide range of Ma and Kn as the reference, the increasing rate of the lift-to-drag ratio after controlling the Gurney flap mounting angle and AOA is shown in Fig. 14. Table IV shows the control law of the Gurney flap mounting angle and angle of attacks at Kn = 0.024, along with the corresponding increasing rates in the lift-to-drag ratio. These findings are highlighted in Fig. 14. The diminishing Mach number entails a corresponding decrease in the AOA control law. As the AOA in the Case5 airfoil was increased from 13.4° to 22.07°, and the Gurney flap mounting angle was reduced from 150° to 131.15°, the lift-drag ratio increased by 29.25%. In addition, case5 compares with the condition of Ma = 2, Kn = 0.024 before and after control optimization; the relative error of the lift drag ratio is reduced from 29.1% to 9.6%. It shows that the effect of the Ma number on the lift-to-drag ratio can be significantly reduced by controlling the Gurney flap mounting angle and AOA.

. | Ma
. | Kn
. | Optimal mounting angle (°) . | Optimal AOA (°)
. | $CLCD$ before optimization . | $CLCD$ after optimization . | Increasing rate (%) . |
---|---|---|---|---|---|---|---|

1 | 2 | 0.0240 | 90 | 13.08 | 1.0058 | 1.1033 | 9.69 |

2 | 2.5 | 0.0240 | 90 | 16.04 | 0.9301 | 1.1028 | 18.57 |

3 | 3 | 0.0240 | 107.52 | 18.63 | 0.8685 | 1.0742 | 23.68 |

4 | 3.5 | 0.0240 | 121.04 | 20.41 | 0.8190 | 1.0306 | 25.84 |

5 | 4 | 0.0240 | 131.05 | 22.07 | 0.7792 | 1.0071 | 29.25 |

. | Ma
. | Kn
. | Optimal mounting angle (°) . | Optimal AOA (°)
. | $CLCD$ before optimization . | $CLCD$ after optimization . | Increasing rate (%) . |
---|---|---|---|---|---|---|---|

1 | 2 | 0.0240 | 90 | 13.08 | 1.0058 | 1.1033 | 9.69 |

2 | 2.5 | 0.0240 | 90 | 16.04 | 0.9301 | 1.1028 | 18.57 |

3 | 3 | 0.0240 | 107.52 | 18.63 | 0.8685 | 1.0742 | 23.68 |

4 | 3.5 | 0.0240 | 121.04 | 20.41 | 0.8190 | 1.0306 | 25.84 |

5 | 4 | 0.0240 | 131.05 | 22.07 | 0.7792 | 1.0071 | 29.25 |

## VI. CONCLUSION

This research investigates the aerodynamic optimization concerning the NACA 0012 airfoil attached gurney at rarefied conditions. Based on the DSMC method, the effects of three geometric parameters (Gurney flap height, mounting angle, mounting location) and a wide range of AOA on aerodynamic characteristics of the NACA 0012 airfoil attached Gurney flap at rarefied gas flow are studied. The results show that the lift-to-drag ratio will increase when increasing the Gurney flap height and mounting angle, and decreasing the distance between the mounting location and trailing edge. Through multi-objective optimization at a single design point, we select a configuration from the Pareto front as the benchmark configuration, which can harmonize the performance of the lift-to-drag ratio and pitching moment coefficient. Moreover, on the basis of the benchmark airfoil, the control of the mounting angle and AOA of the Gurney flap at a wide range of Ma (2–4) and Kn (0.024–0.24) is studied by using the ANN model. The lift-to-drag ratio can be increased by 29.25%, which significantly improves the aerodynamic characteristics of the airfoil, meanwhile reducing the influence of the Ma on the lift-to-drag ratio.

The results of the effect of Gurney flaps on the aerodynamic characteristics of airfoils under rarefied conditions obtained in this paper can provide a reference for future research. The main contribution of this paper is to investigate the aerodynamic optimization by controlling the Gurney flap mounting angle and angle of attack under rarefied and wide-range flight conditions based on the ANN surrogate model. In the future, our work can be extended to the optimization of aerodynamic control of 3D wings and even aircraft under rarefied and wide-range flight conditions.

## ACKNOWLEDGMENTS

The authors acknowledged support from the Sichuan Science and Technology Program (Grant Nos. 2023YFQ0017 and 2023YFG0337) and the Shanghai Academy of Spaceflight Technology (Grant No. USCAST2022-07).

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Keren Lin**: Data curation (lead); Formal analysis (equal); Investigation (lead); Methodology (equal); Software (equal); Writing – original draft (equal). **Songqin Zhang**: Formal analysis (equal); Methodology (equal); Software (equal); Writing – original draft (equal). **Chenfan Liu**: Formal analysis (equal); Methodology (equal); Software (equal); Validation (equal); Writing – original draft (equal). **Haiwei Yang**: Funding acquisition (equal); Writing – review & editing (equal). **Bin Zhang**: Conceptualization (equal); Funding acquisition (equal); Writing – review & editing (equal).

## DATA AVAILABILITY STATEMENT

The data that support the findings of this study are available within the article.

### APPENDIX: INDEPENDENCE STUDY

*N*

_{s}is the number of samples,

*c*

_{p}is the average number of particles in a cell, and

*S*

_{s}is the sampling steps. When the sampling times are 10 000, the statistical error is less than 0.03%, which meets the accuracy requirements of DSMC simulation. According to the research of Hadjiconstantinou and Alexander,

^{48,49}the time step and mesh size are also important factors affecting the accuracy of DSMC simulation.

#### 1. Time step independence study

In this study, the time step independence of DSMC simulation was investigated, considering four different time steps. The computation results of an example with the number of nodes on the airfoil surface as 720, the total number of grids as 304 000, the time step as 1 × 10^{−7}, and the sampling times as 20 000 were taken as the standard. Table V gives the values of lift, drag, and pitching moment coefficients at different time steps when Ma = 2, Kn = 0.0245, and AOA = 20°, as well as the errors relative to the standard example. As can be seen from the table, when the time step decreases from 5 × 10^{−7} to 2 × 10^{−7}, the relative errors of *C*_{L}, *C*_{D}, and *C*_{L} are all less than 0.2%. Therefore, this paper chooses 5 × 10^{−7} as the time step for calculation.

Time step(s) . | 5 × 10^{−6}
. | 1 × 10^{−6}
. | 5 × 10^{−7}
. | 2 × 10^{−7}
. |
---|---|---|---|---|

C_{L} | 0.6429 | 0.6572 | 0.6586 | 0.6587 |

Relative error of C_{L}% | 3.74 | 1.60 | 1.39 | 1.37 |

C_{D} | 0.6417 | 0.5730 | 0.5692 | 0.5682 |

Relative error of C_{D}% | 15.54 | 3.17 | 2.49 | 2.30 |

C_{m} | 0.1404 | 0.1392 | 0.1394 | 0.1394 |

Relative error of C_{m}% | 1.74 | 0.86 | 0.98 | 0.98 |

Time step(s) . | 5 × 10^{−6}
. | 1 × 10^{−6}
. | 5 × 10^{−7}
. | 2 × 10^{−7}
. |
---|---|---|---|---|

C_{L} | 0.6429 | 0.6572 | 0.6586 | 0.6587 |

Relative error of C_{L}% | 3.74 | 1.60 | 1.39 | 1.37 |

C_{D} | 0.6417 | 0.5730 | 0.5692 | 0.5682 |

Relative error of C_{D}% | 15.54 | 3.17 | 2.49 | 2.30 |

C_{m} | 0.1404 | 0.1392 | 0.1394 | 0.1394 |

Relative error of C_{m}% | 1.74 | 0.86 | 0.98 | 0.98 |

#### 2. Mesh independence study

For the same Ma, Kn, and angle of attack, we consider four mesh scales to study the grid independence. The values of lift, drag, and pitching moment coefficient with different mesh numbers are given in Table VI. As can be seen from the table, when the number of nodes on the airfoil surface is 360, its relative errors are all less than 2.5%. Therefore, in order to obtain accurate results and affordable calculation costs, a grid with the number of nodes on the airfoil surface of 360 was selected for simulation.

Number of nodes on airfoil surface . | 180 . | 360 . | 540 . | 720 . |
---|---|---|---|---|

Number of cells | 78 000 | 138 000 | 220 000 | 304 000 |

Computation time (h) | 1.5 | 2.7 | 5 | 6.5 |

C_{L} | 0.6574 | 0.6586 | 0.6629 | 0.6678 |

Relative error of C_{L}% | 1.57 | 1.39 | 0.74 | 0.01 |

C_{D} | 0.5222 | 0.5692 | 0.5688 | 0.5569 |

Relative error of C_{D}% | 5.97 | 2.49 | 2.42 | 0.27 |

C_{m} | 0.1436 | 0.1394 | 0.1382 | 0.1380 |

Relative error of C_{m}% | 4.05 | 0.99 | 0.17 | 0.01 |

Number of nodes on airfoil surface . | 180 . | 360 . | 540 . | 720 . |
---|---|---|---|---|

Number of cells | 78 000 | 138 000 | 220 000 | 304 000 |

Computation time (h) | 1.5 | 2.7 | 5 | 6.5 |

C_{L} | 0.6574 | 0.6586 | 0.6629 | 0.6678 |

Relative error of C_{L}% | 1.57 | 1.39 | 0.74 | 0.01 |

C_{D} | 0.5222 | 0.5692 | 0.5688 | 0.5569 |

Relative error of C_{D}% | 5.97 | 2.49 | 2.42 | 0.27 |

C_{m} | 0.1436 | 0.1394 | 0.1382 | 0.1380 |

Relative error of C_{m}% | 4.05 | 0.99 | 0.17 | 0.01 |

## REFERENCES

*Key Engineering Materials*

*Experimental Flowfields Around NACA 0012 Airfoils Located in Subsonic and Supersonic Rarefied Air Streams*

*The Levenberg-Marquardt Algorithm: Implementation and Theory*

*A Fast Elitist Non-Dominated Sorting Genetic Algorithm for Multi-Objective Optimization: NSGA-II*