Aerodynamic shape optimization using computational fluid dynamics (CFD) is widely used for creating efficient designs in engineering applications such as the automotive or the aerospace domain . Global optimization techniques enable unbiased search for improved shapes and are capable of yielding conceptional new designs. However, they usually require a significant number of CFD simulations for assessing the quality of the proposed designs. This constitutes a substantial drawback for extensive application of global optimization techniques in actual design cycles since realistic CFD simulations are usually rather complex and have a large time and resource demand. Difficulties arise in that context due to the many options available to the designer for setting up such a numerical design optimization. These include the choice of the optimizer algorithm, specific configuration of the parameters for the optimizer algorithm, representation of the design and design changes to be used, and many more.
In practice, one usually employs a specific setup with some more or less well motivated choices. Then optimization is performed and depending on the results some refinements are done and a new optimization run is performed. This is iterated until a satisfactory result is obtained. It is not possible to assess whether the results obtained in this manner are representative. In particular, it is not clear if the results could have been much better, if a slightly different configuration was chosen. Or, if a rather different design could have been obtained which still has a similarly good fitness value. The answers to such questions are of course strongly linked to the qualitative behavior of the fitness function.
In this work, we explore the fitness landscape of a realistic turbofan simulation case. We investigate the outcome of different optimization runs in which we vary the configuration settings. We compare the achievable fitness values obtained with a covariance matrix adaptation evolutionary strategy (CMA-ES)  to a particle swarm optimization (PSO) .
Another key ingredient to shape optimization
approaches is the choice how to represent a specific design and in particular how to generate
sets of different designs to be evaluated during the optimization. In a typical engineering
application, a valid baseline design is readily available and the design process focuses on
improving that design. Deformation methods are especially suited and well-established
for such situations. A crucial decision to be taken by the engineer is the choice of the
deformation method and the details of its realization. In the context of evolutionary
aerodynamic shape optimization, shape morphing methods are popular such as free-form deformation (FFD)
and radial basis functions (RBF) based deformation methods, or shape-function methods like
Hicks-Henne (HH) functions
Another key ingredient to shape optimization approaches is the choice how to represent a specific design and in particular how to generate sets of different designs to be evaluated during the optimization. In a typical engineering application, a valid baseline design is readily available and the design process focuses on improving that design. Deformation methods are especially suited and well-established for such situations. A crucial decision to be taken by the engineer is the choice of the deformation method and the details of its realization. In the context of evolutionary aerodynamic shape optimization, shape morphing methods are popular such as free-form deformation (FFD) and radial basis functions (RBF) based deformation methods, or shape-function methods like Hicks-Henne (HH) functions.
2 Optimization and simulation setup
2.1 Blade representation
During the optimization, the blade geometry has to be changed automatically based on the parameters
proposed by the optimization algorithm. Thus, a blade representation has to be defined. In this
work, we use a deformation approach, where we start from a given baseline blade design and
encode changes to it via the optimization parameters.
Specifically, a blade is represented as a collection of sections stacked on top of each other as shown in
Figure 1 (left). Each section (which can be thought of as a 2D airfoil profile)
is defined by a cylindrical cut of the blade geometry
where the same number of sample points is used for each section.
The 3D geometry is reconstructed from those sections by linearly interpolating between the sections.
(left). Each section (which can be thought of as a 2D airfoil profile) is defined by a cylindrical cut of the blade geometry where the same number of sample points is used for each section. The 3D geometry is reconstructed from those sections by linearly interpolating between the sections.
In the most general situation, all sections would be allowed to change independently to generate a modified blade design. In order to keep the optimization parameters manageable, we select only three sections, one close to the hub, one in the mid-span region, and one close to the shroud, which we deform based on the optimization parameters. The in-between sections are linearly interpolated.
For each of those independently deformed sections, we allow the following changes: (i) rotation of the section around the leading edge point, (ii) movement of the section in the axial-meridional plane, and (iii) deformation of the section profile by adding HH shape functions (see Figure 1 right), which is a well-known representation from 2D airfoil design, see e.g. .The explicit functional form of the HH shape functions we use is
where parametrizes the chord length of each section and controls the location of the maximum of each shape function. The number of shape functions per section will be varied in this work and for a given number of we use equally spaced locations for the maxima, i.e. where .
In order to fully characterize the deformed blade geometry, we need to specify for each of the three independently deformable sections the amplitude for each shape function, the rotation angle, and the shift in axial and rotational direction. Therefore, the total number of free parameters, which is also the number of search variables for the optimization problem, is given by .
2.2 Simulation setup
The turbo-fan simulation setup was inspired by the GE Honda HF120  small business jet engine operated at cruise condition. We used steadyCompressibleMRFFoam, the compressible flow solver from the OpenFOAM CFD suite (specifically, the foam-extend-3.2 version) which was adapted to be more robust . The blade rotation speed was set to a fixed value of 15600rpm and constant mass-flow and radial equilibrium boundary conditions were employed. We used the SST turbulence model. The meshes for all simulations consisted of about 6 million cells and were created by an in-house software based on the OpenFOAM utility blockMesh. The flow equations were solved on 32 CPU cores in parallel, iterated for 12000 solver iterations in order to arrive at a converged solution, and the run-time was typically in the range of 2-4 hours.
The flow-domain of the simulation setup and an example output flow field (pressure) at a circumferential section at a fixed radius is shown in Figure 2. In order to stabilize the CFD solutions we did not include a clearance between rotor blade and shroud geometry. This is somewhat unrealistic as it qualitatively changes the tip vortex caused by such a gap and thus overestimates the blade efficiency. But for the current study it seems a valid simplification which should still provide realistic results except for the blade geometry in the tip region.
2.3 Fitness function and optimization algorithms
The optimization maximizes the isentropic efficiency of the rotor which is calculated from the OpenFOAM CFD solutions as
where and are the mass-flow averaged total pressure and the total temperature at the specified location and is the heat capacity ratio.
The fitness function which is minimized during the optimization is given by
where the efficiency is averaged over the last 1000 flow solver iterations and represents penalty terms which increase (i.e. worsens) the fitness in case the CFD solution does not converge or the generated blade geometry is not feasible.
We employ two different optimization algorithms, the particle-swarm optimization (PSO)  and the covariance-matrix adaptation strategies (CMA-ES) . Unless otherwise stated, for the CMA-ES we always use a populations size of with parents and an initial step size of in relative units where the maximal allowed variation is normalized to one (i.e., amounts to 5% initial variation). For the PSO we use 12 particles (which amounts to 12 evaluations per generation) and the parameters , , and . (see Eq.(6) in ).
3.1 Comparison between CMA-ES and PSO optimization methods
In order to evaluate the impact of the optimization algorithm, we compared the results of the two optimization methods, specifically a () CMA-ES and a PSO with 12 particles. For the blade deformation, we chose =9 which lead to a total of parameters for the optimizer algorithms.
Figure 3 shows the progress of the optimization where the efficiency of each of the proposed individual geometries is shown as function of the optimization generations. The efficiency is normalized to the efficiency of the baseline geometry, i.e. .
As can be seen from the Figure 3, the variation within the population in each generation is much larger for the PSO than for the CMA-ES. This is because we used a rather small initial step size for the CMA-ES, but initialized the particles and velocities for the PSO randomly within the complete search domain. Thus the PSO performs a more global search, whereas the CMA-ES is more local. However, when looking at the progress of the best solution found up to a certain generation, both approaches produce very similar results. Also, the fitness values during later stages of the optimization are approaching very similar results. For example, the best normalized efficiencies were found to be and .
Comparing the actual final optimized geometries, we observe that the two optimization runs converged to different geometries. As can be seen in Figure 4 where the three independently deformed sections are shown, all three sections are different among the baseline, the PSO, and the CMA-ES optimized designs. This leads to substantial differences across the entire blade, which are nicely visualized in Figure 5 where the difference between the geometries is colors-coded. The PSO and CMA-ES designs show qualitatively different changes to the baseline. For example, the PSO has a substantially modified trailing edge region of the fan-rotor whereas the CMA-ES design is mostly changed in the mid-chord upper region of the blade.
3.2 Comparison for different random seed initialization
The stochastic optimization algorithms used in this work need an initialization of the pseudo-random number generator which is done by specifying a random seed. The smallest possible change between two optimization runs is thus given just by using different values for the random seed and leaving everything else unchanged.
We used three different values for the random seed to initialize the CMA-ES optimization run with an otherwise identical setup with =7 which lead to a total of search parameters. The progress of the best solutions produced by the optimization runs is shown in Figure 6. It can be seen that the initial progress is quite similar but after some initial generations, the three runs diverge. During later stages the runs converge again to similar efficiency values, , , .
The comparison of the actual best geometries reveals the same insights as in the previous sections as the final geometries differ strongly between the runs (see Figure 7).
3.3 Variation of the number of evaluations per generation
Another design choice for the optimizer setup is how many evaluations per generation are made. Here we compare two CMA-ES runs with () and () and otherwise identical setup. For the blade deformation, we chose two setups with =9 and =12 which lead to a total of and search parameters, respectively.
The expectation was that the optimization progress would be faster for higher number of evaluations per generation, if measured in terms of generations since a larger region of the search space could be sampled in each generation and the information gathered could be used more efficiently. In practice, with a usually large resource demand for the CFD simulations, there is trade-off between the parallel resources (or licenses) utilized and the overall run-time. This leads to a limit on the number of evaluations running in parallel. Therefore, the more relevant measure for optimization performance is in terms of actual evaluations, i.e. number of CFD simulations.
As can be seen in Figure 8, the optimization progress in terms of actual evaluations was very similar between the and runs for both setups. The runs seem to have slightly larger spread in each generation. It also seems that for larger search dimension, i.e. for =12, a larger population size is beneficial in the initial phase of the optimization. Both tendencies are somewhat instructive, but in order to arrive at a really conclusive evidence, single optimization runs as performed in this work are not sufficient. Instead, one would need to run each optimization multiple times and compute statistics over the results. However, this was not done due to the already high resource demand and long run-times of the CFD simulations.
Comparing the achieved efficiency and the actual optimized blade geometries, the same observations as in the previous sections can be made. Even though the fitness values are comparable, (=9: , ; =12: , ), the resulting geometries differ substantially between the runs over the complete blade for =9 whereas the geometries are rather similar for =12 (see Figure 9).
3.4 Variation of the number of Hicks-Henne shape functions
A very important choice in the course of blade optimization is the specification of the number of HH shape functions per section. It is expected that this parameter has a strong influence on the outcome of the optimization. Choosing a rather small number leads to a very rigid blade representation where only a limited number of designs can be realized which might not include very efficient blades. On the other hand, a very large number might allow to represent a plethora of shapes in principle, including very efficient designs, but the resulting search space has a very high dimension and thus the optimization problem becomes very hard. This usually leads to a much slower progress of the optimization and threatens that a sufficiently improved design will not be found within a reasonable time and the limited resources for the search.
We explored this trade-off by running several optimization runs and varying the number of shape functions ranging from =3 up to =12. This amounted to search dimensions ranging from to , and we used the CMA-ES with a populations size of either or .
The best results as a function of the optimization generations are shown in Figure 10. First of all it can be observed that there is no clear trend for the achievable efficiency as function of . For the smallest number of shape functions per section =3 indeed gives the lowest improvement, but for example, for =4 a rather large improvement is realized whereas for =5 the improvement is as small as for =3. In case of , the best efficiency is actually produced with the lowest =5.
Figure 10 suggests that the optimization runs are not yet converged after 120 generations but instead keep improving with more generations. This is actually expected for such population sizes and search dimensions. Therefore, it still could be that there is a clear relation between the achievable efficiency and , but for all practical purposes, with limited computing resources, this is not relevant. The increased flexibility from a larger number of design parameters is compensated by the increased difficulty of the search problem.
The preliminary insight from this is that already with a rather small number of shape functions per section a considerable improvement in the efficiency can be realized. So the choice of is not determined by the expected efficiency increase but rather based on the geometrical aspects implied by the blade representation and the available optimization run time.
We again observe that all the resulting geometries are rather different over large parts of the blade which can be seen in Figure 11. In particular, the shapes are usually much more curvy and wavy during the optimization for larger values of , but with increasing generations the shapes seem to be smoothed out again even for larger values of (not shown).
We studied multiple shape optimization runs for a realistic turbo-fan blade geometry setup. In order to explore the fitness landscape, we investigated many setups, varying the optimization algorithms, the random seed for the initialization, the population size and the number of Hicks-Henne shape functions per independently deformable blade section.
The number of Hicks-Henne shape functions does not seem to have a pronounced influence on the achievable efficiency improvement, at least with the realistic constraint that only a small number of evaluations (i.e. optimization generations) can be done. This most likely directly reflects the trade-off between the increasing flexibility with a larger number of design parameters and the increasing difficulty of the search problem. For the practical application this implies that a small number of HH shape functions can already be sufficient and geometrical aspects for the location and number of shape functions might be equally or maybe even more important.
The central finding of this work is that for all the tested variants, comparable improvements in the rotor efficiency could be achieved but the actual optimized geometries showed substantial variation over the complete blade geometry. This leads to the conclusion that the fitness landscape of such a realistic turbo-fan optimization is highly multi-modal with many local minima. Even minor changes (e.g. the random seed for the initialization of the optimization algorithm) may lead to very different geometries with comparable efficiencies.
This has strong implications for practical blade design optimization approaches. On the one hand, the designer can easily generate a set of design variations which all have comparable efficiencies by minor changes to the setup. On the other hand, this makes surrogate-assisted optimization approaches  more difficult. First of all, training a surrogate model on a highly multi-modal function requires large amounts of data which are most likely not available thereby reducing the model accuracy drastically. Second, data from finished optimization runs will be dense only around one or a few local minima which makes models trained such data not very helpful for optimizations runs which end up converging to another local minimum.
In the current work we only investigated a single operating condition (fixed mass-flow and rotational speed
at cruise condition) and only the rotor passage. A more realistic scenario should include the downstream
stator blade which is know to have a strong influence on the fan blade design. Also, additional
operating conditions with different mass-flow and rotational speed (e.g. take-off condition) need to be
The expectation of the authors is that the fitness landscape of such a more realistic setup will still be
multi-modal, probably to a lesser extent than for the current case, but the general findings of this
work will still be valid. Such investigations are left for future work.
In the current work we only investigated a single operating condition (fixed mass-flow and rotational speed at cruise condition) and only the rotor passage. A more realistic scenario should include the downstream stator blade which is know to have a strong influence on the fan blade design. Also, additional operating conditions with different mass-flow and rotational speed (e.g. take-off condition) need to be incorporated. The expectation of the authors is that the fitness landscape of such a more realistic setup will still be multi-modal, probably to a lesser extent than for the current case, but the general findings of this work will still be valid. Such investigations are left for future work.
The authors thank Hisato Tanaka, Kunio Nakakita, Jiří Šimonek, Radek Máca and Tomáš Krátký for valuable discussion.
-  See, for eaxmple, Endicott, G., Olhofer, M., Arima, T., Sonoda, T. (2014). A Novel Transonic Fan Swept Outlet Guide Vane Using 3D Design Optimization. Proceedings of the ASME Turbo Expo; Skinner, S., Zare-Behtash, H. (2017). State-of-the-art in aerodynamic shape optimisation methods. Applied Soft Computing 62, pp 933-962.
-  Sieger, D., Menzel, S., Botsch, M., (2014) Constrained Space Deformation for Design Optimization. Procedia Engineering 82, pp. 114-126; Sonoda, T., Schnell, R., Arima, T., Endicott, G., Nicke, E. (2013) A Study of a Modern Transonic Fan Rotor in a Low Reynolds Number Regime for a Small Turbofan Engine Proceedings of ASME Turbo Expo 2013; Kaandorp, M., Menzel, S, Schmitt, S. (2017) An Aerodynamic Perspective on Shape Deformation Methods. AIAA AVIATION Forum (AIAA 2017-3145);
For a recent review see, Bonyadi, M. R., Michalewicz, Z. (2017). Particle swarm optimization for single objective continuous space problems: a review
. Evolutionary Computation 25. pp 1-54.
-  Hansen, N. (2006). The CMA evolution strategy: a comparing review. Towards a new evolutionary computation. Springer, pp 1769-1776.
-  Hicks, R., Henne, P. (1978). Wing Design by Numerical Optimization. Journal of Aircraft 15, pp 407–413; Wu, H., Yang, S., and Liu, F. (2003). Comparison of Three Geometric Representations of Airfoils for Aerodynamic Optimization. Proceedings of the 16th AIAA Computational Fluid Dynamics Conference: 4095
-  See, for example, Honda R&D Technical Review Vol.26 No.2 (2014).
-  Rusche, H., Schmitt, S. (2016). Stability Improvements of Pressure-Based Compressible Solver and Validation for Industrial Turbo Machinery Applications. 4th Annual OpenFOAM User Conference 2016.
-  See, for example, Mengistu, T., Ghaly, W. (2008) Aerodynamic optimization of turbomachinery blades using evolutionary methods and ANN-based surrogate models. Optimization and Engineering 9, pp 239-255.