An effective approach for robust design optimization of wind turbine airfoils with random aerodynamic variables

Abstract

The robust design optimization of an airfoil needs to continuously realize the probability-based aerodynamic simulation for various combinations of geometry and wind climate parameters. The simulation time is lengthy when a full aerodynamic model is embedded for the numerical iteration. To this end, a second-order polynomial-based response surface model is first presented to relate the airfoil performance indicator with geometry and random aerodynamic variables. This allows to quickly evaluate the response moments and optimization constraints. Then, the robust design optimization is formulated to simultaneously maximize the mean aerodynamic performance and minimize the variance of design results due to the variation of geometry and aerodynamic parameters. The robust design optimization based on the NACA63418 and the DU93-W-210 airfoils with random Mach and Reynolds numbers is presented to demonstrate potential applications of this proposed model. Results have shown that the mean-value aerodynamic indicator is generally improved, whereas the variance is minimized to archive the robust design objective. The proposed approach is simple and accurate, suggesting an attractive tool for robust design optimization of airfoils with random aerodynamic variables.

Keywords

Robust design optimization airfoils aerodynamic performance response surface model wind turbines

Introduction

As a key facility to convert wind energy into electricity, the design and analysis of wind turbines have received considerable attention in recent years.¹ Along with probabilistic methods to quantify fatigue and reliability problems of wind turbine structures,^2–4 the consideration of randomly aerodynamic parameters becomes crucial for the design optimization of wind turbine airfoils.

In reality, numerical optimization of an airfoil is typically divided into several sub-problems, that is, geometry modelling, aerodynamic simulation and robust design optimization. In this regard, the geometry modelling of the airfoil can be realized using a group of basis functions to translate the discrete geometry data into a continuous and smooth profile. The coefficients attached to the basis functions can be treated as design variables for design optimizations of the airfoil.^5,6 Methods for the geometry modelling mainly include the Hicks–Henne (HH) method, the class function/shape function transformation (CST) method and the method based on the conformal transform theory (CTT).^7,8 Even though the CTT method is flexible to represent several airfoils, the resultant trailing edge tends to be a straight line.⁹ In this regard, the CST and the HH methods will be used in this article for geometry modelling of airfoils.

The aerodynamic performance of an airfoil is directly related to the power-conversion efficiency of blades. In this regard, design optimization of an airfoil is usually formulated based on aerodynamic indicators, for example, the lift coefficient $C_{l}$ , the drag coefficient $C_{d}$ and the lift-to-drag ratio $C_{ld}$ , whereas results for the pitching moment coefficient $C_{m}$ are determined in numerical examples to check the aerodynamic stability of optimization results.

The aerodynamic simulation of airfoils can be realized based on the compressible or incompressible Reynolds-averaged Navier–Stokes equation and various turbulence models. However, the resultant computational cost will be rather critical, if a large number of aerodynamic simulations are embedded in the design optimization. In this regard, the design optimization of airfoils is alternatively implemented with the XFOIL package.¹⁰ The design-optimization routine presented in this article, however, is potentially applicable to other aerodynamic simulation packages, once a computational procedure for the airfoil aerodynamic characteristic prediction is numerically available.

To further alleviate the computational cost, the design optimization is usually realized based on a surrogate model that links design variables (i.e. geometry parameters) of an airfoil with the aerodynamic performance indicator. Even though the aerodynamic simulation can be directly carried out by means of the computational fluid dynamics (CFD)/XFOIL model, one has to run the full-scale simulation model several times. A surrogate model based either on the Kriging or the artificial neutral network was alternatively used in the literature.^11–15 The article proposes to utilize the second-order response surface model (RSM) to mimic the aerodynamic performance, and its applicability is validated using the aerodynamic optimization of the NACA63418 and the DU93-W-210 airfoils.

Rather than assuming that the aerodynamic parameters are deterministic, measurements for $C_{l}$ and $C_{d}$ have shown uncertainties in the wind tunnel test.^16,17 This motivates the use of a probabilistic model to quantify the design optimization of airfoils in this article. The influence of uncertain geometry parameters on the lift and drag coefficients was investigated by Ernst et al.,¹⁸ whereas Sørensen and Toft¹⁹ presented a risk-informed framework for probability-based design optimization of wind turbine structures. Zhao et al.²⁰ considered the design optimization of the NASA0412 airfoil based on the aerodynamic characteristic that is simulated with the random Mach number $M_{a}$ . As one of the important factors that affects the aerodynamic performance of an airfoil, observational data have justified that $R_{e}$ is usually fluctuated within a range,²¹ whereas the average wind speed has been modelled as the Weibull or the Rayleigh distributed random variables in the IEC 61400-1 Standard.²² Therefore, the predicted aerodynamic indicators, for example, $C_{l}$ , $C_{d}$ and $C_{ld}$ will become random variables due to the uncertain aerodynamic parameters $R_{e}$ and $M_{a}$ .^23,24 To summarize, the design optimization of an airfoil needs to be improved by emphasizing on many uncertain factors.^12,25,26 An approach for probability-based robust design optimization of airfoils is presented in this article.

The objective of this article is to present an effective approach for the robust design optimization of wind-turbine airfoils with uncertain Mach and Reynolds numbers. To implement, the geometry modelling of an airfoil is realized based on the improved Hicks-Henne (IHH) method, whereas a second-order RSM is used to mimic the implicit relation between design variables and aerodynamic indicators. A robust optimization model that maximizes the lift-to-drag ratio and minimizes the uncertain variation of the aerodynamic characteristic $C_{ld}$ is presented. Numerical examples based on the NACA63418 and the DU93-W-210 airfoils are presented to demonstrate potential applications of the proposed approach. Note that the standard coordinate data documented in the literature²⁷ are used to simulate initial result of an aerodynamic indicator for the design optimization of the airfoil.

The rest of the manuscript is organized as follows. Section ‘Geometry modelling of airfoils’ presents a brief summary on the geometry modelling of airfoils. The utility of the second-order RSM to mimic the true aerodynamic characteristic of an airfoil is presented in section ‘A predictive model for the aerodynamic indicator of an airfoil’. With the proposed robust design approach in section ‘Robust design optimization of airfoils’, potential applications of this approach are demonstrated by the design optimization of the NACA63418 and the DU93-W-210 airfoils. Conclusions are summarized in section ‘Conclusion’.

Geometry modelling of airfoils

The geometry modelling is a crucial step towards the design optimization of wind turbine airfoils. In this regard, the original discrete airfoil data can be represented by continuous and smooth curves to facilitate subsequently aerodynamic simulations. Basically, the CST method is commonly used for the airfoil modelling due to its strong approximation ability, whereas the IHH method is suitable for a refinement design of airfoils.²⁸ Specially, an S-shaped trailing edge is desired to improve aerodynamic performance of an airfoil. Note that numerical accuracy of a geometry modelling method is highly related to the investigated airfoil series, and the approximation based on a large number of basis functions are always necessary. To this end, with a brief summary on the geometry modelling technique, the performance of the CST and the IHH methods are examined by a global measure error $ε$ for approximation results for the NACA63418 and the DU93-W-210 airfoils as follows.

The CST method

The origin of the CST method can be traced back to the work by Kulfan.²⁹ It represents the geometry of an airfoil by means of the class and the shape functions. In this regard, the discrete airfoil data on upper and lower surfaces can be generally represented as

{\hat{y}}_{i} (x) = C (x) S_{i} (x) (i = u, l)

(1)

where the subscript $i$ $(i = u, l)$ represents an $i th$ surface of the airfoil.

The class function $C (x)$ here is used to define the geometry of the airfoil

C (x) = x^{n_{1}} (1 - x)^{n_{2}}

(2)

where $n_{i}$ $(i = 1, 2)$ are geometry parameters.

As discussed in the literature,²⁹ the parameter $n_{1} = n_{2} = 1$ realizes either the biconvex airfoil or the ogive body, whereas $n_{1} = n_{2} = 0.75$ is used to define the radius distribution of a Sears–Haack body. In numerical examples, the NACA63418 and the DU93-W-210 airfoils have a semi-circle head and a pointed tail. This determines that the geometry parameters are given as $n_{1} = 0.5$ and $n_{2} = 1$ , respectively. Therefore, based on a class function to control the general profile of an airfoil, the detailed geometry realization is mainly controlled by shape functions $S_{i} (x)$

S_{i} (x) = \sum_{k = 0}^{n} ψ_{ik} P_{k} (x) (i = u, l)

(3)

Herein, parameters $ψ_{ik}$ denotes the expansion coefficients, and the Bernstein polynomials $P_{k} (x)$ are defined as

P_{k} (x) = \frac{n!}{k! (n - k)!} x^{k} (1 - x)^{n - k} (k = 0, 1, \dots, n)

(4)

Substituting for true coordinates of an airfoil, results for unknown parameters $ψ_{ik}$ can be determined based on the least-squares regression analysis. The accuracy of the CST method is highly related to the truncation order $n$ in equation (3), and approximation results based on several realizations of the truncation parameter $n$ are presented as follows.

The IHH method

Rather than directly utilizing all polynomials to represent the airfoil data, the IHH method realizes the airfoil approximation using a perturbation function around the true airfoil data. In this regard, the represented coordinates of an airfoil are expressed as

{\hat{y}}_{i} (x) = y_{i} (x) + \sum_{k = 1}^{n} ψ_{ik} P_{k} (x) (i = u, l)

(5)

where $y_{i} (x)$ represents the true coordinates of the airfoil at the position $x$ , whereas $ψ_{ik}$ denotes the coefficient associated with a $k th$ -order perturbation polynomial

P_{k} (x) = {\begin{matrix} x^{0.25} (1 - x) \exp (- 20 x) & if k = 1 \\ \sin^{3} (π x^{e (k)}) & if k \geq 2 \end{matrix}

(6)

Herein, $e (k) = \log (0.5) / \log (x_{k})$ and $x_{k}$ are equidistant points within the interval. Results in the numeric example have shown that coefficients $ψ_{ik}$ in equation (5) are better to be within the interval $ψ_{ik} \in (- 0.005, 0.005)$ . Otherwise, the crossed trailing edge and non-smooth geometries of the airfoil are possibly generated.

In addition, the perturbation polynomial and its first-order derivative are zero-valued at $x = 1$ , that is, $P_{k} (1) = P'_{k} (1) = 0$ as $k \geq 2$ . This results in the geometry curve being almost a straight line near the position $x = 1$ . To overcome this shortcoming, an additional perturbation polynomial as $k = n$ is introduced²⁰

P_{k} (x) = {\begin{matrix} x^{0.25} (1 - x) \exp (- 20 x) & if k = 1 \\ \sin^{3} (π x^{e (k)}) & if 2 \leq k \leq n - 1 \\ α x (1 - x) \exp [- β (1 - x)] & if k = n \end{matrix}

(7)

In the final polynomial term, the coefficient $α$ controls the slope of the polynomial, whereas $\exp [- β (1 - x)]$ makes the component function attenuate rapidly as approaching the trailing edge $x = 1$ . The attenuation velocity is modelled by the parameter $β$ . In this article, numerical parameters are assumed as $α = 8$ and $β = 10$ , respectively.

Several realizations of the trailing edge based on the HH and the IHH basis functions are pictured in Figure 1. It is observed that an introduction of the last component function in equation (7) is effective to model the variation of the trailing edge to some extents. The validation for the geometry model is further presented as follows by considering the NACA63418 and the DU93-W-210 airfoils.

Figure 1.

An illustration of the trailing edge based on the HH and the IHH perturbation functions: (a) the HH method and (b) the IHH method.

Numerical verification

The design optimization expects a small truncation order $n$ for geometry modelling of the airfoil, in order to reduce the total number of design variables. To this end, the parameter $n$ is assumed to vary from $6$ to $12$ with an incremental step $2$ in numerical simulations. The accuracy measure is generally defined as

ε = \int_{x} | {\hat{y}}_{i} (x) - y_{i} (x) | dx (i = u, l)

(8)

which is expressed as an integral of the absolute error over an entire range of the position parameter $x$ .

Numerical results for the NACA63418 and the DU93-W-210 airfoils are separatively developed based on the CST and the IHH methods in conjunction with several realizations of the truncation order parameter $n$ . As reference to the benchmark result, results for the global error defined in equation (8) are depicted in Figures 2 and 3. It is observed that the global error generally decreases with an increase of the truncation parameter $n$ , and a small accuracy error $(\leq 3 \times 10^{- 3})$ has demonstrated the validation of the CST and the IHH methods for the geometry modelling of airfoils. Note that the performance of the IHH method is generally better than that of the CST method for both two investigated airfoils. In this regard, the IHH method with the truncation parameter $n = 8$ is utilized in subsequent examples for the robust design optimization of airfoils with random aerodynamic variables $M_{a}$ and $R_{e}$ .

Figure 2.

Global error for geometry modelling result of the NACA63418 airfoil: (a) the upper surface and (b) the lower surface.

Figure 3.

Global errors for geometry modelling result of the DU93-W-210 airfoil: (a) the upper surface and (b) the lower surface.

A predictive model for the aerodynamic indicator of an airfoil

The design optimization of an airfoil requires a model to link the design parameters with the aerodynamic characteristic predicted by simulation tools.³⁰ In this section, a second-order RSM with cross-terms is used to approximate the true but computationally demanding aerodynamical model due to its simplicity in numerical realizations.³¹ Its performance is examined against the linear and ordinary second-order polynomial model by representing the aerodynamic coefficients $C_{l}$ and $C_{d}$ of the NACA63418 and the DU93-W-210 airfoils.

The RSM

The RSM has been widely applied to represent complex input–output relations in many engineering fields. The number of training samples for a reliable estimation results needs to be positively proportional to the highest polynomial order of the surrogate model.³² This motivates to use low-order polynomial models based on the efficiency concern in subsequent design optimizations. In this regard, the quadratic RSM including mixed terms is presented as follows

z = a + \sum_{i = 1}^{m} b_{i} u_{i} + \sum_{i = 1}^{m} c_{i} u_{i}^{2} + \sum_{i = 1}^{m - 1} \sum_{j > i}^{m} d_{ij} u_{i} u_{j} + ϵ

(9)

Herein, the variable $z$ represents an aerodynamic indicator, that is, $C_{l}$ , $C_{d}$ and $C_{ld}$ of the airfoil, whereas the vector $u$ contains all design parameters $ψ_{ik}$ (as $k = 1, \dots, n$ , $i = u, l$ ) and random variables in the design optimization. Therefore, totally $l = 2 m + 2$ $(m = 8)$ variables are considered in the surrogate model.

Or in a matrix form, the second-order RSM is expressed as

z = u^{T} β + ϵ

(10)

Following the general procedure of the design of experiment, $p$ realizations of the aerodynamic indicator $z = η (\cdot)$ are determined as a vector $z = [z_{1}, \dots, z_{p}]^{T}$ , together with results for explanatory variables

ξ = {[\begin{matrix} 1 & u_{11} & \dots & u_{1 l}^{2} & \dots & u_{1 (l - 1)} u_{1 l} \\ 1 & u_{21} & \dots & u_{21}^{2} & \dots & u_{2 (l - 1)} u_{2 l} \\ ⋮ & ⋱ & ⋱ & ⋮ & ⋱ & ⋮ \\ 1 & u_{p 1} & \dots & u_{p 1}^{2} & \dots & u_{p (l - 1)} u_{pl} \end{matrix}]}_{p \times q}

(11)

where, the parameter $q$ is defined as $(l + 1) (l + 2) / 2$ .

Unknown coefficients in equation (10) can be calculated with the training data $ξ$ and the model response vector $z$ as

\hat{β} = (ξ^{T} ξ)^{- 1} ξ^{T} z

(12)

as well as the covariance matrix

Cov (\hat{β}) = σ^{2} (ξ^{T} ξ)^{- 1}

(13)

Here, the symbol $σ^{2}$ denotes the global variance of the residual error $ϵ$ . Its unbiased estimator is given as

{\hat{σ}}^{2} = \frac{z^{T} (I_{p} - H) z}{q - p}

(14)

where $I_{p}$ denotes a $p \times p$ identity matrix, and $H = ξ (ξ^{T} ξ)^{- 1} ξ^{T}$ .

Once the training matrix $ξ$ defined in equation (11) and the corresponding aerodynamic response samples $z$ are available, the minimization for the Euclidean-norm of residual errors $ϵ^{(i)}$ (with $i = 1, \dots, p$ ) allows deriving the polynomial-based surrogate model as

\hat{z} = u^{T} \hat{β}

(15)

which is used to mimic the originally true but computationally intensive aerodynamic model response $z = η (u)$ in subsequent design optimizations.

Numerical verification

To verify the effectiveness of the quadratic RSM for the airfoil’s aerodynamic simulation, results for the NACA63418 and the DU93-W-210 airfoils are presented in this section. To implement, the IHH method with the truncation parameter $n = 8$ is used to fit the geometry data. Assume that the angle of attack $(AOA)$ is $7 °$ , whereas Mach and Reynolds numbers follow the normal distribution, that is, $M_{a} ~ N (0.225, 0.025)$ and $R_{e} ~ N (2.75 \times 10^{6}, 0.25 \times 10^{6})$ . Results for the surrogate model of $C_{l}$ and $C_{d}$ are separatively developed for the NACA63418 and the DU93-W-210 airfoils. As reference to the mechanistic model $z = η (u)$ provided by the XFOIL package, the accuracy of the surrogate model $\hat{z}$ can be evaluated by the root-of-mean-square error (RMSE)

RMSE : = \frac{1}{m \bar{z}} \sqrt{\sum_{i = 1}^{m} {(z_{i} - {\hat{z}}_{i})}^{2}}

(16)

where, the mean value of the prediction result is defined as $\bar{z} = (1 / m) \sum_{i = 1}^{m} z_{i}$ .

Results for the RMSE of the $C_{l}$ and $C_{d}$ surrogate models are summarized in Table 1. It is rather remarkable that the small statistical error $(RMSE < 4.0 \times 10^{- 4})$ for both the NACA63418 and the DU93-W-210 airfoils validates the high accuracy of the quadratic surrogate model in approximating the originally aerodynamic characteristic function of airfoils.

Table 1.

Results for the RMSE of the proposed surrogate model.

	NACA63418	DU93-W-210
Lift coefficient $(C_{l})$	$1.0688 \times 10^{- 4}$	$7.5414 \times 10^{- 5}$
Drag coefficient $(C_{d})$	$3.1821 \times 10^{- 4}$	$3.6509 \times 10^{- 4}$

RMSE: root-of-mean-square error.

Numerical validation of the surrogate model for $C_{l}$ and $C_{d}$ is further presented in Figures 4 and 5. It is seen that the predicted and the true model responses are all in close agreements for both the NACA63418 and the DU93-W-210 airfoils. The large value of the coefficient of determination $(\geq 0.9863)$ for all investigated cases has verified the applicability of the proposed RSM in subsequent design optimizations.

Figure 4.

Numerical validation of the RSM for the lift coefficient $(C_{l})$ of the NACA63418 and the DU93-W-210 airfoils: (a) the NACA63418 airfoil and (b) the DU93-W-210 airfoil.

Figure 5.

Numeric validation of the RSM for the drag coefficient $(C_{d})$ of the NACA63418 and the DU93-W-210 airfoils: (a) the NACA63418 airfoil (b) the DU93-W-210 airfoil.

Numerical results for $C_{l}$ and $C_{d}$ predicted based on the linear and quadratic models with/without cross-terms are further depicted in Figure 6. Note that only results for the DU93-W-210 airfoil are presented at here for the sake of brevity. The coefficients of determination for the linear model are determined as $0.988$ and $0.842$ , respectively, for the lift and the drag coefficients. Even though the quadratic model without cross-terms can improve results as $0.992$ and $0.892$ , respectively, the relatively small value of $R^{2}$ for the drag coefficient motivates the use of cross-terms in the quadratic RSM in this article.

Figure 6.

Results for the lift and drag coefficients of the DU93-W-210 airfoil predicted by the linear, the quadratic models with and without cross-terms: (a) results for the lift coefficient $C_{l}$ and (b) results for the drag coefficient $C_{d}$ .

Robust design optimization of airfoils

The aerodynamic performance, for example, the lift-to-drag ratio $C_{ld}$ of an airfoil directly affects the aerodynamic performance of wind turbine blades in engineering reality. However, the wind turbine structure always suffers from uncertain loading information, for example, the random turbulence and the wind speed variables. To account for the input uncertainty, the Reynolds number $R_{e}$ and the Mach number $M_{a}$ are further assumed as normal distribution random variables in this article. The robust design optimization of airfoils is presented as follows. Note that XFOIL package is used to simulate the performance of an airfoil with various geometry and aerodynamic parameters. The proposed approach can be extended to other simulation packages for robust design optimization of airfoils.

The optimization model

Due to the input uncertainty, an aerodynamic performance indicator, for example, the lift or the drag coefficient and the lift-to-drag ratio $C_{ld}$ becomes a stochastic variable

Z = η (ψ, R_{e}, M_{a})

(17)

where, the function $η (\cdot)$ represents the model used to realize the aerodynamic simulation of an airfoil, and the response quantity $Z$ is an implicit function of geometry parameters $ψ$ and the aerodynamic variables $R_{e}$ and $M_{a}$ .

Due to the variation of aerodynamic variables, the performance function $η (ψ, R_{e}, M_{a})$ also has variation. Thus, the robustness of a design objective can be achieved by simultaneously ‘optimizing the mean performance $μ_{Z}$ and minimizing the performance variance $σ_{Z}^{2}$ ’. In other words, the goal of the robust design is to find the most insensitive design to the variation of design variables $ψ$ and random aerodynamic parameters $M_{a}$ and $R_{e}$

{\begin{matrix} Minimize : f (μ_{Z}, σ_{Z}^{2}) \\ Subject to : g_{i} (ψ) \leq 0 \end{matrix}

(18)

where $f (μ_{Z}, σ_{Z}^{2})$ is the robust robust design function of the airfoil, and $g_{i} (x)$ denotes the $i th$ design constraint.

Specially, the robust design function $f (μ_{Z}, σ_{Z}^{2})$ depends on the mean and the variance of the airfoil’s aerodynamic function $Z = η (ψ, R_{e}, M_{a})$ . In this regard, the robust design function can be formulated in various ways based on engineering application types with respect to bi-objectives $μ_{Z}$ and $σ_{Z}^{2}$ .

The following are three important robust function types in reality:³³

Nominal-the-best type

f (μ_{Z}, σ_{Z}^{2}) = w_{1} {(\frac{μ_{Z} - Z_{t}}{μ_{0} - Z_{0}})}^{2} + w_{2} {(\frac{σ_{Z}}{σ_{0}})}^{2}

(19)

Herein, $Z_{t}$ and $Z_{0}$ are the target and the initial values of the aerodynamic quantity $Z$ , respectively, and $w_{1}$ and $w_{2}$ are weight to be determined by the designer. To reduce the dimensionality problem of two objectives, each term is normalized by the initial value $μ_{0}$ and $σ_{0}$ .

Smaller-the-best type

f (μ_{Z}, σ_{Z}^{2}) = w_{1} {(\frac{μ_{Z}}{μ_{0}})}^{2} + w_{2} {(\frac{σ_{Z}}{σ_{0}})}^{2}

(20)

Largest-the-best type

f (μ_{Z}, σ_{Z}^{2}) = w_{1} {(\frac{μ_{0}}{μ_{Z}})}^{2} + w_{2} {(\frac{σ_{Z}}{σ_{0}})}^{2}

(21)

The main goal of this article is to maximize the lift-to-drag ratio, yet to minimize the variance of the uncertain aerodynamic response quantity. Therefore, a mathematical model for the probability-based robust design optimization of an airfoil is generally formulated as

{\begin{matrix} Find : ψ^{*} = [ψ_{i 1}^{*}, \dots, ψ_{i 8}^{*}] (i = u, l) \\ Minimize : w_{1} {(\frac{μ_{0} (ψ, R_{e}, M_{a})}{μ_{Z} (ψ, R_{e}, M_{a})})}^{2} + w_{2} {(\frac{σ_{Z} (ψ, R_{e}, M_{a})}{σ_{0} (ψ, R_{e}, M_{a})})}^{2} \\ Subject to : ψ_{ik} \in (- 0.005, 0.005) (k = 1, \dots, 8; i = u, l) \\ t_{l} \leq t_{max} \leq t_{u}; x_{l} \leq x_{t_{max}} \leq x_{u} \end{matrix}

(22)

In this article, $μ_{0} (\cdot)$ and $σ_{0} (\cdot)$ represent the mean and the standard deviation of the lift-to-drag ratio $C_{ld}$ calculated based on the initial values of the design vector $y = [ψ_{i 1}, \dots, ψ_{i 8}]^{T}$ (for $i = u, l$ ), whereas $μ_{Z} (\cdot)$ and $σ_{Z} (\cdot)$ denotes results for an actual realization of the design vector. Besides, parameters $t_{max}$ and $x_{t_{max}}$ represent the maximal relative thickness and the geometry position, respectively. Note that the weight parameter $w_{1}$ is associated with the mean enhancement of the aerodynamic indicator $Z$ , whereas $w_{2}$ is used to control the robustness of design results, and the realization $w_{1} + w_{2} = 1$ is adopted in numerical simulations.³⁴ Once the second-order regression model in equation (9) is available, the robust design optimization can be efficiently realized based on the approximation result of $μ_{Z}$ and $σ_{Z}$ based on the surrogate model, rather than evaluating the computationally demanding model $Z = η (ψ, R_{e}, M_{a})$ .

Results for the NACA63418 airfoil

The robust design optimization of the NACA63418 airfoil is first implemented based on the proposed mathematical model in equation (22). Note that the genetic algorithm (GA) coded in Matlab © is used to locate optimum results.³⁵

With results summarized in Table 2 and Figure 7, it is observed that results for the geometry characteristic $t_{max}$ and the position coordinate $x_{t_{max}}$ are slightly reduced to some extents as reference to the initial design, whereas the maximum camber $y_{max}$ is increased a little from $2.01 \times 10^{- 2}$ to $2.15 \times 10^{- 2}$ , which is usually beneficial to an improvement of the aerodynamic performance of the NACA63418 airfoil as justified in Figure 8.

Table 2.

Results for geometry characteristics of the optimized NACA63418 airfoil.

	$t_{max}$	$x_{t_{max}}$	$y_{max}$	$x_{y_{max}}$
Initial design	0.180	0.343	$2.01 \times 10^{- 2}$	0.552
Optimized design	0.178	0.309	$2.15 \times 10^{- 2}$	0.449

Figure 7.

The robust design optimization result for the NACA63418 airfoil.

Figure 8.

Results for the empirical distribution of the maximal lift-to-drag ratio $max (C_{ld})$ before and after the design optimization.

Results for the empirical probability density function for the maximal lift-to-drag ratio $max (C_{ld})$ are further depicted in Figure 8. The mean value of $max (C_{ld})$ has been increased from $79.735$ to $85.288$ , whereas the coefficient of variation (COV) is decreases by $38.76 %$ (from $1.86 %$ to $1.14 %$ ) as reference to the initial design result. The reduction of the COV value implies an improved robustness of the optimization result against the input uncertainty due to random variables $R_{e}$ and $M_{a}$ .

Figure 9 further summarizes results for the lift-to-drag ratio $C_{ld}$ for the optimized NACA63418 airfoil in conjunction with various realizations of the $R_{e}$ value. It is observed that the improvement of the aerodynamic characteristic has covered a wide range of the angle of attack of the airfoil for many realizations of the Reynolds number. The validity of the proposed method is further demonstrated by the lift and the drag coefficient result but are not presented at here for the sake of brevity.

Figure 9.

Results of the lift-to-drag ratio $C_{ld}$ for the initial and the optimized NACA63418 airfoil $(M_{a} = 0.15)$ : (a) $R_{e} = 2.0 \times 10^{6}$ , (b) $R_{e} = 2.5 \times 10^{6}$ , (c) $R_{e} = 3.0 \times 10^{6}$ , and (d) $R_{e} = 3.5 \times 10^{6}$ .

Results for the DU93-W-210 airfoil

Numerical application of the proposed approach for the robust design of airfoils is to further extended to the DU93-W-210 airfoil, where the constraint values for $t_{max}$ and $x_{t_{max}}$ in equation (22) are given as $20.8 % \leq t_{max} \leq 21.2 %$ and $0.30 \leq x_{t_{max}} \leq 0.40$ , respectively. Based on the GA algorithm to realize the global optimization, the robust design-optimization result for the DU93-W-210 airfoil is presented as shown in Figure 10.

Figure 10.

Result for the robust design optimization of the DU93-W-210 airfoil.

Figure 10 summarizes the initial and the optimization results for the DU93-W-210 airfoil, whereas the geometry characteristics are presented in Table 3. Besides, the empirical distribution of the maximallift-to-drag ratio simulated based on the initial and designed airfoils are pictured in Figure 11. It is observed that the mean-value of the lift-to-drag ratio has been improved by 6.93%, whereas the variance is approximately reduced by $34.5 %$ , as compared to the initial design result. This is the objective of the robust design optimization procedure by maximizing the mean value aerodynamic indicator and minimizing its uncertain variation given random input variables $M_{a}$ and $R_{e}$ .

Table 3.

The optimized geometry characteristics of the DU93-W-210 airfoil.

	$t_{max}$	$x_{t_{max}}$	$y_{max}$	$x_{y_{max}}$
Initial design	0.2069	0.341	$2.736 \times 10^{- 2}$	0.711
Optimized design	0.2075	0.301	$2.904 \times 10^{- 2}$	0.749

Figure 11.

Results for the empirical distribution of the maximal lift-to-drag ratio $max (C_{ld})$ .

Figure 12 presents the relative ratio of the lift and the drag coefficients for the optimized and the initial results. Herein, the subscript $(0)$ denotes the aerodynamic result simulated based on the initial geometry of the airfoil, whereas the symbol with the subscript (*) represents the robust optimization result. It is observed that the optimized lift coefficient results have been uniformly improved, and the drag coefficients for most of investigated AOA values have been further minimized with relative ratios $C_{d}^{*} / C_{d}^{0} \leq 1$ .

Figure 12.

Results for relative coefficients $C_{l}^{*} / C_{l}^{0}$ and $C_{d}^{*} / C_{d}^{0}$ of the DU93-W-210 airfoil: (a) the ratio for coefficient of lift and (b) the ratio for coefficient of drag.

Figure 13 depicts results for the pitching moment coefficient $C_{m}$ of the optimized DU93-W-210 airfoil. The negative valued $C_{m}$ confirms the aerodynamic stability of the robust optimization result. To summarize, the effectiveness of the proposed approach has been further justified by numerical robust design optimization of the DU93-W-210 airfoil.

Figure 13.

Results of the pitching moment coefficient $C_{m}$ for the optimized DU93-W-210 airfoil.

Conclusion

The article presents an effective approach for the robust design optimization of wind turbine airfoils with random Reynolds and Mach numbers. The geometry property of the airfoil is modelled by the IHH method. And, a second-order polynomial model with cross-terms is used to predict aerodynamic characteristics of the airfoil. With a design optimization model by simultaneously maximizing the mean-value and minimizing the variance of an aerodynamic indicator, the robust design optimization of wind turbine airfoils is demonstrated by considering the NACA63418 and the DU93-W-210 airfoils.

Numerical results have shown that the implicit input–output aerodynamic relation of the airfoil can be accurately represented using the quadratic polynomial model with cross-terms. With the optimized airfoil, the mean-value of the lift-to-drag ratio can be improved by $6.96 %$ and $14.28 %$ for the NACA63418 and the DU93-W-210 airfoils, respectively, whereas the variance of $max (C_{ld})$ has been reduced by $38.76 %$ and $47.49 %$ . The small variance value guarantees the high robustness of optimized airfoil compared to the initial design result. To summarize, numerical simulation results have justified the effectiveness of the proposed methodology for the robust design optimization of wind turbine airfoils subjected to random aerodynamic variables.

Footnotes

Handling Editor: James Baldwin

Declaration of Conflicting Interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The authors would like to appreciate the Chinese National Natural Science Foundation (grant nos. 51775095 and 51405069), the Postdoctoral Science Foundation of China (grant no. 2017M610182) and the Fundamental Research Funds for Central Universities (grant no. N170308028) for financially supporting the research.

ORCID iD

Xufang Zhang

References

Sørensen

. Framework for risk-based planning of operation and maintenance for offshore wind turbines. Wind Energy 2009; 12: 493–506.

Zhu

Liu

Peng

, et al. Computational-experimental approaches for fatigue reliability assessment of turbine bladed disks. Int J Mech Sci 2018; 142–143: 502–517.

Qian

Lei

. A statistical model of fatigue failure incorporating effects of specimen size and load amplitude on fatigue life. Philosoph Magaz 2019; 99: 2089–2125.

Qian

Lei

, et al. Statistical size scaling of breakage strength of irregularly-shaped particles. Theor Appl Fract Mech 2019; 102: 51–58.

Selig

McGranahan

. Wind tunnel aerodynamic tests of six airfoils for use on small wind turbines. J Solar Energy Eng 2004; 126: 986–1001.

Fuglsang

Bak

. Development of the RIS

Φ

wind turbine airfoils. Wind Energy 2004; 7: 145–162.

Masters

Taylor

Rendall

, et al. Review of aerofoil parameterisation methods for aerodynamic shape optimisation. In: Proceedings of the 53rd AIAA aerospace sciences meeting, 2015, https://arc.aiaa.org/doi/abs/10.2514/6.2015-0761

Balu

. Recent trends in geometry parameterization strategies for aerodynamic shape optimization. J Astron Aerosp Eng 2017; 3, http://matjournals.in/index.php/JOAAE/article/view/1168

Chen

Wang

Pang

, et al. Improvement of airfoil design using smooth curvature technique. Renew Energy 2013; 51: 426–435.

10.

Drela

. XFOIL: an analysis and design system for low Reynolds number airfoils. Berlin; Heidelberg: Springer, 1989.

11.

Gao

, et al. Aerodynamic optimization of wind turbine airfoils using response surface techniques. Proc IMechE Part A: J Power and Energy 2010; 224: 827–838.

12.

Ribeiro

AFP

Awruch

Gomes

. An airfoil optimization technique for wind turbines. Appl Math Model 2012; 36: 4898–4907.

13.

Chen

Wang

Zhang

, et al. A new direct design method of wind turbine airfoils and wind tunnel experiment. Appl Math Model 2016; 40: 2002–2014.

14.

Zhang

Wang

Sørensen

. REIF: a novel active-learning function toward adaptive Kriging surrogate models for structural reliability analysis. Reliab Eng Syst Saf 2019; 185: 440–454.

15.

Zhang

Wang

Sørensen

. AKOIS: an adaptive Kriging oriented importance sampling method for structural system reliability analysis. Struct Saf 2020; 82: 101876.

16.

Toft

Svenningsen

Moser

, et al. Assessment of wind turbine structural integrity using response surface methodology. Eng Struct 2016; 106: 471–483.

17.

Ambühl

Sørensen

. Sensitivity of risk-based maintenance planning of offshore wind turbine farms. Energies 2017; 10: 505.

18.

Ernst

Schmitt

Seume

. Effect of geometric uncertainties on the aerodynamic characteristic of offshore wind turbine blades. J Phys, https://iopscience.iop.org/article/10.1088/1742-6596/555/1/012033

19.

Sørensen

Toft

. Probabilistic design of wind turbines. Energies 2010; 3: 241–257.

20.

Zhao

Gao

Z-h

Huang

J-t

. Robust design of natural laminar flow supercritical airfoil by multi-objective evolution method. Appl Math Mech 2014; 35: 191–202.

21.

Abdallah

Natarajan

Sørensen

. Impact of uncertainty in airfoil characteristics on wind turbine extreme loads. Renew Energy 2015; 75: 283–300.

22.

Wind turbines-part 1: design requirements. Geneva: International Electrotechnical Commission, 2005.

23.

Chehouri

Younes

Ilinca

, et al. Review of performance optimization techniques applied to wind turbines. Appl Energy 2015; 142: 361–388.

24.

Kostandyan

Sørensen

. Reliability assessment of solder joints in power electronic modules by crack damage model for wind turbine applications. Energies 2011; 4: 2236–2248.

25.

Wang

Zhang

Chen

. Robust airfoil optimization based on improved particle swarm optimization method. Appl Math Mech 2011; 32: 1245–1254.

26.

Zhang

Pandey

. An effective approximation for variance-based global sensitivity analysis. Reliab Eng Syst Saf 2014; 121: 164–174.

27.

Bertagnolio

Sørensen

Johansen

, et al. Wind turbine airfoil catalogue, 2001, https://orbit.dtu.dk/files/7728949/ris_r_1280.pdf

28.

Sun

. Wind turbine airfoil design using response surface method. J Mech Sci Technol 2011; 25: 1335–1340.

29.

Kulfan

. Universal parametric geometry representation method. J Aircraft 2008; 45: 142–158.

30.

Barrett

Ning

. Comparison of airfoil precomputational analysis methods for optimization of wind turbine blades. IEEE Trans Sustain Energy 2016; 7: 1081–1088.

31.

Zhang

Pandey

. Structural reliability analysis based on the concepts of entropy, fractional moment and dimensional reduction method. Struct Saf 2013; 43: 28–40.

32.

Hosder

Walters

Balch

. Efficient sampling for non-intrusive polynomial chaos applications with multiple uncertain input variables. In: Proceedings of the 48th AIAA/ASME/ASCE/AHS/ASC structures, structural dynamics, and materials conference, 2007, https://arc.aiaa.org/doi/abs/10.2514/6.2007-1939

33.

Choi

Lee

Gorsich

. Dimension reduction method for reliability-based robust design optimization. Comput Struct 2008; 86: 1550–1562.

34.

Marler

Arora

. Survey of multi-objective optimization methods for engineering. Struct Multidiscip Optimiz 2004; 26: 369–395.

35.

The MathWorks Inc. Global optimization toolbox user’s guide. Natick, MA, https://www.mathworks.com/help/pdf_doc/gads/gads_tb.pdf (accessed 17 May 2018).