Abstract
Maintaining the laminar flow on surfaces through active control is a significantly promising technique for reducing fuel burn and alleviating environmental concerns in commercial aviation. However, there is a lack of systematic parameter studies for the hybrid laminar flow control (HLFC) together with natural laminar flow (NLF). To address this need, we optimize the infinite swept wings with different sweep angles and at various conditions, including different Mach numbers, Reynolds numbers, and lift coefficients. The Reynolds-averaged Navier-Stokes (RANS) solver coupled with the linear stability theory is applied for the laminar-turbulent transition prediction, and the traditional optimization method based on evolutionary algorithms is applied for laminar flow wing optimization. The optimization results found that HLFC is required when the NLF fails at a larger sweep angle (35°) and Reynolds number (). The lower pressure peak with boundary-layer suction is found to delay the transition of the regional aviation condition. Besides, the pressure distribution of HLFC is similar to NLF results at the lower Reynolds number () or sweep angle (25°), i.e., a gentle negative pressure gradient near the leading edge and a small favorable pressure gradient behind it. Clarifying the characteristics of laminar flow wings will advance the application of the laminar flow technique within its field.
1. Introduction
Facing the permanent goals to improve fuel economy and reduce their impact on the environment, serious novel technologies have been proposed, such as composite material structures [1], laminar flow wings [2], active aeroelastic wings [3], and adaptive structures [4]. For a typical transport aircraft, total friction drag can be up to 50%, so the laminar flow wing is one of the most promising technologies for reducing drag. Research has shown that the friction drag of a transport aircraft is 18%, 4%, 3%, and 3% for the wing, horizontal tail plane, fin, and nacelles, respectively [5]. If the flow were laminar on 40% of the wing, horizontal tail plane, fin, and nacelles, the total drag of the aircraft would be reduced by 16%.
For the laminar techniques, the NLF and HLFC techniques are widely used in maintaining laminar flow [6–10]. The NLF achieves laminar flow on wings via a suitable profile, whereas the HLFC utilizes the suction before the front spar in combination with a suitable shape of the wing box. Also, thermal control on the natural laminar flow configuration has been investigated [11–14], especially for the transonic shock-boundary layer interaction. Before we proceed with the HLFC design, several phenomena that cause the laminar-turbulent transition need to be involved. There are several phenomena to be considered [15, 16], i.e. attachment line transition, TS instabilities, CF instabilities, and bypass-induced transition. In the laminar flow wing design, the attachment line transition can be diverted with the help of a “Gaster bump” [17]. CF instability is the second phenomenon that induces transition. If the CF velocity profile of the boundary layer has an inflection point [18], the inviscid instability will arise, which is amplified close to the leading edge of a swept wing. When the sweep angle and Reynolds number are increased, the CF instability is amplified significantly and might induce transition very close to the leading edge [19]. The third one, TS instability, can occur in two- or three-dimensional flow. TS instability [15] starts with the development of wave-like disturbances. The linear amplification covers about 75–85% of the distance between the leading edge and the transition starting point. Whether it is TS instability or CF instability, the linear stability theory can capture the main disturbance. For the flow with high-level turbulence, the bypass-induced transition [20, 21] is significant and should also be considered. In this paper, since the low-level turbulent flow is considered for aircraft aviation, the bypass-induced transition is not discussed in detail in the following paragraph.
The TS and CF instabilities are significantly affected by the pressure gradient, so a considerable laminar region can be obtained via the art of shape tailoring [22]. As we all know, the negative pressure gradient amplifies the CF instabilities but suppresses TS instabilities, but the positive pressure gradient has the opposite effect [23]. Thus, a careful balance is required for suppressing both TS and CF instabilities when implementing shape optimization. A series of wind-tunnel experiments and flight tests were performed in Europe, and design guidelines were developed [5]. The Reynolds number and sweep angle with different laminar flow techniques are given and shown in Figure 1. We added the corresponding points of business jet aircraft and short-regional and long-regional airlines. It can be seen that when the Reynolds number and sweep angle are increased, the NLF fails to delay the transition and needs to be replaced by the HLFC technique [24] .
The HLFC and NLF designs have been investigated for several decades. The authors of [9] have designed the NLF wing of a business jet (wing-body configuration) using the rule of cosine. The Reynolds number is below , the cruise Mach number is 0.75, and the lift coefficient is 0.5. The optimized result gains a more-than-40% laminar region and is verified by the wind-tunnel test at DNW-HST. Han et al. [22] designed a laminar flow swept wing at a Mach number of 0.75 and the Reynolds number of with . The sweep angle of the leading edge is 19°, and the MAC (Mean Aerodynamic Chord) is 3.75 m. After optimization, a considerable laminar region, almost more than 60%, is obtained on the upper surface. Shi et al. [8] optimized an infinite swept wing with a 25° sweep angle using an adjoint-based optimization approach. The Mach number is 0.78, and the Reynolds number is . The optimized geometry could gain a laminar region of 59.4% on the upper surface. The balance of controlling TS and CF instabilities in NLF optimization is explored, which shows that a lower pressure peak and the following weak negative pressure gradient correspond to a considerable laminar region.
Besides the NLF optimization, the HLFC optimization and its comparison with NLF have been investigated. Risse et al. [25] applied the HLFC technique on the wing and tails of long-range passenger aircraft. The design’s Mach number ranges from 0.85 to 0.80, and the outboard wing leading edge sweep angle varies from 34° to 28°. With the HLFC technique, the outboard wing of long-range aircraft would obtain a laminar region of more than 50%. As a result, a net block fuel benefit of nearly 11% would be gained, which motivates further investigation of the HLFC technique. Yang et al. [26] have studied the HLFC technique on a wing-glove vertical tail, and the HLFC helps to delay the transition. The wing-glove configuration has a 35° sweep angle. The Mach number is 0.75, and the Reynolds number based on the MAC is . The designed configuration could delay the transition of locations to the position around 25% of the chord. Besides, the design wing with HLFC presents good robustness with the variation of sideslip angles and angles of attack. Sudhi et al. [27] researched the NLF and HLFC applications for a swept wing with a 22.5° sweep angle at Mach number 0.78 and Reynolds number of . The HLFC optimization result gains a more-than-43% drag reduction than NLF. From the investigation of the three-dimensional laminar wing optimization using NLF and HLFC techniques, what is missing is a comprehensive study of how the key parameters, i.e., Reynolds number, lift coefficient, and sweep angle, affect the laminar flow design.
This paper contributes to the study of the effect of the key parameters on the laminar flow design. The HLFC design characteristics and their comparison with NLF will be summarized in this paper. This paper contributes to systemically studying the effects of these key parameters on the laminar flow design. The HLFC design characteristics and their comparison with NLF will be summarized in this paper.
In this work, the 𝑒𝑁 method based on the linear stability theory is used to predict the transition induced by TS and CF instabilities, and the infinite-span swept wings under different flow conditions are optimized. This paper is organized as follows. In Section 2, the complete coupled system, which consists of the Reynolds-averaged Navier-Stokes code, a laminar boundary-layer code, and a stability analysis code, is introduced firstly. Meanwhile, the transition prediction code is verified for both NLF and HLFC cases. The optimization method implemented in HLFC and NLF wing designs is described in Section 3. In Section 4, two optimization cases are done, and the aerodynamic characteristics of HLFC wings and their comparison with those of NLF wings are discussed. Finally, we summarize our work in Section 5.
2. Transition Prediction Method and Validation
In this section, the laminar-turbulent transition (LTT) prediction method is introduced, including the RANS solver, laminar boundary-layer (BL) code, linear stability theory, and 𝑒𝑁 method. The transition module part is composed of the last three methods. The RANS solver and LTT module are coupled by external iteration computation using an intermittency function. Note that the surface suction is simulated using a continuous suction boundary-layer condition [28].
2.1. Linear Stability Theory and 𝑒𝑁 Method
The LTT begins when the boundary-layer flow is disturbed by a set of perturbations (noise, external turbulence, surface roughness, and wall suctioning/blowing) transformed through a process called receptivity into a set of small perturbations which grow or decay in the form of normal modes. For the linear stability theory, a small sinusoidal disturbance (velocity, pressure, density, or temperature) is introduced into the Navier-Stokes equations and defined as where is an orthogonal coordinate system, with being normal to the surface wall. Since the fluctuating quantities are very small, the quadratic terms of the perturbations can be neglected in the N-S (Navier-Stokes) equations. It is also assumed that the mean flow quantities do not vary significantly over a wavelength of the disturbances. The result of the parallel flow approximations is that the complex amplitude functions depend on only.
By applying all the assumptions and submitting the mean flow and into N-S equations, we obtain a local linear system of ordinary differential equations for three-dimensional compressible flows: where , , and are complex coefficient matrices that are determined by mean flow values and the wave numbers and . This equation represents a generalized eigenvalue problem for temporal stability analysis, described as where and can be computed by , , and . For the temporal stability analysis, and are the real wavenumbers, while is the complex frequency, i.e. . This generalized eigenvalue equation is solved using the Rayleigh quotient method.
After that, we obtain the temporal growth . The is transformed into a spatial growth ratio through Gaster’s equation [29]. However, it is not clear in advance which amount of wave amplification finally leads to the transition, so the wave amplifications cannot directly yield a transition location. This missing boundary is introduced by the method, and the amplification factor is defined as where is the arc length from the leading edge. is the envelope of all curves at different frequencies or wavelengths at every chordwise station . The transition location is finally obtained at the position at which the local first exceeds the critical amplification factor .
2.2. Suction Boundary-Layer Condition for HLFC
We simulate the HLFC by setting a steady suction velocity normal to the wall surface [25, 27, 28, 30]. The HLFC technique is realized through a large number of discrete suction holes close to the leading edge. For industrial applications, the continuous surface suction condition is used instead of the hole velocity. The surface continuous velocity is expressed as where is the mass flow rate that is obtained using the flow meter [28], is the density, and is the suction surface area. There exists a transformation between the hole velocity and the surface velocity, which is where denotes the hole suction velocity. The porosity is defined as where is the hole diameter and is the area formed by the adjacent four suction holes. The pressure drop determines the hole velocity through the perforated panel, and this relationship can be described as where and are the outer and inner pressure values of the perforated surface, respectively. In Equation (8), is the surface skin thickness, and and are the coefficients that are determined by the properties of the perforated surface. Generally, the experiments are designed to calibrate the undetermined coefficients [30].
In this work, we use the continuous surface suction condition for the HLFC optimization. The suction condition is set for the laminar BL code, and the nondimensional surface velocity coefficient is set as where is the freestream velocity. We consider the influence of wall suctioning on the steady base flow. The distorted base flow is the input for the linear stability theory.
2.3. Transition Prediction Framework
A laminar boundary-layer code and the RANS solver are required for the transition prediction with surface wall suction, except for the linear stability analysis and boundary-layer condition. We use the conical approximation (quasi-3D) in the laminar BL code, which is a simple and efficient method. The quasi-3D code is accurate for the configuration with a high aspect ratio. The laminar BL equations are nonlinear, and Newton’s method is applied to solving the equations. The input of this code is the section geometry coordinates () and pressure coefficient () distribution, which is provided by the converged solution of the RANS solver.
The RANS solver solves the equations with a finite volume method in both steady and unsteady states for the multiblock-structured overset mesh. It uses cell-centered finite volume formulation, and the inviscid fluxes are discretized with artificially dissipated central differencing, while the viscous fluxes use standard central differencing. An approximate Newton-Krylov (ANK) algorithm or a fully coupled Newton-Krylov (FCNK) method is applied to solve the mean flow and turbulence equations. The steady problem can start using time-marching schemes (DDADI or RK4) or ANK and then switch to FNCK or use ANK in the whole convergence process [31]. Note that we choose the Spalart-Allmaras (SA) turbulence model in this paper.
We start the transition prediction from the RANS solver with fixed transition locations or a fully turbulent condition. The converged RANS solver solution () would provide the at different sections along the spanwise direction for the BL code. With the stability analysis and , the predicted transition locations and lengths are obtained; then, the fixed transition results are updated according to the new results and an iteration relaxation factor. The fixed transition results are returned to the RANS solver through an intermittency function [7], which is computed by interpolating the transition locations and lengths for the whole surface. The RANS solver would continue with the updated transition results, and the iteration stops when both transition results and flow field solution converge (, ). The whole transition prediction framework is shown in Figure 2.
2.4. Validation
In this work, three test cases, the NLF(2)-0415 swept wing, TLFTRM01, and SLFTRM01, are selected to validate the CFD code. These three experiments contain TS and CF instability-induced transitions, which are significant in the following optimization work. The TLFTRM01 and SLFTRM01 models are tested with the HLFC control.
2.4.1. NLF(2)-0415
The transition experiment of the NLF(2)-0415 infinite swept wing with 45° sweep angle has been implemented by Saric et al. [32]. This model is designed to explore the transition induced by CF instabilities. The test angle is set as −4° to generate a negative pressure gradient (Figure 3) that amplifies the CF instabilities. The periodic boundary-layer condition is applied in the spanwise direction. The turbulence intensity is 0.0005, and critical amplification factors for TS and CF instabilities are 9.0 and 6.5, respectively.
We have implemented a grid convergence study here. The geometry and the periodic plane for L0.5, L1, and L1.5 grids are shown in Figure 4. The grid size, drag coefficients, and predicted transition location on the upper surface are illustrated in Table 1. The test condition is and .
(a) L0.5 mesh: 1.75 M cells
(b) L1 mesh: 0.80 M cells
(c) L1.5 mesh: 0.22 M cells
We show the pressure coefficient comparison between the simulation and the experimental data at in Figure 3.
The comparison indicates that the pressure coefficient distribution is well simulated for different grid sizes. From the transition prediction framework (Section 2.3), we know that the transition location is determined by the conditions, geometry, and section pressure coefficient.
As the conditions and geometry are unchanged in the prediction iteration, the transition locations would be close if the pressure coefficient distributions almost coincide. The results of transition locations verify this statement, and the transition location value on the upper surface (0.554) of L1 is almost equal to that of L0.5 (0.553) (Table 1). The good match between and the laminar region determines that the drag coefficient difference would be small. From Table 1, the drag difference between L0.5 and L1 is lower than 0.403 counts, which is very close. Thus, L0.5 has enough accuracy as the increment of grid size would not cause an obvious drag coefficient and transition location change. We use the L0.5 mesh in the following research.
The simulation and experimental condition for verification is shown in Table 2. Figure 5 shows the transition locations/regions for the various experimental approaches and simulations. The transition locations via simulation agree well with experimental data. The transition location moves upstream as the Reynolds number increases. As the transition is induced by CF instabilities except for the case of [33], it verifies the accuracy of the established code for the NLF transition dominated by CF instabilities. In the case of , the transition is induced by laminar separation, which is the point where the BL code stops. This agrees with the results by Grabe et al. [33].
2.4.2. TLFTRM01
A wind-tunnel experiment at the transonic regime is implemented for the HLFC study [34]. The experimental model TLFTRM01 is shown in Figure 6. In the middle of the spanwise direction, an interchangeable panel is laid out close to the leading edge for the NLF and HLFC study. When the panel is perforated, it is used for HLFC research; otherwise, the solid panel is applied for the NLF. The test Reynolds number is around , and the turbulent intensity is approximately 0.1%. Various angles of attack and three Mach numbers were tested, and we chose the Mach number of 0.7 and angle of attack of −3.34° for validation.
Under the condition of , , and AoA = −3.34°, three mass flow rates are included in the experiment, i.e., , 6.45, and 7.27 g/s, and they are measured by the flow meter in the experimental process. The inner pressure is measured in the suction chamber. Subsequently, the coefficients and in Equation (8) were calibrated as 1.289 and 0.185, respectively [28]. The hole velocity is then computed using Equation (8), and the surface suction velocity is obtained via Equation (6). Note that is 0.005454. By integrating the calculated in the suction area, the mass flow rates are obtained, and they are 6.11, 6.35, and 7.35 g/s. The simulated mass flow rates are close to the experimental data (, 6.45, and 7.27 g/s).
Before we verified our CFD results with experimental transition data, we performed grid convergence. The different-level (L0.5, L1, and L1.5) grids are shown in Figure 7. We chose two suction mass ratios ( and ) at the condition of AoA = −3.34°, , and . The pressure coefficient comparisons with are compared in Figure 7(d), and the almost coincides with each other. The mesh sizes, drag coefficients, and transition locations are shown in Table 3. We can see that the differences in drag coefficient between L0.5 and L1 are within 0.8 counts for both the case of and the case of . Thus, L0.5 has enough accuracy in predicting the drag forces. The transition locations for those two cases almost coincide at different grid levels (Table 3). Therefore, we used the L0.5 mesh for research in the following part.
(a) L0.5 mesh: 27.4 M cells
(b) L1 mesh: 13.1 M cells
(c) L1.5 mesh: 3.5 M cells
(d) The pressure coefficient comparison with different grid level for the case of
Figure 8 shows the calibrated surface suction velocity. The critical amplification factors and are 8.7 and 7.0, respectively. The is set as the input of the boundary-layer code, and the corresponding CF and TS amplification results are obtained in Figure 9. When the mass flow rate increases, the CF instabilities are suppressed, and the transition is delayed. From the instability analysis (Figure 9), the transition is induced by CF instabilities. The experimental infrared (IR) images of temperature are shown in Figure 10. The temperature gap (the darker to the bright region in the clean region, except for the premature transition caused by the steps, gaps, and hot film sensors) in Figure 10 indicates the laminar-to-turbulent transition line. The bold yellow line represents the transition with suction, and the bold red line on the right side indicates the transition without suction for comparison. The exact experimental transition locations are defined as the middle between the lowest and highest temperatures in the transition region and are extracted from the IR images. The simulated and experimental transition locations are compared in Figure 9(a).
(a) CF instabilities
(b) TS instabilities
The results show that the transition locations agree well with experimental data for cases where the CF instabilities dominate the transition for the NLF and HLFC.
2.4.3. SLFTRM01
A laminar flow wing-glove model, shown in Figure 11, is designed to demonstrate the HLFC technique through flight experiments. The model consists of the test section (marked in red color), two fairings (marked in blue color), and the main wing. The sweep angle of the leading edge is 5°, which indicates that the TS instabilities dominate the transition for this flight experiment. The porous panel close to the leading edge is used for the boundary-layer suction, delaying the transition.
This flight experiment includes several flight conditions; here, we choose the condition of , , , and AoA =2.2° with different suction mass flow rates for validation (Figure 12). For the porous panel, there are four suction chambers. The continuous surface suction coefficient is determined by the approach introduced in Section 2.4.2. The suction coefficient distributions at the “Sec2” position are shown in Figure 12. These two cases differ mainly in the region close to the leading edge. The different cases are named as “case 1” and “case 2.” We choose to use “case 1” as the grid convergence case. The simulation grid convergence study and the pressure distributions are shown in Figure 13. The pressure coefficient distributions with different grid levels almost coincide with each other, and they match well with the experimental data (Figure 13(d)). The drag coefficients and transition locations are shown in Table 4. The results shows that the L0.5 grid has good accuracy for predicting the drag coefficients and transition locations.
(a) L0.5 mesh: 28.6 M cells
(b) L1 mesh: 13.85 M cells
(c) L1.5 mesh: 3.72 M cells
(d) The pressure coefficient comparison with different grid levels at “Sec2”
The IR images from flight tests are shown in Figure 14, and the transition line is the location where the temperature increases sharply. The bright color represents the laminar region, while the darker region is turbulent. With the boundary-layer suction condition, the transition prediction is implemented, and the TS amplification factor is computed, as shown in Figure 15. The critical TS amplification factor is 9.0, calibrated by Yang et al. [35]. As the streamwise pressure taps trigger (dotted white lines in Figure 11) premature transition in the turbulent wedge-shaped area, the outside clean area is focused on the transition phenomenon detection. Therefore, we choose the outside “Sec2” position for the validation comparison. As shown in Figure 14, the white square point represents the simulated transition location, which agrees well with the results in the IR Images. The comparison results demonstrate that the established transaction prediction tool with HLFC influence is able to handle the transition phenomena induced by TS instabilities.
(a) Case 1
(b) Case 2
3. Optimization Framework
The design method for HLFC wings is the traditional optimization method based on evolutionary algorithms. The multiblock free form deformation (FFD) [36] approach is applied for geometric parameterization in this work, providing consistent parameterization across disciplines. Mesh deformation is realized via the radial basis function (RBF) dynamic mesh approach, which is based on the compact support radial basis function [37]. We adopt the differential evolution (DE) algorithm [38] to find the optimal solution.
All optimization cases are implemented for infinite-span wings. Infinite-span wings are classified as one kind of 2.5D wings [39]. The whole wing is formed by a single airfoil, and the wingspan is set as infinite. Besides, the infinite-span wing might have a given sweep angle. With a large enough sweep angle, CF instabilities dominate the boundary-layer flow. Infinite swept wings are simplified models of real wings. Many factors are beyond our focus and are therefore not taken into account, such as twist angles, lifting devices, and the installed nacelles. This simplicity makes infinite swept wings more appropriate for the research and exploration of the HLFC wing design. We use the CFD solver built in Section 2 to predict the transition locations and aerodynamic forces. A periodic boundary condition in the spanwise direction is set to simulate infinite swept wings.
3.1. FFD Parameterization Method
The FFD approach enables smoothly deforming geometric profiles of any type or degree, including curves, surfaces, and three-dimensional geometry. It embeds models into a lattice of control points and builds a mapping relationship between the physical and parametric spaces. The physical space describes the objective models and lattice control points in physical coordinates. The parametric space presents objective models in the local coordinate defined by the lattice. By modifying the shape of the lattice, a deformation is passed to objective models. Here, the Bernstein polynomial is used to set up the mapping relationship. The deformed position of an arbitrary point on object models is defined as where is the physical coordinates of objective models. is the physical coordinates of lattice control points; and are the local coordinates of objective models; and , , and are, respectively, the -th-, -th-, and -th-degree Bernstein polynomials. The deformation of objective models can be obtained by changing the physical coordinates of lattice control points.
Theoretically, the distribution of lattice control points could affect the deformation ability of the FFD parameterization method. The region with densely distributed control points has a better deformation ability than other regions. As the CF instability is sensitive to pressure distributions at the leading edge, the distribution of FFD control points is denser near the leading edge, as shown in Figure 16. The first and last lines of lattice control points are fixed to make sure that the leading edge and trailing edge are not changed. All the optimization cases are implemented based on this infinite-span wing and FFD control points. The infinite-span wing is formed by a single airfoil, so the two side movements of FFD points at the same chordwise position should be consistent. As a result, the number of the total design variables is 14 in the optimization.
3.2. RBF Dynamic Mesh Approach
RBF dynamic mesh approach realizes the mesh movement through the interpolating technique based on the radial basis function. The traditional RBF dynamic mesh method directly sets up the mapping relationship about the displacement increment between all the nonsurface and surface mesh points through radial basis functions. This feature ensures that the RBF interpolation does not require the connectivity information of mesh points. Thus, the calculation process of interpolation is very easy and robust. However, this feature also makes the coefficient matrix of the interpolating equation a dense matrix, which reduces the computational efficiency, especially when the number of mesh points is huge.
This paper uses the RBF dynamic mesh approach based on the radial basis function with compact support. The value of the radial basis function with compact support decreases continuously as the Euclidean distance to the basis point increases. As long as the Euclidean distance is larger than the defined support radius, the function value is set as zero. In such a way, the coefficient matrix of the interpolating equation is a sparse band matrix, which improves computational efficiency to some extent. The details of the RBF dynamic mesh method based on the radial basis function with compact support are shown by Kady and Takahashi [40].
The interpolating function of the RBF dynamic mesh method based on the compact support radial basis function is described as where is the interpolated value of the displacement for a volume mesh point located at , is the number of RBF points, is the radial basis function, is the Euclidean distance between mesh point and RBF points , is the interpolating coefficient that is needed to be solved, and is the support radial. According to the definition, when , the function value of is zero. Once the interpolating coefficient is solved, the displacement increment could be obtained through the interpolating formula
3.3. DE Algorithm
Differential evolution (DE) is one kind of evolutionary algorithms. Similar to the Genetic Algorithm (GA), mutation, crossover, and selection are the three main operators of DE. A detailed description of the DE algorithm can be found in Yang et al.’s work [35]. The mutation operator is defined as
The constant is used to control the mutation ratio. The crossover operator is
The vector is the result of mutation. is a randomly generated variable (), and CR is a parameter that describes the crossover ratio. The selection operator is defined as where is the -th individual of the -th population and is the dimension of the problem.
The vector is the result of the crossover. We define the function value , which uses as the input parameters.
In each generation, we operate each individual from the parent population based on the following: (1)Each time, select three different random individuals from the parent population , and use the mutation operator to generate new individuals (2)Obtain new individuals from by utilizing the crossover operator(3)Compare with and generate the next population through the selection operator
For HLFC wing design problems with 14 design variables, each generation has 65 individuals. The optimization terminates after 75 generations.
4. Optimization and Analysis of Laminar Flow Wings
We continue to research the HLFC technique based on infinite swept wings using the established transition prediction framework in Section 2 and the optimization tool in Section 3. The HLFC wings have two sweep angles and are designed with different Mach numbers, Reynolds numbers, and lift coefficients. Under the same conditions, we optimize the NLF wings and compare them with the results of the HLFC wings.
4.1. Benchmark Wings
For a conventional transonic transport aircraft, such as the A320 and B787, the Mach number range is basically from 0.78 to 0.85, and the sweep angle varies from 25° to 35°. In this work, the sweep angles of infinite-span wings are accordingly selected as 25° and 35°. The Mach number of the wings with 25° sweep angle is 0.78, while the wings with 35° correspond to the Mach number of 0.85.
NLF wings usually have a lift coefficient of around 0.5 [41, 42]. We select two lift coefficients for each kind of wing to discuss the influence of lift coefficients and explore the difference between supercritical wings and transonic laminar flow wings. One is a smaller lift coefficient, close to the lift coefficient of the NLF wings. The other one is a larger lift coefficient, which corresponds to that of supercritical wings.
Referring to UW-5006, a transonic natural laminar flow wing [9], the small lift coefficient is set as 0.5 for the wing with a 25° sweep angle, while the larger lift coefficient is 0.59. For the case with the 25° sweep angle, the relative thickness of the wing profile is set as 0.125. As for the wings at 35°, the Mach number, lift coefficients, and relative thickness of the airfoil are determined through simple sweep theory [43] and the parameters of wings with a 25° sweep angle. The sweep theory establishes relationships about the Mach number, relative thickness of the airfoil, and lift coefficients between two-dimension airfoils and three-dimension wings, especially infinite swept wings. It is described as where , , and are the Mach number, the relative thickness, and the lift coefficient of two-dimension airfoils and 2.5D infinite-span wings, respectively. The coefficient is the sweep angle of the selected reference line on the wing surface. is a correction factor with different values for different wings.
Using the simple sweep theory, we transform the Mach number, lift coefficient, and relative thickness of the airfoil from the wing with the 25° sweep angle to a two-dimensional airfoil. For the case of the 25° sweep angle, the Mach number of the airfoil is 0.706, the relative thickness is 0.138, the smaller design lift coefficient is 0.608, and the larger one is 0.718. In the same way, the parameters of wings with the 35° sweep angle can be obtained inversely with these design parameters of the two-dimensional airfoil. For the wings with the 35° sweep angle, the Mach number is computed as 0.85, the relative thickness is 0.114, and the smaller and larger lift coefficients are 0.45 and 0.53, respectively.
We select two Reynolds numbers for both infinite swept wings. The first one is , which is similar to the flight Reynolds number of NLF aircraft, and the other one is , which is close to the flight Reynolds number for regional aircraft. In this way, the wings with two different sweep angles have similar design parameters to some extent, so the design results can be compared together. Although there are some differences, the discussion about the influence of sweep angles is still meaningful. The design parameters are given in Table 5. Figure 17 shows a comparison of wings with different sweep angles, which are set as the original models for all the following optimization cases.
(a) The comparison of airfoil geometry
(b) The comparison of planforms
4.2. Optimizations for Wings with 25° Sweep Angle
The influences of the Reynolds number, design lift coefficients, and suction coefficients () on HLFC and NLF are studied for wings with a sweep angle of 25°. According to the previous discussion, eight different design conditions are determined for wings with the sweep angle of 25°, as shown in Table 6. The and denote the Reynolds numbers of and , respectively. The “Low” represents the lower lift coefficient, while the “High” corresponds to the higher lift coefficient.
The computational conditions of the original model are listed in Table 7. The optimization objective is to minimize the drag coefficient. We have the aerodynamic constraint that the nose-down pitching moment coefficient is larger than −0.1 and −0.105 () for the lower and higher lift coefficients, respectively. The geometric constraint is no reduction in airfoil relative thickness. Thus, the mathematical model of the optimization problem is defined as
Tables 8 and 9 show the force and transition location information of both original and optimized geometries, and it also provides as comparison the pressure drag coefficient () and friction drag coefficient (). No matter what the lift coefficient and the Reynolds number are, the drag coefficients of HLFC wings are all reduced, as shown in Tables 8 and 9. After optimization, the total drag is reduced at all conditions. Figures 18–21 show the comparison of the airfoil geometry and pressure distributions between the original models and optimization results.
(a)
(b)
(a)
(b)
Figures 22 and 23 plot the envelope curves () of disturbance amplification curves of TS and CF instabilities on upper and lower surfaces.
(a) The comparison of the TS amplification factor on the upper surface
(b) The comparison of the CF amplification factor on the upper surface
(c) The comparison of the TS amplification factor on the lower surface
(d) The comparison of the CF amplification factor on the lower surface
(a) The comparison of the TS amplification factor on the upper surface
(b) The comparison of the CF amplification factor on the upper surface
(c) The comparison of the TS amplification factor on the lower surface
(d) The comparison of the CF amplification factor on the lower surface
We now move on to the results of the lower lift coefficient and lower Reynolds number. When the Reynolds number is and the lift coefficient is 0.5, i.e., “NLF-Opt-Re1-Low” and “HLFC-Opt-Re1-Low,” the NLF and HLFC optimization results of pressure distributions are similar (Figure 19(a)). On the upper surface, after the negative pressure peak at the leading edge, the pressure distribution has a small adverse pressure gradient, followed by a gentle favorable pressure gradient, as shown in Figure 19(a). Such characteristics of pressure distribution could efficiently control the development of CF and TS instabilities. As a result, the optimized results (“NLF-Opt-Re1-Low” and “HLFC-Opt-Re1-Low”) gain lower pressure drag and friction drag. Compared with the original model, the total drag of “NLF-Opt-Re1-Low” and “HLFC-Opt-Re1-Low” is reduced by 35.7% and 34.1%, respectively. The HLFC geometry gains a longer laminar region on both the upper and lower surfaces than the optimized NLF geometry, which results in a lower friction drag (Table 8).
For the original models (“NLF-Ori-Re1-Low” and “HLFC-Ori-Re1-Low”), Figures 19(a), 22(a), and 22(b) show that the adverse pressure gradient on the upper surface increases the TS instabilities sharply around 25% of the chord. As a result, the TS instabilities induce the transition. However, the NLF optimization result (“NLF-Opt-Re1-Low”) indicates that the favorable pressure gradient after the position of adverse pressure gradients suppresses the TS instabilities. Besides, the suitable adverse pressure gradients on the upper surface, which are just behind the negative pressure peak, would reduce the CF instabilities. Thus, we can conclude that shape optimization successfully delays the transition even without suction control. Note that for cases “NLF-Opt-Re1-Low” and “HLFC-Opt-Re1-Low,” the transition on the upper surface is induced by laminar separation [44, 45] since both TS and CF amplification factors do not exceed the critical values before the separation point (Figures 22(a) and 22(b)). The pressure gradient effects on transition have been investigated by Sengupta and Tucker [44, 45]. The results of “NLF-Opt-Re1-Low” demonstrate the suitable adverse pressure gradient just after the leading edge, and then, a gentle favorable pressure gradient could balance the development of TS and CF instabilities. This is an effective approach in delaying the transition to the NLF design.
Subsequently, we focus on the lower surface without suction. Compared with “NLF-Ori-Re1-Low,” “NLF-Opt-Re1-Low” presents a larger negative pressure gradient close to the leading edge, and then, a weak adverse pressure gradient and a smaller favorable pressure gradient follow, which stabilize the CF instabilities. The pressure distribution characteristics of “NLF-Opt-Re1-Low” result in a larger CF amplification close to the leading edge and then smaller values for the remaining region compared with the original case, i.e., “NLF-Ori-Re1-Low” (Figure 22(d)). Compared with “NLF-Ori-Re1-Low,” the smaller negative pressure gradient in the middle region (“NLF-Opt-Re1-Low”) increases TS instabilities (Figure 22(c)), but does not exceed the TS critical value. For the “HLFC-Opt-Re1-Low,” the almost similar negative pressure gradient distribution results in much smaller CF instabilities compared with that of “HLFC-Ori-Re1-Low” (Figure 22(c)). Compared with “NLF-Opt-Re1-Low,” “HLFC-Ori-Re1-Low” has lower pressure recovery due to a smaller negative pressure gradient in the middle region, which results in a longer laminar region. It can be concluded that for the lower Reynolds number and lift coefficient, both the NLF and HLFC optimized results could gain a considerable laminar flow region.
When the Reynolds number increases (“NLF-Ori-Re2-Low”), the transition of the NLF-optimized geometry occurs close to the leading edge of the upper surface (“NLF-Opt-Re2-Low”), as shown in Table 8. The results indicate that the NLF optimization fails to gain more laminar regions on the upper surface, whereas a considerable laminar region is obtained on the lower surface via optimization. Figure 22 shows that when the Reynolds number is high, the CF and TS instabilities are much higher, especially the CF instabilities close to the leading edge. The fact proves that the development of instabilities is indeed susceptible to the Reynolds number.
So if we want to delay the transition when the Reynolds number is high, a suction control technology is needed. By coupling the suction control and shape optimization, the CF instabilities are successfully suppressed (Figure 22(b)). After the leading edge, the favorable pressure gradient is first reduced a lot and then changed to an adverse pressure gradient distribution, which is aimed at suppressing the CF instabilities. Then, a favorable pressure gradient after the position of adverse pressure decreases the TS instabilities. Finally, the optimization of HLFC delays the transition location from to on the upper surface and to on the lower surface. The pressure drag coefficient and friction drag coefficient are reduced significantly, as shown in Table 8.
In general, when nose-down pitching moments are set as design constraints; higher design lift coefficients mean a higher negative pressure peak at the leading edge. This will increase the difficulty of delaying the transition, especially when the Reynolds number is also high. The phenomenon can be observed by comparing Figures 19 and 21. We take the result of “HLFC-Opt-Re2-High” to testify to this conclusion. Although the suction control technology has been used, a good pressure distribution is still required to effectively suppress the CF and TS instability development. It can be seen that on the upper surface, the CF instabilities almost reach the critical value of the factor ( in Figure 22(b)) just after the leading edge. On the lower surface, the development of CF instabilities maintains a growth trend, as shown in Figure 23(d), and the transition is not pushed aft a lot. For “HLFC-Opt-Re2-High,” it is up to the limitation of the pitching moment coefficient to maintain a reasonable laminar region.
We then focus on the lower lift and Reynolds number to explain the effects of the nose-down pitching moment on laminar wing design. Figures 19 and 22 show that both shock waves of the NLF wing (“NLF-Opt-Re1-Low”) and the HLFC wing (“HLFC-Opt-Re1-Low”) dominate the transition on upper surfaces. Just before and close to the shock wave, TS instabilities present a sharp increase and reach the critical value immediately. If the TS amplification factor is not up to the critical value, the laminar separation is assumed to induce the transition. A downstream movement of the shock wave enables more laminar flow, but this would cause a nose-down pitching moment increase. As shown in Tables 10 and 11, almost all the optimized results have larger pitching moment coefficients than those of the original values, except for “NLF-Opt-Re2-High.” Besides this, a larger nose-down pitching moment means a bigger trim drag. Therefore, for transonic laminar flow wings, nose-down pitching moment constraints and the locations of shock waves significantly affect the design results of laminar flow wings. If the nose-down pitching moment constraint is set too small, laminar flow wings need a higher negative pressure peak at the leading edge to trim the nose-down pitching moment. As we have discussed above, this is harmful to maintaining laminar flow. The conclusions mentioned above are also valid for the wings with as shown in Figures 21 and 23 and Table 11.
Finally, we plot the skin friction drag () distribution as shown in Figure 24. For brevity, we plot the comparison of the original case and optimized result at and in Figure 24(a) and the results at different Reynolds numbers in Figure 24(b). The comparison of original and optimized results in Figure 24(a) shows that the friction drag of laminar flows is much lower than that of turbulent flows. Note that the transition start point is defined as the position where the grows sharply. Figure 24(b) is aimed at comparing the friction drag at different Reynolds numbers, such as “NLF-Ori-Re1-Low” and “NLF-Ori-Re2-Low” and “HLFC-Opt-Re1-Low” and “HLFC-Opt-Re2-Low.” The friction drag coefficient at the higher Reynolds number is much smaller than that at the lower Reynolds number with the same flow station (laminar or turbulent), because the is inversely proportional to the Reynolds number. This kind of phenomenon is much more pronounced in comparison to the turbulent region. This explains why the friction drag of “HLFC-Opt-Re2-Low” is much lower than that of “HLFC-Opt-Re1-Low.” Even the transition maintenance of the higher Reynolds number is much more difficult than that of the lower Reynolds number; a higher Reynolds number with a considerable laminar region could gain lower skin friction.
(a) The friction drag on wings at and
(b) The friction drag comparison of HLFC results at with different Reynolds numbers
We also plot the boundary-layer thickness to elaborate on why the pressure drag is also reduced significantly when the transition is delayed. Figure 25(a) indicates that HLFC cases have thinner boundary-layer thickness than NLF results. Besides, Figure 25(b) shows original and optimized HLFC results at different Reynolds numbers. It can be seen that the boundary-layer displacement thickness of the higher Reynolds number (Re2) is much smaller than that of the lower Reynolds number (Re1). Smaller displacement thickness means lower pressure drag, and this can be testified to by the results of “HLFC-Opt-Re1-Low” and “HLFC-Opt-Re2-Low” and “HLFC-Opt-Re1-High” and “HLFC-Opt-Re2-High.”
(a) The NLF and HLFC results at lower lift coefficients and Reynolds numbers
(b) The comparison with different Reynolds numbers
The characteristics of laminar flows could significantly reduce the friction and the pressure drag, as shown in Tables 8 and 9. Meanwhile, we can summarize that when the laminar region is considerable for both lower and higher Reynolds numbers, the higher Reynolds number corresponds to the smaller friction drag and pressure drag. A higher Reynolds number means a smaller boundary-layer thickness, as shown in Figure 25.
4.3. Optimization for Wings with 35° Sweep Angle
In the case of the large sweep angle, we focus on the lift coefficients of and . The Reynolds numbers are and . The corresponding design conditions of 35° infinite-span wings are listed in Table 12. The suction coefficient () is increased to −0.00045. The original cases are similar to those of 25° (Table 12). The optimization design objective is to reduce the drag coefficient. We have the aerodynamic constraint that the nose-down pitching moment coefficient is larger than −0.125 and −0.165 () at the lower and higher lift coefficients, respectively. The geometric constraint is no reduction in airfoil relative thickness. The mathematical model of the optimization problem is
All the optimization processes converge when they reach the given maximal generation.
The force coefficient results and transition locations are shown in Tables 13 and 14. The airfoil and pressure coefficient () distribution comparisons are shown in Figures 26–29 The corresponding amplification factor results are shown in Figures 30 and 31. It can be summarized that the drag coefficients of optimization results are also much lower than those of the original models, whether it is a low Reynolds number or a high Reynolds number. The NLF cases, i.e., “NLF-Opt-Re1-Low,” “NLF-Opt-Re2-Low,” “NLF-Opt-Re1-High,” and “NLF-Opt-Re2-High” fail to maintain a considerable laminar region on the upper surface.
(a)
(b)
(a)
(b)
(a) The comparison of the TS amplification factor on the upper surface
(b) The comparison of the CF amplification factor on the upper surface
(c) The comparison of the TS amplification factor on the lower surface
(d) The comparison of the CF amplification factor on the lower surface
(a) The comparison of the TS amplification factor on the upper surface
(b) The comparison of the CF amplification factor on the upper surface
(c) The comparison of the TS amplification factor on the lower surface
(d) The comparison of the CF amplification factor on the lower surface
The drag coefficient reduction of these optimized NLF cases is owed to the weaker shock wave, which indicates a lower pressure drag and a longer laminar region on the lower surface (Figures 30 and 31). The NLF optimization results are different from the results of the sweep angle of 25°. Even at the lower Reynolds number (), the transition occurs close to the leading edge due to CF instabilities (Figure 30(b)), demonstrating that the larger sweep angle increases the design difficulty of maintaining laminar flow.
Next, an active flow control technique is required. Compared with the NLF results, the HLFC results (both original and optimized cases) gain a lower drag coefficient. In the matter of the HLFC optimization at the lower Reynolds number and lower lift coefficient (“HLFC-Opt-Re1-Low”), the transition on the upper surface is successfully delayed from to , which is up to the pressure recovery region. The drag is reduced by 32.04%, and both the pressure and fraction drag values are reduced. The CF instabilities (Figure 30(b)) on the upper surface are suppressed and do not exceed to the critical value, and the TS factor is also not up to the critical value (Figure 30(a)). As a result, the transition is assumed to be caused by the laminar separation. On the lower surface, 32% laminar region is gained. Similar to the case of “HLFC-Opt-Re1-Low,” the coupling of the boundary-layer suction and shape optimization ensures the considerable laminar region for “HLFC-Opt-Re2-Low” and “HLFC-Opt-Re1-High.” Except for the leading edge, the weak negative pressure gradient distributions (“HLFC-Opt-Re1-Low,” “HLFC-Opt-Re2-Low,” and “HLFC-Opt-Re1-Highlight”) are similar to the cases of the 25° sweep angle. A higher Reynolds number and higher lift coefficient (“HLFC-Opt-Re2-High”) require a larger negative pressure gradient aimed at suppressing the TS instabilities (Figure 31(a)). The lowered pressure peak is trying to lower the CF instabilities close to the leading edge, and the suction control in the first 20% of the chord could suppress the CF instabilities (Figure 31(b)), which are amplified by the negative pressure gradient. Finally, we obtain a 28.6% chord laminar region on the upper surface. However, the transition occurs at the leading edge on the lower surface. For commercial aircraft, there is no considerable laminar region on the lower surface due to the obstacles formed by the high lift system pylon and nacelle pylon. Thus, the main effort should be made to delay the transition on the upper surface.
For the “HLFC-Opt-Re2-High,” we only obtain a 28.6% laminar region even with suction control, which indicates the difficulty in maintaining the laminar flow when the Reynolds number, sweep angle, and lift coefficient are increased.
Note that for “HLFC-Opt-Re2-Low” and “HLFC-Opt-Re2-High,” it almost reaches the limitation of the nose-down pitching moment (Tables 15 and 16). Similar to the results mentioned before (Section 4.2), the larger lift coefficients mean more effort to delay the transition. The nose-down pitching moment coefficient (a constraint condition) is also one of the most important factors that could limit the further extension of laminar flows. As shown in Figure 29 and Table 14, to maintain considerable laminar flow and reduce drag, the pressure peak is reduced a lot to ensure that the CF instabilities are not amplified sharply close to the leading edge (“HLFC-Opt-Re2-High”).
As a result, the pitching moment is greatly increased and up to the limitation. Besides, the larger sweep exacerbates this phenomenon. We can also summarize that the CF and TS instabilities on wings with a 35° sweep angle are much more difficult to control than those on wings with a 25° sweep angle. At the higher lift coefficient and Reynolds number, a large suction strength is required, or the limitation of the nose-down pitching moment needs to be relaxed.
Similar to the sweep angle of 25°, the higher Reynolds number results in a lower pressure and friction drag. Taking “NLF-Opt-Re1-Low” and “NLF-Opt-Re2-Low” as examples, the “NLF-Opt-Re2” has a short laminar region but has a lower friction drag (Table 9, 33.544 counts and 30.982 counts, respectively). For “HLFC-Opt-Re1-Low” and “HLFC-Opt-Re2-Low,” the “HLFC-Opt-Re1-Low” has a larger pressure drag (16.766 counts), and “HLFC-Opt-Re2-Low” has a stronger pressure recovery (14.189 counts).
5. Conclusions
This paper uses optimization and CFD tools to research HLFC and NLF wing designs. The influences of sweep angles, Reynolds number, pressure distributions, and design lift coefficients on drag reduction and transition delay are studied by comparing the design results of HLFC and NLF wings. Two different lift coefficients are chosen to study and explore the difference between supercritical and transonic laminar flow wings. The flight Reynolds number of NLF aircraft is the lower , while the higher corresponds to regional aircraft.
We establish the laminar-to-turbulent transition by coupling the RANS solver and the transition prediction module. The laminar boundary-layer code, linear stability theory, and 𝑒𝑁 method make up the transition module. The established transition prediction method is verified via the wind tunnel and flight tests. The results indicate that the TS and CF instability-induced transition with or without boundary-layer suction can be accurately captured and is appropriate for the laminar flow configuration optimizations of infinite-span swept wings.
The optimization system consists of the FFD parameterization, RBF dynamic mesh approach, and DE algorithm. NLF and HLFC have both been optimized for comparison. The results show that the HLFC optimization could gain considerable laminar regions under different conditions, but the NLF optimization only succeeds under the condition of a lower Reynolds number and lower sweep angle.
The pressure gradient distribution has a significant effect on delaying the laminar-turbulent transition. For the NLF wings, the pressure distributions with the maximum laminar flow length may have a low negative pressure peak and a gently adverse pressure gradient close to the leading edge, aimed at controlling CF instabilities. Then, behind the adverse pressure gradient, a small favorable pressure gradient is used to suppress the development of TS instabilities. This kind of distribution succeeds in controlling both TS and CF instabilities at the lower Reynolds number and sweep angle, which would be used for future business jet aircraft using the laminar flow technique.
For the lower Reynolds number and sweep angle, the HLFC-optimized is similar to that of the NLF distribution. For the higher Reynolds number at different conditions, the HLFC-optimized 𝐶𝑝 has a much lower pressure peak, together with the suction control in order to suppress the CF instabilities close to the leading edge. A negative pressure gradient is then distributed to suppress the TS instabilities. The lower pressure peak increases the pitching moment up to the limit, especially at the higher Reynolds number and lift coefficient. Thus, the pitching moment limit is an essential constraint for laminar flow wing design. The HLFC could maintain a long laminar region at the lower Reynolds number without significantly increasing the pitching moment. Compared with NLF, HLFC could maintain a considerable laminar region on the upper surface in most cases and is able to obtain a 28.6% chord laminar region for the condition of commercial aircraft. After comparison, the larger sweep angle cases are similar to those of the lower sweep angle to some extent. Still, it is much more challenging to control both CF and TS instabilities, which are the main obstacle to applying the laminar flow technique on commercial aircraft.
In summary, this paper emphasizes the importance of applying HLFC for the condition of a high Reynolds number, high sweep angle, and high lift coefficient. The pressure distribution characteristics at different conditions are explored. The CFD prediction approach, laminar flow optimization tool, and optimization analysis in this paper would be beneficial for the drag reduction of business jet aircraft and short-range and long-range airlines.
Data Availability
Any request to the corresponding author should be directed to https://arc.aiaa.org/doi/full/10.2514/1.J058266https://arc.aiaa.org/doi/full/10.2514/1.J060842.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This research was funded by the Innovative Research Group Project of the National Natural 587 Foundation of China, and the grant number is 12002284. This research was also funded by the 588 certificate that certifies that its holder is awarded the China Postdoctoral Science Foundation fellowship, and the 589 grant number is 2021M692569.