Aerodynamic Optimization of Airfoil S809 Using High Fidelity FFD Method Coupling with CFD

1. Introduction

Wind energy has been harnessed for centuries, yet fossil fuels became the dominant source post-industrialization. The International Energy Agency estimates that coal, natural gas, and oil will deplete in 139, 49, and 54 years, respectively[1]. As environmental concerns over fossil fuels grow, transitioning to renewable energy sources like wind, solar, and nuclear is crucial for sustainability. Policy support for reducing greenhouse emissions is boosting wind energy's popularity as shown in Figure 1. In 2022, onshore wind energy costs fell to USD 0.033/kWh, making it much cheaper than fossil fuels, while offshore wind costs slightly increased to USD 0.081/kWh [2]. Advancements in Horizontal Axis Wind Turbines (HAWT) and Vertical Axis Wind Turbines (VAWT) are driving renewable energy progress, with HAWTs enhancing large-scale efficiency and VAWTs improving urban adaptability.

An essential factor in decreasing energy costs in the wind turbine sector is the advancement of airfoils with optimized aerodynamic characteristics. The design of these airfoils is essential for optimizing energy output and overall efficiency, particularly as wind farms mature. Contemporary wind farms exhibit a decline in power generation over time, with capacity factor reductions between 10% and 15% relative to initial manufacturer forecasts [4]. The losses are largely ascribed to airfoil deterioration, heightened surface roughness, and alterations in aerodynamic performance resulting from wear and environmental exposure. Consequently, improving airfoil designs is crucial for sustaining energy output, minimizing operational expenses, and prolonging the lifespan of wind turbines. Optimizing wind turbine blades involves fine-tuning factors like airfoil shape, chord length, thickness, and twist angle distribution [5].

In recent years, the optimization of airfoil shapes has been a central emphasis in enhancing aerodynamic performance, especially regarding renewable energy. Diverse methodologies have been investigated to enhance airfoil designs, encompassing gradient-based methods and sophisticated evolutionary algorithms. For example, the reshaping of airfoils by parameterization approaches like as PARSEC and CST, in conjunction with these algorithms, has produced encouraging outcomes in improving lift and overall aerodynamic efficacy. Notably, a gradient-free method [6] using the Differential Evolution algorithm achieved a 99.8% reduction in drag-to-lift ratio. Additionally, the S809 airfoil was successfully reshaped using PARSEC and a Nash equilibrium genetic algorithm[7], yielding optimal results on the Pareto front. This research [8] series optimized airfoil chord lengths and twist angles using Bezier curves, a MATLAB GA optimizer, and NREL's FAST Framework, significantly enhancing AEP and reducing costs for the NREL 5MW and a custom turbine.Subsequent studies [9] used a GA with CST and PARSEC methods, and XFOIL calculations, to significantly increase lift coefficient and lift-to-drag ratio, enhancing wind turbine blade design. Another study [10] optimized NREL S-821 and RAE-2822 airfoils using CST parametrization and a MATLAB-based GA, alongside the panel method, leading to significant enhancements in lift-to-drag ratios, validated by XFOIL, NREL data, and complemented by CFD analysis. Research [11] utilized CST parametrization with a MATLAB GA, incorporating fitness estimations from Hess-Smith panel method and XFOIL, validated using Ansys-Fluent, to develop an airfoil with superior aerodynamic attributes. In addition of that , this study [12] employs Bézier curves to generate airfoil polar points, utilizes XFOIL as a flow solver, and implements a genetic algorithm in MATLAB for optimization, resulting in a newly optimized airfoil shape that produces greater lift compared to the original design, with additional flow analysis conducted using ANSYS Fluent. Subsequently, a study [13] introduced a machine learning-based algorithm combining CST with XFOIL for improved lift-to-drag ratio optimization, demonstrating superior convergence and performance over conventional GA methods. Furthermore a multifidelity local search algorithm using manifold mapping (MM) efficiently simulated high-fidelity models with low-fidelity ones, drastically reducing computational demands. Tested against transonic flow problems, it performed comparably to the gradient-based SLSQP but proved more cost-effective and scalable, despite a slight increase in drag count [14] . Additionally, a novel approach introduced a deep learning-based surrogate model, combining CST and Latin hypercube sampling, to optimize rotor airfoil aerodynamics and minimize dynamic stall, leading to significant improvements in drag and moment coefficients while enhancing lift [15]. One study[16] applied the CST technique and a modified XFOIL panel code with a derivative-free Covariance Matrix Adaptation Evolution Strategy to design airfoils for megawatt-class wind turbine rotor blades, achieving equal or improved performance compared to Delft University designs. Lastly, the optimization of the NACA 4412 airfoil was conducted using genetic algorithms and particle swarm optimization in MATLAB, with cubic spline, PARSEC, and CST methods for parameterization, while XFOIL served as the flow solver, and CFD simulations using ANSYS FLUENT validated the results, focusing on maximizing lift and minimizing drag [17]. A notable study [18] introduced a fast, interactive design framework for airfoil aerodynamic optimization, employing a B-spline-based generative adversarial network and various neural network models to predict aerodynamic responses, achieving results that closely align with traditional CFD methods while significantly reducing computation time to just seconds .

High-fidelity aerodynamic shape optimization (ASO) models rely on complex, non-linear Partial Differential Equations (PDEs) essential for both 2D and 3D optimizations. One major challenge in this area is the accurate differentiation of these PDEs, which often limits research in gradient-based high-fidelity ASO. Despite their computational demands, high-fidelity models offer highly accurate predictions, making them essential for applications where even small shape changes significantly impact performance, such as in aerospace engineering and airfoil design. Recent advances in ASO research have shifted towards incorporating evolutionary algorithms and machine learning (ML) techniques, diverging from traditional gradient-based methods. For example, one study [19] utilized a multi-objective Genetic Algorithm (MOGA) integrated with FLUENT software’s Reynolds-Averaged Navier-Stokes (RANS) equations for airfoil optimization, which resulted in notable improvements in both lift and lift-to-drag ratios. Another study [20] employed a continuous adjoint method coupled with CFD to explore the influence of turbulence models on NACA 0012 and DU 93-W-210 airfoil shapes. This study successfully optimized the airfoils by increasing lift coefficients by 20% and reducing drag by 3%, achieved through the use of spline control points.

Furthermore, morphing trailing edge optimization using free-form deformation parameterization was examined in [21] , revealing that an optimal number of control points is crucial for avoiding wavy deformations and enhancing optimization accuracy. In my research, I solve RANS equations with a Spalart-Allmaras turbulence model, applying a gradient-based optimization algorithm and adjoint method to calculate necessary derivatives. Similar to previous research [22], this approach achieved an 8.5% reduction in drag through single-point optimization, while considering lift, pitching moment, and geometric constraints. Moreover, multi-point optimization was used to ensure realistic designs, resulting in drag values within 0.1 counts of each other across various conditions.

Additionally, another pivotal study [23] focused on transonic aerodynamic shape optimization of the Common Research Model (CRM) wing benchmark, addressing both computational intensity and design sensitivity. This research employed six benchmark optimizations, ranging from single to nine-point cases with 768 variables, ultimately achieving a 7.5% drag reduction for single-point designs and enhanced robustness for off-design conditions in the nine-point cases. Finally, recent work [24] developed an efficient discrete adjoint framework for shape optimization using OpenFOAM, effectively applying this method to UAV and car designs, reducing drag by 5.6% and 12.1%, respectively.

This research employs an innovative methodology that integrates Free-Form Deformation (FFD) with computational fluid dynamics (CFD) and the adjoint technique, employing the numerical flow field solver and SLSQP optimizer for the optimization of airfoil shapes. The principal benefit of this technique is its capacity to effectively investigate the design space while maintaining high accuracy in aerodynamic forecasts. In contrast to earlier techniques that frequently depend on surrogate models or evolutionary algorithms, our methodology enables direct gradient-based optimization, markedly improving convergence rates and the accuracy of optimum designs. Furthermore, utilizing the adjoint technique enables the rapid computing of sensitivities concerning various design factors, effectively mitigating the scalability and accuracy deficiencies commonly found in conventional optimization frameworks. This feature facilitates a more comprehensive investigation of the parameter space, resulting in enhanced aerodynamic performance and enabling the design of airfoils that satisfy the requirements of contemporary renewable energy applications. Collectively, these research exemplify the incremental amalgamation of evolutionary algorithms, machine learning, and computational fluid dynamics in enhancing aerodynamic optimization, indicating the possibility for substantial advancements in airfoil design efficacy.

This study is organized as follows. Section 2 describes the methodology for airfoil shape optimization, including the use of Free-Form Deformation (FFD) for parameterization, mesh generation and deformation processes, the CFD solver, the discrete adjoint method, and the SLSQP optimization algorithm. Section 3 presents and discusses the results of the numerical flow field CFD solver validation, as well as the findings from both single-point and multi-point optimization. Finally, the conclusion summarizes the key insights of the study.

2. Methodology

This section provides a brief summary of the components within the optimization framework, as illustrated in Figure 2 using an extended design structure matrix (XDSM) diagram[25]. The study employs a sophisticated framework designed for high-fidelity aerodynamics and structural optimization, utilizing a direct design approach to iteratively refine airfoil shapes, as depicted in the workflow represented by the XDSM diagram.

The optimization process initiates with the assignment of initial design variable values, x(0), to the optimizer ,which then updates the surface mesh using the surface deformation module in every iteration except the first iteration and computes its analytical derivatives in relation to the design variables expressed as dxs/dx ,followed by volume deformation module modifying the volume mesh and determining its derivatives in relation to the surface mesh changes, notated as dxv/dx , while the flow solver processes flow states denoted as wrelayed to the adjoint solver to calculate the total derivative, culminating in the determination of the objective function f and its derivatives df/dx by the optimizer, with surface and volume deformation steps being quick and explicit and flow and adjoint solver steps being iterative and time-consuming , typically involving about o(10²) major iterations, setting a fundamental minimum for the count of CFD solutions and mesh updates required, with each major iteration including additional CFD solutions.

2.1. Geometry Parameterization

The geometry of wind turbine airfoils plays a critical role in the overall performance and efficiency of rotor systems. Among various airfoil designs, the S809 airfoil stands out as a key contributor to advancements in wind turbine technology, particularly within the framework of NREL Phase VI Rotor. This airfoil is specifically engineered to optimize laminar flow, achieving a thickness of 21%, which significantly enhances aerodynamic efficiency.

As illustrated in Table 1 and depicted in Figure 3, the S809 airfoil features a carefully designed tapered and twisted geometry, which is instrumental in maximizing the effectiveness of the rotor's two-bladed configuration [26,27]. Its design principles have been applied consistently in NREL's Phase II, III, and VI (HAWTs), demonstrating its pivotal role in enhancing blade performance throughout their operational length [28]. In light of these attributes, this subsection delves into the geometric parameterization of the S809 airfoil, examining its implications for rotor design and performance optimization in contemporary wind energy systems.

In shape optimization, particularly for airfoil design, geometry parametrization is essential and challenging. The free form deformation (FFD) technique, introduced by Sederberg and Perry [30], is employed to translate design variables into actual geometric changes. It involves applying a grid structure over the object, where modifications to the grid's control points result in the object's indirect deformation as shown in Figure 4. (a) which will be our work on this paper. The FFD volume captures alterations in geometry instead of the geometry itself, leading to a streamlined and concise collection of design variables for the geometry, simplifying the manipulation of intricate shapes.It is straightforward, concise, and effective, with readily calculable analytical sensitivity derivatives suitable for gradient-based optimization techniques.

The subsequent step integrates the method to generate geometry parametrization, constraints, and utility functions for geometry manipulation, particularly notable for its capacity to compute derivatives for all parametrization methods and constraints, as depicted in Figure 4.(b) , where the airfoil section is molded by adjusting 20 FFD control points, distributed evenly between the upper and lower surfaces, thereby yielding 40 variables for shaping the airfoil

2.2. Mesh Generation and Deformation

Before proceeding with the optimization process, mesh generation and deformation are key components in modern CFD applications. It is essential to generate the volume mesh for the baseline S809 airfoil which involves growing or extruding the mesh in successive layers from an initial surface or curve until it extends sufficiently away from the original surface, resulting in the meshing of the entire space surrounding the geometry, as shown in Figure 5.

As shown in Figure 5.b, the O-grid hyperbolic structured mesh generator [31] has superior resilience and convergence speed compared to C-grids. These meshes exhibit exceptional accuracy near edges and demonstrate efficacy under stress situations, provided the spacing is constant. Although there are concerns with element orthogonality in wake regions and numerical instabilities at narrow trailing edges, these problems may often be mitigated by using blunter or rounder trailing edges. The generator demonstrates excellent performance due to its rapid operation and dependability, especially while navigating strongly convex bends.

Throughout the optimization iteration, following the determination of a new airfoil shape using FFD variables, I deform the surface mesh and extend it to the entire volume mesh using an unstructured deformation algorithm, which employs a hybrid method combining algebraic and linear elasticity approaches [32], effectively handling both large-scale and small-scale perturbations in the geometry during the optimization process, while an alternative approach utilizes the inverse distance weighting function (IDW)[33], as demonstrated in Figure 6, showcasing the comparison between the original and deformed mesh of the NREL Phase VI S809 airfoil.

2.3. The High-Fidelity CFD Solver

The high-fidelity CFD solver which has been used in my study, is designed for aerodynamics and multidisciplinary optimization (MDO). It specializes in parallel computing and efficiently processes various equations, including Euler, laminar Navier-Stokes, and Reynolds-averaged Navier-Stokes (RANS), for different simulation types. Utilizing a structured finite volume method, the solver is compatible with structured multiblock meshes and employs the Spalart-Allmaras turbulence model alongside the Jameson–Schmidt–Turkel (JST) scheme for its computations.

$${\displaystyle \frac{\partial \tilde{v}}{\partial t}}+u_j{\displaystyle \frac{\partial \tilde{v}}{\partial x_j}}=C_{b1}(1-f_{t2})S\tilde{v}+{\displaystyle \frac{1}{\sigma }}\left[\nabla \bullet \left(\right(v+\tilde{v})\nabla \tilde{v})+{C_{b2}(\nabla \tilde{v})}^2\right]-C_{\omega 1}{f}_{\omega }{\left({\displaystyle \frac{\tilde{v}}{d}}\right)}^2$$

The four terms represent the convective, production ,diffusion, and near-wall destruction for the turbulent flow, respectively. where $\tilde{\mathrm{v}}$ is related to the turbulent viscosity ${\mathrm{v}}_{\mathrm{t}}={\mathsf{\mu }}_{\mathrm{t}}/\mathsf{\rho }$ via

$$v_t=\tilde{v}{\displaystyle \frac{X^3}{X^3+C_{v1}^3}} X={\displaystyle \frac{\tilde{v}}{v}}$$

Where is the $\tilde{\mathrm{v}}$ is the modified turbulent viscosity,${\mathrm{u}}_{\mathrm{j}}$velocity component ,$\mathrm{S}$is the magnitude of the vorticity,$\mathrm{v}$ is the molecular kinematic viscosity,$\mathrm{d}$ is the distance to the nearest wall ${\mathrm{C}}_{\mathrm{b}1}$,${\mathrm{C}}_{\mathrm{b}2}$,${\mathrm{C}}_{\mathsf{\omega }1}$, and $\mathsf{\sigma }$ are model constants.${\mathrm{f}}_{\mathrm{t}2}$ and ${\mathrm{f}}_{\mathsf{\omega }}$ empirical damping functions within the model.

The high-fidelity solver employs several numerical schemes for flux calculations and solves compressible flow equations using different turbulence models. It utilizes algorithms such as Runge-Kutta, D3ADI, ANK, and NK to ensure efficient residual convergence. The algorithm specifically designed for overset meshes is particularly effective, while another ensures efficient convergence across all mesh types. Various integration schemes, such as explicit Runge-Kutta, implicit LU-SGS, and DDADI methods, are used depending on the simulation requirements. For turbulence ulfilng, a segregated DDADI scheme is employed to manage factorization errors, enabling accurate and efficient solutions.

2.4. High-Fidelity Gradient-Based with Discrete Adjoint Method

In aerodynamic shape optimization, the adjoint method is efficient in managing multiple design variables, in contrast to traditional finite-difference methods, due to its scalability and reduced computational load. The optimization function (f) relies on design variables (x) that shape the airfoil and flow state variables (w), determined by CFD solvers through residuals R(x,w(x))=0. This function is directly influenced by design variables and indirectly affected by changes in state variables needed to ulfil governing equations. Calculating gradients crucial for gradient-based optimization requires accounting for both dependencies, employing the chain rule for total derivative computation.

${\displaystyle \frac{df}{dx}}={\displaystyle \frac{\partial f}{\partial x}}+{\displaystyle \frac{\partial f}{\partial w}}*{\displaystyle \frac{dw}{dx}} \left(2.3\right)$

Here, explicit function derivatives represent partial derivatives, which can be computed analytically at a low computational cost, whereas obtaining the total derivative necessitates iterative solution of the governing equations; ensuring that the total derivative of the CFD residuals R(x,w(x))=0 with respect to x remains zero for a feasible solution allows determination of the complete derivative of the flow solution regarding design variables dw/dx, with the equations being nonlinear and involving millions of state variables, thus necessitating specialized iterative solvers.

${\displaystyle \frac{dR}{dx}}={\displaystyle \frac{\partial R}{\partial x}}+{\displaystyle \frac{\partial R}{\partial w}}*{\displaystyle \frac{dw}{dx}}=0 \left(2.4\right)$

The equation provided facilitates the formulation of a linear system, the solution of which yields the total derivative of the flow solution with respect to the design variables.

${\displaystyle \frac{\partial R}{\partial w}}*{\displaystyle \frac{dw}{dx}}=-{\displaystyle \frac{\partial R}{\partial x}} \left(2.5\right)$

I Substituting the Jacobian's solution from the given equation into the total derivative equation Equations (2.3) yields the following result:

${\displaystyle \frac{df}{dx}}={\displaystyle \frac{\partial f}{\partial x}}-{\displaystyle \frac{\partial f}{\partial w}}{\left[{\displaystyle \frac{\partial R}{\partial w}}\right]}^{-1}{\displaystyle \frac{\partial R}{\partial x}} \left(2.6\right)$

The equation solely comprises partial derivative components, computable analytically with minimal computational expense, and the resolution of the linear system can be accomplished either by computing the solution Jacobian dw/dx from Equations (2.4) or by addressing the adjoint system.

${\left[{\displaystyle \frac{\partial R}{\partial w}}\right]}^T\Psi =-{\left[{\displaystyle \frac{\partial f}{\partial w}}\right]}^{T }\left(2.7\right)$

The adjoint vector $\mathsf{\Psi }$ derived from Equation (2.7) can be applied in the total derivative equation, Equation (2.3). This method proves computationally efficient, particularly in 3D aerodynamic shape optimization (ASO), where typically there are fewer than 10 functions of interest, compared to potentially hundreds of design variables.

In summary, a discrete adjoint involves four main steps:

1)

Compute ${\left[{\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{w}}}\right]}^{\mathrm{T}}$and ${\left[{\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{w}}}\right]}^{\mathrm{T}}$

2)

Solve the adjoint equation ${\left[{\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{w}}}\right]}^{\mathrm{T}}\mathsf{\Psi }=-{\left[{\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{w}}}\right]}^{\mathrm{T}}$to obtain $\mathsf{\Psi }$

3)

Compute ${\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{x}}}$and ${\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{x}}}$

4)

Compute the total derivative ${\displaystyle \frac{\mathrm{d}\mathrm{f}}{\mathrm{d}\mathrm{x}}}={\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{x}}}-{\mathsf{\Psi }}^{\mathrm{T}}{\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{x}}}$.

In summey for above mentioned here summeries in this Figure 7 that show optimization loop using the SLSQP optimizer flowchart , where design variables x are used to update the geometry and mesh, feeding into the flow solver that solves the RANS equation R(x,w(x))=0 to determine state variables w. The adjoint solver then computes the total derivative df/dx using the equations $\mathsf{\Psi }=-{\left[{\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{w}}}\right]}^{\mathrm{T}}{\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{w}}}$ and ${\displaystyle \frac{\mathrm{d}\mathrm{f}}{\mathrm{d}\mathrm{x}}}={\displaystyle \frac{\partial \mathrm{f}}{\partial \mathrm{x}}}-{\mathsf{\Psi }}^{\mathrm{T}}{\displaystyle \frac{\partial \mathrm{R}}{\partial \mathrm{x}}}$ , which are fed back into the optimizer to refine the design variables.

2.5. Optimization Algorithm

In this research, I utilize the Sequential Least Squares Programming (SLSQP) algorithm [34] for all optimization tasks. This algorithm is specifically designed to address nonlinear optimization problems with both equality and inequality constraints, aiming to identify the minimum or maximum of a function subject to these constraints. It achieves this by approximating the objective function and constraints, integrating methods from least squares and quadratic programming. The SLSQP algorithm iteratively refines the solution using the Han-Powell quasi-Newton method, updating the B-matrix according to the BFGS formula and incorporating an L1-test function within its step-length algorithm, which is associated with the reciprocal of the Hessian matrix. The stopping criterion for this algorithm is based solely on the objective function, terminating when the absolute difference between the objective function value at the current iteration (M^(k)) and the previous iteration (M^(k-1)) is less than the product of a predetermined tolerance (e^M) and the current value. For the purposes of our research, this tolerance has been set to 10^-7.

3. Results

3.1. Flow Field Numerical Simulation Validation

The validation of flow field numerical simulation is conducted through a rigorous comparative analysis of its results against empirical data and established computational fluid dynamics (CFD) methodologies. This validation process is essential for ensuring the accuracy and reliability of aerodynamic simulations. The study utilizes an O-type computational domain, with the mesh settings specifying an outer sub-domain radius of 20c, as illustrated in Figure 7. The initial mesh point is positioned at 1.236 × 10^-5 meters from the airfoil, achieving a y+ value of less than 1, which results in a total mesh comprising approximately 67,304 cells. The flow field analysis is performed on the S809 airfoil across angles of attack (AOA) ranging from 0° to 16° at a Reynolds number of 2 × 10⁶.

As depicted in Figure 9.a, the numerical flow simulation demonstrates comparable lift coefficient results that closely align with experimental data up to an AOA of 9 degrees. Beyond this point, the experimental data exhibits a divergence, indicating flow separation that is not accurately captured by the simulations. Figure 9.b illustrates the corresponding trends in drag coefficients, where both solvers exhibit agreement until 9 degrees AOA. Subsequently, the numerical flow simulation indicates a marked increase in drag relative to the experimental observations, suggesting potential overestimation of drag or inadequacies in capturing post-stall flow dynamics.

Both solvers employ numerical simulations to solve the Navier-Stokes equations and the energy equation, utilizing the Spalart-Allmaras turbulence model to predict airflow behavior around the S809 airfoil. The numerical results are validated by comparing the outputs from both solvers against experimental data. Figure 9 through 12 provide visual representations of pressure and velocity contours, streamlines, and pressure variations on the airfoil surface at attack angles of 0.008°, 4.99153°, 9.0275°, and 16.0127°. The comparative analysis indicates that both numerical flow simulations yield closely aligned predictions across the examined angles, thereby corroborating the efficacy of the numerical methods as reliable tools for aerodynamic simulations.

Figure 8. Structured O-Grid for NREL S809 Airfoil.

Figure 8. Structured O-Grid for NREL S809 Airfoil.

Figure 9. Comparative Analysis of S809 Airfoil Using Fluent and Flow Numerical Solvers at Re=2e6 with Experiment Data.

Figure 9. Comparative Analysis of S809 Airfoil Using Fluent and Flow Numerical Solvers at Re=2e6 with Experiment Data.

Figure 10. The Pressure Contour at Different Angle of Attack by (a-d) Flow Numerical Solver,(e-h) Fluent Solver.

Figure 10. The Pressure Contour at Different Angle of Attack by (a-d) Flow Numerical Solver,(e-h) Fluent Solver.

Figure 11. Velocity Contour at different angle of attack by (a-d) Flow Numerical Solver ,(e-h) Ansys Fluent Solver.

Figure 11. Velocity Contour at different angle of attack by (a-d) Flow Numerical Solver ,(e-h) Ansys Fluent Solver.

Figure 12. Velocity Streamlines at Different Angles of Attack by (a-d) Flow Numerical Solver,(e-h) Ansys Fluent Solver.

Figure 12. Velocity Streamlines at Different Angles of Attack by (a-d) Flow Numerical Solver,(e-h) Ansys Fluent Solver.

Figure 13. Pressure Coefficient Distribution: Flow Numerical Solver vs Ansys Fluent Solver at Different Angle of Attack.

Figure 13. Pressure Coefficient Distribution: Flow Numerical Solver vs Ansys Fluent Solver at Different Angle of Attack.

3.2. Airfoil Shape Optimization Design Result

This study presents solutions to the identified optimization challenges and critically evaluates their effectiveness. Our methodology encompasses a series of increasingly complex design iterations aimed at refining the airfoil contour. This refinement considers both the angle of attack and shape parameters across a range of wind velocities. Through this comparative analysis, we seek to enhance our understanding of airfoil shape optimization techniques, specifically focusing on the S809 Airfoil. For this purpose, we have developed a Free Form Deformation (FFD) model, as illustrated in Figure 4.b, and constructed the corresponding volume mesh, depicted in Figure 5.b. Detailed insights into this process are provided in Figure 2, which outlines SLSQP optimization methodology employed in our analysis.

As indicated in Table 2, I present the operational condition parameters selected for both single and three cases of multipoint optimization. Our discussion will begin by examining the results obtained from a single flow condition and will then progress to explore four distinct flow scenarios. In defining thermodynamic conditions, I have employed two distinct methodologies. The first approach involves parameters such as velocity, Reynolds number, Reynolds length, and temperature, while the second approach incorporates velocity, ambient pressure, and ambient temperature. It is worth noting that both approaches are applicable for low-speed flows. Additionally, supplementary quantities are determined using a foundational methodology

Despite the utilization of differing methodologies, our findings consistently demonstrate alignment. In the comprehensive study presented in the accompanying Table 3, the configuration parameters utilized in the optimization of the S809 airfoil shape are delineated, employing SLSQP method for effective optimization. The designated objective function was the minimization of the drag coefficient, meticulously modulated through FFD shape control variables alongside angle of attack adjustments. The optimization was subjected to an array of constraints, encompassing lift, geometric volume, and thickness specifications, Additionally, constraints are placed on the leading edge (LE) and trailing edge (TE) within the FFD box. These constraints are applied exclusively during the single-point and multipoint airfoil shape optimization processes. They ensure that corresponding points on the LE and TE are moved symmetrically, with each pair shifting by an identical measure but in reverse directions, to ascertain the airfoil’s aerodynamic viability and structural feasibility. A systematic augmentation in the complexity of design parameters was methodically executed, extending from a single-point optimization scenario to multiple-point considerations. This incremental complexity aimed to address a spectrum of operational flow conditions, each quantitatively characterized by distinct variables and constraints. Such an approach ensures a holistic optimization strategy, wherein the resulting airfoil configurations are anticipated to demonstrate enhanced aerodynamic efficiency within the defined operational envelope.

3.2.1. Single Point Optimization

This study examines a single point flow condition outlined in Table 2 and in Table 3 which summarizing the optimization problem aimed at minimizing the drag coefficient, involving 41 FFD design variables along with the angle of attack to control the airfoil shape, as depicted in Figure 4. There are 224 constraints in total including constraining the lift coefficient to be 0.4,keeping the local airfoil thickness between 0.85% and 1.25% of the baseline thickness, and ensuring the total volume of the optimized wing is at least equal to that of the baseline wing. Additionally, the leading and trailing edges of the wing are fixed. Where the Free-Form Deformation (FFD) displacement along the y-coordinate is confined within a range of ±0.5, and the angle of attack variable spans from 0 to 16 degrees.

Figure 14 shows the angle of attack design variable and drag coefficient convergence changes history for the optimization iteration .The outcomes demonstrate a decrease in thickness of about -8% and an increase in camber of approximately 18%, as detailed in Table 4, which also indicates a lift-to-drag ratio improvement of 22.25%. Figure 14.(a) demonstrates that the angle of attack changes only slightly during the optimization, suggesting its minimal impact on the results compared to optimizations that exclude it as a variable. Meanwhile, Figure 14.(b) shows that the performance measure (objective function) improves as the optimization iterations progress. In total, I conducted 29 adjoint computations and 59 flow solutions for the optimization

3.2.2 Multiple point optimization

In the previous scenario, it was demonstrated that single-point optimization fails to yield a design that performs satisfactorily across various operating conditions. Consequently, the motivation for employing multipoint optimization arises from the need to encompass a broader spectrum of flow conditions and achieve a more resilient design. In all optimization cases, a high-fidelity gradient-based optimization framework was employed. The design surface was parameterized using the Free-Form Deformation (FFD) method, while volume mesh deformation was managed through an analytic inverse-distance method. The computational fluid dynamics (CFD) simulations were conducted using an adjoint solver with the Spalart-Allmaras (SA) turbulence model.

For local airfoil shape control, 40 FFD control points were used, as shown in Figure 4. Additionally, the angle of attack was incorporated into the optimization for two, three, and four multipoint flight conditions, resulting in 42, 43, and 44 design variables, respectively. Lift coefficients were constrained for each flight condition, and the local airfoil thickness was maintained between 0.85% and 1.25% of the baseline thickness. The total volume of the optimized wing was required to be at least equal to that of the baseline airfoil, with fixed leading and trailing edges. This setup led to 225, 226, and 227 constraints being considered for each case of the multipoint optimization problem, as shown in Table 3. The operating flow condition optimization problem is presented in Table 2. SLSQP is used as the optimizer with accuracy tolerances set to 10^-4. The optimizations run until the tolerances are met or the process aborts due to numerical difficulties.

The design optimization problem parameters and model remain consistent with those of single-point optimizations, albeit with adjustments made to the lift constraint and angle of attack design variable to accommodate the additional flow condition points. The convergence of the objective function history during the optimization iterations for the three flow conditions is apparent, as depicted in Figure 15.a, where the value of the objective function steadily decreases with progressing iterations, indicating effective modification of the airfoil’s shape and resulting in solutions that enhance the performance metric. Additionally, Figure 15.b illustrates the changes in the first angle of attack throughout iterations for two, three, and four multipoint scenarios, revealing a tendency for the angle to approach a minimum value as the number of points increases. Furthermore, with an increase in iteration count, the angle adjustment becomes more refined, signifying a convergence towards optimal solutions. From Figure 15, I can observe that increasing the number of points in the multipoint optimization also increases the number of iterations required to reach the optimal airfoil configuration.

4. Discussion

This study focuses on optimizing the 2D geometry of the S809 airfoil model, considering both single-point and multipoint configurations. To generate 3D meshes in CGNS format, I utilized a foundational methodology, producing meshes with a single cell width and applying symmetry boundary conditions on opposing faces to maintain the quasi-2D nature of the problem . As outlined in Table 4, the optimized airfoil profiles show variations depending on the operational parameters, leading to differences in thickness, maximum thickness location, camber, and maximum camber position. Notably, the optimized airfoils demonstrate a decrease in maximum thickness while shifting towards the leading edge, resulting in reduced drag force. Additionally, improvements in maximum camber and its relocation towards the leading edge enhance lift characteristics, leading to an increase in the lift-to-drag ratio.The aerodynamic coefficients ensuing from the single and multipoint optimization processes are comprehensively delineated in Table 5 I the reference airfoil configuration. These analyses, conducted at a Reynolds number of 3.48e6, illustrate discernibly augmented lift coefficients and cl/cd ratios for the optimized airfoils compared to the baseline S809 airfoil. Noteworthy enhancements in the lift-to-drag ratio, amounting to 22.5%, 22.2%, 21.4478%, and 21.544% for the single-point, two-point multipoint, three-point multipoint, and four-point multipoint optimized airfoils, respectively, further underscore the efficacy of the optimization process.

Next, I examine the shape and pressure (Cp) distributions for both the baseline and optimized geometries in Figure 16, focusing on both single and multipoint optimizations. The optimization process shows a reduction in thickness and a slight increase in camber from single to four-point optimization. This reduction in thickness is expected when prioritizing aerodynamics without considering structural strength constraints. The differences, particularly in thickness and camber parameters, align with Table 4, which indicates a blunt leading edge across all optimization cases. This feature makes the airfoil less sensitive to changes in wind speed and direction, ensuring more consistent performance under varying wind conditions. Additionally, it enhances durability and resistance to damage from debris, rain, and other environmental factors, reducing the risk of leading-edge erosion that can degrade aerodynamic performance over time. Furthermore, a blunt leading edge improves performance at lower wind speeds by delaying flow separation and maintaining smoother airflow over the surface, leading to a higher lift coefficient, which is beneficial for wind turbines in areas with variable wind conditions. A less sharp leading edge also helps reduce noise generated by the wind turbine blades, as sharper edges create more turbulence and higher noise levels. Furthermore, Figure 16.b illustrates the pressure distribution around the airfoil, indicating noticeable improvements resulting from the optimization process. Our simulations were conducted using the Intel i5-8500 CPU (3 GHz), revealing that increased design complexity leads to longer computational times, as evidenced by the significantly higher computational time for the four-point multipoint scenario.Consequent to CFD analyses performed via the Numerical flow field Solver across the entirety of the optimized airfoil configurations, including the baseline, within an angle of attack range spanning from 0 to 20 degrees, flow conditions were established with temperature and static pressure set at standard sea level conditions, specifically 288.15 K and 101325 Pa, respectively. These conditions align with typical operational parameters for wind turbines operating at altitudes not surpassing a few hundred meters. A Reynolds number of 3.48e6 was prescribed at a velocity of 51.48 m/s for all optimized airfoil configurations.

Figure 17 delineates improvements in lift coefficient, lift-to-drag ratio coefficient, drag coefficient, and moment coefficient. Particularly noteworthy is the discernible enhancement in lift-to-drag coefficient observed with the four-point multipoint optimization, especially evident up to an angle of attack of 12 degrees. Conversely, the single-point optimization exhibits superiority in terms of improvement in lift coefficient. Furthermore, Figure 17.c highlights the superior improvement afforded by single-point optimization in the context of lift coefficient.

In Figure 18, the lift-to-drag ratio is displayed for various wind turbine operating conditions, spanning velocities from 5 to 25 m/s and an angle of attack of 5 degrees. Utilizing the flow Numerica solver to demonstrate the lift to drag ratio for various optimization scenarios, including baseline, single-point, two-point, three-point, and four-point optimizations. It is evident that the lift-to-drag ratio improves as the number of optimization points increases, suggesting that multi-point optimization results in superior aerodynamic performance. The results indicate that the four-point optimized airfoil exhibits the highest lift-to-drag ratio across all velocities, showcasing superior overall performance compared to other optimization strategies. This demonstrates the advantages of utilizing advanced optimization techniques to enhance the aerodynamic efficiency of airfoils, a critical factor for the optimal performance of wind turbines within the designated velocity range.

5. Conclusions

This study successfully developed an efficient high-fidelity, RANS-based discrete adjoint method integrated within a constrained shape-design optimization framework, demonstrating its effectiveness through the optimization of the aerodynamic design of the NREL PHASE VI S809 Airfoil Rotor. The optimization framework comprises several key modules, including a geometric parameterization module utilizing Free-Form Deformations (FFDs), an inverse-distance method for volume mesh deformation, a flow simulation solver, a discrete adjoint solver for derivative computation, and a Sequential Quadratic Programming (SLSQP) method for nonlinear constrained optimization.

Through a comprehensive literature review on high-fidelity shape optimization in wind turbine design, this research identified a significant opportunity for advancement in this field. A comparative analysis of two advanced CFD solvers, Fluent and flow field numerical solver, demonstrated robust agreement with experimental data across a variety of angles of attack, validating the reliability of the simulation methods employed.

The study achieved notable results in airfoil optimization. A single-point optimization concerning angle of attack and FFD design variables led to an impressive enhancement in the lift-to-drag ratio by approximately 22.51%. Furthermore, multiple-point optimizations under the same constraints yielded significant improvements in lift-to-drag ratios of approximately 22.2%, 21.45%, and 21.54% for two, three, and four multipoint flow conditions, respectively. These advancements were primarily attributed to a strategic reduction in relative thickness and an increase in camber.

CFD analyses conducted across a range of angles of attack from 0 to 20 degrees for the optimized airfoil shape revealed that the four-point multipoint optimization produced the most substantial improvements in lift-to-drag ratio, lift coefficient, drag coefficient, and moment coefficient. These findings underscore the efficacy of the proposed optimization framework and highlight the potential for further advancements in aerodynamic design within the wind turbine sector.

In conclusion, this research contributes valuable insights to the field of aerodynamic shape optimization, presenting a robust methodology that can be applied to future studies and practical applications in wind turbine design. The demonstrated improvements in aerodynamic performance not only enhance the efficiency of wind turbine rotors but also pave the way for innovative approaches to optimizing complex aerodynamic systems.