See discussions, stats, and author profiles for this publication at: https://www.researchgate.

net/publication/269249043

Metric-Based Mathematical Derivation of Efficient Airfoil Design Variables

Article  in  AIAA Journal · May 2015

DOI: 10.2514/1.J053427

CITATIONS READS

55 428

3 authors:

Daniel J. Poole C. B. Allen

University of Bristol University of Bristol
37 PUBLICATIONS   402 CITATIONS    241 PUBLICATIONS   2,753 CITATIONS   

SEE PROFILE SEE PROFILE

Thomas Rendall
University of Bristol
90 PUBLICATIONS   1,455 CITATIONS   

SEE PROFILE

Some of the authors of this publication are also working on these related projects:

Corinin-1C RNA microinjections View project

SLOWD - Sloshing Wing Dynamics View project

All content following this page was uploaded by Daniel J. Poole on 09 June 2015.

The user has requested enhancement of the downloaded file.

 Metric-Based Mathematical Derivation

of Efficient Airfoil Design Variables∗

Daniel J. Poole† , Christian B. Allen‡ and Thomas C. S. Rendall§

University of Bristol, Bristol, BS8 1TR, U.K.

Within an aerodynamic shape optimization framework, an efficient shape parame-

terization and deformation scheme is critical to allow flexible deformation of the surface

with maximum possible design space coverage. Numerous approaches have been de-

veloped previously for geometric representation of airfoils. A fundamental approach is

considered here from the geometric perspective, and a method presented to allow the

derivation of efficient, generic and orthogonal airfoil geometric design variables. This

is achieved by the mathematical decomposition of a training library. The resulting ge-

ometric modes are independent of parameterization scheme, surface and volume mesh,

and flow solver, thus are generally applicable. However, these modes are dependent on

the training library, and so a benchmark performance measure, the airfoil technology

factor, has also been incorporated into the scheme to allow intelligent, metric-based

filtering, or design space reduction, of the training library to ensure efficient airfoil

deformation modes are extracted. Results are presented for several geometric shape

recovery problems, using two optimization approaches, and it is shown that these math-

ematically extracted degrees of freedom perform particularly well in all cases, showing

excellent design space coverage. These design variables are also shown to outperform

those based on other widely-used approaches; the Hicks-Henne ‘bump’ functions and a

linear (deformative) approximation to the PARSEC parameterization are considered.

∗ A version of this paper was originally presented as AIAA Paper 2014–0114, 10th AIAA Multidisciplinary Design
Optimization Conference, National Harbor, MD, January 2014
† Graduate Student, Department of Aerospace Engineering, University of Bristol, AIAA Student Member
‡ Professor of Computational Aerodynamics, Department of Aerospace Engineering, University of Bristol, AIAA
Senior Member
§ Lecturer, Department of Aerospace Engineering, University of Bristol, AIAA Member

1
 Notation

ai Upper surface PARSEC weighting

bi Lower surface PARSEC weighting

B Number of Hicks-Henne bumps

clx PARSEC lower crest horizontal position

cux PARSEC upper crest horizontal position

clxx PARSEC lower crest curvature

cuxx PARSEC upper crest curvature

cly PARSEC lower crest vertical position

cuy PARSEC upper crest vertical position

CL Lift coefficient

m Number of surfaces in training library

mdef Total number of perturbation vectors

M Mach number

Mdd Drag divergent Mach number

n Number of selected design variables

N Number of surface points

Rle PARSEC leading edge radius

t/c Airfoil thickness to chord ratio

x, y, z Inertial coordinates

αi Design variable scaling factor

ζi Hicks-Henne bump horizontal position factor

θ PARSEC trailing-edge angles

κ Technology factor

ξi Hicks-Henne bump width factor

2
U Modal deformation matrix

V Modal weighting matrix

X Surface coordinate vector [(x1 , . . . , xN )T , (y1 , . . . , yN )T ]T

∆X Training data deformation matrix (2N × mdef )

Σ Modal energy diagonal matrix

I. Background and Introduction

Numerical simulation methods are now used routinely in industrial design, and increasing com-

puter power has resulted in their integration into the optimization process becoming more common

in academic research and industrial design. Integrating an effective geometry control method with

an aerodynamic model and a numerical optimization scheme results in an aerodynamic shape op-

timization approach, which is invaluable during design. The quality and applicability of results

gained by numerical optimization are inherently dependent on the fidelity of the analysis tool used,

and so developing a modular approach means that the complexity and cost of each module required

at various stages of the design process can be changed, in terms of aerodynamic fidelity, number of

surface shape design degrees of freedom, and the complexity of the optimization scheme.

The physics included in the aerodynamic model, and the optimization approach adopted, are

important aspects of any method, and numerous advanced optimizations have been performed using

compressible CFD as the aerodynamic model, for airfoil sections [1, 2], full aircraft [3–5], aeroelastic

aircraft [6], and rotor blades [7–11]. However, the ability for the optimizer to fully interogate the

design space is driven by the ability for the degrees of freedom adopted to represent any shape

within the design space , and so this is a critical aspect of any optimization scheme.

Numerous geometry parameterization and surface control techniques have been developed and

implemented, including constructive methods, such as CST [12] and PARSEC [13], and deformative

methods, such as analytic [1] and discrete [14]. The authors have also presented work in this area,

having developed a domain element-based generic optimization tool, that is flow-solver and mesh

3
type independent, and applicable to any aerodynamic problem [15–18]. The effectiveness of these

various surface control techniques can be measured in terms of being: i) flexible enough to allow

sufficient design space investigation; ii) robust enough to be applicable to any geometry or design

surface; and iii) efficient enough to cover the design space with a small number of design parameters.

The latter is a metric which relates to the computational resources required to perform numerical

aerodynamic shape optimization, particularly if using gradient-based optimizers which require the

computation of the sensitivity of the objective function (which if done using finite-difference is

proportional to the number of design parameters). A popular alternative to finite difference methods

to compute gradients is by solving an adjoint equation [14] such that all gradients are found in one

solution. While this makes the gradient evaluation independent of the number of design parameters

– often an advantage for high-fidelity optimization – large numbers of design parameters can result

in a design space that is highly multimodal and therefore make it difficult for a gradient-based

optimizer to locate a competitive solution. Moreover, if global search optimizers are implemented –

which require large numbers of objective function evaluations – the convergence characteristics are

also related to the number of design parameters. Efficiency of the aerodynamic shape optimization

process therefore necessitates the requirement for a minimum number of design parameters to cover

a maximised design space.

The objective of the work presented here is to introduce a method that can derive a mathe-

matically efficient set of airfoil design parameters for suitable coupling to the aerodynamic shape

optimization process. Using a mathematical technique introduces a more strict and rigorous per-

spective to this process, removing the requirement for user input and intuition. Specifically, the

approach introduced here uses a proper orthogonal decomposition (POD) of a set of airfoil shapes

to obtain design variables, which are then tested for efficiency – which is the extent to which a

surface can be represented with as few design variables as possible – and design space coverage –

which is the extent to which the parameterization can represent any airfoil surface. Furthermore, a

metric-based filtering process, a type of design space reduction via the POD, is also implemented

into the method to ensure maximum efficiency of the extracted shape deformation modes.

4
 II. Surface Representation and Perturbation

Within the aerodynamic shape optimization process, the problems of geometry representation

and perturbation, and the required surface and volume mesh deformations if using a CFD code,

can be considered separate problems, though unified approaches are becoming commonplace, see

for example [19–22].

Methods that consider the geometry creation and perturbation separately are derived from

those that allow a flexible and efficient parameterization of the surface. Perturbation of the base

geometry in the optimization process therefore requires that the new surface be reconstructed which

subsequently requires automatic mesh generation tools for production of a new surface and volume

mesh. For structured multiblock meshes in three dimensions when considering complex geometries

such as a wing with slats and flaps, or a small flow control device on a wing, mesh generation can

require considerable user input and therefore make automatic mesh generation difficult. Though

not impossible, this extra difficulty can make it advantageous to consider approaches that manip-

ulate an existing mesh. An alternative to constructive are deformative methods which unify the

geometry creation and perturbation, which in turn tends to make them simpler to integrate with

mesh deformation tools and allows the use of previously generated meshes – a considerably cheaper

alternative to regeneration – though the mesh deformation is a separate algorithm. A further re-

finement of unifying geometry creation and perturbation is the integration with a mesh deformation

algorithm. Methods of this type typically have some interpolation that describes a link between

the surface and volume, often via a set of control points that are independent of both, such that

deformation of the control points results in deformation of the surface and CFD mesh. Reviews of a

range of parameterization methods have also been presented, and the reader is guided towards the

work of Samareh [23, 24] and Nadarajah [25, 26]. A summary of all of the approaches to geometry

parameterization and control is given in table 1 detailing relative merits of each approach and the

category in which they fit.

To be effective, any parameterization and perturbation technique must be flexible, robust and

efficient such that any shape in the design space can be represented using as few design variables

as possible. These requirements are linked to the suitable selection of design variables for a given

5
Table 1 Comparison of approaches to surface parameterization and control in aerodynamic

shape optimization (D-deformative; C-constructive; U-unified). Efficient dictates number of

design variables required to obtain good design space coverage, coverage relates to variety of

shapes obtainable for small number of design variables (L-low; M-medium; H-high), intuitive

relates to physical identifiability (Y-yes; N-no).

Method Type Efficient/Coverage Intuitive

H/M Y
Analytic[1, 27] D
User set variables for efficiency Define deformation location

M/M N
CST [12] C
Dependent on order, practically often not large coverage Defined by coefficients

L/H Y
Discrete[14] D
Full coverage though at expense of design variables Directly define movement

H/H Y
Domain element[15] U
Hierarchical variables set by user Simple interpolation link

H/H Y
FFD[28] U
Depends on FFD grid though user definable Simple interpolation link

H/L Y
PARSEC [13] C
11 predefined variables, not detailed surface changes Airfoil variables

H/M N
PDE [29] C
Only function BCs, function flexible BCs are non-physical

H/H N
Splines[30] C
User defines efficiency by control points Deformations not related to surface

parameterization technique, though the important characteristic to consider is selecting a set of

design variables that can represent a maximised design space from small surface changes up to

entire planform changes. The discrete method gives full design space coverage by allowing exact

control of every surface point, but leads to large numbers of design variables and makes global

planform changes difficult for an optimization algorithm. The PARSEC method is a specific airfoil

parameterization method defined by a few design variables and can allow large scale surface changes

and therefore alter the global family of the airfoil, though detailed minor local surface changes prove

6
difficult. This problem of various parameterization methods, which have been highlighted in these

two examples, can partly be solved by considering hybridizations or other modifications to such

methods. For example, Zhu and Qin [31] recently presented a method that has the intuitiveness

of PARSEC with the flexibility of the popular CST method and compared their method to both

CST and PARSEC. The issue with hybridized methods often surrounds the large numbers of design

variables, and while in optimization frameworks that use the adjoint approach to compute gradients,

the large number of parameters can be handled, often large numbers of parameters can manipulate

the design space to cause a highly multimodal problem leading to difficulties in locating good optima.

An alternative, however, to surface geometry control is to parameterize based on surface curvature,

such as the approach of Young et al. [32]. This approach may have a more direct effect on the

pressure coefficients as the pressure distribution is related to the curvature hence methods of this

type may have merit. Curvature based surface approximations were tested during this work but

were found not to perform as well as purely geometric surface representations.

Conversely, the more general control point-based methods have the advantage of being able to

use a hierarchical set of design variables, the number of which can be set by the user depending

on the problem, and vary from minor control point position perturbations up to entire groups of

control points moving globally, and the authors have published previous work in this area [15–18].

Further implementations of these types of methods have also been adopted by Zingg and colleagues

[21, 33, 34] and Anderson et al. [35, 36]. However, decoupling the degrees of freedom from the

surface mesh raises the fundamental question of how the design variables are defined. Hence, the

work presented here considers the derivation of fundamental generic airfoil design variables from

a geometric perspective, that can allow both minor surface changes and large scale airfoil family

changes using a small set of design variables. Working in the geometric domain has the advantage

of producing design variables that are independent of the mesh and deformation technique and so

are generally applicable to all necessary parameterization methods.

7
 III. Airfoil Mode Shape Extraction

The derivation of airfoil shape modes, which could lead to shape deformations, has only received

a small amount of attention previously. It is generally accepted that airfoil modes should be orthogo-

nal [37], meaning that each airfoil shape corresponds to a unique set of input parameters; the absence

of orthogonality can lead to the modes having poor design space coverage [25]. Furthermore, the

orthogonality of design variables simplifies the design space considerably and allows for more robust

and faster optimization. Non-orthogonal design variables do still have use, such as the aerofunctions

of Aidala et al. [38] which were derived by inverse design to match certain pressure coefficients, or

the Hicks-Henne ‘bump’ functions [1], though orthogonal modes are preferred. To derive orthogo-

nal modes, either an approach that orthogonalizes non-orthogonal modes, or direct derivation that

produces orthogonal modes, can be used. The post-processing of non-orthogonal modes by Gram-

Schmidt Orthogonalisation has been demonstrated for the derivation of airfoil modes that represent

the NACA 4-series airfoils [39] and NACA supercritical sections [40]. Direct derivation to produce

orthogonal modes could be performed using proper orthogonal decomposition (POD), which was

first introduced by Toal et al. [41] who showed that the design space of the two-dimensional airfoil

shape could be reduced to a few principal modes. Ghoman et al. [42] then showed that POD

could be used to derive airfoil design variables that could represent a specific airfoil family and were

suitable for coupling to the aerodynamic shape optimization problem. However, while the work of

Ghoman et al. is relevant to the current work it is limited to a proof-of-concept study that uses a

set of 10 transonic airfoils to show that transonic shape modes can be used to match a transonic

pressure distribution to within reasonable accuracy.

The aim of the work presented here is to show that a mathematical decomposition approach

using POD can derive generic airfoil design variables based on a training library containing a variety

of airfoil families, and how effective these are at covering the entire airfoil design space. The work is

further considering the effect that suitable selection of data for the training library has by considering

the parameterization of the library based on a performance metric. The suitable selection of data

for the training library is in effect a sampling of the design space for a design space reduction via

the POD, so the investigation of this is paramount in the understanding of this technique.

8
 POD is a method to obtain a low dimensional approximation to a high dimensional space by the

derivation of dominant components, or modes, and can be done by principal component analysis,

Karhunen–Loeve decomposition or singular value decomposition [43], though these methods can be

shown to be equivalent [44]. To derive airfoil design variables, singular value decomposition (SVD)

is the POD method used as it is the simplest to implement for discrete vector data sets. Consider an

airfoil surface parameterized by N surface points, where the i-th surface point has a position in the

space (xi , yi ). A training library of m airfoils has been collated from which the airfoil deformation

modes are extracted. To ensure consistency of the surface description of the training data all airfoils

are parameterized with the same parametric distribution, followed by each airfoil having a rigid

body translation, scaling and then rotation applied to it to map the geometry into a consistent

form where the leading edge is located at the origin and the trailing edge at unit chord along the

horizontal axis. To build the matrix from which SVD is performed, the vector difference of the

i-th surface point between all airfoils needs to be computed, producing mdef = m(m − 1)/2 airfoil

deformations. The x and y deformations are stacked into a single vector of length 2N , for each

airfoil deformation, so a matrix is built of the airfoil deformations which has 2N rows and mdef

columns:

 
 ∆x1,1 · · · ∆x1,mdef 
 
 . .. ..
 ..

 . . 

 
 
∆x · · · ∆xN,mdef 
 N,1 
∆X = 



 ∆y · · · ∆y1,mdef 
 1,1 
 
 . .. ..
 ..

 . . 

 
 
∆yN,1 · · · ∆yN,mdef

Performing a SVD decomposes the matrix into three constituent matrices:

∆X = UΣVT (1)

where U is a matrix of vectors, each of length 2N . The structure is analogous to the decomposed

matrix, so the columns of this matrix are the airfoil mode shapes. Σ is a diagonal matrix of the

9
singular values, arranged in descending order. These can be considered the ‘relative energy’ of the

modes, and represent the ‘importance’ of the mode shapes in the original library. The total number

of possible mode shapes is governed by the number of singular values, which is the minimum of the

number of columns and rows of the decomposed matrix, though the number of modes required to

effectively reconstruct a target airfoil is expected to be considerably fewer. The training library is

parameterized based on deformations and this is an important choice such that design variables that

result from the decomposition are also deformations, ensuring they are independent of the topology

of the airfoils that are used. This allows direct insertion into an aerodynamic shape optimization

framework where deformation of the surface and mesh is important.

IV. Selection of Training Data

Design variables are directly derived from a large library of airfoils, termed the training data,

so it is necessary for the training data to reflect real airfoil design. Little literature is available that

appears to derive design variables from multiple airfoil families, instead simply using a single family

of airfoils – such as the NACA-4 series [39], NACA supercritical airfoils [40], and RAE family [42]

– or using random perturbations in design variables to form a library [41, 45]. The first of these

methods results in design variables that are highly restricted in their overall coverage of the airfoil

design space, and the second has no physical meaning to realistic airfoils.

The derivation of design variables here by POD results in comprehensive modes with charac-

teristics from the training library mapped onto these modes, therefore instigating the demand to

have a training library that is a collection of the most suitable airfoils available. A requirement for

the training library is to have a wide number of airfoils with good performance, and have a link

to physically designed airfoils, which has led to the use of previously designed airfoils that have a

benchmark performance. By including previously designed airfoils, the wealth of knowledge and

experience that has led engineers to design such airfoils will be captured, while also allowing the

flexibility for the derived design variables to cover a large proportion of the physically attainable

design space; it is important that the airfoils that build the training data are a representative sam-

ple of the airfoil design space such that the necessary characteristics are captured by the derived

10
deformation modes. Focussing on physically and intuitively relevant design variables (which might

be thickness and camber, but can also be more general) is an approach that should result in a

comprehensive collection of design variables that are more useful to a designer.

Suitable selection of the training data is imperative to obtain a collection of high performing

design variables. Furthermore, as the POD operator leads to a set of design variables that have

the characteristics of the original training data mapped on to them, then it seems logical to filter

the training data to suit the type of optimization occuring. For example, if the designer wishes to

perform transonic airfoil optimization then a training library can be selected that contains transonic

airfoils and will result in transonic modes. To demonstrate this idea, training library filtering

is performed and subsequent design variables are tested to investigate whether a training library

biased towards a certain performance of airfoils can represent those types of airfoils more than say

no filtering. An unfiltered training library may be more suited to blue-sky airfoil design where the

designer does not know a priori the type of required design.

The performance of the derived design variables are related to the performance of the training

data and so for a filtering process to occur there must be some measurement to act as a training pa-

rameterization. To demonstrate whether a library containing airfoils that perform well at transonic

speeds can better represent transonic type sections a metric approach is proposed to parameterize

the training library in some robust fashion. As the airfoils in the training library are to be a collec-

tion of previously designed airfoils, it seems logical to use a measure that has been used in airfoil

design to determine airfoil performance. Moreover, it is likely, owing to the typical performance of

cruise aircraft, that the optimization of airfoil performance will be driven by the characteristics of

the airfoil in the transonic regime, so a performance measure is sought that describes the transonic

performance of airfoils. A measure that meets this criteria is the airfoil technology factor, κ, which

is described by the Korn equation [46]:

CL t
κ = Mdd + + (2)
10 c

where Mdd is the drag divergence Mach number, CL is the section lift coefficient at the drag diver-
t
gence Mach number, c is the airfoil’s thickness to chord ratio, and typical values of κ are 0.89 for the

11
 NACA-6 series airfoils, and 0.95 for a more advanced supercritical airfoil. A value of κ is therefore

required for as many airfoils as possible leading to a parameterization of the training library, from

which a suitable selection can be chosen for deriving the design variables. The tested airfoils were

all taken from an online database[57].

A. Calculating Drag Divergence Mach Number

The accuracy of the calculation of κ is inherently connected to the accuracy of the estimation of

the drag divergence Mach number as this will also affect the lift coefficient. The calculation of drag

divergence is, however, not a trivial task as it requires many evaluations of a flow solver to obtain an

accurate value. The drag divergence Mach number is the first point, with increasing Mach number,

at which the gradient of the drag coefficient with respect to the Mach number is 0.1[47]:


∂CD 
= 0.1 (3)
∂M M=Mdd

The gradient is found using a central finite-difference formulation. An approach based on

bisection has been implemented using the definition given by equation 3. Overall, the number

of computations required by bisection for a high accuracy of drag divergence Mach number is

considerably smaller than using a sweep approach which would bisect the domain into uniform

intervals and calculate the gradient for all of those intervals.

This method for estimating drag divergence Mach number was validated against known test

data. The first case is a NACA 0015 at CL = 0.0, with validation data from NACA Report 832

[47]. The second is the RAE 2822 at CL = 0.74, with validation from Chapter 8 of Agard Report

256 [46], and Drela [48]. The flow solver used to estimate drag divergence is a structured multiblock

finite-volume, unsteady, inviscid multigrid, upwind code[49, 50] using the flux vector splitting of van

Leer[51].

The results of the validation are given in table 2, showing excellent agreement between the

method used here, and the validation data.

[57] http://www.ae.illinois.edu/m-selig/ads/coord_database.html

12
Table 2 Summary of validation of estimation of drag divergence Mach number by bisection

Validation Bisection

Test case Mdd CL |M =Mdd Mdd CL |M =Mdd

NACA 0015 0.74 0.0 0.74 0.0

RAE 2822 0.70-0.72 0.74 0.71 0.74

B. Calculating Technology Factor

The calculation of the airfoil technology factor can be considered in a few steps. First, the angle

of incidence was found to match the prescribed lift coefficient at low speed, which is stated where

necessary such as in table 3. From this, the drag divergence Mach number was calculated using the

bisection approach described above. The lift coefficient at this Mach number was then found.

C. Training Data Definition Results

The calculation of airfoil performance, as benchmarked by the technology factor, was performed

for 977 sharp trailing edge airfoils from the online database, at three different lift coefficients at

low speed to eliminate any compressibility effects (M=0.2): 0.0, 0.2, and 0.4. The goal of the

design variable derivation is to produce design variables that are highly robust to allow fully flexible

aerodynamic shape optimization to occur. The technology factor is an empirical constant that

describes the trade-off between transonic drag characteristics, lift, and shape, so for robustness

should be roughly constant for lift coefficients within the linear lift region. Hence, any airfoils that

had computed technology factors with a standard deviation greater than 2.5% of the mean value

over the three lift coefficients were eliminated, due to non-robustness. Furthermore, extremely thin

airfoils were also eliminated to ensure that no airfoils that could not realistically be used for design

were included in the library, again for robustness reasons. Results for the full collection of airfoils

before and after the elimination process are given in table 3.

The initial 977 airfoils were thus filtered down to 380 airfoils by this performance analysis.

However, further classification of these airfoils can also be performed by considering the selection of

airfoils based on their performance for the training data. Four sets of training data have therefore

13
 Table 3 Results of airfoil technology factor for all tested airfoils before and after elimination

CL =0.0 at M=0.2 CL =0.2 at M=0.2 CL =0.4 at M=0.2

κ Mdd κ Mdd κ Mdd

Max 0.937 0.869 0.930 0.822 0.934 0.784

Min 0.480 0.409 0.617 0.466 0.675 0.503

Before elimination
StDev 0.0791 0.0821 0.0479 0.0549 0.0362 0.0416

Mean 0.812 0.691 0.846 0.695 0.843 0.666

Max 0.920 0.844 0.930 0.822 0.934 0.784

Min 0.715 0.559 0.719 0.541 0.724 0.522

After elimination
StDev 0.0319 0.0480 0.0322 0.0461 0.0322 0.0419

Mean 0.855 0.727 0.861 0.703 0.861 0.666

been collated comparing various densities of design space sampling by empirical techniques, and are

compared for design space coverage:

• Library A: Full set of 977 airfoils[58];

• Library B: Reduced set of 380 airfoils by considering constant κ;

• Library C: Concentrated set of 100 airfoils with highest technology factor;

• Library D: Concentrated set of 100 airfoils uniformly distributed throughout the range of

technology factors.

The first of these cases will test the effect that having a very large number of airfoils in the

training data has, with a theoretically very large design space coverage; this represents a dense

sampling of the airfoil design space (library A). The second (library B) will also test a large sampling

of the design space, though in a more comprehensive manner than the first case. The third case

(library C) will test the effect that having entirely high performance airfoils in the training data

has on the modes extracted, and this demonstrates the effect of the metric-based filtering approach.

[58] airfoils that were used in the shape recovery tests below were individually omitted from this database to ensure
they were not contained within the original training data

14
The majority of these high performance airfoils are, however, supercritical sections so the ability to

represent supercritical sections should be high at the expense of representing other types of section.

This type of section has been used for decades as a basis for transonic design of airfoils, however,

optimization is interested in selecting the designs that a designer has yet to consider, which may be

similar to current sections, but may also be different. It is for this reason that airfoils with varied

performance were also tested for the final case (library D).

The first six mode shapes extracted using the SVD approach for all libraries are shown in figure

1, with the modal energy variations shown in figure 2. It should be noted that some of the trailing

edges of the surfaces in figure 1 cross over, however this is due to scaling up of the design variables

for visualization purposes. While appearing slightly different, the first four modes of these libraries

parameterized based on performance are representative of similar deformations. Modes 1 and 2

of the libraries represent thickness and camber respectively, though it is interesting that the full

airfoil library has a highly dominant thickness mode. Modes 3 and 4 of all libraries are similar and

represent trailing edge thickness distributions, and the addition of some form of trailing edge cusp

similar to the trailing edge angles of the PARSEC approach (see later section). Later mode shapes

begin to represent higher frequency surface oscillations, which are intuitively difficult to understand,

though allow detailed surface changes that are required for full scoping of the design space.

V. Geometric Shape Recovery

The fundamental geometric design variables extracted through the SVD are independent of flow

solver, optimizer, and surface parameterization/control/deformation methods. As such, the quality

of these design variables is also independent of these modules so analysis need not take these into

account. The analysis performed needs to answer the question of whether these modes provide

sufficient design space coverage, and to what extent does using a reduced set of modes affect this.

These issues can be explored by considering geometric inverse problems, so demonstrating the extent

to which the design variables can recover arbitrary airfoil shapes.

To compare the effect that the training data has on the derived design variables, inverse geo-

metrical design problems were considered, to analyse how effectively and efficiently the extracted

15
 (a) Library A mode 1 (b) Library B mode 1 (c) Library C mode 1 (d) Library D mode 1

(e) Library A mode 2 (f) Library B mode 2 (g) Library C mode 2 (h) Library D mode 2

(i) Library A mode 3 (j) Library B mode 3 (k) Library C mode 3 (l) Library D mode 3

(m) Library A mode 4 (n) Library B mode 4 (o) Library C mode 4 (p) Library D mode 4

(q) Library A mode 5 (r) Library B mode 5 (s) Library C mode 5 (t) Library D mode 5

(u) Library A mode 6 (v) Library B mode 6 (w) Library C mode 6 (x) Library D mode 6

Fig. 1 First six geometric mode shapes (design variables) of the dense sampling (library A),

medium sampling (library B) and sparse sampling filtered by high and varied performance

(libraries C and D respectively) airfoil libraries.

102
A-Full Library
B-380 Reduced
C-Best 100
D-Vary 100
1
10
Energy (%)

100

-1
10

0 10 20
Mode

Fig. 2 Modal energy of libraries.

variables can be used to recover an arbitrary airfoil shape from another, i.e. reconstruct a target

16
airfoil using a weighted combination of the modes added onto a base airfoil.

A. Optimization Problem

Using the extracted geometric design variables, which are surface deformations, the reconstruc-

tion is a solution of the problem:

n
X
Xrecon = Xbase + αj Uj (4)
j=1

where n is the selected number of modes used, and αj are the weightings of each extracted mode.

This is the description of a deformative parameterization, which is simpler to implement within the

ASO process. Hence, to test design space coverage, the values of the weighting coefficients α need

to be found that minimise the error between the target airfoil and the reconstructed airfoil. An

optimization approach is used to evaluate the weighting coefficients.

The evaluation of the error needs careful consideration. One evaluation of error could simply

be the vector difference between the surface points of the reconstructed and target airfoil, however,

this may constrain the problem too much by asking for an exact point to point mapping, whereas an

initial point need only move to the target surface, thus a representation of the whole surface is needed.

To simplify evaluation of the error an implicit surface [52] was constructed of the target airfoil that

gives a function value throughout the entire domain, where that function has a value of zero on

the surface. The zero-th level set of the function, which is the airfoil surface, is constrained to pass

through the points on the surface. The implicit surface is constructed by minimising the bending

energy around the defined surface points and so does not allow the modelling of sharp corners. The

sharp trailing edge of the airfoil was accounted for by constructing two implicit surfaces, one each for

the upper and lower surfaces, where curvature at the leading edge was maintained by overlapping the

upper and lower surfaces, see figure 3. The error is evaluated by marching in the normal direction of

the reconstructed surface until the function becomes zero – the norm of this is then the individual

point error.

17
 0.6 0.6

0.4 0.4

0.2 0.2

0 0

-0.2 -0.2

-0.4 -0.4

-0.6 -0.6
0 0.5 1 0 0.5 1
X X

(a) Upper surface (b) Lower surface

Fig. 3 Implicit surface construction.

B. Optimization Algorithm

An optimization approach was taken to consider the recovery problem, though the correct choice

of optimizer is a paramount issue. Two types of optimizer typically found in optimization theory

are gradient-based and agent-based. Gradient-based methods use the local gradient as a basis upon

which to construct a search procedure and terminate when the gradient becomes zero, which occurs

in any minima. On the other hand, agent-based methods avoid the computation of the gradient

by having a number of agents that investigate the search space and use the position and fitness of

those agents to progress the optimization in time. Agent-based methods are more suited for global

optimization as they have less chance of terminating in local minima, though the cost of using these

algorithms is more than a gradient-based approach.

The recovery problem investigated here was tested initially using both algorithms stated above.

It is unknown whether the problem is unimodal (one local optimum which is the global optimum) or

multimodal (many local optima, one or more of which is the global optimum) so a simple gradient

approach was tested along with an agent-based search system. To ensure the correct computation

of the recovery problem the implicit surface approach was used as a target surface representation,

however, to guide the points onto the surface the point-to-point error was also introduced. In the

gradient approach this was done by first running the point-to-point problem to convergence followed

by running the implicit surface problem from this point until that had converged. In the agent-

based approach a compound objective was used that was initially point-to-point, but then biased the

18
implicit surface more through the optimization by a partition of unity approach until the objective

was purely the implicit surface near to the end of the process.

As the objective is cheap to compute, optimizing the algorithms to ensure good convergence

characteristics is not an issue though the accurate location of the optimum solution is important.

To minimise the error, as measured by the root-mean-squared (RMS) of the difference between the

target and recovered airfoils, the gradient approach used is a conjugate gradient approach with the

Polak-Ribiere step direction where the step size is found by an Armijo line search procedure [53],

and the agent-based method used is the gravitational search algorithm (GSA) [54]. The conjugate-

gradient was run until either convergence or 1000 evolutions using a central finite difference for

gradient computation, and GSA with 500 particles for 1500 evolutions. To compare the two op-

timizers, both were run on a shape recovery problem of a NACA0012 to KC135c using library D

introduced above for various numbers of dimensions, and the results are given in table 4.

Table 4 Comparison of optimizers for shape recovery for NACA0012 to KC135c using library

D. Final RMS error of reconstructed section against target section is given as well as average

number of objective evaluations per timestep (cost).

RMS error Cost per evolution

Modes Gradient GSA Gradient GSA

8 4.41 × 10−4 4.41 × 10−4 21 500

12 9.18 × 10−5 9.22 × 10−4 29 500

16 9.02 × 10−5 9.15 × 10−4 37 500

20 6.91 × 10−5 7.44 × 10−5 45 500

From the results, the problem appeared well behaved. The difference in results between the two

optimizers is zero to within the tolerance shown in the table when considering only a few modes,

though when considering higher dimensional search spaces the agent-based algorithm has difficulty

fully converging to the true optimum solution. The agent system has still produced a high quality

solution. This result indicates that the design space of the recovery problem is convex, possibly as

a result of the orthogonality of the design variables that are due to the SVD approach implemented

19
here. As the gradient approach is always as good, and sometimes slightly better than the agent-based

method, this algorithm was used to test recoverability of the design space. Furthermore, the cost of

implementing the gradient-based method is considerably less than the global search algorithm.

C. Modal Convergence

This work deals with the ability for design variables to cover the design space by considering the

ability that a set of design variables has to recover an arbitrary airfoil. When considering the design

variables derived by the SVD approach, this involves testing the shape reconstruction of an arbitrary

airfoil not contained within the training data that was used to construct the design variables. If the

target and initial airfoils were contained within the original library then this type of optimization

would be described by the recovery of the column in the original deformation matrix that related

to the deformation between the two airfoils. This problem has a known solution, which is exact

recovery, and as such the developed optimization algorithm has been tested on this for convergence

of the modes. The sparsely sampled airfoil training data sets have been chosen, and the problem

was the recovery of the YS900 to the RAE5214 – though this is an arbitrary selection and could be

any two airfoils contained within both sets – using an increasing number of design variables. Figure

4 shows the results of this investigation and also shows the energy of the modes from each library.

The results show that the convergence of the error and energy follow roughly the same trend, and

using all the modes (defined by the energy to within machine tolerance) results in a machine zero

error, as expected.

D. Alternative design variables

The four sets of design variables derived by the SVD approach (those from 977, 380 and 100

airfoils selected by various reduction methods, termed libraries A, B, C and D) have also been

compared to two other deformative approaches that are common in aerodynamic shape optimization,

and are often accepted to be good approaches to the parameterization and surface control problem.

The first is the Hicks-Henne method [1], which is characterised by a set of deformation design

variables. The second is the PARSEC [13] method which is a constructive method common in

optimization. The Hicks-Henne approach was selected to test the effect of using a method that

20
 101 101 101 101

-1 -1 -1 -1
10 10 10 10

-3 -3 -3 -3
10 10 10 10

Energy (%)

Energy (%)
-5 -5 -5 -5
10 10 10 10

RMS

RMS
-7 -7 -7 -7
10 10 10 10

10-9 10-9 10-9 10-9

10
-11
Energy 10
-11
10
-11
Energy 10
-11

RMS RMS
-13 -13 -13 -13
10 10 10 10

-15 -15 -15 -15

10 10 10 10
20 40 60 80 100 20 40 60 80 100
Mode Mode

(a) Library C (b) Library D

Fig. 4 Modal energy of varied κ library and reconstruction error of airfoils within library.

has large design space coverage but which requires large numbers of design variables (Sripawadkul

et al. [37] tested 250 airfoil shapes and Hicks-Henne was able to give very good representation

using 32 design variables). This is compared to the PARSEC approach, which is found in many

aerodynamic shape optimization problems [45, 55], and performed almost perfectly in many of the

tests in the review paper of Sripawadkul et al.. This gives an indication of the effect that using

few design variables has, though at the expense of not accounting for detailed surface changes [25].

Furthermore, the PARSEC approach tends to be commonly used as it allows a designer to have

an intuitive link between the design variables and specific airfoil shape parameters, such as leading

edge radius and upper surface curvature.

1. Hicks-Henne

The Hicks-Henne approach describes the surface perturbation as a combination of a set of

smooth ‘bump’ functions. The perturbation of the jth function from a total of B functions is

described by:

  log 0.5 ξj

Uj = sin πx log ζj (5)

where ζ is the horizontal position of the bump and ξ controls the width of the bump. The prescribed

horizontal position of the bumps allows the user to constrain or relax control over certain areas of

21
the airfoil, however the approach used here is interested in global shape matching so for flexibility

a uniform position of the bumps across the chord is implemented. ξ is set as a design variable to

allow a bump global or local control to recover the target surface. For an airfoil, which must have

upper and lower surface perturbations, each bump has a weighting and a width so the total number

of design variables is given as (2Blower + 2Bupper ), and in this work (B = Blower = Bupper ).

2. PARSEC

The parametric section (PARSEC) airfoil parameterization, developed by Sobieczky [13], is

considered one of the most intuitive methods for parameterizing the surface of an airfoil, by defining

two sixth order polynomials using 11 design variables, which all relate to airfoil features, as shown

in figure 5.
Y

cxxu

cuy
u
Rle clx yy
l te
0
0 cux 1
X
cly
cxxl

Fig. 5 PARSEC airfoil geometric allocation

The upper and lower surfaces are given by:

6
X
yu (x) = an xn−0.5 (6)
n=1
6
X
yl (x) = bn xn−0.5 (7)
n=1

where a and b are the coefficients of the polynomial, which relate to the airfoil design variables by

the following system of equations [25]:

22
     
p
 1 0 0 0 0 0  a1   2/Rle 
    
    
1 1 1 1 1 1  a2  yte + ∆yte 
     

    
    
u

 (cux )1/2 (cux )3/2 (cux )5/2 (cux )7/2 (cux )9/2 u 11/2   
(cx )  a3  
 cy 

   =   (8)
    

 1/2 3/2 5/2 7/2 9/2 11/2  a4   tan(θu ) 
   

    
    
1 u −1/2 3 u 1/2 5 u 3/2 7 u 5/2 9 u 7/2 11 u 9/2    0
2 (cx ) 2 (cx ) 2 (cx ) 2 (cx ) 2 (cx ) (c ) a
  
 2 x   5   
    
     2  
d y
− 41 (cux )−3/2 3 u −1/2
4 (cx )
15 u 1/2
4 (cx )
35 u 3/2
4 (cx )
53 u 5/2
4 (cx )
99 u 7/2
4 (c x ) a 6 2
dx x=cu

x
    
p
 1 0 0 0 0 0  b1   − 2/Rle 
    
    
1 1 1 1 1 1  b2  yte − ∆yte 
     

    
    
l

 (clx )1/2 (clx )3/2 (clx )5/2 (clx )7/2 (clx )9/2 (clx )11/2  b3  
    cy 

   =   (9)
    

 1/2 3/2 5/2 7/2 9/2 11/2  b4   tan(θl ) 
   

    
    
1 l −1/2 3 l 1/2 5 l 3/2 7 l 5/2 9 l 7/2 11 l 9/2    0
2 (cx ) 2 (cx ) 2 (cx ) 2 (cx ) 2 (cx ) 2 (cx )  b5  
  
 
    
     2  
d y
− 41 (clx )−3/2 3 l −1/2
4 (cx )
15 l 1/2
4 (cx )
35 l 3/2
4 (cx )
53 l 5/2
4 (cx )
99 l 7/2
4 (cx ) b6 dx2  l
x=cx

The formulation is simple to implement to produce airfoil shapes though the issue with PAR-

SEC, and all constructive based methods, is that they describe a surface by a function instead of

perturbations away from a geometry. This makes their integration more difficult, as automatic mesh

generation tools are required. The PARSEC method can be slightly modified to act as the basis of a

deformative technique by approximating the perturbations as linear. For this work, which considers

surface perturbations as design variables, the PARSEC parameters are perturbed individually and

then remeshed to provide a set of surfaces from which the deformations can be directly obtained and

are shown in figure 6. The airfoils tested all have sharp trailing edge with zero vertical perturbation

at the trailing edge, so the PARSEC description reduces to nine parameters.

VI. Design Space Coverage Analysis

For the demonstration of full design space coverage it should be possible for the derived modes

to recover an arbitrary airfoil not in the original training library. To analyse this, inverse shape

matching was considered for a number of different airfoil shapes that represent typical profiles used

in aircraft design. Furthermore, the base airfoil, NACA0012 used here, was also not in the training

23
 (a) Leading edge radius, Rle (b) Upper crest position, cu
x (c) Lower crest position, clx

d2 y 

(d) Upper crest magnitude, cu
y (e) Lower crest magnitude, cly (f) Upper crest curvature, dx2 x=cu

x

d2 y 

(g) Lower crest curvature, dx2 x=cl
 (h) Trailing edge upper surface (i) Trailing edge lower surface
x

direction, θu direction, θl

Fig. 6 Nine linearly approximate PARSEC deformation parameters for sharp trailing edge

airfoils

library.

The effect of the metric-based filtering on the derived design variables is investigated by con-

sidering separately the ability for the various sets of modes to recover airfoils of a certain family.

The analysis can be split into high transonic performance shapes with high technology factors, and

other general shapes with varied transonic performance. Six different sections are tested to analyse

these criteria. The first three are typical high transonic performance and supercritical section air-

foils. The first is an airfoil from the Boeing KC135, termed the KC-135c from the online database

that the airfoils were taken from. This airfoil is taken directly from a current in-service aircraft

which has been designed to optimize some criteria and should give a good indication as to how

the derived design variables can match this type of airfoil. The second is from the NACA-6 series

and is the 66(2)-415, which was designed for high speed performance at a specific set of operating

conditions, though these airfoils have relatively high drag divergent Mach numbers and as such have

high technology factors also. The final airfoil of the transonic sections tested is the NPL airfoil from

ARC CP-1372[56]. This is a typical supercritical section – flat upper surface, rounded leading edge

and high trailing edge camber – that is typical in transonic airfoil design and so the derived modes

should be able to represent this common section type. The effect that high technology factors has on

24
the recoverability using the design variables is investigated here, with the first three airfoils having

average technology factors of 0.88, 0.87 and 0.90 respectively.

The final three airfoils considered are from three different airfoil families not considered to have

necessarily good transonic performance and, even if they do, they are considerably different from

airfoils contained within the training libraries. The first of these is the DAE-11 designed for the wing

of the Daedalus HPA project, and was optimized for low Reynolds number flows to produce high lift.

Subsequently, the section is has a high maximum thickness, so a large variation in thickness along

the chord, and a high camber, representing a particularly difficult problem. The second section if the

NLF-415, which like the NACA-6 series airfoils also to be tested was designed to maximise laminar

flow over the chord, however this section has a much later maximum thickness position with high

curvature at the trailing edge so demonstrates a considerable difference in design between sections to

achieve the same goal, and its testing is therefore important. The final section considered is a typical

wind-turbine section which is designed for low Reynolds number flows like the DAE-11, though the

section shape is thinner and less cambered. The final three airfoils have average technology factors

of 0.75, 0.84 and 0.80 respectively, so a considerable range of airfoil shapes with varied performance

levels are tested.

The initial airfoil in all cases was taken as the NACA 0012. All the airfoils are shown in figure

7, and it should be stressed that none of the airfoils used appear in any of the libraries tested, and

so this is a wide-ranging test of design parameter space coverage.

NACA 0012

Boeing KC-135c
DAE-11

NACA 66(2)-415
NLF-415

NPL from CP-1372[56] SG6042

Fig. 7 Airfoils used for testing

A 97 point surface was used for each airfoil where the distribution of points is as per the original

25
airfoil library. The full set of recovery results is presented in figure 8 and the surface errors of a high

performance airfoil (NACA 66(2)-415) and low performance airfoil (SG6042) are given in figures 9

and 10 respectively and the lowest error result is highlighted in bold. The results for dense (library

A), medium (library B) and sparse (libraries C and D) sampling libraries are presented as well as

those for Hicks-Henne (HH) and linearly approximate PARSEC (L-PARSEC) results. The surface

error plots are computed using only 12 design variables in each case (except PARSEC which has a

constant nine variables). The wind tunnel tolerance value of 5 × 10−4 was taken from Kulfan [12]

to be the tolerance which the results should be within to be able to state that the airfoil shape had

been recovered sufficiently.

Library A
Library B
Library C
Library D
HH
L-PARSEC

-2 -2 -2
10 10 10

-3 -3 -3
RMS

RMS

RMS
10 10 10

-4 -4 -4
10 10 10

5 10 15 20 5 10 15 20 5 10 15 20
Design Variables Design Variables Design Variables

(b) Boeing KC-135c (c) NACA 66(2)-415 (d) NPL

-2 -2 -2
10 10 10

10-3 10-3 10-3

RMS

RMS

RMS

10-4 10-4 10-4

5 10 15 20 5 10 15 20 5 10 15 20
Design Variables Design Variables Design Variables

(e) DAE-11 (f) NLF-415 (g) SG6042

Fig. 8 Root-mean-squared (RMS) error of surface recovery tests for varied numbers of design

variables.

The results when using an increasing number of modes for recovery exhibit the expected trend of

26
 Library A
Library B
Library C
Library D
HH
L-PARSEC

0.002 0.002

0.0015 0.0015
|Error|

|Error|
0.001 0.001

0.0005 0.0005

0 0
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/c x/c

(b) Upper - SVD modes (c) Upper - other

0.002 0.002

0.0015 0.0015
|Error|

|Error|

0.001 0.001

0.0005 0.0005

0 0
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/c x/c

(d) Lower - SVD modes (e) Lower - other

Fig. 9 Surface errors from recovery of NACA 66(2)-415 for 12 design variables and 9 for

PARSEC.

producing lower error. The increased number of modes adds flexibility to the system and therefore

allows for greater surface recovery, however the caveat to this is that the low energy modes represent

high frequency surface oscillations and therefore when used in an optimization scheme should often

be avoided. All the airfoil recovery problems have produced good results using at least one of the

libraries tested, demonstrating the efficient manner in which design variables are produced using

the SVD approach.

The results of the high density sampling library (library A) are good, often outperforming both

Hicks-Henne and PARSEC for an equivalent number of design variables, though the medium sam-

27
 Library A
Library B
Library C
Library D
HH
L-PARSEC

0.002 0.002

0.0015 0.0015
|Error|

|Error|
0.001 0.001

0.0005 0.0005

0 0
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/c x/c

(b) Upper - SVD modes (c) Upper - other

0.002 0.002

0.0015 0.0015
|Error|

|Error|

0.001 0.001

0.0005 0.0005

0 0
0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1
x/c x/c

(d) Lower - SVD modes (e) Lower - other

Fig. 10 Surface errors from recovery of SG6042 for 12 design variables and 9 for PARSEC.

pling library (library B) outperforms this. This indicates that while the high density sampling covers

the design space, the modes produced are characterised by large noise variations resulting in the di-

lution of minor surface changes. The reduction from this full set to a more manageable set of airfoils

(libraries C and D) by the consideration of the metric is therefore required to produce a generous

trade-off between wide sampling of the design space and concentrated, usable design variables that

capture minor surface changes and large shape changes. The sparsely sampled data sets perform

better than both the more densely sampled sets further indicating that a trade-off between fidelity

of surface changes in the modes needs consideration. The sparsely sampled library containing en-

tirely high performance design variables (library C) is characterised by primarily supercritical and

28
NACA-6 series sections, and as a result can represent these type of sections better than all the other

libraries. The varied performance library (library D), which contains a wider selection of airfoil

section families can subsequently represent more varied sections other than supercritical compared

to the high performance library. The advantage of the method developed in this paper is that a

user can tailor a set of design variables to a specific problem where detailed surface change is ideal,

or where large topology changes may be useful. The comparison of all four sets of design variables

from the SVD compared to the other two parameterization methods (Hicks-Henne and PARSEC)

further demonstrates that this method is effective at producing design variables that can efficiently

cover a wide portion of the airfoil design space.

It is promising to note that the errors around the leading edge are small for the best results from

each of the surface error plots which demonstrates that the SVD methods for deriving the design

variables allows control of the leading edge, as opposed to PARSEC which has large errors; this was

also a conclusion of Castonguay and Nadarajah [25]. Furthermore, the errors have few oscillations

over the surface and are well behaved, and given the use of only 12 modes represents the use of a

highly efficient set of design variables.

VII. Conclusions

The definition of efficient design degrees of freedom for use in aerodynamic shape optimization

has been considered. A proper orthogonal decomposition approach has been developed and pre-

sented, to allow mathematical extraction of effective and orthogonal airfoil design variables from

a library of training data. The resulting geometric modes are independent of parameterization

scheme, surface and volume mesh, and flow solver, thus are generally applicable.

An intelligent, metric-based approach has also been developed, using the airfoil technology

factor, to perform a type of design space reduction to the training data library. Geometric inverse, or

shape matching, problems have been considered to analyse the effectiveness of these extracted modes,

and results presented for six airfoil shape recovery problems, using a gradient-based optimizer,

with modes extracted from four different airfoil libraries with various sampling densities: high

density sampling containing 977 airfoils; medium sampling of 380 airfoils reduced by considering

29
the performance metric; 100 airfoils concentrated by considering high performance airfoils; 100

airfoils concentrated by considering a variety of performances.

It is shown that the mathematically extracted mode shapes perform very well in all cases and are

able to represent high curvature regions accurately, such as around the leading edge. It is possible

to represent a wide variety of airfoil shapes to within a small tolerance using as few as eight SVD

design variables, and this is a particularly good result. Furthermore, including a wide range of airfoil

shapes is essential to ensure the best design space coverage of a variety of sections, though a user

can tailor design variables to a specific problem. It is, however, important that a sampling density

is chosen to ensure that minor surface changes and large topology changes can both be captured

effectively.

References

[1] Hicks, R. M. and Henne, P. A., “Wing Design by Numerical Optimization,” Journal of Aircraft, Vol. 15,

No. 7, 1978, pp. 407–412.

doi:10.2514/6.1977-1247.

[2] Li, W., Huyse, L. W., and Padula, S., “Robust Aerofoil Optimization to Achieve Consistent Drag

Reduction Over a Mach Range,” Tech. rep., 2001, NASA/CR–2001–211042. ICASE report number

2001–22.

doi:10.1007/s00158-002-0212-4.

[3] Chung, H. S. and Alonso, J. J., “Multiobjective Optimization Using Approximation Model-Based Ge-

netic Algorithms,” 10th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference, Al-

bany, New York, 2004, AIAA Paper 2004–4325.

doi:10.2514/6.2004-4325.

[4] Wong, W. S., Le Moigne, A., and Qin, N., “Parallel adjoint–based optimisation of a blended wing body

aircraft with shock control bumps,” The Aeronautical Journal, Vol. 111, No. 1117, 2006, pp. 165–174.

[5] Qin, N., Vavalle, A., Le Moigne, A., Laban, M., Hackett, K., and Weinerfelt, P., “Aerodynamic Con-

siderations of Blended Wing Body Aircraft,” Progress in Aerospace Sciences, Vol. 40, No. 6, 2004, pp.

321–343.

doi:10.1016/j.paerosci.2004.08.001.

[6] Jameson, A., Vassberg, J. C., and Shankaran, S., “Aerodynamic-Structural Design Studies of Low-

Sweep Transonic Wings,” Journal of Aircraft, Vol. 47, No. 2, 2010, pp. 505–514.

30
 doi:10.2514/1.42775.

[7] Le Pape, A. and Beaumier, P., “Numerical optimization of helicopter rotor aerodynamic performance

in hover,” Aerospace Science and Technology, Vol. 9, No. 3, 2005, pp. 191–201.

doi:10.1016/j.ast.2004.09.004.

[8] Dumont, A., Le Pape, A., Peter, J., and Huberson, S., “Aerodynamic Shape Optimization of Hovering

Rotors Using a Discrete Adjoint of the Reynolds Averaged Navier–Stokes Equations,” Journal of the

American Helicopter Society, Vol. 56, No. 3, 2011, pp. 1–11.

doi:10.4050/JAHS.56.032002.

[9] Tatossian, C., Nadarajah, S. K., and Castonguay, P., “Aerodynamic Shape Optimization of Hovering

Rotor Blades Using a Non-Linear Frequency Domain Approach,” Computers and Fluids, Vol. 51, No. 1,

2011, pp. 1–15.

doi:10.1016/j.compfluid.2011.06.014.

[10] Nielsen, E. J., Lee-Rausch, E. M., and Jones, W. T., “Adjoint Based Design of Rotors in a Noninertial

Frame,” Journal of Aircraft, Vol. 47, No. 2, 2010, pp. 638–646.

doi:10.2514/1.46044.

[11] Choi, S., Lee, K. H., Potsdam, M., and Alonso, J. J., “Helicopter Rotor Design Using a Time–Spectral

and Adjoint Based Method,” Journal of Aircraft, Vol. 51, No. 2, 2014, pp. 412–423.

doi:10.2514/1.C031975.

[12] Kulfan, B. M., “Universal Parametric Geometry Representation Method,” Journal of Aircraft, Vol. 45,

No. 1, 2008, pp. 142–158.

doi:10.2514/1.29958.

[13] Sobieczky, H., “Parametric Airfoils and Wings,” Notes on Numerical Fluid Mechanics, Vol. 68, 1998,

pp. 71–88.

[14] Jameson, A., “Aerodynamic design via control theory,” Journal of Scientific Computing, Vol. 3, No. 3,

1988, pp. 233–260.

doi:10.1007/BF01061285.

[15] Morris, A. M., Allen, C. B., and Rendall, T. C. S., “CFD–based Optimization of Aerofoils Using Radial

Basis Functions for Domain Element Parameterization and Mesh Deformation,” International Journal

for Numerical Methods in Fluids, Vol. 58, No. 8, 2008, pp. 827–860.

doi:10.1002/fld.1769.

[16] Morris, A. M., Allen, C. B., and Rendall, T. C. S., “Domain–Element Method for Aerodynamic Shape

Optimization Applied to a Modern Transport Wing,” AIAA Journal, Vol. 47, No. 7, 2009, pp. 1647–

31
 1659.

doi:10.2514/1.39382.

[17] Allen, C. B., Rendall, T. C. S., and Morris, A. M., “Computational-Fluid-Dynamics-Based Twist Op-

timisation of Hovering Rotors,” Journal of Aircraft, Vol. 47, No. 6, 2010, pp. 2075–2085.

doi:10.2514/1.C000316.

[18] Allen, C. B. and Rendall, T. C. S., “Computational-Fluid-Dynamics-Based Optimisation of Hovering

Rotors Using Radial Basis Functions for Shape Parameterisation and Mesh Deformation,” Optimization

and Engineering, Vol. 14, 2013, pp. 97–118.

doi:10.1007/s11081-011-9179-6.

[19] Samareh, J. A., “Aerodynamic Shape Optimization Based on Free-Form Deformation,” 10th

AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference, Albany, New York, 2004, AIAA

Paper 2004–4630.

doi:10.2514/6.2004-4630.

[20] Lassila, T. and Rozza, G., “Parametric free-form shape design with PDE models and reduced basis

method,” Computer Methods in Applied Mechanics and Engineering, Vol. 199, No. 23–24, 2010, pp.

1583–1592.

doi:10.1016/j.cma.2010.01.007.

[21] Leung, T. M. and Zingg, D. W., “Aerodynamic Shape Optimization of Wings Using a Parallel Newton-

Krylov Approach,” AIAA Journal, Vol. 50, No. 3, 2012, pp. 540–550.

doi:10.2514/1.J051192.

[22] Gagnon, H. and Zingg, D. W., “Two-Level Free-Form Deformation for High-Fidelity Aerodynamic Shape

Optimization,” 12th AIAA Aviation Technology, Integration and Operations (ATIO) Conference and

14th AIAA/ISSMO Multidisciplinary Analysis Optimization Conference, Indianapolis, Indiana, 2012,

AIAA Paper 2012–5447.

doi:10.2514/6.2012-5447.

[23] Samareh, J. A., “Status and Future of Geometric Modelling and Grid Generation for Design and

Optimization,” Journal of Aircraft, Vol. 36, No. 1, 1999, pp. 97–104.

doi:10.2514/2.2417.

[24] Samareh, J. A., “Survey of Shape Parameterization Techniques for High Fidelity Multidisciplinary

Shape Optimization,” AIAA Journal, Vol. 39, No. 5, 2001, pp. 877–884.

doi:10.2514/2.1391.

[25] Castonguay, P. and Nadarajah, S. K., “Effect of Shape Parameterization on Aerodynamic Shape Op-

32
 timization,” 45th AIAA Aerospace Sciences Meeting and Exhibit, Reno, Nevada, 2007, AIAA Paper

2007–59.

doi:10.2514/6.2007-59.

[26] Mousavi, A., Castonguay, P., and Nadarajah, S. K., “Survey of Shape Parameterization Techniques and

its effect on Three Dimensional Aerodynamic Shape Optimization,” 18th AIAA Computational Fluid

Dynamics Conference, Miami, Florida, 2007, AIAA Paper 2007–3837.

doi:10.2514/6.2007-3837.

[27] Pickett, R. M., Rubinstein, M. F., and Nelson, R. B., “Automated Structural Synthesis with a Reduced

Number of Design Coordinates,” AIAA Journal, Vol. 11, No. 4, 1973, pp. 494–498.

doi:10.2514/6.1973-336.

[28] Samareh, J. A., “Novel Multidisciplinary Shape Parameterization Approach,” Journal of Aircraft,

Vol. 38, No. 6, 2001, pp. 1015–1024.

doi:10.2514/2.2888.

[29] Bloor, M. I. G. and Wilson, M. J., “Generating Parameterizations of Wing Geometries Using Partial

Differential Equations,” Computer Methods in Applied Mechanics and Engineering, Vol. 148, 1997, pp.

125–138.

doi:10.1016/S0045-7825(97)00033-9.

[30] Braibant, V. and Fleury, C., “Shape Optimal Design Using B–Splines,” Computer Methods in Applied

Mechanics and Engineering, Vol. 44, No. 3, 1984, pp. 247–267.

doi:10.1016/0045-7825(84)90132-4.

[31] Zhu, F. and Qin, N., “Intuitive Class/Shape Function Parameterization for Airfoils,” AIAA Journal,

Vol. 52, No. 1, 2014, pp. 17–25.

doi:10.2514/1.J052610.

[32] Young, D. P., Huffman, W. P., Melvin, R. G., Hilmes, C. L., and Johnson, F. T., Nonlinear Elimination

in Aerodynamic Analysis and Design Optimization, Springer Berlin Heidelberg, 2003.

doi:10.1007/978-3-642-55508-4_2.

[33] Hicken, J. E. and Zingg, D. W., “Aerodynamic Optimization Algorithm with Integrated Geometry

Parameterization and Mesh Movement,” AIAA Journal, Vol. 48, No. 2, 2010, pp. 400–413.

doi:10.2514/1.44033.

[34] Han, X. and Zingg, D. W., “An Adaptive Geometry Parameterization for Aerodynamic Shape Opti-

mization,” Optimization and Engineering, Vol. 15, No. 1, 2014, pp. 69–91.

doi:10.1007/s11081-013-9213-y.

33
[35] Anderson, G. R., Aftosmis, M. J., and Nemec, M., “Constraint-based Shape Parameterization for

Aerodynamic Design,” 7th International Conference on Computational Fluid Dynamics, Big Island,

Hawaii, 2012, ICCFD7–2001.

[36] Anderson, G. R., Aftosmis, M. J., and Nemec, M., “Parametric Deformation of Discrete Geometry

for Aerodynamic Shape Design,” 50th AIAA Aerospace Sciences Meeting including the New Horizons

Forum and Aerospace Exposition, Nashville, Tennessee, 2012, AIAA Paper 2012–0965.

doi:10.2514/6.2012-965.

[37] Sripawadkul, V., Padulo, M., and Guenov, M., “A Comparison of Airfoil Shape Parameterization Tech-

niques for Early Design Optimization,” 13th AIAA/ISSMO Multidisciplinary Analysis and Optimization

Conference, Fort Worth, Texas, 2010, AIAA Paper 2010–9050.

doi:10.2514/6.2010-9050.

[38] Aidala, P. V., Davis, W. H., and Mason, W. H., “Smart Aerodynamic Optimization,” AIAA Applied

Aerodynamics Conference, Danvers, Massachusetts, 1983, AIAA Paper 1983–1863.

doi:10.2514/6.1983-1863.

[39] Chang, I. C., Torres, F. J., and Tung, C., “Geometric Analysis of Wing Sections,” Tech. rep., NASA

Ames Research Centre, Moffett Field, California, 1995, NASA Technical Memorandum 110346.

[40] Robinson, G. M. and Keane, A. J., “Concise Orthogonal Representation of Supercritical Airfoils,”

Journal of Aircraft, Vol. 38, No. 3, 2001, pp. 580–583.

doi:10.2514/2.2803.

[41] Toal, D. J. J., Bressloff, N. W., Keane, A. J., and Holden, C. M. E., “Geometric Filtration Using Proper

Orthogonal Decomposition for Aerodynamic Design Optimization,” AIAA Journal, Vol. 48, No. 5, 2010,

pp. 916–928.

doi:10.2514/1.41420.

[42] Ghoman, S. S., Wang, Z., Chen, P. C., and Kapania, R. K., “A POD-based Reduced Order Design

Scheme for Shape Optimization of Air Vehicles,” 53rd AIAA/ASME/ASCE/AHS/ASC Structures,

Structural Dynamics and Materials Conference and Co-located Events, Honolulu, Hawaii, 2012, AIAA

Paper 2012–1808.

doi:10.2514/6.2012-1808.

[43] Chatterjee, A., “An Introduction to the Proper Orthogonal Decomposition,” Current science, Vol. 78,

No. 7, 2000, pp. 808–817.

[44] Liang, Y. C., Lee, H. P., Lim, S. P., Lin, W. Z., Lee, K. H., and Wu, C. G., “Proper Orthogonal

Decomposition and its Applications–Part I: Theory,” Journal of Sound and Vibration, Vol. 252, No. 3,

34
 2002, pp. 527–544.

doi:10.1006/jsvi.2001.4041.

[45] Andres, E., Salcedo-Sanz, S., Monge, F., and Perez-Bellido, A. M., “Efficient Aerodynamic Design

Through Evolutionary Programming and Support Vector Regression Algorithms,” Expert Systems with

Applications, Vol. 39, 2012, pp. 10700–10708.

doi:10.1016/j.eswa.2012.02.197.

[46] Slooff, J. W., Thibert, J. J., Schmidt, W., Sacher, P., Papailiou, K. D., Bucciantini, G., Borsi, M.,

van der Vooren, J., Ashill, P. R., Holst, T. L., and Boppe, C. W., “Technical Status Review on Drag

Prediction and Analysis from Computational Fluid Dynamics: State of the Art,” Tech. rep., Advisory

Group for Aerospace Research and Development (AGARD), 1988, AGARD Advisory Report Number

256.

[47] Graham, D. J., Nitzberg, G. E., and Olsen, R. N., “A Systematic Investigation of Pressure Distributions

at High Speed Over Five Representative NACA Low-Drag and Conventional Airfoil Sections,” Tech.

rep., National Advisory Committee for Aeronautics, 1945, NACA TR–832.

[48] Drela, M., Pros and Cons of Airfoil Optimization, World Scientific, 1998.

[49] Allen, C. B., “Parallel Universal Approach to Mesh Motion and Application to Rotors in Forward

Flight,” International Journal for Numerical Methods in Engineering, Vol. 69, No. 10, 2007, pp. 2126–

2149.

doi:10.1002/nme.1846.

[50] Allen, C. B., “Multigrid Convergence of Inviscid Fixed- and Rotary-Wing Flows,” International Journal

for Numerical Methods in Fluids, Vol. 39, No. 2, 2002, pp. 121–140.

doi:10.1002/fld.282.

[51] Parpia, I. H., “Van-Leer Flux Vector Splitting in Moving Coordinates,” AIAA Journal, Vol. 26, 1988,

pp. 113–115.

doi:10.2514/3.9858.

[52] Turk, G. and O‘Brien, J. F., “Modelling with Implicit Surfaces that Interpolate,” ACM Transactions

on Graphics, Vol. 21, No. 4, 2002, pp. 855–873.

doi:10.1145/571647.571650.

[53] Gilbert, J. C. and Nocedal, J., “Global Convergence Properties of Conjugate Gradient Methods,” SIAM

Journal on Optimization, Vol. 2, No. 1, 1992, pp. 21–42.

doi:10.1137/0802003.

[54] Rashedi, E., Nezamabadi-pour, H., and Saryazdi, S., “GSA: A Gravitational Search Algorithm,” Infor-

35
 mation Sciences, Vol. 179, 2009, pp. 2232–2248.

doi:10.1016/j.ins.2009.03.004.

[55] Jahangirian, A. and Shahrokhi, A., “Aerodynamic Shape Optimization Using Efficient Evolutionary

Algorithms and Unstructured CFD Solver,” Computers and Fluids, Vol. 46, 2011, pp. 270–276.

doi:10.1016/j.compfluid.2011.02.010.

[56] Lovell, D. A., “A wind-tunnel investigation of the effects of flap span and deflection angle, wing planform

and a body on the high-lift performance of a 28 degrees swept wing,” Tech. rep., Aeronautical Research

Council, 1977, ARC/CP-1372.

36

View publication stats