Copyright © IFAC Large Scale Systems, Rio Patras, Greece, 1998
STOCHASTIC DESIGN OPTIMIZATION FOR A SIMPLE 2-D AEROFOIL
W. A. Wrigbt, CM.E. Holden
SO\l'erby Research Centre British Aerospace Plc, Bristol Andy. [email protected]
uk, [email protected]
Abstract: This paper describes the results of a study. which considers the aerodynamic optimization of a 2-D aerofoil, given a number of realistic but simple constraints. It is shown that this simple problem. with a realistic objective function, has a number of local optima. It is demonstrated that stochastic optimisation methods provide a distinct advantage over the more conventional methods and it is further shown that the optimal results are obtained when these methods are combined with conventional methods. Copyright © 1998 IFAC Keywords: Aerodynamic Optimisation. Wing Design. Genetic Algorithms. Simulated Annealing.
the relative performance of the Genetic Algorithm and Simulated Annealing methods on this problem.
Although there are a large number of investigations into the use of stochastic optimisation methods (e.g. Genetic Algorithms (4) and Simulated Annealing (7» in design (for example see (5) and (2)) few of these look at the performance on real aerospace problems. So far only limited studies have been reported in the open literature (1) and most of these concentrate on the use of more conventional optimisation methods.
2. PROBLEM DESCRIPTION 2.1 B-sp/ine representation The complexity of a shape optimisation problem is often governed by the representation used for the shape of the object to optimised. Where possible it is desirable to incorporate any s:mmetry or constraints inherent in the shape of the object In the aerospace industry standards (STEP) have been laid down which outline appropriate representations for aerodynamic geometries. However. what constitutes (within this standard) the most appropriate representation is still an open question.
To test further the performance of these more advanced optimisation strategies a problem is considered here which looks at the optimisation of the drag performance of this 2-D aerofoil. Although relatively simple the 2-D aerofoil provides a flexible and real shape optimisation problem which. not only requires the uses of real objective, or cost, function, but also a number of constraints (on both the general shape and SUUcture of the aerofoil together with its aerodynamic characteristics) to ensure that a realistic and plausible aerodynamic shape is produced. This paper reports the findings of an investigation into
Here the aerodynamic 2-D cross-sectional shape of a NACA 0012 aerofoil was defined by a set of B-spline poles (see figure 1). These poles were fixed in x and varied in y to allow the shape of the aerofoil to be adjusted. The B-spline representation provides a number of advantages. In particular the method
on the fuel usage and this is related to the aircraft drag. Here the drag coefficient. Co. of the aerofoil for a given lift coefficient. CL. is used as the objective function . This provides a convenient nondimensional. measure of tltis aerofoil'S perfonnance.
allows a smooth but flexible shape to be generated from a relatively small set of similar parameters. The similar nature of all the parameters also allows for the easy encoding of the optimisation problem since no one parameter requires a different treatment. However. the flexibility of the representation. as will be shown. does lead to difficulties. To optimize this representation an inviscid Euler method was used to compare the pressure distribution induced by the flow round the aerofoil for various B-spline representations with that of the true geometry. From this it \vas detennined that a Bspline representation with 25 8th order poles provided a geometry with sufficient accuracy for the purpose of this study. It should be noted that during the optimisation. it was found possible to reduce the dimensionality of the representation, with the resulting sa\oing in computational costs. without unduly effecting the quality of the results. For the complete aerofoil optimisations. it was found that 15 3rd order and 23 5th order B-spline poles were sufficient for the two cases presented. In all cases, it was also found that minor modification to the shape representation (i.e .. adjustment of the number of poles and their order) had little effect on the quality of these results. For the GA optimisation each pole position was encoded using 12 bits.
9 S" -
For the purpose of this study t\\·o objective functions were constructed. • A single point objective function . The drag coefficient at a lift coefficient of 0.5 (i.e. (C D) C, =1) )' This drag coefficient is typical for an aerofoil in cruise. •
A multi-point objective function. The sum of the drag coefficients at t\vo specified lift coefficients of 0.2 and 0.6 (i.e. (Co) C,-O_ _ , +(C o ) C, -_0 6 ).
In addition to these objective functions a number of simple, but real, constraints were also chosen. These were that: •
the cross sectional area of the aerofoil ~4) bet\veen tlle leading and trailing edge spars (i.e. the wing box area) nondimensionalised \\ith respect to chord (the thickness of the \\ing) was greater than 0.052. This restriction ensures that the wing is not too heavy, as a thinner wing section incurs a weight penalty and that there is sufficient capacity for fuel Vtithin the \ving-box.
the value of the zero lift pitching moment coefficient (C",o) was restricted such that -0.1
The thickness at 15% and 65% chord (corresponding notionally to the positions of the front and rear spar respectively) were also limited such that they did not go below 9.66% chord and 7.34% chord respectively. This also limits the thinness of the wing. Note the weight of a wing is inversely proportional to its thickness.
. 2 1125 ~
Fig. 1: B-Spline representation of the NACA 0012 Aerofoil. The B-spline poles are shown here numbered from 1 to 25 .
2.2 Objective Function For this study it was important to obtain a realistic objective function (cost function) against which to undertake the optimisation. In general most aircraft optimisation problems use complex objective functions such as Direct Operating Cost. This measures at the "in service" cost of tlle resulting aircraft to the purchaser. Not surprisingly a significant element of this complex function depends
The objective functions were calculated using the panel method of Eppler and Somers (3) . This technique is an inviscid potential method wltich coupled to an engineering boundary layer allows for the viscous effects in the flow around the aerofoil.
(c) Between Middle and Hi9hest
----= - '
. . - --
~- . -.
_ _ _"'-
Figure 2: Four isometric plots of the objective function plotted as a function of the y position of pole 2 (x a.xis) and pole 24 (y axis) for pole 1 fixed at a val ue of (a) -0 .016079, (b) -0.23743 x 10- 4 , (c) 0.0080039, (d) 0.016032. Transition was fixed in the calculation of the penaltv function.
The method was used in incompressible flow mode (i .e zero Mach number M=OO). For the purpose of this study a Reynolds number. based on aerofoil chord. of Re.=6.0 x 10° (where chord. c=l.O) was used. It should be noted that in cruise the typical modern jet passenger aircraft will fly at compressible speeds (i.e. greater than .\/=0.3) any "wave drag" contributions to the total aerofoil drag have. therefore. been ignored here. Codes .vhich are able to cope \vith the shocks produced by the compressible flow together with appropriate viscous effects are more expensive and so were not used here ' . Ho\vever. the study is still relevant since the effect on profile drag. which is not a function of compressible flow. are still of interest and can be examined with this model. Furthermore. this work is directly applicable to the design of high lift devices and UA Vs (Wl-manned airborne vehicles). Both of which operate in a regime where the flO\v is incompressible. For both the simulated annealing and genetic algorithm methods the constraints \vere combined with the objective function using the modified Fiacco and McCormick penalty function
As an alternative t\yo conventional optimisation methods were also considered. These were: • PoweIrs Direct Search Method (10) . \vhich uses a direct search in conjugate directions with quadratic convergence. •
To obtain an Wlderstanding of the complexity of tllis simple problem the objective function as a function of the y position of poles 2 and 2.+ for four fixed values of pole 1 (the pole at end of the trailing edge) \Vas plotted (see figure 2). From these plots it can be seen that even for this problem the objective function behaves in a complex manner. TIllS is due to the separation of the bOWldary layer. For most fixed values of pole 1 the objective function clearly has more than one minima when ploned as a function of the other variables (poles 2 and 2'+). Furthermore, as pole 1 is varied the shape and position of these minima evolve. Because of the presents of this local minima the non-stochastic methods did not perform well on the full aerofoil problems. In fact it was found that for both objective functions the NACA 0012 aerofoil sits in a local minima. The simplex and gradient based methods were ineffective even returning the aerofoil to the original design when the shape was perturbed by moving the poles by 3% of their original position.
Three test cases were examined. In the first the three poles at the trailing edge were allowed to move and transition was fixed in the boundary layer calculation for simplicity. In the second test case, a 15 3rd order B-spline poles representation was used with fixed transition in the boundary layer and the performance of the different optimisation methods investigated. The final case was a mUlti-point optimisation problem, here 23 5th order B-spline poles were used and transition in the bOWldary layer was not fixed such that natural laminar flow could occur.
Initially other problems were encountered. On using the stochastic optimisation methods non-physical designs were obtained. Most aerodynamic methods only provide an accurate approximation of the drag within certain strict assumptions. The Eppler method is unable to cope with separation of the bOWldary layer caused by ''bumps'' on the leading edge of the aerofoil. Such "bump" geometries were often produced by the stochastic algorithms. For such geometries the Eppler method breaks down and returns an unrepresentative low estimate of drag. The generation of such shapes. therefore, had to be eliminated to prevent the stochastic methods from finding a minimum wllich was not physical . TIllS was achieved by applying an additional constraint on the shape of the leading edge. The constraint simply restricted the shape of the leading edge of the aerofoil ensuring that the y coordinate of the upper surface of the leading edge increased towards the front spar and conversely the lower surface y coordinate decreased to\vards the front spar thus preventing "bumps" from occurring.
To undertake this study a number of optimisation methods were utilized from the OPTIONS suite of prograrrunes (6) produced by Dynamics Consulting Ltd. Two stochastic methods were used. •
Simulated Annealing (SA). This method used a temperature dependent Metropolis algorithm (8) \\ith a geometric temperature schedule.
Genetic Algorithm (GA). TillS was an elitest method wllich uses both niche forming and clustering.
Nelder and Meade' s Simplex Method (9).
I Other \vork using a transonic methodology is being undertaken. but cannot be presented due to commercial reasons.
Table I . Table describing the first annealing schedule used for bv the Simulated Annealing algorithm for the oQtimisation of the comQlete aerofoil No of Iterations TemE No TemE 50. 200 1 10. 200 2 1. lOO 3 OA lOO -+ lOO 0.08 5 0.016 200 6 .., I 0.0032 200 200 0.00064 S 0.000128 200 9 0.0000256 200 10 0.0000000001 200 LI
I ,:,.; .... 0.04 1 I.~: 0.02
---~~. : .. •. •.....
Fig. 3: Shape of the complete aerofoil produced by tile Genetic Algorithm .. compared with that produced by the Simulated Annealing optimisation. The result \\'as produced using a geometry represented by 15 3rd order B-spline poles. The penalty function was calculated with transition fixed. The uneven shape is smoothed through further optimisation using the PDS method. Note that the aspect ratio of the aerofoil has been altered artificially to excentuate the differences between the shapes.
3.1 Single Point Optimisation (transition fixed) For the single point objecti\'e function both the simulated annealing and genetic algorithm methods were used. 15 3rd order B-spline poles were found sufficient for tilis problem. The calculation of the penalty function used fixed transition in the boundary layer method. Both optimisation methods were able to give physically realistic results. Of the two the simulated annealing method (with a relatively quick annealing schedule, see table 1) produced the largest reduction in the objective function. The method reduced (CD)
- - "ACAOO12
...... .. ..
after 2200 object function evaluations compared to the genetic algorithm (with a population of 100) wllich achieved a reduction of 13% after 8000 evaluations. The resulting GA and SA aerofoils (see figure 3) have some unevenness around the leading edge. This was smoothed by further local optimisation using the PDS method reducing tile objective function to 30% of the original NACA 0012 value. 3.2 Multi-point Optimization (transition free)
C L ",0. ~5
aerodynamic designs resulting from the optimisation using the genetic algorithm. These shapes are characterized by the aerofoil maximum tllickness having moved aft to maintain large tracts of laminar flow. Figure 5 shows tile impact of large tracts of laminar flow on the aerodynamic characteristics of the aerofoil. The single point design result has a narrow "bucket" giving a low drag at just one point
(a two point design) and,
for comparison, tile optimum of (CD)
Fig. 4: Initial and final aerofoil shapes produced after optimisation of the complete aerofoil by the Genetic Algorithm, for the single and two point opuIIUsations. The geometry was represented using 23 5 th order B-spline poles. The penalty function was calculated with free transition in the boundary layer. Again not the artificially increased aspect ratio.
Here only the GA was used. It was found necessary to use 23 5th order B-spline poles because a more course representation gave less smooth geometries, which trip the boundary layer, with a resultant drag penalty. The population size was increased to 300 to allow for the high dimensionality of the optimization task. To ensure a smooth geometry on the leading edge and pennit natural laminar flow the transition position was calculated by the boundary layer method when detennining the penalty function . Tv,:o scenarios were computed: The optimum of tile sum
(CD )CL=O ~ +(C D )C L=06
single point design) . Figure 4 shows the final
in the graph. The two point design has a wider "bucket" giving a better result across a range of CL values. Below a CL of about 0.1 and above a CL of 0.8, the Co is better for the NACA 0012 aerofoil.
believable answers. Here the maximum thickness of the aerofoil moved to the rear in order to maintain large tracts of laminar flow. Indeed the drag polars contained the sort of laminar flow "buckets" that may be expected from a laminar flow aerofoil.
of . CONCLUSION
It is considered that results obtained using these methods could provide the designer with a useful starting point in those situations where the velocity profile of a similar aerofoil is not available. In reality of course an actual design study would be more complex than the test case investigated here. This example uses a relatively simple objective function with only a limited number of constraints. Clearly further investigation on more complex (possibly multi-disciplinary) problems needs to be undertaken to detennine the true potential of these methods.
This im'estigation has looked at the optimisation of a simple but realistic problem, the shape optimisation of a 2-D aerofoil. The use of a gradient based method gave no change to the NACA 0012 aerofoil this suggests that the geometry sits in a local minima. Because of the presents of local minima it was found necessary to use stochastic optimisation method to search the design space. For the first full aerofoil problem two methods were used; Simulated Annealing and a Genetic Algoritlun. Here it was found that the Simulated Annealing gave a, marginally, better solution with fewer objective function evaluations.
ACKNOWLEDGIvlENTS The authors would like to thank Professor Keane of the University of Southampton for both his help in using the OPTIONS suite of programs and advice on how best to approach this project. CH would also like Clyde Warsop of SAe Sowerby Research Centre without whose help this work would not have been possible .
To allow physical drag calculations to be produced it was necessary constrain the shape of the '10'
--""ACAOO'2 . - .• ••••• . ''''01'' PO'"t
- - - - - - . !'IIIOOOIt'If t e,1I't
AGARD. Optimum design methods for aerodynamics, volume AGARD-R-803 , 1994.  D. Dasgupta and Z. Michalewicz. Evolutionary Algorithms in Engineering Applications. SpringerVerlag, 1997.  R. Epple;- and D. M. Somers. A computer program for the design and analysis of low speed airfoils. Technical Report NASA-Tlvf-B0210, NASA Langley Research Center, 1980 . .  J. H. Holland. Adaption in natural artificial systems. University ofMichgan Press, 1975.  A. J. Keane. Experiences with optimizers in strutural design. In Conference on Adaptive Computing in Engineering Design and Control, pages 14-27, 1994.  A. J. Keane. The OPTIONS Design Exploration System: Reference Manual and User Guide. 1995 .  C. D. Kirkpatrick. C. Dd. Gellatt Jr, and M. P. Vecchi. Optimisation by simulated annealing. SCience. 220 :671-680, 1983.  N. Metropolis, A Rosenbluth, Rosenbluth M., A Teller, and E. Teller. Equation of state calculations by fast computing machines. 1. Chem. Phys., 6:1087, 1953 .  J. A. Nelder and R. Meade. A simplex method for nmction minimisation. Computer Journal, 7 :308313 , 1965 [IOJ M. J D. Powel!. An efficient method for fmding the minimum of a function of several variables wi thout calculating derivatives. Computer Journal, 7(4)3 03-307. 1964. [ I]
., f r
'0, , • \
'- . ... . . :
Fig. 5: Graph of the drag coefficient versus lift coefficient for the NACA 0012 aerofoil and the single and two point designs from the Genetic Algorithm. Note: Re c =6 x 10 6 trial aerofoils by restricting the geometry of the leading edge. The need for such constraints highlights a potential difficulty when using these methods. If the space of designs is not appropriately constrained then the optimisation method may produce shapes which fall outside the working assumptions of the aerodynamic code used to calculate its viability. This was seen in these experiments where the aerodynamic code produced very low drag values for unphysical aerofoil shapes. In the second full aerofoil experiment both the single and two point optimisation results gave physically