Turbine Aerothermal Optimization Using Evolution Strategies
by
Drew A. Curriston
A thesis submitted to the Graduate Faculty of
Auburn University
in partial ful llment of the
requirements for the Degree of
Master of Science
Auburn, Alabama
May 4th, 2014
Keywords: Turbine, Optimization, Aerothermal, B ezier
Copyright 2014 by Drew A. Curriston
Approved by
Roy Hart eld, Chair, Walt and Virginia Woltosz Professor of Aerospace Engineering
Brain Thurow, Associate Professor of Aerospace Engineering
Andrew Shelton, Assistant Professor of Aerospace Engineering
Abstract
The use of adaptive optimizers in turbine blade geometry design is discussed. Many
di erent adaptive optimizers have been applied to turbine blade design in the past, but
often the optimizer does not t the problem and unnecessary constraints are implemented.
This work analyzes the popular choices for turbine blade optimization, select the most t-
ting optimizer that will lead to the lowest number of constraints, and apply it to an airfoil
optimization problem for validation. After validation, the optimizer is used in a turbine
blade optimization problem considering aerodynamics and convective heat transfer in the
objective function. The algorithm uses a 2D Navier-Stokes blade-to-blade analysis to evalu-
ate the e ectiveness of the solution encoding and optimizer performance. By using a highly
e ective optimizer for the speci ed problem, the computation time is drastically reduced
while producing comparable or better results, and these bene ts will also allow the use of a
higher delity analysis without excessive computational costs.
ii
Acknowledgments
I would like to thank my wife Meg, and my children Cally, Joshua, Jacob, and Asher for
their continued support in completing this paper. I would also like to thank Dr. Hart eld
for the time and e ort have gave while advising me on this research, and Dr. Shelton for
lending his expertise in helping me understand and implement CFD.
iii
Table of Contents
Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ii
Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . iii
List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi
List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x
1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.1 Objectives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.1.1 Objective 1: Produce and Validate Optimizer . . . . . . . . . . . . . 2
1.1.2 Objective 2: Turbine Aerodynamic Optimization . . . . . . . . . . . 2
1.1.3 Objective 3: Turbine Aerothermal Optimization . . . . . . . . . . . . 2
1.2 Literature Review . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3 Problem Background . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2 Evolution Strategies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.1 Algorithm Selection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.2 Solution Encoding and Move Operator . . . . . . . . . . . . . . . . . . . . . 12
2.3 Optimizer Validation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.3.1 Validation Case 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.3.2 Validation Case 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
2.3.3 Conclusions From Optimizer Validation . . . . . . . . . . . . . . . . . 26
3 Turbine Aerodynamic Optimization . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.1 Turbine Fundamentals Review . . . . . . . . . . . . . . . . . . . . . . . . . . 28
3.2 Preliminary Inviscid Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.2.1 Model and Assumptions . . . . . . . . . . . . . . . . . . . . . . . . . 30
3.2.2 Grid Generation and CFD Code . . . . . . . . . . . . . . . . . . . . . 31
iv
3.2.3 Objective Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3.2.4 Results and Optimizer Performance . . . . . . . . . . . . . . . . . . . 33
3.3 Viscous Turbine Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
3.3.1 Model and Assumptions . . . . . . . . . . . . . . . . . . . . . . . . . 38
3.3.2 Optimizer Updates and Objective Function . . . . . . . . . . . . . . . 39
3.3.3 Structural Considerations . . . . . . . . . . . . . . . . . . . . . . . . 42
3.3.4 Results and Optimizer Performance . . . . . . . . . . . . . . . . . . . 43
4 Turbine Aerothermal Optimization . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.1 Heat Transfer Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.2 Adaptations to Optimizer . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.2.1 Objective Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.3 Preliminary Aerothermal Optimization Results . . . . . . . . . . . . . . . . . 51
4.3.1 Issues With Preliminary Results . . . . . . . . . . . . . . . . . . . . . 53
4.3.2 Updates to Objective Function . . . . . . . . . . . . . . . . . . . . . 54
4.3.3 Stage Aerothermal Results . . . . . . . . . . . . . . . . . . . . . . . . 55
5 Conclusions and Future Applications . . . . . . . . . . . . . . . . . . . . . . . . 62
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
Appendices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
A Pseudo Code for ESTurb . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68
B Example GRAPE and RVCQ3D input les . . . . . . . . . . . . . . . . . . . . . 69
C Additional Plots for Optimization Runs . . . . . . . . . . . . . . . . . . . . . . 71
D O -Design E ects on Aerothermal Blade . . . . . . . . . . . . . . . . . . . . . . 74
E Brief Grid Re nement Study . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
v
List of Figures
1.1 Example Velocity Diagrams for a Stator[21] . . . . . . . . . . . . . . . . . . . . 6
1.2 Example Velocity Diagrams for a Rotor[21] . . . . . . . . . . . . . . . . . . . . 6
2.1 Search space comparison for minimization problem . . . . . . . . . . . . . . . . 12
2.2 PARSEC encoding . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
2.3 Example of using B ezier curves to encode blade[23] . . . . . . . . . . . . . . . . 14
2.4 Spline connected polynomial encoding angles and radii[21] . . . . . . . . . . . . 14
2.5 Spline connected polynomial encoding de ning points[21] . . . . . . . . . . . . . 14
2.6 Blade surface for example solution encoding . . . . . . . . . . . . . . . . . . . . 15
2.7 Move operator example for a single control point . . . . . . . . . . . . . . . . . 17
2.8 Second control point mutation option 1 . . . . . . . . . . . . . . . . . . . . . . . 18
2.9 Second control point mutation option 2 . . . . . . . . . . . . . . . . . . . . . . . 18
2.10 Flow chart for ESTurb algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.11 Test case one objective function results (sought to minimize) . . . . . . . . . . . 21
2.12 Test case one computation time comparison . . . . . . . . . . . . . . . . . . . . 21
2.13 XFoil output for converged airfoil . . . . . . . . . . . . . . . . . . . . . . . . . . 22
vi
2.14 History of objective function value over generations for initial eight runs . . . . 23
2.15 Example of variation in nal population from Figure 2.13 . . . . . . . . . . . . . 23
2.16 Pareto front results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.17 Pareto front progression with computation time . . . . . . . . . . . . . . . . . . 27
3.1 Elliptical H-grid . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.2 Leading edge without vertical clustering . . . . . . . . . . . . . . . . . . . . . . 31
3.3 Leading edge with vertical clustering . . . . . . . . . . . . . . . . . . . . . . . . 31
3.4 Euler code validation results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 32
3.5 Objective function values over 50 generations . . . . . . . . . . . . . . . . . . . 35
3.6 Final blade geometry using spline-connected polynomial encoding . . . . . . . . 35
3.7 Final blade geometry using B ezier curve encoding . . . . . . . . . . . . . . . . . 36
3.8 Example 169x45 grid for VK-LS89 turbine blade . . . . . . . . . . . . . . . . . 40
3.9 Close-up of trailing edge for VK-LS89 turbine blade grid . . . . . . . . . . . . . 41
3.10 Blade geometry and pressure contour results for single stage . . . . . . . . . . . 45
3.11 Mach Number contours for stator . . . . . . . . . . . . . . . . . . . . . . . . . . 46
3.12 Absolute Mach Number contours for rotor . . . . . . . . . . . . . . . . . . . . . 47
4.1 Ideal engine power vs. pressure ratio . . . . . . . . . . . . . . . . . . . . . . . . 51
4.2 Relationship for normalized power ratio vs turbine inlet temperature . . . . . . 51
vii
4.3 Aerothermal run 1 pressure contours . . . . . . . . . . . . . . . . . . . . . . . . 52
4.4 Aerothermal run 2 pressure contours . . . . . . . . . . . . . . . . . . . . . . . . 52
4.5 Aerothermal run 3 pressure contours . . . . . . . . . . . . . . . . . . . . . . . . 52
4.6 Preliminary aerothermal results: Heat transfer coe cient . . . . . . . . . . . . . 52
4.7 Aerothermal run 1 Mach number distribution . . . . . . . . . . . . . . . . . . . 54
4.8 \Bump" in velocity distribution for turboprop blade[3] . . . . . . . . . . . . . . 54
4.9 Stage aerothermal Run 1 stator Mach number distribution . . . . . . . . . . . . 56
4.10 Stage aerothermal Run 1 rotor Mach number distribution . . . . . . . . . . . . 56
4.11 Comparison of nal stator geometries . . . . . . . . . . . . . . . . . . . . . . . . 57
4.12 Comparison of nal rotor geometries . . . . . . . . . . . . . . . . . . . . . . . . 57
4.13 Stage pressure contours for Run 1 . . . . . . . . . . . . . . . . . . . . . . . . . . 58
4.14 Types of leading edges and incidence angle e ects[3] . . . . . . . . . . . . . . . 59
4.15 Heat transfer coe cient for optimized stator . . . . . . . . . . . . . . . . . . . . 60
4.16 Convergence history for both stator and rotor . . . . . . . . . . . . . . . . . . . 61
C.1 Pressure contours for Run 2 stage results . . . . . . . . . . . . . . . . . . . . . . 72
C.2 Isentropic Mach number distribution for Run 2 stator . . . . . . . . . . . . . . . 72
C.3 Isentropic Mach number distribution for Run 2 rotor . . . . . . . . . . . . . . . 72
C.4 Pressure contours for Run 3 stage results . . . . . . . . . . . . . . . . . . . . . . 73
viii
C.5 Isentropic Mach number distribution for Run 3 stator . . . . . . . . . . . . . . . 73
C.6 Isentropic Mach number distribution for Run 3 rotor . . . . . . . . . . . . . . . 73
D.1 Mach number surface distribution over a range of incidence angles . . . . . . . . 75
D.2 Close-up of Figure D.1 suction side . . . . . . . . . . . . . . . . . . . . . . . . . 75
E.1 RRMS plot for 169x45 grid . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 77
E.2 Mach number distribution for multiple grid sizes . . . . . . . . . . . . . . . . . 78
E.3 Heat transfer coe cient distribution for multiple grid sizes . . . . . . . . . . . . 78
ix
List of Tables
2.1 Example solution encoding . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
4.1 Preliminary Aerothermal Results Summary . . . . . . . . . . . . . . . . . . . . 51
4.2 Stage Aerothermal Results Summary . . . . . . . . . . . . . . . . . . . . . . . . 56
x
List of Symbols
C Absolute Velocity
CD Coe cient of drag
CL Coe cient of lift
Cm Moment coe cient
Cp Speci c Heat Capacity
gc Gravitational Constant
h Heat transfer Coe cient
h Speci c Enthalpy
k Thermal Conductivity
_m Mass ow rate
P Pressure
q Heat transfer
r radius from axis of rotation
Rn Stage Reaction
Re Reynold?s number
St Stanton Number
u Tangential Blade Speed
V Absolute Velocity
W Relative Velocity
xtr Axial distance from leading edge to
ow separation
Greek
Angle of attack
Speci c Heat Ratio
E ciency
O spring population size
Parent population size
c Compressor Pressure Ratio
Density
Flow Coe cient
Work or Loading Coe cient
! Shaft angular velocity
! Total Pressure Loss Coe cient
Subscripts
0 Stagnation Condition
1 Freestream Condition
Tangential direction
ex Exit Condition
in Inlet Condition
m Meridional direction
s Static Condition
t s Total-to-static
t t Total-to-total
xi
Chapter 1
Introduction
Gas turbine engines are used widely across the world for power generation, commercial
and military aircraft, and even ground vehicles. Their high power to weight ratio prove
especially tting for aviation purposes, and the demand for e cient and powerful engines is
extremely high. Optimization of these systems in the preliminary design stage is generally
well understood and several books have been written on the subject, but the optimization
process for the detailed design of speci c components is still an area of ongoing research. A
deterministic optimization approach is often not possible due to the high number of design
variables and complexity of the ow analysis, as well as the high possibility of the search
stopping in a local optimum. Global optimization techniques are often e ective for problems
such as turbine design, but are ine cient in comparison, often requiring thousands of function
evaluations to produce a worthwhile result. Although this may prove possible with a simple
through- ow analysis code, using a full 3D Navier-Stokes analysis could require months or
years to complete, which is an unreasonable amount of computation time. Finding a more
e cient optimizer that nds a global optimum in only a few hundred function evaluations
or less would prove extremely useful in allowing a higher delity analysis and optimization
without excessive computational costs.
1.1 Objectives
The objectives of this paper include producing and validating an optimizer that is
both e ective and e cient at evaluating turbine cascade problems, conduct an aerodynamic
optimization of a multi-stage turbine, and conduct an aerothermal optimization of the same
multi-stage turbine. These objectives are discussed in detail below.
1
1.1.1 Objective 1: Produce and Validate Optimizer
This objective includes evaluating a number of global optimizers, noting their strengths
and weaknesses, and selecting the appropriate optimizer for a turbine design problem that
capitalizes on it?s canonical strengths. This objective focuses on using a canonical, or close
to canonical, optimizer and limit the use of unreliable and problem speci c optimizer en-
hancements. The optimizer has been written by the author and validated using a number
of test cases where the results are compared against published articles to evaluate optimizer
performance.
1.1.2 Objective 2: Turbine Aerodynamic Optimization
The optimizer has been applied to a turbine aerodynamic optimization problem. This
analysis includes a preliminary 2D inviscid analysis as well as a 2D Navier Stokes analysis. In
order to provide a comparison for results, the optimized blade geometry was evaluated against
the blade geometry of the T55-GA-714A engine currently used on the CH-47D Chinook
helicopter.
1.1.3 Objective 3: Turbine Aerothermal Optimization
The third objective improves on the aerodynamic optimization by investigating the
e ects of blade geometry on the aerodynamic e ciency and overall convective heat transfer
to the blade. This optimization seeks to produce turbine blade geometries that allow a
higher turbine inlet temperatures due to reduced hot spots, while maintaining a high level
of aerodynamic performance. Higher turbine inlet temperatures may be directly correlated
to increased total power output and decreased speci c fuel consumption[1]. This objective
includes an investigation into whether a blade geometry that allows higher turbine inlet
temperatures and slightly lower aerodynamic performance may prove to be a more optimal
blade.
2
1.2 Literature Review
The literature reviewed for this research included a variety of books on turbine design
and a number of published articles researching optimization, turbine design, and CFD spe-
ci c topics. Although a large number of books were referenced in this paper, the books
most heavily referenced for detailed turbine blade design were Turbine Aerodynamics by
Ronald Aungier[2] and Principles of Turbomachinery in Air-Breathing Engines by Erian
Baskharone[3]. Both of these books provided a thorough walk-through of detailed blade
design that was well written and easy to understand. Baskharone?s text also provided case
studies that were helpful in understanding the fundamentals of blade design. Because both
of these references were written within the past 10 years, they represent a fairly up to date
review of current design techniques and processes, including detailed blade geometry design.
Wilson?s text[4] was also referenced heavily for preliminary and stage design. The methods
presented in these books were viewed by the author as the most current methods for blade
parameterization and design, and are therefore used as a basis for comparison.
CFD is a quickly growing eld of study and is a very integral part of current turbo-
machinery research. Fluid Dynamics and Heat Transfer of Turbomachinery[5] by Lakshmi-
narayana along with the user?s manual for the CFD code were the main references used by
the author to ensure responsible CFD utilization. All CFD validation studies were conducted
in reference to Arts[6] due to the fact that it is commonly cited and used for CFD validation
in other publications. Articles by Horlock & Denton[7] as well as Denton & Dawes[8] provide
a background of CFD use in turbomachinery design compared to traditional methods over
the last 50 years. Because the focus of this thesis is on the performance of the optimizer
and a 2D blade-to-blade optimization, a well known and validated CFD code will be used
rather than conducting a literature study to nd the most current CFD techniques. Several
articles have been written by Dr. Chima on the use and validation of the primary CFD code
used in this thesis, Rotor Viscous Code Quasi-3D (RVCQ3D), and all have been reviewed in
detail by the author.[9][10][11][12]
3
Optimization of turbomachinery using metaheuristics is a much more recent topic and
is not yet preferred over the traditional design methods. Research in this eld is much more
sporadic and to the author?s knowledge no consolidated reviews of this topic exist. A number
of articles using stochastic algorithms to optimize turbine blades were reviewed by the author
and a summary of the more relevant articles is given in this literature review. The earliest
articles using metaheuristics to optimize blade shape appeared in the mid-1990?s. Turbine
Airfoil Design Optimization by Goel[13] was published in 1996 and used both circular arcs
and B ezier curves to represent the blade geometry. Several articles followed in 1998 and
1999 including two by Giannakoglou[14] [15], and articles by Pierret[16], Trigg[17], and
Dennis[18]. Each conducted an aerodynamic optimization using either a Genetic Algorithm
or Arti cial Neural Network and used a very similar blade parameterization technique with
large restrictions on the movement of the control points. More recent work includes 3D
blade optimization using GAs as in articles by Chen[19] and Oksuz[20], which use Non-
Uniform Rational B-Splines (NURBS) to parameterize the blades but maintain circular
arc segments at the leading and trailing edges. In addition to those listed, several other
articles were reviewed and all had two common factors. The rst is the use of constraints
to restrict the blade parameterization, either through limiting the movement of the control
points, using prescribed circular arc segments, or using a lower delity blade parameterization
method altogether. The second is conducting a purely aerodynamic optimization without
considering heat transfer and sometimes structural requirements as well. It was also apparent
that Genetic Algorithms were the optimizer of choice, with only one article using a di erent
optimizer.
The main reason metaheuristic optimizers are not preferred over traditional blade de-
sign methods is the high computational cost. However, analyzing several objective functions
at once, such as aerodynamic, thermal, and mechanical considerations, may make the com-
putational time more reasonable compared to doing each step manually. Less restrictions on
4
the search space may also nd better solutions, especially when multiple objectives are con-
sidered. Therefore, more research is needed in this area, where there are very few published
articles available.
1.3 Problem Background
Turbine design is a very complex subject, but is also well documented by several books
that break down the steps of preliminary design into a manageable problem[4][1][2][3]. De-
tailed stage or blade design is a much more ambiguous topic with numerous methods available
in literature. These methods of detailed design usually rely on a large amount of practical
experience on the part of the engineer and knowledge of current best practices. The direct
blade design is usually preferred over the inverse design method, and one of a number of blade
curvature procedures is used along with manual iteration until a desirable blade is achieved.
This process provides ample opportunity for improvement upon using recent advanced in 2D
and 3D turbine modeling programs.
In order to understand the applicability of the detailed design process discussed in this
paper, it is important to rst understand the preliminary design process. The preliminary
design process outlined here is a summary of the discussion in The Design of High-E ciency
Turbomachinery and Gas Turbines by Wilson and Korakianitis[4]. The design process starts
with a number of given turbine speci cations, which are listed below.
1. the uid to be used
2. the uid stagnation temperature and pressure at inlet
3. the uid stagnation or static pressure at outlet
4. uid mass ow or turbine power output required
5. shaft speed
5
Figure 1.1: Example Velocity Diagrams
for a Stator[21]
Figure 1.2: Example Velocity Diagrams
for a Rotor[21]
An initial investigation is usually conducted to determine if a single stage con guration
will su ce for the design. If the stage loading is too high, a multi-stage turbine must be
used and a stage work split is speci ed. With this given information and an assumption
of the inlet and exit ow angles, the velocity triangles for the turbine may be speci ed.
Velocity triangles depicting the in ow and out ow velocities and angles for each blade row
are a fundamental concept in turbine design, and an example is shown in Figures 1.1 and
1.2 . The velocity triangles are determined from the design work output requirements for
the engine and then the blade pro les are designed to match the velocity triangles, not the
other way around. The concept of velocity triangles is based on the Euler turbine equation.
A modi ed version of this equation is shown as Eqn. 1.1, which states that the enthalpy
change across a blade row is a function of the shaft speed, mean blade radius, and the change
in tangential velocity.
H2 H1 = !(r2C 2 r1C 1) (1.1)
In order do de ne the velocity triangles three parameters are required: the work (or
loading) coe cient , the ow coe cient , and the stage reaction Rn. These parameters
are discussed in detail in Wilson[4] and calculated using Eqns. 1.2-1.4. Velocity triangles are
often calculated from initial guesses for ow angles and then the stage reaction is checked. If
6
the stage reaction is unsatisfactory, the ow angles are adjusted and iterations are completed
until the desired stage reaction is obtained.
= gc
ex
inh0
u2 =
ex
in(uC )
u2
(1.2)
= Cmu (1.3)
Rn = hs;rotor h
0;stage
(1.4)
After velocity triangles for the mean line and other radial locations are obtained, the
number of blades per row, blade shape, and blade size may be determined. This detailed
design of the blade rows will be the focus of this paper. In order to focus on the detailed
blade design, it will be assumed that the prior preliminary design steps have already been
completed. The data from the T55-GA-714 engine will be used as a starting point for this
analysis, including the velocity triangles, annulus geometry, shaft speed, and thermodynamic
properties at each stage location. It is important to note that the information used in this
analysis from the T55 engine design report is proprietary to Honeywell International Inc.
and will not be disclosed in this report. All comparisons between optimized results and the
original T55 geometry will be presented as ratios or percentages to protect this information.
7
Chapter 2
Evolution Strategies
2.1 Algorithm Selection
A variety of meta-heuristic optimizers are available to apply to this problem, and choos-
ing the appropriate one will have a great in uence on the nal results. Although many
di erent types of optimizers have been applied to turbine optimization problems in the past,
not all of them t the problem at hand. Because this is a continuous problem, optimizers
that are strong at solving combinatorial problems but are weak in a continuous solution space
will not be considered. These include Tabu Search, Ant Colony optimizers, and others. To
narrow the choice of optimizers further it is desirable to obtain a variety of near-optimal
results at the end of an optimization process, and therefore non-population based optimiz-
ers such as Simulated Annealing would not be as advantageous. Three primary population
based optimizers that often perform well for aerodynamic problems are Genetic Algorithms
(GA), Evolution Strategies (ES), and Particle Swarm Optimizers (PSO). A short descrip-
tion of each is given below, and a much more detailed description may be found in the book
Metaheuristics for Hard Optimization[22].
Genetic Algorithms are based on the premise of natural selection and small evolutionary
changes. In the algorithm, members of an initial population of solutions are combined in
some sort of fashion to produce \o spring" that bear some of the characteristics of the
parents. These solutions are canonically represented as a bit string, but may also be a series
of real numbers that de ne a solution for the problem. Several o spring are produced and
then evaluated based on some function that determines tness, an this objective function
may either be maximized or minimized for the particular problem. At this point the group of
parents and o spring are ordered by tness and the most t solutions become the parents for
8
a new generation while the remainder are discarded. This process is repeated until some pre-
determined stopping criteria is met, such as a maximum number of generations. A mutation
is usually included as well, which selects and perturbs a solution at random to increase
diversity in the population. This mutation rate is often very low, such as a few percent of
the o spring or less. After a large number of generations, the solutions will continue to get
better until they converge on a near-optimal solution. GAs have proven to be very successful
in nding the global optimum for a variety of problems. They do not easily get stuck in local
minimums, provide multiple solutions, and are easy to implement. However, GAs can be
very ine cient for certain problems, especially if the solution space is not well de ned or the
initial solutions are di cult to randomly place across the breadth of the solution space. GAs
generally require a larger population size as well to ensure enough variation in the o spring
to make an e ective search.
Particle Swarm Optimizers are inspired by birds ocking together or sh swimming
in a school. In this algorithm a population of solutions, or \particles", move through the
solution space using a velocity that is updated after each move. Each member?s velocity
is updated based on it?s own inertia, the location of the best solution seen by the particle,
and the location of the best solution seen by anyone in communication with the particle.
The algorithm uses this collective knowledge of the swarm to locate the global optimum in
the solution space. There are a number of variations and enhancements to the PSO algo-
rithm, such as limiting the communication to neighborhoods in the particles and limiting the
maximum velocity. In general, the PSO algorithm is very easy to implement and converges
fairly fast. Disadvantages of this algorithm include being stuck in local optimum easier than
other meta-heuristics and sensitivity to input parameters. The ability of a PSO to nd the
global optimum is highly dependent on the location of the initial population in relation to
the optimum and the starting inertia of each particle.
Evolution Strategies algorithms are very similar to GAs, and in fact were developed
in Europe around the same time as GAs in the 1970s. As in GAs, ES algorithms use a
9
population of solutions to produce o spring, sort the o spring based on a tness value, and
select a certain number of o spring to proceed as the parents of the following generation.
The di erence between GAs and ES algorithms is the move operator, which is the method of
producing o spring. In an ES algorithm 100% of the o spring are mutated from the parent
solutions. This is accomplished through a random perturbation using a normal distribution
and a standard deviation value . After the designated number of o spring have been
mutated from the parent solutions, the solutions are sorted and the most t will proceed to
the next generation. This move operator will be discussed in more detail in the following
sections. ES algorithms are easy to setup, can easily handle problem constraints and nd
global optimum, and have the ability to self-adjust parameters during the optimization
process to adapt to the search space. They also perform well with micro-populations of
only a few members because the value controls the diversity of the search instead of the
diversity in the parent population. Disadvantages of ES algorithms include the possibility
of getting stuck in a local optimum if the parent population is far from the global optimum
and the values are not su cient enough to move the solutions the required distance.
Although any of these three algorithms would be capable of nding a optimal or near
optimal solution, performance may vary greatly between the three based on their canonical
form and how they are chosen to be implemented. In an aerodynamic optimization problem,
an algorithm that searches fast and requires relatively few function evaluations is desired
because each function evaluation may be extremely costly, especially if a Navier-Stokes
solver is used for the evaluation. The limits of the search space are also very di cult to
determine in aerodynamic problems, if they can be determined at all. Sometimes random
initial solutions are used, but more often a known good solution is used as the starting point.
Therefore, an algorithm that can search outward as e ciently as it can search inward would
be desirable. A simple, robust algorithm is also very bene cial due to the long optimization
times. If the optimization process only required a few seconds or minutes to complete it
10
would be simple to restart if a problem occurred, but for an optimization process that takes
several hours or days a robust algorithm becomes extremely important.
Considering the characteristics desired in an algorithm for turbine optimization, Evo-
lution Strategies is the most appropriate meta-heuristic for this problem. The larger pop-
ulation sizes required to make GAs and PSOs e ective would greatly increase computation
time while using a micro-population with an ES algorithm would drastically reduce the num-
ber of required function evaluations. Airfoil and blade geometry poses an inherent problem
with the move operators in genetic algorithms and PSOs because these algorithms work best
when the initial solutions are randomly placed across the breadth of the search space. This
typically does not happen in airfoil optimization because the limits of the search space are
often unknown, and the author often chooses a set of known airfoil geometries as the starting
points. Move operators in GAs and PSOs are more inclined to move inward towards the
local and global optima and do so very quickly and e ciently, while outward movement past
the initial starting points is slower and more di cult. In order to move outside the area of
the search space bounded by the parent solutions, GAs must wait for a random mutation
to occur that moves the solution in the right direction. Canonically the mutation rate in a
GA is very low, usually less than 2%, and the mutation is randomly generated such as a bit
ip or a random number generation. This is often overcome by raising the mutation rate to
much higher levels, but this only results in less e ective version of an ES algorithm and an
ine cient search. PSOs also have di culty moving outside the bounds of the initial parent
population if the inertia from the initial velocities cannot carry the particles far enough to
get close to the global optimum. An example of a minimization problem is shown in Figure
2.1 with white stars depicting the initial solutions. The example on the left is well suited
for a genetic algorithm or a PSO because the limits of the search space are well de ned and
random initial solutions can be placed across the breadth of the search space. The example
on the right is more typical for an airfoil optimization problem, where the initial solutions
11
Figure 2.1: Search space comparison for minimization problem
are known well-performing airfoils and the limits of the search space are not known. Well-
constructed genetic algorithms and PSOs should nd the optimum to this problem, but a
better suited optimizer such as an ES algorithm may be able to nd a comparable or better
solution at a much lower computation cost. ES algorithms are also very simple and do not
have as many issues with robustness as PSOs, which can fail due to parameter sensitivities.
2.2 Solution Encoding and Move Operator
An essential part of the algorithm used in this research, which will be referred to as
ESTurb, is determining how the optimizer will move the solutions through the search space
and how those solutions will be represented. The de nition of the move operator and solution
encoding in an adaptive optimizer has a tremendous e ect on the algorithm?s performance.
Usually, the type of move is determined by the type of optimizer, such as crossover for a
genetic algorithm or velocity updates for a PSO. The user still has a great degree of freedom
to customize the move operator to t the problem. It is important to ensure the move is
stochastic in nature for an adaptive optimizer to nd the global optimum rather than getting
stuck in local optima. Developing a move operator comes with many considerations, such
as using a discrete or continuous move, searching only feasible solution space or allowing
infeasible solutions, and ensuring scaling is not an issue. The solution encoding is always
12
problem dependent, but must work well with the move operator to produce an e ective
algorithm.
When selecting a solution encoding there are several factors to be considered. The
encoding must fully represent the solution, which is the upper and lower surface curves of
the turbine blade, in as few variables as possible. Choosing two few variables would not
fully represent the solution and adding too many variables, such as a series of coordinate
pairs to describe the blade surface, would add so many dimensions to the problem that the
search space would be unreasonably large. The solution encoding must also limit unnec-
essary restrictions to the search space which could not allow an optimum solution to be
represented. Obviously, if the optimum solution cannot be encoded, then it cannot be found
by the algorithm. A good solution encoding for an aerodynamic problem should be able
to represent almost any feasible airfoil using a small number of variables, and maintain a
smooth continuous curve for the blade surfaces. There are a number of options available,
such as the 11 variable PARSEC encoding shown in Figure 2.2, spline connected polyno-
mial segments shown in Figures 2.4 and 2.5, or a series of polynomial curves such as B ezier
Curves as shown in Figure 2.3. Most parameterization methods have inherent restrictions
on representing certain geometries, such as the problem PARSEC encoding has with fully
describing the trailing edge of an airfoil. B ezier Curves are used in ESTurb due to it?s ability
to represent almost any curve, ease of use, and small number of variables. B ezier Curves
are a form of parametric curve where the polynomial coe cients are controlled through a
number of control points. When B ezier Curves are used with other optimizers the control
points are almost always limited to certain regions in order to de ne the limits of the search
space and allow a more e ective search for a GA or PSO. An example of this is the boxes
depicting the limits for each control point shown in Figure 2.3. However, these restrictions
prevent the algorithm from representing a large portion of the search space, which may
include the optimum solution. By using an ES algorithm in ESTurb these restrictions are
not necessary because the algorithm searches outward very e ectively. This may allow the
13
Figure 2.2: PARSEC encoding
Figure 2.3: Example of using B ezier
curves to encode blade[23]
Figure 2.4: Spline connected polynomial
encoding angles and radii[21]
Figure 2.5: Spline connected polynomial
encoding de ning points[21]
algorithm to nd better solutions that simply could not be represented by other algorithms.
Bernstein polynomials may also be a good option for airfoil representation because they
are the mathematical basis for B ezier curves, but the binomial coe cients would be more
di cult for the algorithm to move through the search space when compared to the B ezier
control points. For this reason, B ezier curves were seen as the better option for use with an
ES algorithm. A generalization of B ezier splines also published by French engineer Pierre
B ezier known as Non-Uniform Rational B-Splines (NURBS), would be an excellent option
for representing 3D shapes and is commonly used in the industry.
14
Table 2.1: Example solution encoding
Figure 2.6: Blade surface for example so-
lution encoding
Because a B ezier curve requires a set of control points to de ne a curve, the solution
encoding for ESTurb consists of a set of (x;y) coordinate pairs for both the upper and lower
surface of the blade. The number of control points used for each surface is important in
allowing the curve to easily represent any shape while not over-expanding the search space
by adding too many points. After completing a series of runs in the validation cases discussed
below, using six points for each surface proved to be e ective in representing almost any shape
while not adding too many variables to the problem. An example of a solution encoding is
given in Table 2.1, where each row is a speci c control point and the columns represent the
x and y coordinates for each surface. Figure 2.6 shows the curves generated by the given
solution as well as the location of the control points.
The canonical move operator for an Evolution Strategies algorithm is used in ESTurb.
The move applies a random perturbation to the solution based on a normal distribution,
where is the desired deviation. The value is adjusted based on the results found according
to the one fth rule, which is integral part of a canonical Evolution Strategies algorithm. The
one fth rule looks at the tness of the o spring compared to the parents that generated each
o spring. If greater than one fth of the o spring are more t than the parents, then there
are lots of areas in the search space that provide good solutions and the value is increased
to diversify the search and prevent the algorithm from being stuck in a local optima. If
less than one fth of the o spring are more t than the parents, then the algorithm is not
15
nding a lot of good solutions and the value is decreased to intensify the search. This
typically happens when the algorithm nds a solution near the global minimum. To decrease
the value of it is typically multiplied by 0:85, and is divided by 0:85 for an increase.
A method of preventing stagnation, which will be referred to as \ restart," is also used in
ESTurb. Stagnation may occur if the value of becomes so small that it has no e ect on
the solution, and if this occurs the subsequent generations will not improve results. If the
value of drops below some predetermined level, \ restart" will reset the value of to the
initial value or some other designated number. This method works well in problems where
the solution space is not a smooth gradient and there may be many local optima that are
close in vicinity and tness to the global optimum.
Because the solution in ESTurb involves a matrix of numbers, every time an o spring is
generated the move is applied in both x and y directions for each control point as shown in
Figure 2.7. In an evolutionary algorithm the terms \move" and \mutate" are synonymous,
and are used interchangeably in this paper. The rst control point is omitted in the move in
order to keep the leading edge located at the origin and the x direction for the nal control
point is also omitted to maintain a unit axial chord. However, the move operator is applied
to the trailing edge in the y direction for the upper surface and the trailing edge for the
lower surface is set equal to the upper trailing edge if a cusp is desired.
The mutation of the second control point is of particular interest because it determines
the angle that the curve leaves the rst control point. For each B ezier curve, the function is
tangent to the line drawn between the edge control point and the rst interior control point.
Therefore, the location of the second control point will determine if the leading edge of the
blade is blunt or has an angle where the rst control points for each surface meet. By not
restricting the movement of the second control point a sharp leading edge may be produced,
which may be applicable for a supersonic turbine but is not desirable in this case. Ensuring
the rst and second control point for each surface lies on the same line ensures a blunt leading
edge, which is necessary for a subsonic turbine. Two di erent methods were used to ensure
16
Figure 2.7: Move operator example for a single control point
a blunt leading edge on the blade, the rst ensuring the ow is perpendicular to the ow
at the location of the leading control points, and the other method simply ensures the rst
two control points for each surface lie on the same line. The rst method is accomplished by
setting the x value of the second control point on each surface (CP2;lower &CP2;upper) to zero
and not allowing the second control point to mutate in the x direction. This ensures the
curves going into the rst control points (CP1;lower & CP1;upper) are vertical, and therefore
a smooth leading edge is obtained. The second method allows CP2;upper to mutate freely.
After the mutation is complete, the distance between (CP2;lower and CP1;lower) is recorded
and mutated based on the current value. Then CP2;lower is relocated to a position on the
line including CP1;lower, CP1;upper, and CP2;upper the new mutated distance away from the
leading edge. The second method is desirable whenever the ow is not horizontal, such as
entering a turbine rotor after passing through a nozzle guide vane. The rst method was
used in ESTurb for the optimizer validation cases and the inviscid turbine analysis, but was
updated to the second method prior to the viscous turbine analysis. Examples of an airfoil
leading edge using these two mutation methods are shown in Figures 2.8 and 2.9.
17
?0.15 ?0.1 ?0.05 0 0.05 0.1 0.15
?0.1
?0.05
0
0.05
0.1
0.15 CP2,upper
CP1,lower and CP1,upper
CP2,lower
Figure 2.8: Second control point mutation
option 1
?0.06 ?0.04 ?0.02 0 0.02 0.04 0.06
?0.04
?0.02
0
0.02
0.04
0.06CP2,upper
CP1,lower and CP1,upper
CP2,lower
Figure 2.9: Second control point mutation
option 2
2.3 Optimizer Validation
After determining the type of optimizer to be used, along with the solution encoding
and move operator, the initial version of ESTurb was written following the canonical ES
optimizer form. In order to validate the optimizer and determine performance, two separate
test cases were performed. Each test case evaluated the optimizer performance compared to
the results in a published article. Considering the reasons Evolution Strategies was picked in
section 2.1, the algorithm is expected to conduct an e cient search outward from the seeded
solutions and should converge quickly. The algorithm should also nd improved solutions
that the other solution encodings in the reference papers may not be able to represent due
to the restrictions placed on geometry due to the parameterization method. A owchart for
ESTurb is shown in Figure 2.10.
2.3.1 Validation Case 1
The rst validation case compared optimizer performance against the results published
in \Designing Airfoils using a Reference Point based Evolutionary Many-objective Particle
18
Figure 2.10: Flow chart for ESTurb algorithm
19
Swarm Optimization Algorithm," written by Wickramasinghe[24]. Wickramasinghe?s paper
attempted to optimize the geometry of the NLF0416 airfoil used on a UAV with respect to
certain criteria. Six di erent objectives were evaluated and a Particle Swarm Optimizer was
used for his algorithm. The ow solver used for both the reference paper and this analysis
was XFoil, a well known subsonic airfoil analysis code written by Dr. Mark Drela of MIT.
The objectives for this problem are listed below, which are described by Wickramasinghe as
desirable characteristics during several di erent ight modes for a UAV.
1. Minimize Coe cient of Drag (at CL = 0:5, Re = 4:0 106, and Mach= 0:3)
2. Minimize CD
C
32
L
, (at oating , Re = 4:0 106, and Mach= 0:3)
3. Minimize Cm (at CL = 0, Re = 4:0 106, and Mach= 0:3)
4. Minimize 1C2
L;max
(at oating , Re = 4:0 106, and Mach= 0:3)
5. Minimize 1C2
L
(at = 5o, Re = 2:0 106, and Mach= 0:15)
6. Minimize 1 xtr (at = 5o, Re = 2:0 106, and Mach= 0:15)
The ES algorithm performed well in nding near-optimal solutions and performed ex-
tremely well in computational time compared to Wickramasinghes PSO algorithm. Eight
runs were initially performed to evaluate the algorithm while changing certain parameters,
such as maximum number of generations, initial sigma values, and toggling the one fth
rule on or o . It was found that 50 generations were su cient with the original ESTurb
algorithm to obtain good results. This algorithm included the canonical one- fth rule with
restart, but no other enhancements such as recombination were necessary. Three di erent
runs were completed with the baseline optimizer and compared against the results published
in Wickramasinghe?s paper. Each run began with a seeded initial solution consisting of six
NACA series airfoils (0012, 2412, 4415, LS(1)-0417, 64-210, 63-015), a NLF0416 airfoil, and
a RAE100 airfoil. The NLF0416 airfoil was used as the seeded initial solution in the refer-
enced paper. The algorithm used a + replacement scheme, where both the parents and
20
1 2 3 4 5 610?3
10?2
10?1
100 Best Airfoil Objective Values
Objective Number
Obj Function Value
ES AirfoilReference
PSO Airfoil
Figure 2.11: Test case one objective function
results (sought to minimize)
Figure 2.12: Test case one computation
time comparison
o spring were considered for selection into the next generation. Five o spring were gener-
ated per parent in each generation, which resulted in 40 function evaluations per generation
and 2000 function evaluations across 50 generations.
The results for both the ES algorithm and the benchmark algorithm are summarized
in Figure 2.11, where all objective functions are sought to be minimized. Figure 2.13 shows
and example of the converged XFoil results for the best solution. Each of the three runs pro-
duced very similar results in terms of objective function values, so only one airfoil is shown.
The algorithm converged on a solution by about the 50th generation, so any subsequent
generations generally added computation time without improving results. Figure 2.14 shows
the sum of the objective functions over generations. By approximately generation 50 the
algorithm converged to a near-optimal solution. \ restart" improved the consistency of the
algorithm by ensuring it did not easily get stuck in local optima, but did not signi cantly
improve the nal results or computation time. Random initial solutions were also tried in
addition to seeded initial populations. This did not signi cantly impact the nal results, but
did slightly increase the computation time because the algorithm required approximately 10
generations to move from the random shapes to an airfoil shape equivalent to the NLF0416
airfoil. Even though the objective function results were very similar between the runs, the
21
Figure 2.13: XFoil output for converged airfoil
airfoil geometry was noticeably di erent in some of the runs. An example airfoil geometry
from another run is shown in Figure 2.15.
The computation time required for the algorithm was very reasonable considering the
large computational cost of the ow solver. The fastest run consisting of 50 generations,
which also produced the best objective function values, had a run time of 1.11 hours. All
computations were performed in Matlab on a 2.0 GHz dual-core processor. The PSO algo-
rithm required approximately 36 hours to run on a 2.3 GHz dual-core processor. Because
the XFoil ow solver accounted for over 80% of the computational time per run, the main
di erence in computation time between the ES and PSO algorithms came from the number
of function evaluations required. The PSO algorithm required 10,000 function evaluations
22
Figure 2.14: History of objective function value over generations for initial eight runs
Figure 2.15: Example of variation in nal population from Figure 2.13
23
while the ES algorithm required about 2,000. The use of a Hyper-Volume method for eval-
uating multiple objectives rather than using a linear combination of the objective functions
also increased time for the PSO algorithm.
2.3.2 Validation Case 2
The second paper selected for comparison is \Multi-Objective Evolutionary Optimiza-
tion of Subsonic Airfoils by Meta-Modeling and Evolution Control, written by D?Angelo[25].
In this paper a multi-objective optimization is performed using a Pareto Front to minimize
the coe cient of drag and maximize the coe cient of lift. This paper also includes work
in meta-modeling, arti cial neural networks, and evolution control, but the comparison will
focus on the solution encoding and results. For problem two, the ES algorithm was generally
unchanged with the exception of adding a Pareto archive instead of a parent population, and
the parents were selected randomly from the rst rank of this Pareto archive to create o -
spring. The algorithm still used a micro-population of eight members and made no changes
to the solution encoding, move operator, or any other fundamental aspects of the algorithm.
Fifteen di erent runs were completed in order to gain a su cient amount of data to evaluate
the performance. Three runs for each of the ve scenarios below were completed:
1. 75 generations, 300 member archive
2. 150 generations, 300 member archive
3. 250 generations, 300 member archive
4. 150 generations, 100 member archive
5. 150 generations, 500 member archive
The data from these runs is too much to display in one gure, but in general the
larger archives performed better with no change in computation time. Obviously, more
generations led to better results, but in general returns diminished substantially after about
24
150 generations. There was no single run that dominated the other runs when comparing
Pareto fronts, but the run with 150 generations and a 500 member archive performed well
and will be used as the primary comparison against D?Angelo?s results.
Figure 2.16 shows the results from D?Angelo, both using a ve parameter airfoil rep-
resentation and a seven parameter representation. Each representation used two B ezier
curves, but instead of de ning the airfoil geometry the curves de ned the mean camber line
and airfoil thickness. The number of parameters is the number of B ezier control points used.
D?Angelo?s best result using the ve parameter representation was selected for comparison,
as well as the seven parameter result utilizing Bayesian learning, which performed signi -
cantly better. The ES results for 150 generations and a 300 member archive achieved much
better results in terms of the Pareto Front. Although not shown, each of the remaining 15
runs dominated at least a portion of the seven parameter Pareto front, and all but three of
the runs completely dominated the ve parameter Pareto front from D?Angelo. The three
remaining runs dominated a majority of the ve parameter front, but did not dominate the
whole front.
The ve parameter model from D?Angelo required approximately two hours to complete
using just XFoil for a ow solver on a Pentium IV 2.8 GHz processor. The computational time
for the seven parameter model was not stated, but it can be assumed with the higher degrees
of freedom that the algorithm required at least as long as the ve parameter model. The
ES algorithm required approximately one hour to complete 50 generations, which resulted
in a total computation time of 2.98 hours for the 150 generation, 300 member archive run.
However, even though the entire algorithm required more time to complete, the Pareto front
surpassed the ve parameter front in approximately 18 minutes and dominated the seven
parameter front in approximately 25 minutes, as shown in Figure 2.17.
25
0.004 0.006 0.008 0.01 0.012 0.014 0.016 0.018 0.020.5
1
1.5
2
2.5
3
Pareto Front
Cd
Cl
DAngelo 5 Par
DAngelo 7 Par
ES Front
Figure 2.16: Pareto front results
2.3.3 Conclusions From Optimizer Validation
Both test cases proved that ESTurb is both e cient and e ective at solving aerodynamic
optimization problems. By evaluating the strengths and weaknesses of di erent adaptive op-
timization methods, an ES algorithm was expected to converge on the optimal solution faster
than GA or PSO algorithms due to the move operator. The use of B ezier curves for the
solution encoding was expected to allow the algorithm to better represent an aerodynamic
shape, and possibly nd better results that other parameterization methods could not rep-
resent. Although using B ezier curves would increase the size of the search space, the move
operator in the ES algorithm was expected to be e ective enough to still allow an e cient
search with decreased computation time.
As seen in both test cases, the computation time was drastically reduced. This was
mainly due to a reduction in the number of function evaluations since the cost of each
was extremely high. The increased e ciency in the search allowed the number of function
evaluations to be greatly reduced while still maintaining comparable or better results. In
test case two especially, the ES algorithm well outperformed D?Angelo?s algorithm, nding
26
0.004 0.006 0.008 0.01 0.012 0.014 0.016 0.018 0.02
1
1.5
2
2.5
3
Pareto Front vs. Computation Time
Cd
Cl
DAngelo 5 Par.
DAngelo 7 Par.
Initial Seeds
8.73 min
17.6 min
24.7 min
1.06 hrs
1.95 hrs
2.96 hrs
Figure 2.17: Pareto front progression with computation time
a Pareto front that was well in front of his best results. This validation shows that ESTurb
has the potential to be very e ective for turbine optimization, and the drastic reduction in
computation time will be even more prominent as the cost of function evaluations go up,
especially when using an inviscid or viscous CFD analysis code. Additional details on each
validation case may be found in a previous article by the author[26], including algorithm
parameters used, a more detailed discussion on evolutionary computation, and algorithm
sensitivity to parameter values.
27
Chapter 3
Turbine Aerodynamic Optimization
3.1 Turbine Fundamentals Review
This paper will not spend an extensive amount of time on turbine fundamentals since
a number of great books have been written that discuss this topic in a much more detailed
and eloquent way[4][2][3]. This section will focus on reviewing the relevant aspects of turbine
fundamentals that are applicable to the optimization process used in this paper. It is assumed
that the reader is familiar with the conservation of mass, momentum and energy equations,
as well as the equation of state and the second law of thermodynamics.
The basic principle behind the turbine section in a gas turbine engine is to extract
useful power from the thermal energy in the gas leaving the combustor. A turboshaft engine
normally has a gas generator turbine (GGT) and a power turbine (PT), and each may
consist of one or more stages. The GGT is almost always the high pressure turbine directly
downstream of the combustor and it extracts power in order to drive the compressor. The
PT extracts power from the remaining thermal energy to drive the main transmission or
other components. Turbofan and turboprop engines operate in the same manner, except the
PT shaft is connected to the fan or propeller either directly or through a gearbox.
Each stage in a turbine section consists of a stator, sometimes referred to as a Nozzle
Guide Vane (NGV), and a rotor. The stator is stationary with respect to the frame and
imparts a angular momentum in the uid with minimal losses. The rotor is rotating with an
angular velocity equal to the shaft speed and the product of the pressure force on the blade
and the angular momentum of the rotor disk result in shaft power. By applying Newton?s
2nd law the pressure force on the blade equals the net momentum change in the uid, and
this premise leads to the well known Euler Turbine Equation (Eqn. 3.1)
28
Power = ! _m(rinC ;in rexC ;ex) (3.1)
The Euler Turbine Equation will be the basis for the rotor objective functions to increase
shaft power. The other major component that will be used is e ciency. There are several
di erent means of assessing turbine performance, including total-to-total e ciency, total-
to-static e ciency, polytropic e ciency, kinetic energy loss coe cient, total pressure loss
coe cient, entropy based e ciency, and others. Polytropic e ciency is commonly used
because it is a means of relating engines with di erent compressor ratios, but because we
are only concerned with one engine this is not necessary. The total-to-total (t t) e ciency
is su cient for all but the last stage, where the thermal energy at exit is expected to be
used. For the last stage, where any exit swirl or other thermal energy will be wasted, the
total-to-static (t s) e ciency must be used. Note that the only di erence between the
total-to-total and total-to-static e ciencies is using the static exit pressure rather than the
total exit pressure in the denominator. The total pressure loss coe cient provides a good
measure of the e ciency of a stator, where there is no work being done, and is shown in
Equation 3.4.
t t =
1 T0;exT0;in
1
P0;ex
P0;in
1
(3.2)
t s =
1 T0;exT0;in
1
Ps;ex
P0;in
1
(3.3)
! = P0P
0;in
= 1 P0;exP
0;in
(3.4)
29
3.2 Preliminary Inviscid Analysis
The purpose of the inviscid analysis is to investigate the potential for ESTurb to conduct
optimization for turbomachinery. ESTurb performed very well for simple airfoil optimiza-
tions, but it important to apply the concepts of turbine thermodynamics to an objective
function in as simple a manner as possible. Therefore, an inviscid analysis was completed
to adapt the algorithm to using a CFD code and test the new objective functions prior
to adding complications such as turbulence and heat transfer. Another main focus of the
inviscid analysis is to determine the e ectiveness of using B ezier curves for turbine solution
encoding compared to other conventional blade parameterization methods. This analysis
is much easier to complete using the less costly inviscid CFD code rather than the full
Navier-Stokes CFD analysis.
3.2.1 Model and Assumptions
The application of ESTurb is very modular, so applying it to a turbine problem is as
simple as switching out the evaluation subroutine with a new one. The evaluation subroutine
can be viewed as a \black box", receiving the solution encoding for a particular o spring as
input and producing an objective function value as output. The remainder of the ESTurb
algorithm remains the same. Because this will only involve a single objective function, the
ESTurb version used for validation case one will be used for all turbine problems. Because a
CFD code is used as the ow solver, the evaluation subroutine consists of receiving the input,
generating a grid, running the CFD code, extracting the results, and calculating the objective
function. An in-house Euler CFD code written by the author was used for this portion of
the analysis. The ow is assumed to be inviscid and the blades evaluated are assumed to
be stationary. A distance of a half chord was used both upstream and downstream of the
blade limits and the cascade was linear. There were also no structural requirements for this
analysis, only aerodynamics were considered in the objective function.
30
Figure 3.1: Elliptical H-grid
Figure 3.2: Leading edge
without vertical clustering
Figure 3.3: Leading edge
with vertical clustering
3.2.2 Grid Generation and CFD Code
The details on grid generation and the Euler code are expounded upon in a previous
article by the author[27] and discussed at length in the Handbook of Grid Generation by
Thompson[28], but will be paraphrased here for continuity. The grid generation technique
utilized was a di erential equation method using an elliptic scheme. A H-grid control volume
was used that placed the mesh in the blade passage while using the blade surfaces as the grid
boundary, as well as the area extending in front of the leading edge and behind the trailing
edge. Hyperbolic Tangent Stretching Functions (HTSF) were used to help cluster the grid
around the leading and trailing edges, while Trans nite Interpolation (TFI) was used to
generate an initial solution to the di erential equation. Blended normalized arc lengths were
also used to help the grid lines not cross the high curve of the leading edge for the initial
solution. The Poisson equation was solved numerically to generate a smooth elliptical grid.
In order to cluster the grid vertically at the leading and trailing edge, the smoothed grid was
split into three blocks and the process was repeated on each block to produce a nal grid as
seen in Figure 3.1. Figures 3.2 and 3.3 show a close-up of the leading edge before and after
the multi-block technique was used. An 80x40 grid was used for evaluations, but a 60x25
grid is shown in Figure 3.1 for clarity.
The Euler code used the AUSM+[29] ux vector splitting method and a Min-Mod
limiter was applied. Characteristic subsonic in ow and out ow boundary conditions were
used on the inlet and exit of the control volume. On the north and south face of the control
31
0 0.2 0.4 0.6 0.8 1 1.2 1.40
0.2
0.4
0.6
0.8
1
1.2
1.4
Blade Velocity Distributions
Surface Arc?length (s/c)
Mach (isentropic)
Euler (for validation)
Experimental Data
Arts Euler
Figure 3.4: Euler code validation results
volume a periodic boundary condition was used prior to the leading edge and after the
trailing edge while a slip-wall boundary condition was applied to the surfaces of the turbine
blade. The CFD code was validated by comparing to experimental data and results from an
Euler CFD code presented in an well known article by Arts[6]. The validation results are
shown in Figure 3.4 for a VK-89 turbine blade at an isentropic exit Mach number of 0:875.
3.2.3 Objective Function
The Euler code provided ow data at the inlet and exit of the control volume, which
was used to compute the net change in momentum of the ow. Because the focus of this
article is to evaluate the e ectiveness of the solution encoding and move operator, an ideal
analysis is su cient for this application. Equation 3.5 gives the force (F) per unit span for
the turbine blade due to change in momentum and is used as the objective function, which is
sought to be maximized. The summation is across all j indices, which span vertically across
the inlet and exit of the control volume. For Eqn. 3.5, nj is the total number of j indices, u
32
and v represent the horizontal and vertical uid velocities respectively, and l represents the
vertical cell height.
F =
njX
j=1
in;j(uin;jlin;j)vin;j
njX
j=1
ex;j(uex;jlex;j)vex;j (3.5)
Because grid generation and CFD codes can sometimes fail to converge properly, a few
safeguards were put in place to ensure ESTurb did not fail altogether. First, the solution
was checked after mutation to ensure the upper and lower surfaces did not cross and the
o spring was feasible. If either the grid generation code or the CFD code failed or did not
converge properly an output le would not be created. A simple timer was started when the
evaluation began, and if the timer expired without an output le present the solution would
be assigned an objective function value of zero and the optimizer would move on to the next
o spring. Convergence failures were fairly uncommon and did not impact the algorithm?s
ability to converge on a solution.
3.2.4 Results and Optimizer Performance
In order to evaluate the ability of the algorithm to optimize turbine blade geometry
and the e ectiveness of the solution encoding, runs were completed using the B ezier curve
encoding as well as the spline-connected polynomial (SCP) encoding shown in Figures 2.4
and 2.5. Each run consisted of a micro-population of 4 members and a = (o spring to
parent) ratio of 4, for a total of 16 o spring per generation. As in the validation cases, 50
generations were evaluated and the standard deviation, , was initially set to 0.03 with the
one- fth rule and \ restart" included. Initial solutions were seeded, the in ow angle was
set to 30o, and the exit isentropic Mach number was set to 0:875. In order to set the throat
distance a desired throat was entered and a bisection method was used to adjust the pitch
until the correct throat was obtained. The throat in turbine design is the minimum distance,
or area for three-dimensions, between blades. The throat is not necessarily choked, but this
is often the case for one or more stages in a turbine near design operating conditions. The
33
trailing edge of the turbine blades is cusped, and a maximum tangential chord of 1 was
allowed for feasible solutions.
Three runs were completed using each solution encoding method and the results were
then compared against each other. The objective function values over the range of genera-
tions are shown in Figure 3.5 for both the B ezier and spline-connected polynomial encodings.
Both encoding options began improving very quickly over the rst few generations, which
was expected based on the performance of the algorithm in the validation cases. As dis-
cussed in section 2.2, B ezier solution encoding is advantageous because it is able to represent
a lot of solutions that SCP and other parameterization method cannot. This can be seen in
the results for the turbine optimization, where both encoding method start improving very
quickly, but the SCP encoding method stalls shortly after. Once the search space, which is
signi cantly smaller for SCP encoding, has been su ciently searched, the optimization will
stagnate because it cannot nd any better solutions. The B ezier solution encoding can rep-
resent almost any shape, and the optimizer continues to nd solutions that the SCP simply
cannot represent. The nal geometries for each method and the pressure contours are shown
in Figures 3.6 and 3.7.
The computation time for 50 generations was approximately 12.2 hours using the spline
encoding and approximately 13.7 hours using the B ezier curve encoding. The slightly longer
time for the B ezier curve encoding was due to the increased number of iterations required to
converge the mesh on certain geometries. Over 96% of the computation time was spent on
mesh generation and evaluating the blade, which means the total computation time depends
mostly on the mesh generation and ow solver algorithm used. This analysis successfully
showed the bene ts of using a B ezier curve solution encoding and validated the objective
function. However, several factors must be addressed to provide compelling results for the
objectives in this paper. Because boundary layers are very in uential in the throat location
and ow in the blade passage, a viscous analysis must be performed. Blunt trailing edges
34
0 10 20 30 40 500.35
0.4
0.45
0.5
0.55
0.6
0.65
0.7
Generation
Obj Function Value
Obj Function vs. Generation
Spline Run 1
Spline Run 2
Spline Run 3
Bezier Run 1
Bezier Run 2
Bezier Run 3
Figure 3.5: Objective function values over 50 generations
?1.5 ?1 ?0.5 0 0.5 1 1.5 2 2.5
0
0.5
1
1.5
2
2.5
3
0.4
0.5
0.6
0.7
0.8
0.9
1
1.1
1.2
Figure 3.6: Final blade geometry using spline-connected polynomial encoding
35
?2 ?1 0 1 2 3
?1
?0.5
0
0.5
1
1.5
2
2.5
3
3.5
4
0.4
0.5
0.6
0.7
0.8
0.9
1
1.1
1.2
1.3
Figure 3.7: Final blade geometry using B ezier curve encoding
36
and simple structural requirements must also be included. These will be discussed in the
next section.
3.3 Viscous Turbine Analysis
The viscous analysis will provide a much more thorough evaluation of the turbine blades,
including boundary layer e ects, ow separation, and frictional losses. This higher delity
analysis also comes with a computational cost, especially when a large number of iterations
are required. Because the application of a Navier-Stokes analysis requires a much greater
amount of expertise and is outside the scope of this paper, two codes generated by NASA were
used. The grid generation code was replaced with GRAPE, a well known grid generation
software also based on solving Poisson?s equation that was written by Reece Sorenson at
NASA Ames Research Center and was modi ed by Dr. Rodrick Chima. This code follows
the same methodology as the code written by the author, but is much more user friendly
and robust. GRAPE allows the use of a \C", \O", or \H" grid topology, but the C-grid will
be used for this application because it works well with the Navier-Stokes code.
The CFD code used is Rotor Viscous Code Quasi-Three Dimensional (RVCQ3D), which
was written by Dr. Rodrick Chima at NASA Glenn Research Center. Both GRAPE and
RVCQ3D are available on the NASA Glenn Research Center Software Repository. RVCQ3D
is a code speci cally developed for turbomachinery that solves the thin-layer Navier-Stokes
equation on a blade-to-blade surface of revolution. It uses central or upwind nite dif-
ferencing techniques and includes three turbulence models and the ability to predict heat
transfer. The AUSM+ upwind di erencing scheme was used for this analysis along with the
low Reynold?s number Wilcox k ! turbulence model. An explicit multistage Runge-Kutta
scheme is used with local time-stepping and the ow variables are non-dimensionalized us-
ing freestream conditions. Multiple articles have been written by Dr. Chima[9][10][11][12]
validating the results of RVCQ3D with experimental data for each of the turbulence models
as well as heat transfer predictions. These articles also show the sensitivity of RVCQ3D
37
to changing variables such as Reynold?s Number, eddy viscosity ratio, freestream turbu-
lence, grid re nement, and a comparison between the turbulence models. Full details on the
methodology of RVCQ3D are provided in the User?s Manual[30] and are elaborated on in
the other articles by Dr. Chima. Care was taken to check and ensure the wall grid spacing
(y+) was approximately less than 1 for each of the nal solutions and the eddy viscosity
ratio was in the range used by Dr. Chima in other published articles[10]. The nal grid for
each solution was also visually checked to ensure the grid spacing was adequate for the nal
blade shape and spacing.
3.3.1 Model and Assumptions
The viscous model accounts for frictional losses through the use of the k ! turbulence
model, and uses an updated non-re ecting boundary condition for inlet and exit as well as
an updated no-slip boundary condition on the turbine blade wall. The analysis is Quasi-3-D,
taking into account the stream tube geometry from the T55 mean radial path, but does not
include any losses speci c to full 3D analysis such as tip leakage. RVCQ3D also allowed
the rotational velocity of the turbine blades to be entered for rotors, while the stator blades
were kept stationary. The numerical parameters recommended in the User?s Manual for the
most robust scheme were used, including a Courant number of 2.5. Example input les for
RVCQ3D are given in Appendix B.
The grid used for the analysis was 169x45, which has been shown by several articles
to be su cient for this type of analysis[31] and veri ed by testing a variety of grids up to
383x55 with no signi cant change in the results. Grid parameters recommended by the
RVCQ3D User?s Manual for heat transfer applications were used, including a grid spacing of
a hundredth of a chord around the trailing edge and four hundredths of a chord around the
leading edge. For a unit axial chord, the grid height at the blade surface was 1:0 10 5. The
remainder of the input parameters are listed in the example GRAPE input le in Appendix
38
B. An example 169x45 grid for a VK-LS89 turbine blade as described in Arts[6] is shown in
Figure 3.8.
Stator/rotor interactions are very di cult to model, and an unsteady solution is neces-
sary to capture the e ects of the pressure uctuations as the rotor blades interact with the
stator wakes. RVCQ3D is a single-block steady state code, which limits the accuracy of the
stator/rotor interaction considerably. The analysis was completed in the direction of ow
using the pressure ratios from the T55 design report, and the energy averaged ow properties
at the stator exit were used for the inlet conditions of the rotor. This analysis is rudimen-
tary at best, but will be su cient in determining the e ectiveness of ESTurb for turbine
optimization. If ESTurb is shown e ective in turbine optimization, a higher delity analysis
using a multi-block code such as SWIFT would allow a more accurate model of a stage or
multi-stage analysis. A multi-block code would evaluate all stages at once by transferring
properties between blade rows at grid interfaces as described in Calculation of Multistage
Turbomachinery Using Steady Characteristic Boundary Conditions by Dr. Chima[32].
3.3.2 Optimizer Updates and Objective Function
In order to incorporate GRAPE and RVCQ3D the evaluation subroutine was updated
again, but the remainder of the ESTurb code remained the same. The evaluation subroutine
consisted of accepting inputs from the main program, writing input les and executing
GRAPE and RVCQ3D, extracting data from the output les, and calculating the objective
function. As in the inviscid case, if either of the executables failed the objective function
would be set to zero and the algorithm would move on to the next o spring.
The objective function was re ned to included the losses calculated by the turbulence
model. The objective function used for stators was as combination of the uid tangential
momentum change and the ratio of P0;2=P0;1, which is also equal to (1 !). Because the
uid mass ow rate varied slightly due to stream tube thickness and other factors, only
the di erence in uid tangential velocity was considered for the momentum change. The
39
X
Y
0 1 2
-3.5
Figure 3.8: Example 169x45 grid for VK-LS89 turbine blade
40
X
Y
2.03 2.04 2.05 2.06 2.07 2.08 2.09
-0.73
-0.72
-0.71
-0.7
-0.69
-0.68
Figure 3.9: Close-up of trailing edge for VK-LS89 turbine blade grid
41
tangential momentum change, which the rotor will later convert into power, is multiplied by
the square of the total pressure ratio to yield the stator objective function equation, shown
below as Equation 3.6.
OBJstator =jC ;ex C ;inj
P
0;ex
P0;in
2
(3.6)
The objective function used for the rotor analysis was a function of the tangential
momentum change and the total-to-total e ciency. The tangential momentum change is
calculated in the same way as the stator objective function, using the di erence in tangential
velocities between inlet and exit. This momentum change in the uid is equal to the force
exerted on the blade, and when combined with the shaft speed and mean blade radius the
work output may be calculated. The total-to-total e ciency discussed in section 3.1 will be
used, which is described in Equation 3.2. The rotor objective function is shown in Equation
3.7, which holds the same form as the stator objective function.
OBJrotor =jC ;ex C ;inj 2t t (3.7)
3.3.3 Structural Considerations
Although an in depth structural analysis is outside the scope of this research, taking
structural considerations into account during the optimization process is important to ensure
feasible results are obtained. For this reason, a simple structural check was completed as
part of the feasibility check shown in the algorithm ow chart. If the o spring does not pass
the structural portion of the feasibility check, the o spring is discarded and a new o spring
is generated. The structural check consists of three parts, a minimum blade thickness check,
blade area check, and surface distance to blade area ratio check. The minimum thickness
check nds the minimum thickness of the blade and ensures it is above a minimum value.
The minimum thickness allowed for this research was equal to the diameter of the trailing
edge, which was equal to that of the T55 blade. The minimum blade area check was a simple
42
way to ensure the blade could withstand the stress of the pressure and inertial forces on the
blade and ensured the optimized blade area was no less than the area of the T55 blade minus
a 10% bu er. The last check measured the ratio of the surface distance to the blade area.
This check was not only for structural purposes, but during the aerothermal optimization
this was a simple way of keeping the ratio from increasing to a point where the blade could
not be cooled properly. The surface to area ratio was not allowed to be more than that of the
T55 blade plus a 20% bu er. Although this check does not replace a full structural analysis,
it does provide a simple estimate of the structural feasibility of the blade and is su cient
for the purposes of this research.
3.3.4 Results and Optimizer Performance
The single stage viscous analysis was completed using the stream tube and thermody-
namic data for the rst stage of the gas generator turbine (GGT) in the T55 engine. The T55
stator and rotor blades were used for the seeded solutions in a micro-population of 4. The
= ratio was set to 2 for a total of 8 o spring per generation and 400 function evaluations
over 50 generations. The analysis required approximately 14 hours to complete per blade
row on a 2.0 GHz dual core processor, for a total of 28 hours for the entire stage. The stator
was optimized rst over 50 generations using the inlet conditions and exit static pressure in
the T55. The energy averaged ow variables at stator exit were used as the inlet conditions
for the rotor analysis along with the rotor exit static pressure and the blade angular velocity
from the T55 design report. The rotor optimization process was also 50 generations. As
discussed earlier, all results compared to the original T55 engine will be given in percentages
to protect the proprietary information.
The optimized geometry and pressure contours for the rst stage are given in Figure
3.10. Blanking variables are used to omit the overlapped region between the blades in
Figure 3.10 even though the grid extended half a chord length axially both before and after
the blade. The Mach contours for both the stator and rotor are shown separately with the
43
full computational grid in Figures 3.11 and 3.12. Note that several computational domains
were stacked in a cascade for viewing purposes, but the control volume analyzed consisted
of only one blade passage. ESTurb performed well at optimizing both the stator and rotor
in reasonable amounts of time. As expected the optimizer increased the ow angle exiting
the stator to impart more momentum into the uid without separating the ow over the
suction surface of the stator. The rotor calculation also performed well at providing a large
momentum change with little or no turbulence. After completing two runs, the average
optimized blade geometry resulted in a 1.1% increase in isentropic e ciency and a 7.0%
increase in power output for the stage when compared to the T55. The results for the
two runs were very similar, di ering by less than 0.1%, so only the geometry for the rst
run is shown. The increased performance is mainly due to a greater change in tangential
momentum by the stator and the rotor, but it is important to note that the optimizer was
able to nd a blade geometry that could signi cantly increase the momentum change while
also increasing the isentropic e ciency by over one percent. Example GRAPE and RVCQ3D
input les are included in Appendix B, but the stream tube geometry for the T55 has been
removed and must be replaced or the \nblade" and \nmn" variables must be updated for it
to run. All remaining data in the input les is non-proprietary.
44
m (in)
Th
eta
*ra
diu
s
1 1.5 2 2.5 3
-1.5
-1
-0.5
0
0.5
1
1.5
Pressure
0.680.65
0.620.59
0.550.52
0.490.46
0.420.39
Figure 3.10: Blade geometry and pressure contour results for single stage
45
X
Y
1 1.5 2 2.5 3-1.5
-1
-0.5
0
0.5
Mach Number
0.80.75
0.70.65
0.60.55
0.50.45
0.40.35
0.30.25
0.20.15
0.10.05
Figure 3.11: Mach Number contours for stator
46
X
Y
2 2.5 3 3.5
-1.5
-1
-0.5
0 Mach Number
0.70.65
0.60.55
0.50.45
0.40.35
0.30.25
0.20.15
0.10.05
Figure 3.12: Absolute Mach Number contours for rotor
47
Chapter 4
Turbine Aerothermal Optimization
An aerothermal analysis is completed in this report in order to investigate the trade-
o between aerodynamic and thermal performance during optimization. This analysis is
completed for two primary reasons, the rst being the lack of research in this area using
an adaptive optimizer. Of all the articles reviewed by the author, none could be found
that used an adaptive optimizer and evaluated both aerodynamic performance and heat
transfer. The second reason is the potential for B ezier curve blade representation to provide
optimal leading edge geometries for this problem. Of the 7 turbine design books and many
more journal articles reviewed by the author, including those listed in the literature review,
every one suggests de ning the leading edge geometry of a blade using a circle or an ellipse.
Although this provides a great aerodynamic advantage and therefore is commonly used,
intuition would suggest that having a more exible parameterization method might lead to
better designs, especially when multiple objectives are considered. This greater exibility in
the geometry may lead to a more thermally resistant blade which would allow higher turbine
inlet temperatures. Because turbine inlet temperature is directly related to maximum power
output, this could be advantageous. The analysis conducted in this paper is assuming the
turbine blade is not cooled for simplicity purposes.
4.1 Heat Transfer Model
The heat transfer model used for this analysis is incorporated into the RVCQ3D analysis
code. The code uses the k ! turbulence model with transition e ects to predict the Stanton
Number at locations around the blade. As seen in Eqn. 4.1, a slightly adapted equation for
the Stanton Number de nition from the RVCQ3D User?s Manual[30], the Stanton Number
48
may be used to calculate the heat transfer coe cient around the blade. A thorough validation
of this heat transfer model has been conducted by Dr. Chima in comparison to experimental
data[10], and the author replicated the results within reasonable accuracy. The details of
the heat transfer model are outside the scope of this paper, but are described in the article
by Dr. Chima[10]. It is important to note that turbulence and heat transfer in turbine
blades is an extremely complex subject that is the focus of much current research, and an
accurate and robust heat transfer model that works in a variety of situations is still sought
after. Every heat transfer model has limitations and is an estimate at best, therefore the
research being conducted in this paper is primarily to identify the potential for conducting
further experimental and high delity research in this area, rather than providing a speci c
optimized blade shape.
St = k
T
njwall
inV0inCp(T0;in Twall) =
h
inV0inCp (4.1)
4.2 Adaptations to Optimizer
As with the previous sections, the main portion of the ESTurb algorithm remains the
same while the objective function is updated to match the new optimization criteria. Because
the evaluation code remains the same, the aerodynamic portion of the objective function will
not change. The previous objective function from section 3.3.2 will simply be multiplied by a
heat transfer objective function to produce the overall tness function for the optimization.
Because the main goal of the aerothermal optimization is to investigate the trade-o between
thermal and aerodynamic performance, the heat transfer objective function relate the heat
transfer coe cient on the blade to the ideal engine power output
4.2.1 Objective Function
The equation of convective heat transfer is given below in Eqn. 4.2. Assuming the blade
is at the maximum allowable blade temperature and the turbine inlet temperature as well as
49
the heat transfer coe cient are known, the heat transfer into the blade may be calculated.
For evaluating an optimized blade it can be assumed that the heat transfer into the blade
and the blade temperature is held the same, and therefore a change in the heat transfer
coe cient will result in a change in allowable turbine inlet temperature. The new turbine
inlet temperature may then be moved to the left hand side of the equation as shown in Eqn.
4.3
q = h(T0 Twall) (4.2)
T0;new = hnewh
old
(T0;old Twall;old) +Twall;old (4.3)
The ideal engine analysis discussed in Mattingly[1] provides a relation between a turbine
inlet temperature increase and power increase. Equation 4.4 is a simpli ed version of the
expression in Mattingly for shaft power from an ideal turboshaft engine. Plotting this equa-
tion over a range of compressor ratios and turbine inlet temperatures results in Figure 4.1.
If a compressor ratio is selected a linear relationship between the turbine inlet temperature
ratio and shaft power ratio may be obtained through a simple curve t. The relationship
used for this analysis is described in Eqn. 4.5 and shown in Figure 4.2. After running the
evaluation on a blade and the peak heat transfer coe cient is obtained, Eqn. 4.3 is used to
calculate the maximum turbine inlet temperature and Eqn. 4.5 is used to obtain a ratio for
_Wnew= _Wold. This ratio is then multiplied with the aerodynamic objective function to obtain
an overall objective function.
_Wshaft = _m0CpT0
"
T0;turbine
T0;1
1 1
1
c
!
1
c 1
#
(4.4)
_Wnew
_Wold = 1:5935
T0;new
T0;old 0:5935 (4.5)
50
2 4 6 8 10 12 14 16 18 2060
80
100
120
140
160
180
200
220 Ideal Engine Power for Varying Turbine Inlet Temperatures
Compressor Pressure Ratio
Shaft Power to Air Mass Flow Rate Ratio (BTU/(lbm/s))
2800K2700K
2600K2500K
2400K2300K
2200K
Figure 4.1: Ideal engine power vs. pres-
sure ratio
0.85 0.9 0.95 1 1.05 1.10.75
0.8
0.85
0.9
0.95
1
1.05
1.1
1.15 Ideal Engine Power at Specified Pressure Ratio
Normalized Turbine Inlet Temperature (T0,new/T0,old)
Normalized Shaft Power to Air Mass Flow Rate Ratio (Power
new/Power
old)
Figure 4.2: Relationship for normalized
power ratio vs turbine inlet temperature
4.3 Preliminary Aerothermal Optimization Results
The preliminary aerothermal results are shown in Table 4.1. Three separate runs were
completed on the stator using the data from the T55 engine for the thermodynamic variables
and stream tube geometry. The blade geometry from the T55 was also used for an initial
seed for the rst population. As expected, the algorithm produced a blade with a more
blunt leading edge to lessen the e ects of heat transfer in that region. Of the three runs, the
results were fairly similar, as shown in Figures 4.3-4.5.
Table 4.1: Preliminary Aerothermal Results Summary
Run # OBJ P0;ex=P0;in T0;new=T0;old _Wnew= _Wold
1 107.8 65:4o 0.9872 1.128 1.204
2 107.3 64:2o 0.9867 1.129 1.206
3 105.3 63:9o 0.9793 1.123 1.196
The adjusted shape of the leading edge had a signi cant impact on the maximum heat
transfer coe cient on the blade. Figure 4.6 shows the heat transfer coe cient across the
blade surface for each of the optimizer runs as well as the numerical results for the T55
stator geometry. The optimized blade shape predicts a signi cant allowable rise in turbine
51
m
Theta
*r
1 1.2 1.4 1.6 1.8 2 2.2
-1
-0.5
0
Pressure0.68
0.640.60.56
0.520.480.44
0.4
Figure 4.3: Aerothermal run
1 pressure contours
m
Theta
*r
1 1.2 1.4 1.6 1.8 2 2.2
-1
-0.5
0
Pressure0.68
0.640.60.56
0.520.480.44
0.4
Figure 4.4: Aerothermal run
2 pressure contours
m
Theta
*r
1 1.2 1.4 1.6 1.8 2 2.2
-1
-0.5
0
Pressure0.68
0.640.60.56
0.520.480.44
0.4
Figure 4.5: Aerothermal run
3 pressure contours
?1 ?0.5 0 0.5 1
Heat Transfer Coefficient for ES Blade
Arc Length Distance from Leading Edge (Pressure/Suction)
Heat Transfer Coefficient "H" (BTU/lbm.hr.ft
2 )
Run 1
Run 2Run 3
T55
Pressure Side Suction Side
Figure 4.6: Preliminary aerothermal results: Heat transfer coe cient
52
inlet temperature of over 12%, which corresponds to an almost 20% power increase in an
ideal engine, but this blade shape also introduces some problems.
4.3.1 Issues With Preliminary Results
The CFD results are promising in a thermal respect, but fail to capture some important
aerodynamic aspects. Flow separation and transition e ects are very hard to model and are
often sensitive to input parameters such as the eddy viscosity ratio. One very important
aspect in detailed turbine blade design is ensuring a monotonic increase in velocity along
the suction surface. \Apart from mechanical considerations, it is the primary objective, in
designing a turbine cascade, to produce a monotonically accelerating ow along the suction
side, in particular. Such a velocity distribution will give rise to a continually decreasing
pressure that substantially slows the boundary layer buildup and minimizes the possibility
of ow separation."[3] Figure 4.7 shows the isentropic Mach number distribution across the
blade surface for Run 1, which clearly shows a large decrease in velocity directly downstream
of the sharp leading edge curve on the suction side. A much smaller \bump" in the velocity
distribution of a turboprop blade, as shown in Figure 4.8, separated ow without reattach-
ment and decreased e ciency by over three points. Changing the blade geometry to produce
a \continually accelerating nozzle-like ow eld [provides] an environment where boundary
layer buildup is signi cantly suppressed, and the likelihood of boundary layer separation is
practically nonexistent, at least under design point operation."[3] Maintaining a monotoni-
cally increasing velocity along the pressure side of the blade is not as critical, since there is
always a strong monotonic increase in velocity as the channel approaches the throat which
will reattach the ow.
Another important aspect of blade design not considered in the preliminary analysis is
the uncovered turning angle. This angle is often incorporated into the blade parameterization
method, but this is not the case when using B ezier curves. The uncovered turning angle is
the angle between the tangents to the suction surface at the throat and the trailing edge
53
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Adverse Pressure Gradient"Bump" in Velocity Distribution
Figure 4.7: Aerothermal run 1 Mach num-
ber distribution
Figure 4.8: \Bump" in velocity distribu-
tion for turboprop blade[3]
merging point. It is referred to as \uncovered" because the region of the blade suction
side after the throat is not con ned within the blade passage. A small uncovered turning
angle facilitates accelerating ow over the suction side after the throat, which increases
performance. However, too large of an angle may have detrimental e ects.
4.3.2 Updates to Objective Function
Updating the objective function to account for these two issues is easily implemented
in the form of a penalty function. After evaluating each blade in RVCQ3D, the velocity
distributions for each surface are checked. If the suction surface decreases velocity prior to
the throat location the gradient of the velocity decrease is used to produce a penalty factor
while is less than one. The penalty factor is then multiplied with the objective function
value to decrease it proportional to the severity of the di usion on the blade. Because the
pressure side is more forgiving, a similar process is used including a bu er. If the velocity
decreases more than 15% of the exit Mach number along the suction surface, the penalty is
applied to the objective function in the same manner.
A limit is applied to the suction side di usion as suggested by Baskharone[3] to manage
the uncovered turning angle. The suction side di usion is calculated by Eqn. 4.6 for a
stator and 4.7 for a rotor and it is unpenalized under 0.25. If the di usion exceeds 0.25, the
54
objective function is penalized based on the value of the angle over the limit in a similar
manner to the pressure side velocity.
stator = (
P
P0 )ex (
P
P0 )maxV
1 ( PP0 )maxV (4.6)
rotor = (
P
P0 )ex (
P
P0 )maxW
1 ( PP0 )maxW (4.7)
4.3.3 Stage Aerothermal Results
Three additional optimization runs were completed for an entire stage using the updated
objective function incorporating the penalty factors. The results for these three runs are
summarized in Table 4.2. As in the preliminary runs, the T55 geometry for the annulus and
blades were used for initial conditions, as well as the T55 thermodynamic variables. The
momentum portion of the objective function sought to match the T55 stage work output
while the e ciency and allowable inlet T0 were maximized. This is to highlight the e ects of
the thermal improvements while maintaining the same velocity triangle angles as the T55.
Note that the t t and t s are arti cially high when compared to a typical gas turbine engine
due to the 2D analysis which does not include tip losses or end wall losses. The referenced
T55 e ciencies are those produced by the CFD code using the T55 blade geometry for
comparison purposes, not the actual e ciencies of the engine stage including all losses.
Because RVCQ3D is a single block code, the energy-averaged ow variables calculated for
the stator exit were used as the inlet conditions for the rotor. A multi-block code such
as SWIFT would provide a more accurate model of the entire stage at once, but with the
pressure ratios already de ned by the T55 data this approach will be su cient for the
objectives of this analysis.
The results for each of the runs were fairly similar, so the results for Run 1 will be
presented in detail in the main text. Additional plots for the other runs are included in
Appendix C. The stage geometry and pressure contours for Run 1 are shown in Figure 4.13,
55
Table 4.2: Stage Aerothermal Results Summary
Run # OBJstator OBJrotor stator;ex rotor;ex t t t s T0;new=T0;old _Wnew= _Wold
1 102.9 101.6 64:2o 12:3o 94:6% 83:8% 1.063 1.101
2 102.9 100.8 63:9o 10:8o 93:9% 82:1% 1.050 1.079
3 102.4 100.9 64:0o 12:9o 94:3% 83:7% 1.055 1.088
T55 - - - - 93:7% 80:5% 1 1
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9 Run 1 Stator Isentropic Mach Number
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure 4.9: Stage aerothermal Run 1 sta-
tor Mach number distribution
0 0.2 0.4 0.6 0.8 1 1.2 1.40
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8 Run 1 Rotor Isentropic Mach Number
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure 4.10: Stage aerothermal Run 1 ro-
tor Mach number distribution
and an example of the input les for GRAPE and RVCQ3D are shown in Appendix B. As
stated earlier, the stream tube geometry for the T55 has been removed and must be replaced
or the \nblade" and \nmn" variables must be updated in order to run the input les. The
di erences in the nal solutions for the three runs are shown in Figures 4.11 and 4.12, as
well as the radius of curvature used for the T55 stator blade leading edges for comparison.
When compared against the preliminary aerothermal results, the updated objective function
drastically reduced the ow acceleration around the corners of the leading edge, and Figures
4.9-4.10 show a monotonic ow acceleration along the suction side from the leading edge
to the throat. This updated blade geometry will be very unlikely to have ow separation
on the suction side during design operating conditions. The decreasing velocity along the
pressure surface is reasonable and less than the velocity decrease for the original T55 blade.
The suction side di usion was also kept under the limit of 0.25 for each blade.
56
0 0.2 0.4 0.6 0.8 1?0.9
?0.8
?0.7
?0.6
?0.5
?0.4
?0.3
?0.2
?0.1
0
0.1
Optimized Stator Geometry
m
theta*r
Run 1
Run 2Run 3
Approximate T55 leading edge radius
Figure 4.11: Comparison of nal stator
geometries
0 0.2 0.4 0.6 0.8 1
?0.6
?0.5
?0.4
?0.3
?0.2
?0.1
0
0.1
0.2
Optimized Rotor Geometry
m
theta*r
Run 1
Run 2Run 3
Figure 4.12: Comparison of nal rotor ge-
ometries
Figure 4.11 shows the di erence in the leading edge between the optimized blades and
the original radius of curvature. The optimized blades all developed a small at area on
the leading edge which reduced the maximum heat transfer coe cient. The distribution
of the heat transfer coe cient across the surface of the blade and the comparison with
the T55 blade is shown in Figure 4.15. Traditionally, thicker leading edges performed well
with changes in incidence angles but produced larger kinetic energy losses, while thinner
leading edges sometimes developed signi cant issues with positive incidence angles. Usually
a turbine blade is designed to provide a compromise between the two based on aerodynamic
performance. However, if the optimized turbine blades provide an increase in turbine inlet
temperature that outweighs the kinetic energy loss from the thicker leading edge, the overall
performance would be better. It is important to note that of all the blade parameterization
methods reviewed by the author in 7 turbine design books and numerous articles, every
one used either a radius of curvature or ellipse to de ne the leading edge. The three typical
types of leading edges are shown in Figure 4.14 from Baskharone[3]. Therefore, the optimized
blade produced by B ezier curve encoding through ESTurb could not be represented by any
of those parameterization methods. Using B ezier curves to encode the solution did not only
prove bene cial in nding better solutions, but in this case it was required to produce the
nal results presented in this paper.
57
m (in)
Th
eta
*ra
diu
s
1 1.5 2 2.5 3
-1.5
-1
-0.5
0
0.5
1
1.5
Pressure
0.680.65
0.620.59
0.550.52
0.490.46
0.420.39
Figure 4.13: Stage pressure contours for Run 1
58
Figure 4.14: Types of leading edges and incidence angle e ects[3]
Although a very large increase in performance is predicted for design conditions, very
little of an engines life is spent at or near design conditions. Therefore a brief study of
o design conditions for the optimized stator from Run 1 is conducted in Appendix D by
altering the in ow conditions. A brief grid re nement study is also conducted in Appendix
E. A variety of numerical studies could be performed on each of the blades over a large
range of conditions, but due to the large amount of computation time per evaluation only a
brief o -design analysis was conducted to provide enough compelling information to warrant
future experimental research in this area.
The computation time to conduct an optimization using ESTurb was very reasonable.
Because over 96% of the computation time was spent during objective function evaluation,
the number of objective function calls was extremely important. As shown in the validation
cases, a much smaller number of function evaluations was required by ESTurb. A popu-
lation size of 4 and a ratio of 2 was used for the aerothermal optimization, which
resulted in 16 function evaluations per generation and only 800 function evaluations over 50
generations. The convergence history is shown in Figure 4.16 for both the stator and rotor,
although the objective functions for the stator and rotor di ered. The optimization required
59
?1 ?0.5 0 0.5 1Arc Length Distance from Leading Edge
Heat Transfer Coefficient "H"
Run 1
Run 2Run 3
T55
Pressure Side Suction Side
Figure 4.15: Heat transfer coe cient for optimized stator
approximately 14 hours to complete on a 2.2 GHz dual core processor and only 5.5 hours to
complete on a 3.4 GHz 8-Core processor per blade row. This resulted in a total of approx-
imately 28 hours or 9 hours of computation time respectively for the entire stage. This is
extremely reasonable given a Navier-Stokes CFD code was used for the function evaluation,
and would allow much higher delity codes to be used without making the computation time
excessively high. Parallel processing was not implemented for this analysis, but could easily
be incorporated into the algorithm in the future.
60
0 5 10 15 20 25 30 35 40 45 5060
65
70
75
80
85
90
95
100
105
Generation
Objective Function Value
Stator Run 1Stator Run 2
Rotor Run 1Rotor Run 2
Stator Run 3Rotor Run 3
Figure 4.16: Convergence history for both stator and rotor
61
Chapter 5
Conclusions and Future Applications
The three objectives for this thesis were to produce and validate an optimizer, conduct
a turbine aerodynamic optimization, and conduct a turbine aerothermal optimization. The
ESTurb optimizer produced performed remarkably well during the validation cases on airfoil
geometries. The ES optimization algorithm coupled with B ezier curve solution encoding was
well suited for aerodynamic problems and performed a very e cient and e ective search to
produce optimal or near optimal solutions. The inviscid turbine aerodynamic optimization
successfully created an objective function to evaluate turbine blade geometry, and validated
the use of B ezier curves in representing turbine blades. The viscous analysis showed the
potential of nding optimal turbine blades that increased e ciency by over 1% and increased
shaft power by 7% over the original T55 blade geometry.
The aerothermal analysis provided the most compelling results of the paper, showing
that using B ezier curves to represent the blade geometry may produce turbine blades that
have better heat transfer characteristics at the leading edge. This improved geometry may
allow up to a 6% increase in turbine inlet temperature, which would correlate to over a 10%
increase in power for an ideal engine. The algorithm also performed an e cient search that
only required 512 hours to complete on a 3.4 GHz desktop computer. Based on these ndings,
the following major conclusions may be made:
1. Using an Evolution Strategies algorithm coupled with a B ezier curve encoding is a very
e cient and e ective technique for aerodynamic optimization
2. The reduced number of function calls required for an ES algorithm may allow higher
delity evaluation methods, such as 3D Navier-Stokes CFD codes, to be used for opti-
mization without unreasonable computation times.
62
3. The aerothermal optimization results show compelling evidence that further research
is needed in this area. Speci cally, conducting experimental research to evaluate if the
optimized blades produced in this paper show improved performance over the existing
turbine blades.
Future work in this area should include application of ESTurb or similar algorithms
to other problems in order to determine e ectiveness and exibility across a wide range of
aerospace or other engineering disciplines. Continued research using ESTurb to optimize
turbine blades with higher delity 3D models will be pursued by the author, speci cally
using SWIFT for a 3D ow solver. Continued research into the validity of the optimized
aerothermal blades by those experienced in turbine design would greatly bene t this area
of research, especially experimental tests to con rm or deny the predictions made by CFD.
Research into the e ects of blade geometry on lm cooled turbine blades would also prove
very bene cial.
63
Bibliography
[1] Mattingly, Jack D. Elements of Gas Turbine Propulsion. New York: McGraw-Hill, 1996.
Print.
[2] Aungier, Ronald H. Turbine Aerodynamics: Axial- ow and Radial-in ow Turbine De-
sign and Analysis. New York: ASME, 2006. Print.
[3] Baskharone, Erian A. Principles of Turbomachinery in Air-breathing Engines. Cam-
bridge: Cambridge UP, 2006. Print.
[4] Wilson, David Gordon., and Theodosios Korakianitis. The Design of High-e ciency
Turbomachinery and Gas Turbine. Upper Saddle River, NJ: Prentice Hall, 1998. Print.
[5] Lakshminarayana, B. Fluid Dynamics and Heat Transfer of Turbomachinery. New York:
Wiley, 1996. Print.
[6] Arts, T., and M. Lambert De Rouvroit. "Aero-Thermal Performance of a Two-
Dimensional Highly Loaded Transonic Turbine Nozzle Guide Vane: A Test Case for
Inviscid and Viscous Flow Computations." Journal of Turbomachinery 114.1 (1992):
147. Print.
[7] Horlock, J. H., and J. D. Denton. "A Review of Some Early Design Practice Using
Computational Fluid Dynamics and a Current Perspective." Journal of Turbomachinery
127.1 (2005): 5. Print.
[8] Denton, J. D., and W. N. Dawes. "Computational uid dynamics for turbomachinery
design." Proceedings of the Institution of Mechanical Engineers, Part C: Journal of
Mechanical Engineering Science 213.2 (1998): 107-124.
[9] Chima, Rodrick V. Analysis of Inviscid and Viscous Flows in Cascades with an Explicit
Multiple-grid Algorithm. Washington, D.C.: National Aeronautics and Space Adminis-
tration, 1984. Print.
[10] Chima, Rodrick V. A K-! Turbulence Model for Quasi-three-dimensional Turbomachin-
ery Flows. Washington, D.C.: National Aeronautics and Space Administration, 1995.
Print.
[11] Chima, Rodrick V., P. W. Giel, and Robert J. Boyle. An Algebraic Turbulence Model
for Three-dimensional Viscous Flows. Washington, DC: National Aeronautics and Space
Administration, 1993. Print.
64
[12] Chima, Rodrick V. "Application of the K-omega Turbulence Model to Quasi-three-
dimensional Turbomachinary Flows." Journal of Propulsion and Power 12.6 (1996):
1176-179. Print.
[13] Goel, S., Cofer, J.I., & Singh H. (June 10-13 1996). Turbine Airfoil Design Optimization,
International Gas Turbine and Aeroengine Congress and Exposition, Birmingham, UK.
[14] Giannakoglou, C.K. \A Design Method for Turbine Blades Using Genetic Algorithms
on Parallel Computers", ECCOMAS 98 John Wiley & Sons, 1998.
[15] Giannakoglou, K. C. "Designing turbomachinery blades using evolutionary methods."
ASME Paper (1999).
[16] Pierret, Stphane, and R. A. Van den Braembussche. "Turbomachinery blade design
using a NavierStokes solver and arti cial neural network." Journal of Turbomachinery
121.2 (1999): 326-332.
[17] Trigg, M. A., G. R. Tubby, and A. G. Sheard. "Automatic genetic optimization approach
to two-dimensional blade pro le design for steam turbines." Journal of turbomachinery
121.1 (1999): 11-17.
[18] Dennis, Brian H., George S. Dulikravich, and Zhen-Xue Han. "Constrained shape opti-
mization of airfoil cascades using a Navier-Stokes solver and a genetic/SQP algorithm."
ASME International Gas Turbine and Aeroengine Congress and Exhibition, Indianapo-
lis, IN. Vol. 7. No. 10. 1999.
[19] Chen, Bo, and Xin Yuan. "Advanced aerodynamic optimization system for turboma-
chinery." Journal of turbomachinery 130.2 (2008): 021005.
[20] ksz, zhan, and Ibrahim Sinan Akmandor. "Multi-Objective Aerodynamic Optimization
of Axial Turbine Blades Using a Novel Multilevel Genetic Algorithm." Journal of Tur-
bomachinery 132.4 (2010): 041009.
[21] Aungier, Ronald H. Axial-Flow Compressors: A Strategy for Aerodynamic Design and
Analysis. New York: ASME, 2003. Print.
[22] Dreo, Johann. Metaheuristics for Hard Optimization: Methods and Case Studies. Berlin:
Springer, 2006. Print.
[23] Giannakoglou, K. C. "Designing turbomachinery blades using evolutionary methods."
ASME Paper (1999).
[24] Wickramasinghe, U., Carrese, R., Li, X. Designing Airfoils using a Reference Point based
Evolutionary Many-objective Particle Swarm Optimization Algorithm. Paper presented
at the IEEE World Congress on Computational Intelligence, Barcelona Spain, July
18-23, 2010.
[25] DAngelo, S., Minisci, E. Multi-objective evolutionary optimization of subsonic airfoils
by meta-modeling and evolution control. Proceedings of the Institution of Mechanical
Engineers, Part G: Journal of Aerospace Engineering Vol. 221, 2007, pp. 805-814.
65
[26] Curriston, Drew A., and Alice E. Smith. "Airfoil Optimization by Evolution Strategies."
Evolutionary Computation (CEC), 2013 IEEE Congress on. IEEE, 2013.
[27] Curriston, Drew A., and Roy J. Hart eld. "Power Turbine Blade Aerodynamic Op-
timization Using Non-Restrictive Evolution Strategies." Joint Propulsion Conference,
2013.
[28] Thompson, Joe F., B. K. Soni, and N. P. Weatherill. Handbook of Grid Generation.
Boca Raton, FL: CRC, 1999. Print.
[29] Liou, Meng-Sing. "A Sequel to AUSM: AUSM+." Journal of computational Physics
129.2 (1996): 364-382.
[30] Chima, Rodrick V. "RVCQ3D Rotor Viscous Code Quasi-3-D User?s Manual and Doc-
umentation." (1999).
[31] Boyle, R. J. "NavierStokes Analysis of Turbine Blade Heat Transfer." Journal of Tur-
bomachinery 113.3 (1991): 392. Print.
[32] Chima, Rodrick V. Calculation of multistage turbomachinery using steady characteris-
tic boundary conditions. Vol. 206613. National Aeronautics and Space Administration
Langley Research Center, 1998.
[33] Boyce, Meherwan P., Gas Turbine Engineering Handbook. Houston: Gulf Publishing
Co., 1982.
[34] Cohen, Henry., Rogers, G.F.C., Saravanamuttoo, H.I.H. Gas Turbine Theory 3rd edi-
tion. New York: Wiley & Sons, 1987.
[35] Fielding, L. Turbine Design, The E ect on Axial Turbine Performance of Parameter
Variation. New York: ASME Press, 2000.
[36] Pritchard, L. J. An Eleven Parameter Axial Turbine Airfoil Geometry Model. ASME
Paper No. 85-GT-219.
[37] Ye, Z.-Q. A Systematic Computational Design System for Turbine Cascades, Airfoil
Geometry and Blade-to-Blade Analysis. Trans. ASME, Journal of Engineering for Gas
Turbines and Power (July 1984): 598-605.
66
Appendices
67
Appendix A
Pseudo Code for ESTurb
Randomly select or seed initial population of size
Evaluate tness of initial population
Perform generation loop until stopping criteria is met f
for i = 1 to f
randomly select parent based on uniform probability
mutate to generate = o spring
evaluate tness of each o spring
g end for
sort o spring and parents based on tness
best solutions from + replace population
implement one fth rule
g End generation loop
Save and display results
68
Appendix B
Example GRAPE and RVCQ3D input les
Example GRAPE input le
&grid1 jmax=169 kmax=45 ntetyp=3 nairf=5 nibdst=7 nobshp=7 jairf=229 jtebot=41 jtetop=129 norda=0 3 maxita=0 300
nout=4 xle=2.335500 xte=3.058800 xleft=1.973850 xright=3.420450 rcorn=0 dsi=1.000e-05 &end
&grid2 nobcas=0 nle=35 nte=24 dsra=0.500000 dsle=0.006000 dste=0.001000 pitch=0.608508 yscl=1.000000 xtfrac=1.000000
dsobi=0.013522 aaai=0.550000 bbbi=0.550000 ccci=0.700000 dddi=0.700000 csmoo=0 jcap=18 dswex=0.008000 jwakex=1 kwakex=1
joble=9 fswake=0.600000 &end
&grid3 airfx=
7.20223359E-01 7.19743354E-01 7.19250123E-01 7.18748279E-01 7.18242519E-01 7.17737573E-01 7.17238166E-01 7.16748968E-01
7.16274558E-01 7.15819372E-01 7.15387670E-01 7.14983490E-01 7.14610614E-01 7.14272529E-01 7.13972398E-01 7.13713029E-01
7.11596086E-01 7.09150004E-01 7.06377036E-01 7.03279756E-01 6.99861057E-01 6.96124144E-01 6.92072523E-01 6.87710001E-01
6.83040675E-01 6.78068925E-01 6.72799408E-01 6.67237054E-01 6.61387053E-01 6.55254857E-01 6.48846164E-01 6.42166917E-01
6.35223298E-01 6.28021715E-01 6.20568804E-01 6.12871414E-01 6.04936605E-01 5.96771641E-01 5.88383982E-01 5.79781278E-01
5.70971361E-01 5.61962241E-01 5.52762094E-01 5.43379263E-01 5.33822245E-01 5.24099685E-01 5.14220373E-01 5.04193234E-01
4.94027321E-01 4.83731811E-01 4.73315994E-01 4.62789272E-01 4.52161147E-01 4.41441217E-01 4.30639169E-01 4.19764771E-01
4.08827866E-01 3.97838368E-01 3.86806249E-01 3.75741539E-01 3.64654313E-01 3.53554691E-01 3.42452826E-01 3.31358897E-01
3.20283109E-01 3.09235677E-01 2.98226826E-01 2.87266781E-01 2.76365762E-01 2.65533977E-01 2.54781614E-01 2.44118835E-01
2.33555770E-01 2.23102509E-01 2.12769095E-01 2.02565520E-01 1.92501716E-01 1.82587547E-01 1.72832806E-01 1.63247205E-01
1.53840369E-01 1.44621832E-01 1.35601024E-01 1.26787272E-01 1.18189787E-01 1.09817662E-01 1.01679860E-01 9.37852135E-02
8.61424117E-02 7.87599983E-02 7.16463624E-02 6.48097326E-02 5.82581699E-02 5.19995613E-02 4.60416127E-02 4.03918425E-02
3.50575746E-02 3.00459318E-02 2.53638289E-02 2.10179664E-02 1.70148231E-02 1.33606499E-02 1.00614629E-02 7.12303646E-03
4.55089684E-03 2.35031512E-03 5.26300639E-04 -9.16405794E-04 -1.97333863E-03 -2.64031430E-03 -2.91343801E-03 -2.78911040E-03
-2.26403437E-03 -1.33522172E-03 0.00000000E+00 1.37412562E-03 3.57375250E-03 6.54814423E-03 1.02483270E-02 1.46270603E-02
1.96388076E-02 2.52397073E-02 3.13875429E-02 3.80417145E-02 4.51632088E-02 5.27145702E-02 6.06598716E-02 6.89646847E-02
7.75960512E-02 8.65224532E-02 9.57137840E-02 1.05141319E-01 1.14777686E-01 1.24596836E-01 1.34574015E-01 1.44685733E-01
1.54909736E-01 1.65224975E-01 1.75611580E-01 1.86050828E-01 1.96525113E-01 2.07017920E-01 2.17513793E-01 2.27998307E-01
2.38458037E-01 2.48880532E-01 2.59254282E-01 2.69568692E-01 2.79814049E-01 2.89981498E-01 3.00063006E-01 3.10051339E-01
3.19940030E-01 3.29723348E-01 3.39396272E-01 3.48954460E-01 3.58394218E-01 3.67712477E-01 3.76906754E-01 3.85975132E-01
3.94916225E-01 4.03729152E-01 4.12413505E-01 4.20969321E-01 4.29397052E-01 4.37697540E-01 4.45871979E-01 4.53921894E-01
4.61849109E-01 4.69655715E-01 4.77344045E-01 4.84916641E-01 4.92376228E-01 4.99725682E-01 5.06968004E-01 5.14106285E-01
5.21143683E-01 5.28083391E-01 5.34928606E-01 5.41682504E-01 5.48348206E-01 5.54928752E-01 5.61427070E-01 5.67845948E-01
5.74188003E-01 5.80455655E-01 5.86651094E-01 5.92776251E-01 5.98832772E-01 6.04821987E-01 6.10744877E-01 6.16602053E-01
6.22393719E-01 6.28119644E-01 6.33779138E-01 6.39371016E-01 6.44893574E-01 6.50344555E-01 6.55721125E-01 6.61019838E-01
6.66236611E-01 6.71366694E-01 6.76404639E-01 6.81344271E-01 6.86178662E-01 6.90900097E-01 6.95500047E-01 6.99969141E-01
7.04297134E-01 7.08472880E-01 7.12484301E-01 7.16318360E-01 7.19961027E-01 7.23397257E-01 7.23499088E-01 7.23552557E-01
7.23557164E-01 7.23512865E-01 7.23420076E-01 7.23279664E-01 7.23092943E-01 7.22861660E-01 7.22587978E-01 7.22274458E-01
7.21924032E-01 7.21539980E-01 7.21125893E-01 7.20685647E-01 7.20223359E-01
airfy=
-4.05503505E-01 -4.05663650E-01 -4.05776673E-01 -4.05841518E-01 -4.05857576E-01 -4.05824698E-01 -4.05743191E-01 -4.05613819E-01
-4.05437790E-01 -4.05216753E-01 -4.04952775E-01 -4.04648325E-01 -4.04306251E-01 -4.03929754E-01 -4.03522357E-01 -4.03087869E-01
-3.99086672E-01 -3.94436305E-01 -3.89183038E-01 -3.83371881E-01 -3.77046599E-01 -3.70249713E-01 -3.63022513E-01 -3.55405064E-01
-3.47436215E-01 -3.39153605E-01 -3.30593676E-01 -3.21791674E-01 -3.12781666E-01 -3.03596539E-01 -2.94268016E-01 -2.84826659E-01
-2.75301878E-01 -2.65721943E-01 -2.56113986E-01 -2.46504015E-01 -2.36916916E-01 -2.27376470E-01 -2.17905350E-01 -2.08525139E-01
69
-1.99256332E-01 -1.90118348E-01 -1.81129536E-01 -1.72307182E-01 -1.63667521E-01 -1.55225742E-01 -1.46995995E-01 -1.38991406E-01
-1.31224076E-01 -1.23705096E-01 -1.16444551E-01 -1.09451531E-01 -1.02734139E-01 -9.62994957E-02 -9.01537522E-02 -8.43020952E-02
-7.87487562E-02 -7.34970196E-02 -6.85492308E-02 -6.39068043E-02 -5.95702319E-02 -5.55390907E-02 -5.18120517E-02 -4.83868876E-02
-4.52604807E-02 -4.24288319E-02 -3.98870679E-02 -3.76294501E-02 -3.56493824E-02 -3.39394191E-02 -3.24912739E-02 -3.12958272E-02
-3.03431345E-02 -2.96224350E-02 -2.91221592E-02 -2.88299372E-02 -2.87326070E-02 -2.88162226E-02 -2.90660621E-02 -2.94666359E-02
-3.00016949E-02 -3.06542385E-02 -3.14065230E-02 -3.22400696E-02 -3.31356725E-02 -3.40734072E-02 -3.50326385E-02 -3.59920290E-02
-3.69295468E-02 -3.78224739E-02 -3.86474144E-02 -3.93803026E-02 -3.99964111E-02 -4.04703589E-02 -4.07761198E-02 -4.08870305E-02
-4.07757984E-02 -4.04145104E-02 -3.97746404E-02 -3.88270577E-02 -3.75420356E-02 -3.58892588E-02 -3.38378321E-02 -3.13562883E-02
-2.84125964E-02 -2.49741700E-02 -2.10078750E-02 -1.64800383E-02 -1.13564553E-02 -5.60239882E-03 8.17373397E-04 7.93861017E-03
1.57975588E-02 2.44309570E-02 3.38760246E-02 4.00273850E-02 4.63543678E-02 5.28155309E-02 5.93713228E-02 6.59840470E-02
7.26178273E-02 7.92385724E-02 8.58139411E-02 9.23133070E-02 9.87077234E-02 1.04969888E-01 1.11074109E-01 1.16996269E-01
1.22713788E-01 1.28205593E-01 1.33452078E-01 1.38435074E-01 1.43137807E-01 1.47544870E-01 1.51642185E-01 1.55416965E-01
1.58857685E-01 1.61954041E-01 1.64696918E-01 1.67078356E-01 1.69091512E-01 1.70730626E-01 1.71990987E-01 1.72868896E-01
1.73361632E-01 1.73467419E-01 1.73185385E-01 1.72515533E-01 1.71458705E-01 1.70016542E-01 1.68191455E-01 1.65986585E-01
1.63405774E-01 1.60453521E-01 1.57134957E-01 1.53455802E-01 1.49422332E-01 1.45041348E-01 1.40320136E-01 1.35266432E-01
1.29888391E-01 1.24194547E-01 1.18193783E-01 1.11895291E-01 1.05308540E-01 9.84432393E-02 9.13093054E-02 8.39168245E-02
7.62760190E-02 6.83972117E-02 6.02907912E-02 5.19671764E-02 4.34367815E-02 3.47099813E-02 2.57970753E-02 1.67082535E-02
7.45356066E-03 -1.95713846E-03 -1.15141943E-02 -2.12082076E-02 -3.10300645E-02 -4.09709715E-02 -5.10224907E-02 -6.11765749E-02
-7.14256026E-02 -8.17624132E-02 -9.21803418E-02 -1.02673255E-01 -1.13235584E-01 -1.23862364E-01 -1.34549264E-01 -1.45292626E-01
-1.56089496E-01 -1.66937665E-01 -1.77835698E-01 -1.88782973E-01 -1.99779715E-01 -2.10827029E-01 -2.21926940E-01 -2.33082423E-01
-2.44297440E-01 -2.55576977E-01 -2.66927077E-01 -2.78354874E-01 -2.89868630E-01 -3.01477772E-01 -3.13192921E-01 -3.25025934E-01
-3.36989934E-01 -3.49099349E-01 -3.61369943E-01 -3.73818853E-01 -3.86464628E-01 -3.99327256E-01 -3.99822919E-01 -4.00326101E-01
-4.00832095E-01 -4.01336168E-01 -4.01833602E-01 -4.02319746E-01 -4.02790051E-01 -4.03240117E-01 -4.03665733E-01 -4.04062919E-01
-4.04427958E-01 -4.04757435E-01 -4.05048267E-01 -4.05297735E-01 -4.05503505E-01
&end
Example RVCQ3D input le (with stream tube data omitted)
?Title?
&nl1 m=169 n=45 mtl=41 mil=76 &end
&nl2 nstg=2 ivdt=1 irs=1 eps=1.5000 cfl=2.5000 avisc2=0.0000 avisc4=0.5000 ndis=1 ipc=0 pck=0.3000 refm=0.5000 icdup=1 ausmk=0.6000 &end
&nl3 ibcin=1 ibcex=3 itmax=6500 iresti=0 iresto=1 ires=10 icrnt=50 ixrm=0 &end
&nl4 amle=0.6679 alle=64.1600 bete=-56.2000 prat=0.5923 p0in=0.9904 t0in=0.9984 g=1.4000 &end
&nl5 ilt=5 jedge=15 renr=2.267e+06 prnr=0.7000 tw=0.6804 vispwr=0.6670 cmutm=14.0000 itur=2 &end
&nl6 omega=0.0697 nblade=1 nmn=1 &end
&nl7 tintens=0.0100 tmuinf=1.000e-01 hrough=0.0000 &end
70
Appendix C
Additional Plots for Optimization Runs
71
m (in)
Theta
*radi
us
1 1.5 2 2.5 3
-1.5
-1
-0.5
0
0.5
1
1.5 Pressure0.680.65
0.620.590.55
0.520.490.46
0.420.39
Figure C.1: Pressure contours for Run 2 stage results
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure C.2: Isentropic Mach number dis-
tribution for Run 2 stator
0 0.2 0.4 0.6 0.8 1
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure C.3: Isentropic Mach number dis-
tribution for Run 2 rotor
72
m (in)
Theta
*radi
us
1 1.5 2 2.5 3
-1.5
-1
-0.5
0
0.5
1
1.5 Pressure0.680.65
0.620.590.55
0.520.490.46
0.420.39
Figure C.4: Pressure contours for Run 3 stage results
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure C.5: Isentropic Mach number dis-
tribution for Run 3 stator
0 0.2 0.4 0.6 0.8 10
0.1
0.2
0.3
0.4
0.5
0.6
0.7
Surface Arc Length
Isentropic Mach Number
PressureSuction
Throat Loc
Figure C.6: Isentropic Mach number dis-
tribution for Run 3 rotor
73
Appendix D
O -Design E ects on Aerothermal Blade
A very brief o -design analysis was conducted on the optimized stator for Run 1. The
Mach number distribution along the blade while varying the incidence angle from 15o to
+15o is shown in Figure D.1. These results show the blade performs well at negative incidence
angles and has a monotonic velocity increase on the suction side up to almost +10o. This
performance is comparable to that of the T55 and most turbine blades, which generally have
a drastic reduction in performance at incidence angles over +10o. A close-up of the suction
surface is shown in Figure D.2
74
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
Isentropic Mach Number for Varying Incidence Angles
Surface Arc Length
Isentropic Mach Number
?15 deg
?10 deg
?5 deg
0 deg
5 deg
10 deg
15 deg
Throat Loc
Figure D.1: Mach number surface distribution over a range of incidence angles
0.15 0.2 0.25 0.3 0.35
0.64
0.66
0.68
0.7
0.72
0.74
0.76
Isentropic Mach Number for Varying Incidence Angles
Surface Arc Length
Isentropic Mach Number
?15 deg?10 deg
?5 deg0 deg
5 deg10 deg
15 degThroat Loc
Figure D.2: Close-up of Figure D.1 suction side
75
Appendix E
Brief Grid Re nement Study
Along with the referenced articles[31] showing a grid of 169x45 is su cient for turbine
analysis and using the suggested grid in the RVCQ3D User?s Manual[30], a brief grid re ne-
ment study is included to verify the mesh used. The Fortran program Verify available on
the NPARC Alliance CFD Veri cation and Validation Web Site was used to supplement this
analysis. In addition to the Verify program output, plots of the heat transfer coe cient and
Mach number distribution over a variety of grids are provided to show the convergence of the
grid. All RVCQ3D evaluations were allowed to converge at least four orders of magnitude
for the speci ed grid. The program output shows that the 169x45 grid used was within
asymptotic range.
Verify Program Output
--- VERIFY: Performs verification calculations ---
Number of data sets read = 3
Grid Size Quantity
1.000000 0.992430
1.200000 0.992160
1.440000 0.991520
Order of convergence using first three finest grid
and assuming constant grid refinement (Eqn. 5.10.6.1)
Order of Convergence, p = 4.734972
Richardson Extrapolation: Use above order of convergence
and first and second finest grids (Eqn. 5.4.1)
Estimate to zero grid value, f_exact = 0.9926269
76
Grid Convergence Index on fine grids. Uses p from above.
Factor of Safety = 1.250000
Grid Refinement
Step Ratio, r GCI(%)
1 2 1.200000 0.024801
2 3 1.200000 0.058818
Checking for asymptotic range using Eqn. 5.10.5.2. A ratio of 1.0 indicates asymptotic range.
Grid Range Ratio
12 23 0.999728
--- End of VERIFY ---
0 1000 2000 3000 4000 5000 6000 700010
?9
10?8
10?7
10?6
10?5
10?4
Iterations
RRMS
Figure E.1: RRMS plot for 169x45 grid
77
0 0.5 1 1.50
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
Surface Arc Length
Isentropic Mach Number
169x45
250x50
383x70
Figure E.2: Mach number distribution for multiple grid sizes
?1 ?0.5 0 0.5 1
Arc Length Distance from Leading Edge
Heat Transfer Coefficient "H" (BTU/lbm.hr.ft
2 )
85x23, dsi=1e?4
169x45, dsi=1e?5
250x50, dsi=1e?5
360x68, dsi=1e?5
383x70, dsi=5e?6
Pressure Side Suction Side
Figure E.3: Heat transfer coe cient distribution for multiple grid sizes
78