Optimal Propeller Design for Future HALE UAV

The main roles of unmanned air vehicles (UAVs) include: observation, surveillance, transportation, remote sensing and various security tasks. Improved, augmented type of UAVs are high-altitude long-endurance (HALE) aircraft capable and designed, as their name suggests, for lengthy flights at higher altitudes (which also usually implies subsonic cruising velocities). Different variants, in both size and applied technical solutions, have been tried. Common approach incorporates standard wing-fuselage-aft empennage configuration and propelled flight as the most efficient for the required speed range. The paper gives a brief overview of a preliminary aerodynamic analysis of the main lifting surfaces as well as a detailed description of the performed multi-objective optimization of the propeller capable of producing a sufficient amount of thrust at the cruising altitude and speed. Aerodynamic performances of the investigated propellers are estimated by a simple blade element momentum theory (BEMT). The chosen optimizing method, genetic algorithm (GA), is suitable for dealing with a large number of input variables.


Introduction
ODAY, HALE UAVs, capable of flying at altitudes higher than 15000 m for more than 24 hours, present one of the most exciting segments of aerospace industry [1][2]. Operation in deteriorated atmospheric environment as well as quite restrictive demands on efficiency and reliability require: the design of aerodynamically and structurally improved solutions, the use of novel materials and propulsion, etc [2][3][4][5][6]. Interesting examples of ongoing projects can be found in [1,5,6].
Conceptual design of a new aircraft is an iterative process comprising several steps. Everything starts with the clear definition of the initial requirements. This is followed by the weight estimation, aerodynamic simulations, the choice and dimensioning of the propulsion system, flight mechanics, stability and control assessment, structural analysis of aircraft parts, equipment definition, cost-benefit analysis, etc. which is repeated until reaching a viable design. Therefore, the structure of this paper follows the initial phases of the conceptual design of a novel HALE UAV.
The initial aircraft requirements are: -Cruising altitude h = 15000 m, where air density is ρ = 0.193674 kg/m 3  Brief description of the assumed aircraft concept and initial weight estimation are presented in the following section. The chapter "Preliminary wing performance computation" describes the performed numerical simulations of the flow around a clean wing, while the chapter "Propeller design" explains how the optimal geometrical shape of the HALE UAV propeller was attained. Brief conclusions are given in the end.

Aircraft concept
After the analysis and comparison of similar aircraft, in particular Boeing: Phantom Eye [5] and ALTOY: Global Observer [6] presented in Fig.1 (more examples are listed in [1]), a standard wing-fuselage-aft empennage configuration was adopted. The main lifting surfaces are the wings, and they make up most of the aircraft structure. Horizontal flight is made possible by the two propellers installed symmetrically on two halves of the wing.While the choice of actual power plant (modified piston, turbine or hybrid) remains undecided for now, the present study also provides an initial assessment of the required power. T By assuming a simple mission where the aircraft spends most of the time in cruising or loitering, the initial weight estimation was performed. As described in [1], the design take-off weight of the aircraft W o can be considered as a sum of the empty weightW e , fuel weight W f and useful payload weight W pl . For the fixed mass of useful payload W pl = 200 kg and the desired range R = 1300 km the following values can be obtained: -Fuel mass W f = 1960 kg (also roughly corresponding to E> 16 h with the assumed average specific fuel consumption c = 0.3 kg/kWh and available power of 400 kW), Take-off gross mass W o = 4507 kg.By adopting an average value of wing loading from the range W/S = [25 kg/m 2 , 35 kg/m 2 ] typical for this kind of aircraft, it is possible to make a starting guess of the required reference area of the wing S = 139.6 m 2 . Since the cruising air density and aircraft speed are known, cruising lift coefficient can be obtained (C L = 0.346) and a choice of an adequate airfoil can be made.Here, a reliable NACA 4415 airfoil, with well proven aerodynamic characteristics, was chosen in the initial wing design phase. Since its drag is minimal at an approximate angle-of-attack α ≈ −1.5°, and drag force directly influences the required propeller thrust force intensity, the wing built-in angle (and cruising angle-of-attack) are planned to be 0°. The final definition of the wing's geometrical shape requires a more detailed study of the effects of wing aspect ratio AR and taper ratio λ to its aerodynamic characteristics. Since the performed process is described in [1], here, only the adopted values are summarized: -Aspect ratio AR = 22, taper ratio λ = 0.6.

Preliminary wing performance computation
Numerical simulations of the incompressible, turbulent flow around an untwisted wing with the mentioned geometrical characteristics were performed in the commercial software package ANSYS FLUENT [7] where the governing flow equations, Navier-Stokes and additional transfer equations for turbulent quantities, are solved by the finite volume method.
The computational half-domain is shaped like a cuboid stretching 5c before and 20c aft of the wing, where c is the wing mean chord. The domain also extends 5c above and below the wing, and 2b from the wing tip to the lateral boundary, Fig.2a.
Different sizing functions are defined along the wing edges and surfaces. Twenty layers of prismatic cells encompass the wall boundaries and enable the resolution of the boundary layer flow, Fig.2b. The total number of cells is somewhat over 400000. Whilst it might seem rough, it is quite sufficient for an initial computation of a clean, steady subsonic flow (free of any problematic flow phenomena).
The flow is modeled as steady, incompressible and turbulent. k−ω SST turbulence model [8], as a usual choice in aeronautical applications, is used for the closure of the system of basic flow equations. The cruising speed V c = 97.22 m/s is defined at the inlet boundaries, while the atmospheric pressure at the altitude of h = 15000 m is set at the outlet. Wing walls are modeled as rigid and no-slip. For the coupling of pressure and velocity fields a segregated solver and SIMPLE scheme were used. All spatial discretizations are of 2 nd order. The computations were performed until the convergence of wing aerodynamic coefficients of lift and drag, C L and C D respectively, was achieved. Computed wing lift C L and drag C D curves as functions of angle-of-attack (AoA) α are sketched in Fig.3.Two additional aerodynamic characteristics, C L /C D and C L 1.5 /C D , are usually considered when maximum range or endurance is of interest. However, in this case these factors reach their extreme values at α = 4° and α = 6°, respectively. Given that at those anglesof-attack the drag is nearly doubled compared to the value at α = 0° (which also implies that the required propeller thrust force should be proportionally higher) they were not used here for the definition of the cruising regime. While the lift of the wing can be projected to the whole aircraft, a certain increase in drag can be expected with the introduction of the subsequent aircraft components such as fuselage and tail surfaces. In steady, cruising regime, the thrust forces generated by the propellers should be no less than the total drag force. To be on the safe side (by assuming an increased value of the drag of the whole aircraft compared to the wing alone, i.e. by decreasing the finesse of the whole aircraft) the minimal total thrust force for maximal take-off weight should at least be T min = 4340 N, where aircraft lift-to-drag ratio is assumed to be L/D = 10. For lesser aircraft weights, the total thrust generated by propellers can be somewhat reduced.

Propeller design and optimization
Propeller aerodynamic performances primarily depend on the blade number N b and geometrical shape, i.e. span-wise distributions of relative thickness and curvature (airfoils), chord c(r) and pitch θ(r) along the blade. The geometries in question are generally quite complex and curved, and therefore require a large number of parameters for appropriate description. On the other hand, if an adequate selection of input variables is made and their optimal values determined, it is possible to drastically increase propeller efficiency η. Furthermore, reliable operation of the designed propeller in certain, specified, even unconventional working regimes also becomes attainable. Some examples of the research of high altitude propellers can be found in [9][10][11].
Again, to commence the process of the conceptual design of the propeller, the design task in this example can be formulated as: the development of an optimal propeller that provides an approximate thrust force of 1800 N ≤ T ≤ 2200 N at the cruising altitude of h = 15000 m and cruising speed of V c = 97.22 m/s (one of two existing symmetrically located on both left and right sides of the wing). The choice of two engines/power plants is made for redundancy and cost reduction (since larger engines cost more and spend more fuel). Due to the compressibility effects the blade tip speed is limited to V tip = 200 m/s (M tip ≈ 0.678). The cruising advance ratio is then expected to be J = πV/V tip = 1.523 with a certain, required reserve on both sides.
The imposed requirements call for an advanced design approach in the form of multi-objective optimization procedure described below. However, it is first necessary to give a brief overview of the employed computational model used for the estimation of blade aerodynamic performances.

Computational BEMT model
The formulation of BEMT computational approach is derived from the combination of momentum theory (MT) and blade element theory (BET), hence the name BEMT. Rotor disk is divided into ring-shaped stream tubes of small width over which the changes of flow quantities are observed. Thevalues of induced velocities in axial and tangential directions are of particular interest. Aerodynamic contributions (elemental forces and moments) are estimated from the aerodynamic characteristics (lift and drag coefficients) of local airfoils. Fig.4 illustrates the fluid flow around a representative blade element with characteristic velocities, angles and forces clearly denoted. Induced velocities are usually quantified by induction factors, a and a', as: Local increments of the thrust force and torque appearing along the elemental ring (from MT) are: The local induced angle φ and angle-of-attack α are respectively: The relative (total) local velocity is: On the other hand, from BET: By equating expressions (2) and (5), it is possible to compute implicitly defined induction factors as: ,  1  cos  sin  4sin cos  1 ,  ,  1 sin cos for known local values of φ, c L (α) and c D (α) where σ' = N b c/(2yπ) denotes local solidity. In general, eq.6 can be iteratively solved by conducting the following steps: initial assumption of a and a' (e.g. 0), calculation of local φ and α (eq.3), estimation of local c L (α) and c D (α), computation of new a and a' (eq.6). The process is repeated until reaching converged values of induction factors, which is followed by the computation of relative velocity (eq.4) and aerodynamic contributions (eq.5). In the end, global values of thrust and torque are obtained by integration (summing) along the blade. If necessary, dimensionless coefficients of the propeller can also be computed: , .
Description of model parameters As previously mentioned, the blade is fully geometrically (and aerodynamically) described by the three distributions along the blade: airfoil (thickness and curvature), chord and pitch. Theoretically ideal longitudinal changes follow the shape of hyperbolic functions. However, this design is not feasible for a multitude of working regimes, and produces oversized values near the blade root. Instead, standard solutions mostly assume polynomial relations/transformations that can be easily represented by few independent parameters depending on the order of the polynomial.
The work presented in [9] uses a single airfoil (constant function) and 4th order distributions of both chord and pitch. Similarly, [10] employs four predefined airfoils together with 3 rd order functions for c(r) and θ(r). This research describes a possible approach to the optimization of all three relevant blade features without the use of any constant or predefined values. However, for simplicity some assumptions and simplifications are made: -Blade diameter is set to D = 5 m which is quite similar to the diameter of the Boeing: Phantom Eye of 4.9 m. For the adopted V tip = 200 m/s, the angular velocity of the propeller is constant and equal to Ω = 80 rad/s ≈764 rpm. However, it must be mentioned that the initial design objectives can be achieved in more than one way with somewhat smaller (but thicker, more curved or faster) propellers. A different combination of Ω, D and V tip would alter and shift the Pareto frontier in x-and y-directions. -Airfoil distribution is defined by two cross-sections, located at the blade root (y/R = 0.2) and tip (y/R = 1.0). Since airfoils are complex, spline curves originally defined by a large number of points, it is necessary to use a technique that adequately represents them by the smallest possible number of input parameters. One such approach is Class-Shape Transformation (CST) developed by Brenda Kulfan [12]. In short, the method substitutes complex curves by the sums of Bezier polynomials multiplied by the corresponding coefficients that actually represent the input parameters that are going to be optimized. Here, each segment (lower and upper side) of the airfoil is defined by three coefficients, resulting in 12 independent variables used for the description of the blade thickness and curvature. In order to avoid structurally unacceptable solutions (airfoil shapes), lower and upper boundaries of possible parameter values are defined. Maximal relative thickness of the root airfoil must take values from the range [6%, 35%] while the tip airfoil can be twice thinner. -Relative chord distribution is assumed as 2 nd order polynomial c(r)/R = ar 2 + br + c defined by three coefficients a, b and c. -Pitch distribution is also a quadratic function, but with an additional condition θ(0.7) = 0° resulting in the form θ(r) = a(r 2 − 0.7 2 ) + b(r − 0.7) and two input parameters. Overall, the total number of input parameters is 17 (= 12 + 3 + 2).They are all assumed as continuous and allowed to take any value from the defined ranges. Lower and upper boundaries of the domain of each variable are defined in a way they make physical sense (e.g. max(θ) <60°, rotor solidity σ< 0.15).Number of blades is set to N b = 4 (again, similar to Boeing: Phantom Eye aircraft [5]).

Optimization
Apart from the choice and definition of input parameters, in order to reach the optimal solution it is equally important to make a proper selection ofthe cost functions and constraints, i.e. the criteria of evaluation of each considered design. In this case, although simpler, single-criterion optimization does not represent an adequate choice since the objective of maximal efficiency max(η) might result in insufficient intensity of thrust force T. On the other hand, the goal of maximal thrust force max(T) might result in too high values of required power P (and cost). Therefore, in this case, a two-objective optimization of a HALE aircraft propeller in the cruising regime (corresponding to J = 1.527 and Δθ = 40°) was performed: Two, equally important goals of the study are the maximization of efficiency and maximization of thrust force. Listed, purely geometrical constraints are imposed to ascertain the selection of feasible designs whose aerodynamic performances can, with sufficient accuracy, be estimated with the chosen computational model. Aerodynamic computations are performed by BEMT model, a satisfactory compromise between simplicity, ease of implementation and performance and accuracy. In this step, no corrections for blade root and tip were added to the basic model because the authors wanted to keep the problem definition as simple as possible and test the proposed methodology without any additional modifications and effects. However, in future work, a simple Prandtl function can easily be implemented. Each blade was segmented into 50 elements. As a necessary prerequisite for the complete blade simulation, aerodynamic characteristics of each generated root and tip airfoil are obtained by a simple and practical tool − XFOIL [13] that assumes potential flow but incorporates certain viscosity effectswith respect to the imposed constant value of Reynolds number of 1 MRe. Since aerodynamic computations were performed for angles-of-attack −5° ≤ α ≤ 20°, the estimated lift and drag curves, c L (α) and c D (α) respectively, had to be additionally extrapolated. This was done by simple combinations of constant and linear functions in order to simulate the odd parity of the lift function and even parity of the drag function. To further simplify the computational process, the aerodynamic characteristics of inner segments located between the root and the tip were not separately estimated. Instead, it was regarded that blade airfoil characteristics can be applied to segments in the first half of the blade (r< 0.5) while tip characteristics were used for 0.5 ≤ r ≤ 1.
Genetic algorithm (GA), based on the principles of biological evolution and natural selection [14][15], is adopted as the optimization method. It is widely used in problems involving both continuous and discrete variables, non-linear cost functions and constraints as well as insufficiently known factors or physical phenomena.It operates over a population constituted from individuals described by the vectors of input parameters that are randomly formed at the beginning of the procedure.The basic steps (phases) of the method include: matching and reproduction of the individuals − parents (usually, individuals that have better fitness function values have the advantage in mating), crossover (formation of the children − next generation individuals) and mutations that introduce a level of stochasticity into the process and ascertain that the algorithm does not get stuck in local extrema. Additional options include elitism that means that a certain number of the best individuals survives into the next generation without any change. Here, the population numbers 200 individuals since higher values of input parameters and cost functions require larger populations. The number of computational steps, i.e. generations is 150. This value was adopted after repeating the optimization process several times (with different combinations of population size and generation number) and establishing that no significant changes in Pareto frontier appear after 150 computational steps with a population of this size (i.e. 200).

Optimal results
The results of the multi-objective optimization in the form of Pareto frontier, a set containing equally good individuals by both criteria, are illustrated in Fig.5. Scores of the last generation are marked by blue dots. As expected, higher values of thrust imply lower efficiencies. However, for the range of interest T[1800 N, 2200 N], the propeller efficiency can be expected within the bounds η[0.940, 0.946] which can be considered quite satisfactory.These overly optimistic values of efficiency are obtained for several reasons, initial model assumptions and implemented simplifications: -Blade extent, 0.2 ≤ r ≤ 1, is large. In reality, the pitching mechanism and modified blade tips would reduce the operational part of the blade. -Losses due to wake separation at the blade root and tip were not considered. -Mach number (compressibility) effects and wave drag were not considered. -Fluid viscosity was only partially considered through certain corrections that might be sufficient for 2D stationary airfoils, but apparently not so for 3D rotating lifting surfaces. However, even though propeller efficiencies will certainly be lower in real application, this optimization framework still enables the determination of optimal solutions. Since all optimal solutions are equally good, and perform quite similarly, one representative individual is chosen for further analysis and inspection. Attained blade root and tip airfoils are sketched in Fig.6. Both maximal relative thickness of the root airfoil of 19% and its maximal camber of 1.4% areachieved approximately at 35% chord length. Its computed zero-lift angle-of-attack is α n = −0.61°, lift gradient a = 5.77, stall angle α cr = 18.5°, while minimal drag coefficient of c D,min = 0.0065 is achieved at α = 0.5°. Analogously, for the tip airfoil, the following characteristics are obtained: maximal relative thickness is 6% at 35% chord, extreme negative camber is −0.5% at 32%, α n = 0.5°, increased a = 7.22, delayed stall with decreased minimal drag coefficient of c D,min = 0.0035 again achieved at α = 0.5°.Furthermore, since the airfoil in question is thin its behavior near stall is smooth and gradual without a definite value of critical angle-of-attack which can be termed as delayed stall.
Optimal pitch and chord distributions along the blade are given in Fig.7.Global rotor solidity ratio ofσ = 0.115 is achieved with nearly rectangular blades. The explanation for this shape which is not typical for propeller blades lies in the adopted goal functions and disregard of the blade structure. Namely, since one of the imposed goals is maximal thrust and blade tips aerodynamically contribute the most, without taking into account the structural limitation and bending torques that will inevitably appear at the blade root, the greatest possible chords will be obtained in the vicinity of the blade tip. This shape would certainly change if an integral aerodynamic-structural optimization was performed. The most important and interesting propeller feature that requires additional investigation is the complete map of its aerodynamic performances obtained at different combinations of advance ratio J and collective pitch Δθ. Although the optimization was performed solely for the cruising regime, Fig.8 demonstrates that the estimated propeller operation at other working regimes can also be considered satisfactory with efficiency η mostly taking the values in the range [0.75, 0.95]. Small decrease of efficiency compared to the values presented in Fig.5 is due to the interpolation of airfoils and their aerodynamic characteristics along the blade. Still, optimized propeller seems able to produce 2000 N of thrust force at the velocity of approximately 325 km/hat Δθ = 40° and the "cost" of 200 kW. Similarly, it should be able to produce the required 350 km/h at somewhat increased collective pitch Δθ ≈ 42°. With the increase of collective angle, the same amount of thrust can be maintained at higher aircraft speeds and the same required mechanical power. Analogously, at constant Δθ, the aircraft velocity can be increased by decreasing the required thrust and vice versa.

Conclusions
The paper presents the initial phases of the conceptual design of a HALE UAV. The aircraft in question is designed to fly at high altitudes with a long endurance. Several existing similar examples have been mentioned and overall aircraft configuration chosen. Its characteristic weights have been estimated, main lifting surfaces dimensioned and basic aerodynamic performances computed.
The next step is the definition of an adequate propeller. For that purpose, a multi-objective GA has been employed. For the defined number of blades and rotor diameter, blade geometry can be described by 17 independent variables describing the longitudinal distributions of relative thickness, curvature, chord and pitch along the blade. One possible choice of cost functions is the simultaneous maximization of efficiency and thrust at a certain operating regime (e.g. cruise). The optimization study confirmed that it is possible to achieve the requested 2000 N of thrust force at the velocity of 350km/h for the required power of approximately 200 kW per propeller. The performed study also serves as a base for power plant definition and final choice. However, the described, new geometry realized through optimization is a direct consequence of the imposed goal functions and constraints. For any other combination or choice of requirements, different geometries could be obtained. Also, for a different total number of propellers (e.g. 4 instead of the assumed 2), value of rotor diameter D or blade number N b a new set of optimal solutions would be achieved.
These are some of the main advantages of the multiobjective approach: the optimal solution is not unique and many effects can be simultaneously considered! However, this also implies that, in order to define a viable, usable and truly most favorable geometry, deep knowledge of the investigated problem is needed. The presented work can be additionally improved by: upgrading the computational aerodynamic model so that it takes more physical effects into account (in particular drag and fluid viscosity) through various semi-empirical corrections, considering structural aspects of the blade, performing integral aerodynamicstructural optimization of the blade, increasing the number of input parameters, performing an optimization study of the complete aircraft instead of an isolated propeller, etc.