PROPELLER PERFORMANCE ANALYSIS AND MULTIDISCIPLINARY OPTIMIZATION USING A GENETIC ALGORITHM Except where reference is made to the work of others, the work described in this dissertation is my own or was done in collaboration with my advisory committee. This dissertation does not include proprietary or classified information. ____________________________________ Christoph Burger Certificate of Approval: _______________________ _______________________ John E. Burkhalter Roy J. Hartfield, Jr., Chair Professor Emeritus Professor Aerospace Engineering Aerospace Engineering _______________________ _______________________ Robert S. Gross Ronald M. Barrett Associate Professor Associate Professor Aerospace Engineering Aerospace Engineering University of Kansa _______________________ George T. Flowers Interim Dean Graduate School ii PROPELLER PERFORMANCE ANALYSIS AND MULTIDISCIPLINARY OPTIMIZATION USING A GENETIC ALGORITHM Christoph Burger A Dissertation Submitted to the Graduate Faculty of Auburn University in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy Auburn Alabama December 17, 2007 iii PROPELLER PERFORMANCE ANALYSIS AND MULTIDISCIPLINARY OPTIMIZATION USING A GENETIC ALGORITHM Christoph Burger Permission is granted to Auburn University to make copies of this dissertation at its discretion, upon request of individuals or institutions at their expense. The author reserves all publication rights. _______________________ Signature of Author _______________________ Date of Graduation iv VITA Christoph Burger, son of Richard and Baerbel (Kundler) Burger, was born on September 18, 1972 in Munich, Bavaria (Germany). Christoph graduated from High School in St.Georgen, Freiburg (Germany) in June of 1993. From October of 1993 until April of 1995 he served as an EMT in fulfillment of the mandatory military service. He attended the University of Stuttgart in October of 1995 and graduated with a Diploma Engineer (Master of Science) Degree in Aerospace Engineering, in 2001. In the fall of 2001, Christoph started work as a Research Scientist at the Aerospace Engineering Department at Auburn University, Alabama. After about two years of research and development work on unmanned aerial systems, he enrolled as a student at Auburn University to complete his Doctor of Philosophy Degree in Aerospace Engineering. v DISSERTATION ABSTRACT PROPELLER PERFORMANCE ANALYSIS AND MULTIDISCIPLINARY OPTIMIZATION USING A GENETIC ALGORITHM Christoph Burger Doctor of Philosophy December 17, 2007 (Diploma Engineer (M.S.), University of Stuttgart 2001) 181 Typed Pages Directed by Roy J. Hartfield A propeller performance analysis program has been developed and integrated into a Genetic Algorithm for design optimization. The design tool will produce optimal propeller geometries for a given goal, which includes performance and/or acoustic signature. A vortex lattice model is used for the propeller performance analysis and a subsonic compact source model is used for the acoustic signature determination. Compressibility effects are taken into account with the implementation of Prandtl-Glauert domain stretching. Viscous effects are considered with a simple Reynolds number based model to account for the effects of viscosity in the spanwise direction. An empirical flow vi separation model developed from experimental lift and drag coefficient data of a NACA 0012 airfoil is included. The propeller geometry is generated using a recently introduced Class/Shape function methodology to allow for efficient use of a wide design space. Optimizing the angle of attack, the chord, the sweep and the local airfoil sections, produced blades with favorable tradeoffs between single and multiple point optimizations of propeller performance and acoustic noise signatures. Optimizations using a binary encoded IMPROVE ? Genetic Algorithm (GA) and a real encoded GA were obtained after optimization runs with some premature convergence. The newly developed real encoded GA was used to obtain the majority of the results which produced generally better convergence characteristics when compared to the binary encoded GA. The optimization trade-offs show that single point optimized propellers have favorable performance, but circulation distributions were less smooth when compared to dual point or multiobjective optimizations. Some of the single point optimizations generated propellers with proplets which show a loading shift to the blade tip region. When noise is included into the objective functions some propellers indicate a circulation shift to the inboard sections of the propeller as well as a reduction in propeller diameter. In addition the propeller number was increased in some optimizations to reduce the acoustic blade signature. vii ACKNOWLEDGMENTS I would like to thank Dr. Roy J. Hartfield for his exceptional support and guidance during the exertive time. Without his help and especially his patience it would have been impossible to complete the Degree. I also would like to thank Dr. John E. Burkhalter for his advice and help with the programming work. viii Style manual or journal used: The American Institute of Aeronautics and Astronautics Journal Computer software used: Microsoft Word 2003, IMPROVE ? 3.1 Genetic Algorithm, Compaq Visual FORTRAN, Tecplot 10, Microsoft Excel 2003 ix TABLE OF CONTENTS LIST OF FIGURES ........................................................................................................... xi LIST OF TABLES............................................................................................................ xv NOMENCLATURE ........................................................................................................ xvi 1 INTRODUCTION ...................................................................................................... 1 2 BASIC THEORY........................................................................................................ 5 2.1.1 Lifting-Line Theory .................................................................................... 6 2.1.2 Lifting Surface Theory................................................................................ 9 2.2 Propeller Wake.................................................................................................. 13 2.3 Propeller Performance Prediction by VLM ...................................................... 16 2.4 Propeller Geometry........................................................................................... 30 2.4.1 2-D Airfoil Geometry ............................................................................... 30 2.4.2 3-D Propeller Geometry............................................................................ 35 2.5 Compressibility................................................................................................. 39 2.6 Viscous Model .................................................................................................. 46 2.7 Flow Separation Model..................................................................................... 48 2.8 Acoustic Signature Model................................................................................. 53 3 PROPELLER AERODYNAMIC PERFORMANCE PREDICTION PROGRAM. 63 4 CODE VERIFICATION........................................................................................... 66 4.1 Gauss Seidel Solver convergence ..................................................................... 66 4.2 Wake investigation............................................................................................ 70 4.3 Number of Chordwise and Spanwise Panels .................................................... 73 4.4 Validation of the numerical method ................................................................. 75 4.4.1 Comparison to lifting line and lifting surface method.............................. 75 4.4.2 NACA 109622 Propeller........................................................................... 78 5 GENETIC ALGORITHM ........................................................................................ 82 5.1 Optimization methods....................................................................................... 82 5.2 Genetic Algorithm ............................................................................................ 85 x 6 IMPLEMENTATION OF THE PROPELLER PERFORMANCE PROGRAM INTO THE GENETIC ALGORITHM ............................................................................. 88 7 OPTIMIZATION RESULTS.................................................................................... 91 7.1 Comparison of propeller optimizations with and without viscous and compressibility effects taken into account.................................................................... 92 7.1.1 Single point Optimization for Cruise Condition: Binary encoded GA..... 92 7.1.2 Single point Optimization for Cruise Condition: Real encoded GA ...... 100 7.2 Single point Optimization Launch Condition................................................. 107 7.3 Dual point Optimization Cruise and Launch Condition ................................. 111 7.4 Single Point Optimization in Cruise considering Propeller Noise.................. 116 7.5 Single/Dual Point Optimization in Cruise/Launch Condition considering Propeller Noise and Performance ............................................................................... 119 8 CONCLUSION & RECOMMENDATIONS......................................................... 129 REFERENCES ............................................................................................................... 132 APPENDIX A: BIOT-SAVART LAW.......................................................................... 142 APPENDIX B: PROPELLER GA VARIABLES .......................................................... 147 APPENDIX C: PROPELLER WAKE INVESTIGATION............................................ 149 APPENDIX D: PROPELLER PANEL DISCRETIZATION......................................... 150 APPENDIX E: NACA 109622 PROPELLER GEOMETRY ........................................ 151 APPENDIX F: NACA 0009 AIRFOIL DATA .............................................................. 154 APPENDIX G: PROPELLER GEOMETRY OF SKYHAWK 172............................... 156 APPENDIX H: WIND TUNNEL DATA CAM 13x7 PROPELLER ............................ 160 xi LIST OF FIGURES Figure 2.1 Lifting-line model with single horseshoe vortex on finite wing ....................... 6 Figure 2.2 Distribution of horseshoe vortices, lifting-line model....................................... 7 Figure 2.3 Induced drag due to induced downwash of the trailing vortices....................... 8 Figure 2.4 Lifting surface method by horseshoe elements ............................................... 10 Figure 2.5 Vortex ring element arrangement.................................................................... 11 Figure 2.6 Propeller blade grid discretization................................................................... 13 Figure 2.7 Propeller wake geometry................................................................................. 14 Figure 2.8 Propeller wake horseshoe trailer discretization............................................... 16 Figure 2.9 Induced velocity at point P(x,y,z) by a single vortex element ........................ 18 Figure 2.10 Scanning sequence on propeller blade panels ............................................... 21 Figure 2.11 Normal vector on propeller blade panel ........................................................ 22 Figure 2.12 Propeller boundary condition and orientation angles.................................... 23 Figure 2.13 Vortex circulation sequence on a two bladed propeller ................................ 24 Figure 2.14 Trailing vortex segments responsible for induced downwash ...................... 26 Figure 2.15 Total thrust and drag of the propeller............................................................ 28 Figure 2.16 Class function 2D design space..................................................................... 32 Figure 2.17 Bernstein polynomial representation of the unit shape function................... 33 Figure 2.18 Class/Shape function single and dual surface geometries............................. 35 Figure 2.19 Propeller sweep functions based on CST method ......................................... 38 Figure 2.20 Thin propeller blade geometry by Class/Shape function scheme ................. 39 Figure 2.21 Comparison of lift coefficient for incompressible and compressible................ flow with direct application of the Prandtl-Glauert rule (Source ......................................... from Ref.[45] ).................................................................................................................. 40 Figure 2.22 Prandtl-Glauert domain stretching in flow direction..................................... 41 Figure 2.23 Prandtl-Glauert domain stretching of vortex ring elements .......................... 42 Figure 2.24 Spanwise lift coefficient normalized by blade radius.................................... 43 Figure 2.25 Generic two bladed propeller with wake....................................................... 44 Figure 2.26 Compressibility effects vs propeller tip Mach number ................................. 45 Figure 2.27 Spanwise lift coefficient normalized by blade radius [23] ............................... 48 xii Figure 2.28 Lift coefficient of NACA 0009 in post stall.................................................. 50 Figure 2.29 Drag coefficient of NACA 0009 in post stall................................................ 51 Figure 2.30 Post stall thrust behavior of a NACA109622 propeller with empirical ............ stall model......................................................................................................................... 52 Figure 2.31 Post stall power behavior of a NACA109622 propeller with empirical ........... stall model......................................................................................................................... 52 Figure 2.32 Locations of noise sources............................................................................. 54 Figure 2.33 Acoustic signature analysis geometry ........................................................... 57 Figure 2.34 Linearized airfoil thickness distribution of a APC 14x7 propeller................ 59 Figure 2.35 1C160 quarter scale Cessna 172 propeller .................................................... 60 Figure 2.36 Pressure perturbation comparison experiment and theory from........................ Miller [5] ............................................................................................................................. 61 Figure 2.37 Pressure perturbation of 1C160 quarter scale propeller ................................ 62 Figure 3.1 Flow chart propeller aerodynamic performance program............................... 65 Figure 4.1 Summation of LHS and RHS difference after matrix sweeps......................... 68 Figure 4.2 Accuracy of the solution after matrix sweeps by Gauss Seidel solver............ 69 Figure 4.2-4 Wake vortex sheet extensions for ?, 2 and 5 rotations downstream ........... 71 Figure 4.5 Propeller thrust variation due to wake vortex sheet extension........................ 72 Figure 4.6 Propeller Discretization schemes .................................................................... 73 Figure 4.7 Lift force convergence by discretization scheme and panel number .............. 74 Figure 4.8 Wing and wake discretization ......................................................................... 76 Figure 4.9 Lifting line and lifting surface method cl and cd of flat plate wing AR=6 ..... 77 Figure 4.10 NACA109622 propeller geometry with wake............................................... 79 Figure 4.11 Propeller efficiency comparison at different advance ratios ......................... 80 Figure 4.12 Propeller thrust coefficient at various advance ratios.................................... 81 Figure 4.13 Propeller power coefficients at various advance ratios ................................. 81 Figure 6.1.Flow chart propeller performance program implementation into GA............. 89 Figure 7.1 Optimized cruise propeller geometry from binary encoded GA; inviscid, ......... incompressible .................................................................................................................. 95 Figure 7.2 Optimized cruise propeller geometry from binary encoded GA; inviscid, ......... compressible ..................................................................................................................... 96 Figure 7.3 Optimized cruise propeller geometry from binary encoded GA; viscous,.......... incompressible .................................................................................................................. 97 Figure 7.4 Optimized propeller spanwise circulation distributions for cruise.................. 98 Figure 7.5 Objective function convergences binary encoded propeller ........................... 99 xiii Figure 7.6 Objective function convergence inviscid, incompressible, for binary ................ encoded GA .................................................................................................................... 100 Figure 7.7 Optimized cruise propeller geometry from real encoded GA; inviscid, ............. incompressible ................................................................................................................ 102 Figure 7.8 Optimized cruise propeller geometry from real encoded GA; inviscid, ............. compressible ................................................................................................................... 103 Figure 7.9 Optimized cruise propeller geometry from real encoded GA; viscous, .............. incompressible ................................................................................................................ 104 Figure 7.10 Optimized circulation distributions for cruise propeller, real encoded ............. GA ................................................................................................................................. 105 Figure 7.11 Objective function convergences real encoded GA, cruise condition......... 106 Figure 7.12 CAM 13x7 propeller efficiency at 19 [m/s] free stream velocity ............... 107 Figure 7.13 Optimized launch propeller geometry: viscous, incompressible case......... 109 Figure 7.14 Optimized propeller spanwise circulation distribution viscous,........................ incompressible case ........................................................................................................ 110 Figure 7.15 Objective function convergence from real encoded GA: viscous,.................... incompressible case ........................................................................................................ 111 Figure 7.16 Optimized cruise and launch propeller; viscous, incompressible case........ 114 Figure 7.17 Objective function convergences real encoded GA, cruise condition......... 115 Figure 7.18 Objective function convergences real encoded GA, cruise condition......... 116 Figure 7.19 Spanwise circulation distribution of the optimized propeller in cruise............. condition ......................................................................................................................... 118 Figure 7.20 Optimized cruise propeller from real coded GA; viscous, incompressible....... with noise considered...................................................................................................... 119 Figure 7.21 Optimized cruise propeller; viscous, incompressible case with noise and........ performance considered.................................................................................................. 122 Figure 7.22 Optimized launch propeller; viscous, incompressible case with noise and....... performance considered.................................................................................................. 123 Figure 7.23 Optimized launch / cruise propeller; viscous, incompressible case with .......... noise and performance considered.................................................................................. 124 Figure 7.24 Spanwise circulation distribution of the optimized propeller in cruise....... 125 Figure 7.25 Spanwise circulation distribution of the optimized propeller at launch...... 126 Figure 7.26 Spanwise circulation distribution of the optimized propeller in cruise............. and launch condition....................................................................................................... 127 Figure 7.27 Objective function convergences, cruise, launch and combined case......... 128 xiv Figure A.1 Velocity at point P due to a vortex distribution............................................ 144 Figure A.2 Induced velocity at point P by a vortex segment.......................................... 145 xv LIST OF TABLES Table 7.1 Fixed input parameters cruise condition........................................................... 92 Table 7.2 Optimized propellers performance parameter results cruise condition ............ 94 Table 7.3 Optimized propellers performance parameter results real encoded GA......... 101 Table 7.3 Fixed input parameters for launch condition .................................................. 108 Table 7.4 Optimized propellers performance parameter results for launch condition ... 108 Table 7.5 Fixed input parameters for dual point optimization........................................ 112 Table 7.6 Optimized propeller performance parameters: cruise and launch .................. 113 Table 7.7 Optimized propellers performance parameters cruise condition with noise......... considered ....................................................................................................................... 117 Table 7.8 Propellers performance parameter results with noise optimization................ 121 xvi NOMENCLATURE Roman Letters A 1-10 Bernstein polynomial coefficients, 2-D airfoil section [-] A n Shape function coefficients [-] zyx AAA ,, Influence coefficient matrices [1/m] AOA Angle of attack function [?] aoa n AOA function coefficients [-] Au, Al Coefficient functions for shape function [-] au n , al n Shape function coefficients [-] b Wing span [m] b Unit propeller radius [-] B Vector with right hand side of boundary condition [-] zyx BBB ,, Influence coefficients horseshoe trailers [1/m] c Section propeller chord [m] c d Drag coefficient [-] c d_max maximum drag coefficient [-] xvii c f Skin friction coefficient [-] c Maximum chord [m] 1 2 N N C Class function with coefficients [-] Chord Propeller chord function [m] chord n Propeller chord function coefficients [m] c l Section lift coefficient [-] c p Propeller power coefficient [-] c t Propeller thrust coefficient [-] d Propeller diameter [m] D Total induce drag [N] D i Induced rag [N] dl length of the vortex element [m] D visc Viscous drag [N] ?D Panel drag force due to induced downwash [N] ?F Panel force perpendicular to the onset velocity [N] ?Q Panel torque [Nm] ?S Panel width, spanwise [m] ?T Panel thrust [N] ?V Velocity induced by single vortex element [m/s] ?y Panel width in the spanwise direction [m] F ? Propeller section drag force [N] F R Propeller section radial force [N] xviii F T Propeller section thrust force [N] h Rankine body diameter [m] horeshoe_x,y,z x, y,z coordinates of the horseshoe elements [m] i, j, k Program counters [-] inf_coeff x,y,z Influence coefficients in x, y and z direction [1/m] L Lift [N] M Mach number [-] m Wake gradient [-] M r Mach number with respect to observer [-] M s Mach number relative to propeller [-] n Propeller revolutions [rad/sec] n? Surface normal [-] N 1 , N 2 Class function exponents [-] nodal_x,y,z x, y,z coordinates of the propeller nodal elements [-] nu n , nl n Class function coefficients [-] nxbl Number of propeller blades [-] P Propelr power [W] p-p 0 Presure perturbation [Pa] P visc Power due to viscous torque [W] Q Propeller torque [Nm] qm Mas flux [kg/s] Q visc Torque due to viscous drag [Nm] r Distance collocation point to vortex element; observer [m] xix Re c Section chord Reynolds number [-] r i,j Distance of center of vortex ring panel to origin [m] S Surface area [m 2 ] S r,n Bernstein polynomial term [-] Sweep Propeller sweep function [m] sweep n Propeller sweep function coefficients [-] T Propelr thrust [N] t Airfoil thickness [m] V ? Free Stream velocity [m/s] V local Local velocity [m/s] V onset Product of free stream velocity and propeller RPM [m/s] w ind Induced downwash velocity [m/s] x, y, z Coordinate system x, y and z axis [-] x Vector with circulations ? [m 2 /s] Greek Letters ? Flow angle / Angle of source with y axis [?] ? onset Angle of V onset with x-axis [?] ? i Induced downwash angle [?] ? Wake angle [?] ? Circulation [m 2 /s] xx ? Onset flow angle with x axis [?] ? Angle of panel normal with z-axis [?] ? Source/sink distance [m] ? Propeller efficiency [-] ? Angles in the x,y-plane [?] ? Propeller panel dihedral angle [?] ? Viscosity coefficient [-] ? Density of air [kg/m 3 ] ? Retarde time [s] ? Angle of observer with x axis [?] ? Fraction of local chord x/c [-] ? Angular velocity [1/sec] Mathematical Expressions ? Increment [-] t? ? First order partial derivative in time [-] Abbreviations AR Aspect ratio BP Bernstein polynomials xxi CST Class/Shape function transformation GA Genetic Algorithm LHS Left hand side NACA National Advisory Committee for Aeronautics PSO Particle Swarm Optimization RHS Right hand side RPV Remotely piloted vehicles UAS Unmanned aerial system VTL Vortex Lattice Method 1 1 INTRODUCTION Even though aircraft propellers have been designed for over a century, investigations of the performance and the design of propellers are as important as ever. The continued interest in propeller propulsion is due to the fact that propellers are more efficient at low and modest flight speeds when compared to turbofan propulsion. A relatively new field in propeller design is the prediction and reduction of propeller noise, which is driven by the public living close to airports or large wind turbines. Aircraft propeller noise is also of special interest to the military and was first investigated in 1919 by Lynam and Webb. [1] This military interest in minimizing noise continues today with the recent fielding of low flying Unmanned Aerial Systems (UAS) on reconnaissance missions. This work describes an effort to optimize the propeller performance of small thin airfoil propellers, which are found on electric powered UAS?s and Remotely Piloted Vehicles (RPV). Previous work on 2D airfoil optimization using a Genetic Algorithm (GA) has been investigated by Refs. [2]- [4]. Fanjoy and Crossley [3] generated 2D airfoil geometries using a B-spline method with 80 panels distributed over the airfoil surface and used a panel method to do the aerodynamic analysis. Penalty functions were added to provide a minimum thickness thus structural integrity is guaranteed. The two Genetic Algorithm?s used are a binary tournament selection with uniform crossover and a real encoded GA. The shortcoming is following: If there are cases for which the 2 analysis over predicts the desirable features of an airfoil, the population will move toward those cases. A 3D propeller analysis and optimization was done by Miller [5] . In his approach, the Vortex Lattice Method (VLM) using a single panel in the chordwise direction was employed to determine the aerodynamic propeller performance. The propeller geometries were generated using an extended version of the classical method in which the blade shank is the center of rotation. Airfoil camber is not considered due to the limitation of using only a single chordwise panel. The work includes a subsonic noise analysis using sectional forces located at the quarter chord to represent the blade loading and point mass sources and sinks to represent blade thickness. The multivariable optimization method used was developed by Powell [6] and is a conjugate directions technique. The propeller optimization is limited because of constraints of the propeller geometry functions as well as the constraints of the lifting line method to account for chordwise geometric shapes. Olsen [7] described ship propeller optimizations using the VLM where the propeller blade is replaced by a lattice of quadrilateral panels with constant circulation and the horseshoe vortices follow regular helices. The propeller surface geometry is defined by a vector which contains the rake, the skew of the midchord and the horizontal pitch angle. The propeller optimization solves a variational problem to find the radial circulation distribution which provides the lowest torque for a given thrust, thus generating propellers with minimum losses. The propeller performance computation did not include viscous effects and the propeller geometries were generated using simple cosine and sine functions. 3 Hampsey [1] investigated the optimization of small wind turbines employing a three dimensional panel method with distributions of source and doublet singularities over the discretised blade surface. Airfoil geometries were generated by least squares B-spline curve fitting through specified data points. The spanwise geometry was defined by stringing a handful of scaled and rotated airfoils along the blade span and then skinning it with a similar B-spline surface in the radial direction. The optimization scheme used was a differential evolution algorithm which allows for multi-objective optimization. In addition to the aerodynamic performance optimization, Hampsey included a wind turbine blade structure model. The main focus of the present work was: ? Develop a propeller program which accurately computes the propeller aerodynamic performance parameters while being computationally efficient; ? Implementation and development of a new way to generate propeller geometries using the Class/Shape function Transformation (CST) methodology; ? Increase computational accuracy by including viscosity, compressibility and flow separation models; ? Implementation of an acoustic signature model. The developed propeller performance program is attached to both a binary encoded GA and a real coded GA. A fundamentally new method to generate airfoil and propeller geometries based on a Class/Shape function methodology is utilized. The VLM using lifting surfaces with 4 vortex ring elements has been used to determine the aerodynamic performance parameters. In an effort to reduce computation times, a detailed wake extension convergence model has been developed and a panel density investigation has been completed. The focus of this work is the optimization of propellers with multiple objectives, which has not been researched extensively. Finally the implementation of a validated acoustic propeller signature model in combination with the VLM performance analysis and the optimization scheme using a GA has been developed. This investigation can be split up into four major areas: propeller geometry determination, aerodynamic analysis, noise analysis and an optimization scheme. In each Chapter, the background in the field begins with a discussion of the basic issues, followed by a detailed description of the analysis method used, and concluded by the presentation of the work done to verify the method. Optimization results include single and multipoint performance and acoustic signature optimizations with performance trade-off?s. 5 2 BASIC THEORY The first and simplest method to predict propeller performance was developed by Rankine [9] . The basic version is called axial momentum theory, in which the propeller is replaced with a disc and the thrust is considered uniformly distributed. Even though the theory was extended to account for swirl, friction losses [10] and compressibility [11] , it still ignores many real world effects, which makes it useful only for predicting the upper limits of propeller efficiency. The blade element theory was formulated by Froude [12] and Drzewiecki. [13] Anderson [14] said ?Wilbur Wright was the first person to recognize ? that a propeller is nothing more than a twisted wing.? Using information from their wind tunnel tests, Wilbur Wright had invented ?blade element theory,? which is the idea that at each point along the span, a propeller meets the air at a different angle and speed. Blade element theory divides the propeller blade into elements, each of which is analyzed separately by an analogy drawn between propeller and finite wing. The method takes axial and associated rotational velocities at each spanwise section into account and uses 2D airfoil coefficients to compute local and overall propeller performance. Even though this method is more accurate than the axial momentum theory it does not consider the influence of the propeller wake. Glauert [15] combined the momentum theory with the blade element theory and assumed that the propeller is lightly loaded, hence slip stream contractions are small 6 and the radial components of the induced velocities are negligible. The Vortex Theory as discussed above does not include tip losses. An improvement can be made using lifting line theory which is discussed in the following section. 2.1.1 Lifting-Line Theory The lifting-line theory was the first practical theory for incompressible, inviscid flow to predict aerodynamic properties of unswept three dimensional wings. [16] The theory was developed by Ludwig Prandtl and colleagues [17] between 1911 and 1918 and is essentially an implicit potential flow solution. The wings are modeled as a single bound vortex line, which is located at the quarter chord position and has a shed vortex sheet extending to infinity to satisfy the Helmholtz theorem, [16] shown in Figure 2.1 below. Figure 2.1 Lifting-line model with single horseshoe vortex on finite wing ? V ? ? c b ? ? ? . Collocation point 7 To realistically simulate the downwash distribution of a wing, a large number of horseshoe vortices are superimposed, each with a different bound vortex length located at the quarter chord. Solution methods to the lifting line theory are described in numerous texts including Katz & Plotkin [16] , Kuethe and Chow [18] and Anderson. [14] Figure 2.2 Distribution of horseshoe vortices, lifting-line model Figure 2.2 shows the geometric arrangement of the horseshoe elements over the wing span with circulation ? changing along the lifting-line. The unknown circulation ? of the individual horseshoe vortices is obtained using the Biot-Savart Law [16] with the boundary condition, 0? =?nV r (zero flow normal to the surface) at a pre-defined collocation point (also referred to as a control point). The resulting form of the Biot-Savart Law is: 3 4 r rdl dV ? ? ? = ? (2.1) 1 ? 1 ? 2 ? 2 ? 3 ? 3 ? -b/2 b/2 x y z 1 ? 21 ?+? 321 ?+?+? Lifting line 8 The local lift is obtained from the Kutta-Joukowski Theorem [18] , which states that lift, L per unit span, is directly proportional to the circulation ?, and is always perpendicular to the local freestream velocity. ???= ? VL ? (2.2) The total lift is the summation of the local lift, L, which is perpendicular to the free stream velocity. An illustration of the velocities and related forces is shown in Figure 2.3. To determine the induced drag, the resultant velocity, V, which has two components, the free stream and the induced velocity, w i , must be determined. The induced angle, ? i , is the angle between the free stream velocity and the resultant velocity. With the induced angle of attack known, the induced drag is obtained as shown in Figure 2.3. Figure 2.3 Induced drag due to induced downwash of the trailing vortices The induced downwash velocity, w i , at any span location can be determined once the circulation ? is known. This allows for the computation of the local flow angles on V ? x y V ? w i ? i ???= ? VL ? ( ) ii VD ?? sin?????= ? ???= VL ? V 9 the entire wing. A consequence of the downwash flow is that the direction of action of each section?s lift vector is rotated relative to the freestream direction. The local lift vectors are rotated backwards and hence give rise to a lift induced drag. By integrating the component of section lift that acts parallel to the freestream across the span b, the induced drag D i can be found as () ? ? ????= 2/ 2/ sin b b ii dyVD ?? (2.3) 2.1.2 Lifting Surface Theory Due to the limitations of predicting accurately aerodynamic performance for only straight medium to high aspect ratio wings, the Lifting-Line Method was extended by placing a series of lifting vortex lines along different wing chord stations. [2], [14] Each lifting line has two associated trailing vortices which are parallel to the x axis. As the location for lift calculations shifts downstream in the chordwise direction, starting at the leading edge, additional superimposed trailing vortices are added. Thus vortex strength dependency is gained in both the x and y direction. The vortex sheets cover both, the spanwise and the chordwise direction and make up the lifting surface. Figure 2.4, a wing with two chordwise vortex elements is shown. 10 Figure 2.4 Lifting surface method by horseshoe elements As with the lifting-line method, the unknown circulations ? must be determined based on the flow-tangency condition, ( )0? =?nV r , at all points on the wing, so that lift and induced drag can be calculated. The central problem with lifting surface theory is to solve for all the unknown circulations since circulation strengths are mutually dependent. One numerical method is the Vortex Lattice Method (VLM) which was developed to compute aerodynamic performance data of various wing geometries. Numerous papers [19]- [22] have been published which demonstrate the accuracy and applicability of the VLM. In this method, the wing is covered by a net of bound and horseshoe elements which is called a vortex lattice system. When the Biot-Savart law and the flow tangency condition is applied at all collocation points, a system of simultaneous equations is obtained which then can be solved for the unknown circulation strengths, ?. x y z n r 11 Lifting Surface Method by Vortex Ring or Quadrilateral Elements In this method, vortex ring elements are distributed over the wing surface instead of the horseshoe elements as previously shown. This simplifies the vortex model as well as reduces computation times due to the lower number of vortex filaments. The boundary condition on each ring vortex is the same as before: no flow through the wing at a control point. Since all vorticity based solutions are essentially implicit potential flow solutions, flow separation, viscous and compressibility effects cannot be simulated and must be taken into account separately, as discussed in Chapters 2.5 ? 2.7. A single ring vortex consists of four straight bound vortex elements, which are arranged in a geometric closed form. Positive circulation is defined according to the right-hand rotation rule, which is illustrated in Figure 2.5. Figure 2.5 Vortex ring element arrangement In previous efforts [5] the lifting surface solution by horseshoe elements was used, which is computationally more expensive but yields similar results when compared to vortex ring elements. The vortex geometry method used in this effort has been applied to ship propellers [7] , aircraft propellers [23] and yacht sails [24] and is proven to give accurate Element 2 Element 4 Element 1 ? Element 3 12 aerodynamic performance predictions. In the present work, a single propeller blade is divided up into chordwise and spanwise panels. The corners of the panels are defined by nodal points and the geometric center by a collocation point. For thin propeller airfoils, a single layer of panels is placed on the mean camber line. For thick airfoil propeller blades (>10%), the upper and lower blade surface is covered with a single layer of vortex elements. Vortex ring elements are placed on the propeller surface starting at the blade leading edge and extending to the second last element in the chordwise direction. The last chordwise panel elements are horseshoe elements with a spanwise bound vortex and two trailers. The horseshoe trailers extend downstream to infinity following the pitch of the last chordwise panel. Figure 2.6 illustrates the discretization of a propeller blade by vortex elements. 13 Figure 2.6 Propeller blade grid discretization The aerodynamic performance of the propeller is obtained, as described for the lifting surface method with horseshoe elements by applying the boundary condition to each panel that the flow normal to the panel surface is zero. For a detailed description of the numerical computation of the aerodynamic performance of the propeller is given in Section 2.3. 2.2 Propeller Wake Generally there are two types of wake models, the prescribed wake models [23], [25], [26] , which define the wake by a function and the free wake methods [24], [27], [28] , which adjust the wake structure to account for wingtip roll up. Blade trailing edge Horseshoe element Ring element ? 1 Collocation point Nodal point ? 11 ? 2 ? 3 ? 12 ? 5 ? 8 ? 10 ? 7 ? 4 ? 9 ? 6 14 In this effort, the propeller wake geometry follows Goldstein?s [25] helical model, which assumes lightly loaded propellers. Thus the interference flow of the vortex system is small compared to the velocities of the blade. The pitch of the wake is constant and there is no slipstream contraction. In Ref. [23], it is shown that aerodynamic propeller performance can be accurately predicted using constant pitch wakes. To consider wing tip rollup, a free wake model must be chosen. Most free wake models use an iterative process to determine the exact wake geometry. The free wake models are not considered within this optimization effort due to the increased computation time and minimal gains in accuracy. The wake of the propeller is simulated by horseshoe vortices which extend from the last chordwise panel extending two revolutions downstream. More details on propeller wake extensions are discussed in Section 4.2. The wake gradient follows the trailing edge panel slope to closely approximate the Kutta condition. Figure 2.7 Propeller wake geometry z x Horseshoe trailer sare divided up by 10 o rotation elements ? First streamwise horseshoe trailer Collocation Points Nodal Points 15 The trailers are divided up into 10 degree rotational increments which are implemented in Equations 2.4-2. The y and z coordinates of the horseshoe trailers are computed similarly but offset by 90 degrees. The x coordinates are determined by the gradient of the last chordwise vortex elements multiplied by the distance traveled during a 10 degree rotation. Once the gradient, m, of the last chordwise panels is obtained and the angles ? of each last chordwise nodal point with the y-axis are determined, the x, y and z-coordinates of the horseshoe vortex elements can be written as: m i jncynodaljncxnodaljixhorseshoe ? ? ???+= 360 )(10 ),(_2),(_),(_ ? (2.4) ) 180 ),()(10 cos()),(_),(_(),(_ 5.022 jnci jncznodaljncynodaljiyhorseshoe ?+? ?+= (2.5) ) 180 ),()(10 sin()),(_),(_(),(_ 5.022 jnci jncznodaljncynodaljizhorseshoe ?+? ?+= (2.6) with m being the gradient of the last streamwise propeller panel ),1(_),(_ ),1(_),(_ jncznodaljncznodal jncxnodaljncxnodal m ?? ?? = (2.7) and ? the angle in the x-y-plane associated with the last chordwise nodal point. ? 180 )),(_),(_( ),(_ sin),( 5.022 ? ? ? ? ? ? ? =? jncznodaljncznodal jncznodal ajnc (2.8) 16 The counter i runs to 72 to generate two full wake rotations and j from 1 to nr which are the spanwise number of nodal points. Figure 2.8 Propeller wake horseshoe trailer discretization Initially, a single propeller blade including horseshoe trailers is generated. The geometry is then rotated about the origin depending on the required number of propeller blades. 2.3 Propeller Performance Prediction by VLM For the computation of the propeller loading, three different methods are available. The methods differ in accuracy and computation time. One method is the Trefftz plane analysis [29] in which it is assumed that the induced flow far downstream from the propeller does not depend on the streamwise coordinate and thus is considered to be two dimensional in all planes normal to the trailing vortex sheet. The induced velocities of the X Y Z 10 degree rotation increments ? 17 bound vortex elements are neglected since they are determined far downstream in a cross plane normal to the wake. Thus, the analysis is essentially reduced to finding the flow field of a two dimensional vortex sheet to satisfy certain boundary conditions. The second method is the pressure summation method [30] which uses the velocity difference on the upper and lower propeller surface to obtain pressure differences. The total force acting on the propeller is obtained by surface integration. The third method is the Joukowski method which determines the forces on each discrete vortex line in the lattice and through summation obtains the total propeller loading. There are different ways of determining the induced velocity which are described in Refs. [31]- [33]. A comparison of the three computational methods with respect to accuracy and computational effort, applied to yacht sails is found in Ref. [24]. The analysis concludes that the Joukowski method has good accuracy and provides the force distributions in both the spanwise and the chordwise directions. The Joukowski method was also used by Cheung and generated good accuracy in the performance prediction of propellers. Based on this analysis, the Joukowski method was selected to determine propeller performance parameters. Lifting-Surface Solution by Vortex Ring Elements: The Propeller blade is divided into chordwise and spanwise vortex ring elements. The ring elements cover the entire blade except the last chordwise panels, which are horseshoe elements. The horseshoe trailers leave the propeller blade with the same pitch as the last chordwise panels as discussed in Section 2.1. The velocity induced by a vortex element of the strength ? with a length of dl on a collocation point P(x,y,z) is calculated by: 18 ( ) 3 10 10 4 rr rrdl V ? ?? ? ? ? =? ? (2.9) with ( ) 10 rr ? being the vector from the collocation point to the midpoint of the vortex element. The circulation ? and the induced velocities are initially unknown and the last part of Equation 2.9 contains the influence coefficients which are obtained with: ( ) 3 10 10 4 1 inf_ rr rrdl coeff ? ?? ? ? = ? (2.10) The influence coefficients describe the geometric relation between the length of a vortex element and the point at which the induced velocity is to be determined. Figure 2.9 illustrates the geometric relations of a vortex element and the point of velocity evaluations. Figure 2.9 Induced velocity at point P(x,y,z) by a single vortex element dl ds ? P(x,y,z) r 0 r 1 Origin r 0 ? r 1 ? ?V Vortex element 19 For a single quadrilateral vortex element the total influence coefficients in the x, y and z direction are the summations of the influence coefficients in x, y and z direction of the four straight bound vortex legs. Cheung [23] determines them by ()()()()[] ()()()() ()()()()Ixxyyyyxxcoeff Izzxxxxzzcoeff Iyyzzzzyycoeff startendcolocationstartstartendcolocationstart i z startendcolocationstartstartendcolocationstart i y startendcolocationstartstartendcolocationstart i x ????????= ????????= ????????= ? ? ? = = = 4 1 4 1 4 1 inf_ inf_ inf_ (2.11) with () ? ? ? ? ? ? ? ? ? +?+ + ? ?? = c b cba ba bca I 2 1 2 (2.12) and ()()( ) ()( )()( )()( ) ()() 222 222 ncollocatiostartncollocatiostartncollocatiostart ncollocatiostartstartendncollocatiostartstartendncollocatiostartstartend startendstartendstartend zzyyxxc zzzzyyyyxxxxb zzyyxxa ?+?+?= ???+???+???= ?+?+?= (2.13) with I, a, b, c being variables which describing the geometric relations of a vortex element with respect to a collocation point. The same applies to the horseshoe elements, with the difference that the summation of the influence coefficients consists of one bound vortex element and two trailers which extend two revolutions downstream and which are divided into 10? sections. It is of importance to follow the right hand rule when defining start and end 20 point of a single vortex element. The direction of the circulation within the panels cannot change. Also once a direction of the circulation is selected, all following panel circulations must be in the same rotational direction. The sequence of calculating the influence coefficients with respect to one collocation point is as follows: The first collocation point on the first panel of the first propeller blade is selected and the influence coefficients in the x, y and z direction of the first vortex ring element is computed. The obtained x, y and z values of the influence coefficients are then stored in the first location of the matrices A x 1,1 , A y 1,1 , and A z1,1 . The next step is to keep the first collocation point and determine the influence coefficients of the second vortex ring element, which is the second panel in the spanwise direction of the propeller blade. Again the influence coefficients are stored in the matrices at the location A x 1,2 , A y 1,2 , and A z 1,2 . This process is repeated until all influence coefficients in x, y and z direction are determined with respect to all vortex ring elements and horseshoe elements of all propeller blades. Then the process is repeated with the collocation point of the second vortex ring element selected. The individual influence coefficients are determined and stored in the matrices A x 2,n , A y 2,n , and A z 2,n . The computations of all the influence coefficients ends with the last collocation point on the last propeller blade and the last horseshoe panel on the last propeller blade. Figure 2.10 the scanning sequence of the propeller panels is shown. 21 Figure 2.10 Scanning sequence on propeller blade panels All influence coefficients in x, y and z directions are stored in the matrices A x , A y , and A z with the rows representing the dependencies with respect to a single collocation point and the columns representing the dependencies of a single panel with respect to all collocation points of all propeller blades. Thus the Biot-Savart law becomes: ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ???? ? ? ? ? =? MNMNMNMNMN MN MN AAA AAA AAA V 2 1 ,2,1, ,22,21,2 ,12,11,1 4 1 ? (2.14) To calculate the unknown circulations, ?, the boundary condition that the velocity at each collocation point normal to the individual panel is zero, or: ? x z i i+1 i+2 j j+1 j+2 j+3 M N Counter k=M x N Tip Root 22 0? =?nV (2.15) The surface normal of each quadrilateral and horseshoe element is determined as the cross product of the diagonal vectors of each panel and is given as: kk kk k BA BA n ? ? = (2.16) Figure 2.11 Normal vector on propeller blade panel The boundary condition applied to each collocation point is: []()() kkkkonsetzijzijyijyijxijxij M i k N j VnAnAnA ??? cossin 11 ???=?+?+??? ?? == (2.17) with ? k being the angle of attack of the onset velocity V onset_k with the x- axis, ? k the angle of each panel normal with the z-axis and ? the orientation of the panels. Figure 2.12 illustrates the geometric relations in the propeller coordinate system. i i+1 i+2 j j+1 j+2 A k B k n k 23 Figure 2.12 Propeller boundary condition and orientation angles The onset velocity V onset changes with each spanwise panel, since it is a combination of the free stream velocity and the rotational speed of the propeller. ()( ) 2 1 2 , 2 , jionset rVV ji ?+= ? ? (2.18) Equation 2.14 is a system of linearly dependent equations which are solved for the unknown circulations ? with a Gauss Seidel solver. Section 4.1 describes in detail the Gauss Seidel Solver. [ ][ ][ ]BxA =? (2.19) with A being the influence coefficients [ ] zijzijyijyijxijxij M i N j nAnAnAA ?+?+?= ?? ==11 (2.20) ? onset x z x z ? y z ? n n - + V onset 24 and x the unknown circulations: [] k x ????= ,...,, 321 (2.21) and B the boundary condition zero flow normal to the panel surface: ( ) ( ) ()() ()() ? ? ? ? ? ? ? ? ? ? ? ? ??? ??? ??? = kkkkonset onset onset V V V B ??? ??? ??? cossin ..... cossin cossin 2222 1111 (2.22) The returned vector x contains the circulations ? k of each propeller blade panel as shown in Figure 2.13 below. Figure 2.13 Vortex circulation sequence on a two bladed propeller With the returned circulations, the normal force of each quadrilateral and horseshoe element is obtained by the use of the Kutta-Joukowski theorem [16] . Kutta-Joukowski theorem states that the forces are always perpendicular to local onset flow, which is the combination of the free steam velocity and the local rotational speed of the propeller. The forces for the panels on the propeller blade, except for the forces on the leading edge panels are: ? ? 1 ? 2 ? 3 ?4 ?5 ? 6 ? 7 ? 8 25 ( ) jjijijionsetji yVF ???????=? ?? ,1,,_, ? (2.23) For the leading edge panel (i=1) the forces are: ( ) yVF jionsetji ?????=? ? ,, ? (2.24) For the computation of the thrust and required power of the propeller, the induced downwash at each collocation point must be determined. The influence coefficients from only the trailing vortex segments with respect to each collocation point are computed in a manner similar to the approach used for the quadrilateral panels. These coefficients are stored in the matrices B x , B y , and B z . Figure 2.14 shows a propeller blade with related horseshoe elements. 26 Figure 2.14 Trailing vortex segments responsible for induced downwash Since the circulations ? k are known from the previous calculations, the induced downwash on each collocation point is determined by: ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? ? ? kkkkk k k indk ind ind BBB BBB BBB w w w 2 1 ,2,1, ,23,21,2 ,12,11,1 2 1 ... ............ ... ... ... (2.25) The circulations used for the computations of the induced downwash are obtained from the Gauss Seidel solution. Trailing vortex segments Propeller geometry 10? horseshoe wake segments . . . . . . . . . . . . . . . . . . Collocation point 27 With the induced downwash, the induced drag on each panel is determined with: () 1, ,,1,,_,_ >????????=? ? iywD jijijijiindjiind ? (2.26) () 1, ,,,_,_ =??????=? iywD jijijiindjiind ? (2.27) The total induced drag is the summation of all individual panel contributions of one propeller blade multiplied by the number of the propeller blades. nxblDD N j jiind M i ? ? ? ? ? ? ? ?= ?? == 1 ,_ 1 (2.28) To compute the thrust and the required power of the propeller the forces and the induced drag forces on each panel must be rotated in the propeller coordinate system as illustrated in Figure 2.15. 28 Figure 2.15 Total thrust and drag of the propeller The angle of rotation is defined by the onset flow angle ? j which changes in the spanwise direction. Thus the propeller thrust contributions of the individual panels are: ( ) ( ) jjiindjjiji DFT ?? cossin ,_,, ?????=? (2.29) V ? x z ? ind ??r i w ind ?F ? ind V onset V local ?D ind ?T = ?F sin(?) - ?D ind ?cos(?) ?D = ?F cos(?) + ?D ind sin(?) ? ?D ind = -?V local ?sin(? ind ) ?F = ?V local ?cos(? ind ) = ?V onset ? Rotor shaft ? 29 and the drag contributions are: ( ) ( ) jjiindjjiji DFD ?? sincos ,_,, ??+??=? (2.30) The total thrust of the propeller is obtained by the taking the summation of the individual panel contributions multiplied by the number of the propeller blades: nxblTT N j ji M i total ? ? ? ? ? ? ? ?= ?? == 1 , 1 (2.31) To compute the required propeller power, the torque ?Q i,j generated by the individual drag component of the panels, is multiplied by the distance from the collocation point to the axis of rotation. The summation of the torque contributions of the individual panels multiplied by the number of the propeller blades nxbl will determine the total propeller torque Q total . ` nxblQQ rDQ N j ji M i total jijiji ? ? ? ? ? ? ? ?= ??=? ?? == 1 , 1 ,,, (2.32) The required propeller power is the total torque multiplied by the rotational speed ? of the propeller. ??= total QP (2.33) 30 With the propeller thrust and power known, the thrust and power coefficients are: 42 dn T c total t ?? = ? ? (2.34) 53 dn P c p ?? = ? ? (2.35) with n being the revolutions in rad/sec and d the propeller diameter. The propeller section lift coefficient, c l , and the propeller efficiency are other parameters which are commonly used in the performance analysis. () cV rc il ? ?? = ? 2 (2.36) P TV total ? = ? ? (2.37) 2.4 Propeller Geometry In general the geometric definition of propellers can be reduced to the airfoil type, the chord length, the sweep and the angle of attack function. In this effort the CST (Class/Shape Function Transformation) method for defining general aerodynamic shapes has been applied to the generation of random propeller shapes which is discussed as follows. 2.4.1 2-D Airfoil Geometry Several methods have been developed to represent general airfoil shapes for the use in aerodynamic design optimizations. The goal is to define a simple analytic function 31 which efficiently describes the entire design space for airfoils. Some of the methods used in previous efforts can be found in Ref. [34]- [39]. The method chosen is based on the work of Ref. [40] and Ref. [41]. In this approach a simple and well behaved analytic unit shape function, based on Bernstein polynomials (BP) is introduced. This shape function directly controls key airfoil parameters including leading edge radius, thickness and trailing edge angle. A unit class function is added to the unit shape function to allow for a wider variety of general body shapes. The streamwise class function in the design space is defined as ( ) ( ) [ ] 211 2 1 NNN N C ??? ??= (2.38) with ? being the fraction of the local chord. In the physical space the unit class function is 21 1 NN c x c x c x C ? ? ? ? ? ? ?? ? ? ? ? ? ? = ? ? ? ? ? ? (2.39) where c is the chord and x is the local chord fraction having a range from 0-1. The first term of Equation 2.39 defines the shape of the airfoil leading edge and the second term can be used to ensure a sharp trailing edge. If N 1 = 0.5 and N 2 = 1.0 a round airfoil leading edge and a sharp trailing edge, is enforced. Figure 2.16 shows some examples for different values of the class function coefficients N 1 and N 2 . Further geometry class determinations due to N 1 and N 2 variation can be found in Ref. [40]. 32 Figure 2.16 Class function 2D design space The unit shape function is defined by a BP of the order n with the variable x ranging from 0-1.0. The BP?s were chosen due to the mathematical property of ?Partition of Unity? as described in Ref. [41]. The first term of the shape function defines the binomial coefficients with increasing order n of the BP. For each order n of BP there exist n+1 terms which are defined by Equation 2.40. () () ? = ? ? ? ? ? ? ? ?? ? = n r rn nr c x rnr n xS 0 , 1 !! ! (2.40) The individual terms of the BP?s can be illustrated by means of a Pascal?s triangle as shown in 2.17. N 1 =1.0 , N 2 =1.0 N 1 =0.5 , N 2 =1.0 N 1 =0.5 , N 2 =0.5 c x 33 Figure 2.17 Bernstein polynomial representation of the unit shape function The first term of Equation 2.40, when written in polynomial form as shown in Figure 2.17, defines the leading edge radius and the last term the trailing angle. The other terms are shaping terms which do not influence airfoil leading or trailing edge shape. The entire airfoil can be represented by one upper and one lower unit class function (if thickness is included) multiplied by a unit shape function with Bernstein polynomials. In Ref. [41] it was found that a BP of the order n of 6-9 matched the airfoil geometries and aerodynamic forces. It was further suggested that for optimization purposes the order n could be lowered to 4 which reduces the design variables and thus improves computation times. Based on the findings in Ref. [41] the order of the BP is set to n=4 for all further propeller optimizations. Equations 2.41 and 2.42 define the upper and lower airfoil geometry. The first two terms are the class function which sets the leading and trailing shape of the airfoil. The 1 x?1 x 2 )1( x? 3 )1( x? 4 )1( x? 5 )1( x? xx)1(2 ? 2 x xx 2 )1(3 ? 2 )1(3 xx? 3 x xx 3 )1(4 ? 22 )1(6 xx? 2 )1(4 xx? 4 x xx 4 )1(5 ? 23 )1(10 xx? 32 )1(10 xx? 4 )1(5 xx? 5 x n=0 n=1 n=2 n=3 n=4 n=5 Leading Edge Radius Trailing Edge Boattail Angle ?Partition of Unity? () 1 0 , = ? = n r nr xS ?Shaping Terms? which do not affect leading or trailing edge geometry 34 remaining terms are the BP for n=4 with the individual coefficients as shown in Figure 2.17. When thin airfoil propellers are considered only the upper airfoil geometry function is used to define the propeller shape. A single layer of vortex elements is placed onto the geometry function to represent the mean chord of the thin bladed propeller. Upper surface definition: () [ ] 4 5 23 4 22 3 3 2 4 1 114 1614 1 21 ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ? = c x A c x c x A c x c x A c x c x A c x A c x c x xy NN upper (2.41) Lower Surface definition: () [ ] 4 10 23 9 22 8 3 7 4 6 114 1614 1 43 ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ? ?= c x A c x c x A c x c x A c x c x A c x A c x c x xy NN lower (2.42) with c being the maximum chord, x the local variable ranging from 0 ? 1 and A 1-10 the coefficients (with a range of 0-1). This 2-D class/shape function requires 14 variables for a thick airfoil (upper and lower surface) and 7 variables for a thin propeller with a single 35 surface geometry function. Figure 2.18 illustrates two random airfoils from the design space of the class/shape function, one thick and one thin. Figure 2.18 Class/Shape function single and dual surface geometries 2.4.2 3-D Propeller Geometry The class/shape function methodology of representing 2D airfoils is extended here to define general 3D propeller shapes. This approach is based on the work done by Ref. [40] and [41] and then extended further to allow for a more open design space in the creation of general propeller geometries. To accomplish spanwise geometry variation the coefficients A 1-10 and N 1-4 for the top and the bottom airfoil side must be dependent on the local spanwise position. Thus for each of the dependent variables A 1-10 a BP function with n=3 describing the spanwise coefficient variation, is defined. For the coefficients of the Class function, N 1-4 , a BP function with n=1 is introduced, to provide smooth spanwise coefficient transitions. Thus the variables for the upper and lower airfoil geometry are: Class/Shape function single mean chord airfoil geometry Class/Shape function upper lower airfoil geometry 36 [ ] [ ] ? ? ? ? ? ? ?+ ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? ? ? ?+ ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ??? + ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ?= ? ? ? ? ? ? ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ??? + ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ?= ? ? ? ? ? ? + + ++ + ++ + b y nl b y nl b y Nu b y nu b y nu b y Nu b y al b y b y al b y b y al b y al b y Al b y au b y b y au b y b y au b y au b y Au jjj jjj ii iii ii iii 1 1 113 13 113 13 2 2 3 15 2 10 2 5 3 3 15 2 10 2 5 3 (2.43) b is the unit propeller radius, y is the local spanwise station, and i and j are the variable counters from 1-5 and 1-2. The total number of variables to describe the upper and lower airfoil geometry is 44. Twenty for each, the upper and lower airfoil shape function and 4 variables to define the two class functions. The geometry of a single surface is described by the unit class/unit shape function for the upper airfoil side as: () ( ) ( ) [ () () () () () ] 4 5 23 4 22 3 3 2 4 1 114 1614 1, 21 ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ? = c x yA c x c x yA c x c x yA c x c x yA c x yA c x c x yxy yNuyNu upper (2.44) An angle of attack (AOA), chord length and sweep variation function is added to the unit class/shape function to further extend the propeller design space. All three functions are defined in the same manner as the unit shape function variables A 1-5 , with 37 the addition of a multiplication factor to define the max value of the function. Each function requires 5 additional variables, four to define the BP with n=3 and one to set the max value of the function. The angle of attack, chord length and sweep function can be written as: () [ ] 3 4 2 3 2 2 3 15 113 13 ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ??? + ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ??= b y aoa b y b y aoa b y b y aoa b y aoaaoayAOA (2.45) () [ ] 3 4 2 3 2 2 3 15 113 13 ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ??? + ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ??= b y chord b y b y chord b y b y chord b y chordchordyChord (2.46) with aoa 5 and chord 5 being the max angle of attack and max chord length. () [ ] ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ??+ ? ? ? ? ? ? ? ? ? ? ? ? ? ??? + ? ? ? ? ? ? ? ? ? ? ? ? ? ???+ ? ? ? ? ? ? ??= rootatchord b y sweep b y b y sweep b y b y sweep b y sweepsweepySweep 3 4 2 3 2 2 3 15 113 13 (2.47) The variable, sweep 5 , defines the maximum sweep in multiples of the propeller root chord. Figure 2.19 illustrates how the sweep function influences the propeller shape. This brings the total number of variables to 66. Sixty three for the blade geometry, one for the propeller radius, one for number of propeller blades and one for the position, which the individual spanwise airfoil sections are rotated about. 38 Figure 2.19 Propeller sweep functions based on CST method In this effort, propeller optimizations for small electric powered UAS are investigated thus only thin airfoil propellers are considered. This eliminates one surface geometry function and reduces the number of variables to 40. An example of a single surface propeller blade with wake is shown in Figure 2.20. Propeller root Propeller tip Sweep backwards 39 Figure 2.20 Thin propeller blade geometry by Class/Shape function scheme 2.5 Compressibility Even though only small thin bladed propellers are considered in this effort, air density variations, ?, due to Mach numbers higher than M > 0.3 can have noticeable effects on the propeller performance. The general, Prandtl-Glauert compressibility correction uses the Mach number to correct pressure, lift and moment coefficients at Mach numbers up to M=0.7 for thin airfoil and at small angle of attack. The derivation of the compressibility correction is described in Ref. [14]. Gennaretti and Morino [42] showed that at higher Mach numbers the Prandtl-Glauert rule greatly over predicts the pressure, moment and lift coefficients at helicopter rotor tip sections when applied directly to incompressible results from lower Mach numbers as shown in Figure 2.21. Further background information on compressibility is described in Ref. [31] [43] [44]. 40 Figure 2.21 Comparison of lift coefficient for incompressible and compressible flow with direct application of the Prandtl-Glauert rule (Source from Ref. [45] ) Szymendera [45] derived the Prandtl-Glauert rule starting with the continuity equation based on Ref. [31], [43], [44] and then applied it to the Biot-Savart law. This resulted in a domain stretching in the flow direction which is dependent on the Mach number and has a direct effect on the induced velocity computation. Equation 2.48 and Figure 2.22 show the relation between domain stretching and Mach number. 2 ' 1 M r r p p ? = (2.48) As shown in Fig. 2.21 compressible results from Szymendera show better agreement with the coefficients than the results computed directly by applying the Prandtl-Glauert rule. ? Experiment ? ? Incompressible Method ----- Prandtl-Glauert Correction ____ Compressible Method 41 Figure 2.22 Prandtl-Glauert domain stretching in flow direction The implementation of compressibility into the propeller performance program is done after the propeller geometry is determined based on the input parameters. Since the tangential velocity increases along the propeller span, the stretching is more pronounced at the tip where the velocities are the highest, as illustrated in Figure 2.23 below. . r 1 r 2 r 0 r 1 ? r 1 ? r 0 ? Stretched domain analogous to incompressible flow ? Actual domain analogous to compressible flow P r p r p ? 42 Figure 2.23 Prandtl-Glauert domain stretching of vortex ring elements The propeller performance code was run with and without compressibility taken into account. The propeller diameter was 0.3 [m] and the propeller rpm was 10,000, which results in a tip Mach number of 0.47. The section lift coefficient for both the compressible and the incompressible solution is shown in Figure 2.24. This result agrees with the findings of Gennaretti and Morino [42] . The section lift coefficient is higher for the compressible case since the domain shrinking based on higher Mach numbers results in a smaller chord which in return increases c l . Stretched domain, incompressible Actual domain, compressible 43 Radial position r/R [-] S e ct io n lif t c o e f f i c i e n t C l [ - ] 0.2 0.4 0.6 0.8 0.06 0.08 0.1 0.12 0.14 Figure 2.24 Spanwise lift coefficient normalized by blade radius To determine the influence of the compressibility on the overall performance of small thin airfoil propellers at lower Mach numbers considered in this effort, a second investigation has been done. A generic propeller was selected as shown in Figure 2.25, with a diameter of 0.25 [m] operating at 36 [m/s] free stream velocity. The propeller speed was raised in steps starting from 4000 to 10000 revolutions per minute. The parameters of interest are the propeller efficiencies which are shown in Figure 2.26 for the compressible and incompressible case. It is observable that at Incompressible Compressible 44 Mach numbers above 0.25 efficiencies are slightly higher in the compressible case. The higher the tip Mach number the larger the discrepancy between incompressible and compressible solution. X Z Y Figure 2.25 Generic two bladed propeller with wake 45 Mach number [-] P r o p e lle r E f f ic ie n c y ? [% ] 0.1 0.15 0.2 0.25 0.3 0.35 0.4 65 70 75 80 Compressible solution Incompressible solution Figure 2.26 Compressibility effects vs propeller tip Mach number Most of today?s electric propulsion systems used on small UAS?s such as the US Army?s Raven system (manufactured by AeroVironment) use propeller speeds in cruise below 8000 revolutions per minute and propeller diameters of 0.2-0.4 [m]. This generates a maximum propeller tip speed in the range of M=0.3 and small propeller efficiency differences, of about 1% between incompressible and compressible solution, are noticeable. Figure 2.26 illustrates the reduced efficiency when compressibility is considered at different propeller speeds. 46 2.6 Viscous Model Since the VLM is an inviscid solution technique the physically accurate treatment of viscous flows is beyond its capability. Inviscid aerodynamic solvers predict lift and drag well for high Reynolds numbers (Re) and moderate angle of attacks, but these predictions become increasingly less accurate as the Re decreases and the sectional angle of attack increases. Viscous-inviscid interaction methods [46]- [47] increase the range of the prediction accuracy of the boundary element methods while maintaining computational economy, but they are beyond the scope of this work. Most of the viscous-inviscid interaction methods use the solution obtained by a boundary method and apply a numerical boundary layer model to account for viscous effects and flow separation. In this effort, a simple Reynolds number based model is implemented in an attempt to account for the effects of viscosity as it varies in the spanwise direction. Since the propeller blades considered in this effort all have thin airfoils, a flat plate boundary layer method [14] is used to determine the skin friction coefficients c f . It is assumed that laminar flow exists over the first 15% of the chord and that the remaining 85% the boundary layer flow is turbulent. [48] Thus the skin friction coefficients for laminar and turbulent flows are: c arlaf c Re 328.1 min_ = (2.49) 5 1_ Re 074.0 c turbulentf c = (2.50 47 with Re c being the Reynolds number at the spanwise collocation points with respect to the propeller section chord length c. ? ? cV onset c ?? =Re (2.51) Once the Reynolds number is determined, the skin friction coefficient on each spanwise collocation position is computed. The viscous drag is then obtained from: [ ] fonsetvisc cSVD ????= 2 2 1 ? (2.52) with S being the summation of all sectional wetted areas at one spanwise station. With the viscous drag known at each spanwise propeller blade position, the torque and the power increase due to viscous effects are: nxblQP rDQ viscvisc iviscvisc ??= ?= ? (2.53) The viscous torque is multiplied by 2.0 to account for both the top and the bottom side of the propeller surface. Finally the total power due to viscous effects is calculated by multiplying the viscous torque by the number of propeller blades. To gain confidence in the viscous flow modeling, a NACA 109622 straight blade propeller was simulated and compared with results obtained from lifting line method with and without viscous effects considered. To account for viscosity, Cheung [23] used drag data obtained from Ref. [49] of a NACA16-004 airfoil section. Figure 2.27 shows both the results of from Cheung obtained with lifting line method and the flat plate viscous 48 computation method shown above. There is good agreement between the two inviscid (lifting line and lifting surface) methods and small discrepancies in the results of the two viscous methods at lower advance ratios. 0.1 0.15 0.2 0.25 0.3 0.35 1.7 1.8 1.9 2 2.1 2.2 Advance Ratio J [-] Po w er C o efficien t c P [-] Lifting Line Method no Thickness and cd=0 Lifting Line Method no Thickness with Drag Data Lifting Surface Method no Thickness and cd=0 Lifting Surface Method viscous no thickness Figure 2.27 Spanwise lift coefficient normalized by blade radius [23] 2.7 Flow Separation Model The flow field associated with propellers is highly three dimensional and unsteady in nature and, as yet, not well understood. [50]- [53] Of particular difficulty is the prediction of propeller performance at high angles of attack when flow separation occurs. In recent years, new efforts have been made to better predict post stall behavior since flow separation on inboard regions of wind turbines are common during much of their 49 operational times. [54]- [59] Most of the developed methods for predicting flow separation use an empirical or semi empirical stall model due to the complexity of the flow. Tangler and Kocurek [60] developed a global post stall method based on the Viterna?s equations. [61]- [62] The disadvantage is that the slope and the maximum value of the aerodynamic coefficients must be known for the model to be used. Further if the aerodynamic coefficients do not match the behavior of a flat plate, the method becomes unusable up to 20? angle of attack. In Ref. [60], it was demonstrated that even when thick airfoil sections are considered, post stall behavior of c l /c d of a rotor, beyond 20 degrees angle of attack, follows flat plate theory. In this effort, the first part of the Viterna equations is used, since it predicts the maximum section drag coefficient as a function of the aspect ratio of the propeller. The maximum global drag coefficient prediction, shown in Equation 2.54, was used by Ref. [60] for sectional post stall predictions and showed good agreement when compared to experimental data obtained from a rotor with five spanwise pressure sensors: ARc d ?+= 018.011.1 max_ (2.54) For the computation of the stall and post stall behavior of the flow, data obtained from a NACA 0009 symmetrical airfoil [63] were used to generate the empirical post stall model. The NACA 0009 airfoil was selected because experimental data exist for it and since it is still considered a thin airfoil which matches the thin propeller blade simulation efforts of this work. This was done with a simple curve fitting of the lift and drag coefficient slopes starting at ~11? and going up to 90? angle of attack. Figure 2.28 and 2.29 illustrates the curve fitting of the lift and drag coefficients. In Appendix F the figures 50 generated from the experimental data from Ref. [63] are shown. If a single vortex panel exceeds 13? angle of attack, the flow separation model is applied to the corresponding panel of the propeller blade and new lift and drag forces are computed. This allows for performance computations of partially stalled propellers. During most of the optimizations, the panel angles of attack stay below 13?, since at larger Angle of Attack propeller efficiencies degrade. Thus the empirical flow separation model is mainly used when fixed pitch propellers designed for high cruise speeds are analyzed at low free stream velocities where local panel angles of attack exceed 13 degrees. y = -0.0000000108x 5 + 0.0000029014x 4 - 0.0002898088x 3 + 0.0127624897x 2 - 0.2301008464x + 2.1352985456 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 10203040506070809010 Angle of Attack [deg] Lift C o ef fici en t C l [ - ] Data experiment Polynomial trend line Figure 2.28 Lift coefficient of NACA 0009 in post stall 51 y = -0.00000290x 3 + 0.00034614x 2 + 0.00345460x - 0.00587855 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 10203040506070809010 Angle of Attack [deg] Dr ag coef f i ci ent Cd [ - ] Data experiment Polynomial trend line Figure 2.29 Drag coefficient of NACA 0009 in post stall The NACA109622 straight blade propeller was used to quantitatively verify the effects of the flow separation model. Figure 2.30 shows that below an advance ratio of 1.3 certain propeller regions are above 13? angle of attack, which activates the flow separation model, thus reducing the lift and increasing the required power. Figure 2.30 and 2.31 illustrate two runs with the propeller performance program where lift and power with and without the flow separation model activated are investigated. 52 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 0 0.5 1 1.5 2 2.5 Advance Ratio J [-] Thr us t T [N ] Trust of NACA propeller no flow separation Thrust of NACA propeller with flow separation Figure 2.30 Post stall thrust behavior of a NACA109622 propeller with empirical stall model 2.0E+05 2.5E+05 3.0E+05 3.5E+05 4.0E+05 4.5E+05 5.0E+05 00.511.522.5 Advance Ratio J [-] Po wer P [W] Power of NACA propeller no flow separation Power of NACA propeller with flow separation Figure 2.31 Post stall power behavior of a NACA109622 propeller with empirical stall model 53 2.8 Acoustic Signature Model The aerodynamic noise prediction of propellers dates back to 1919 when Lynam and Webb first published their work, which was driven by the requirement to have aircraft flying undetected over enemy territory. Since then, several different empirical, semi empirical and theoretical methods of propeller noise prediction have been developed and continuously improved. In recent years, propeller noise has again become of public interest due to the rapid increase in air traffic and the popularity of wind turbines for electric power generation. A summary of some of the different methods can be found in Reference [64]. More recently, the introduction of UAS?s as reconnaissance platforms in military operations, brought attention to the propulsion and propeller noise, since low flying reconnaissance UAS?s have low detectability requirements. Engine related noise, especially for small unmanned aircraft, is avoided by switching to an electric propulsion system, which leaves the propeller as one of the main sources of noise. [65] In this effort, a subsonic steady analysis using point sources (point forces located at the quarter chord for the elemental forces and point mass sources and sinks to represent propeller thickness) is implemented. The location of the point sources is illustrated in Figure 2.32 54 Mass sources Point forces Mass sinks Figure 2.32 Locations of noise sources The model is based on the work of Lowson [66] and was further developed by Miller [5] to account for noise generated by radial forces. It is assumed that the forward flight velocity is small M<0.1 which seems appropriate because only small low airspeed aircraft are considered in this work, thus a static analysis technique is applicable. For optimization efforts only the pressure perturbation, p-p 0, is considered and implemented into the objective function to reduce propeller noise. If the perceived noise level dBA of a human ear is of interest, the pressure perturbation must be converted to sound pressure levels (SPL) and then adjusted based on the frequency to represent human ear sensitivity. The pressure perturbation after Miller [5] is: 55 () () () ()() () () () ()() () () ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ?? ? ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ?? ? ??? +??? ? ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ?? ? ??? +??? ? ? ? ? ? ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ?? ? ?? ?? ? ? ? ? ? ? ? ? ???? =? r MM t M aM r aq r M t M aM rR r a F r M t M aM r r a F r M t M aM xx F rM pp srr r m sr r obss obsR sr r obs obs sr r s T r 2 0 0 2 00 2 00 2 0 22 0 1 1 11 1 sinsin cossin 11 1 sinsin cossin 11 1 14 1 ?? ?? ?? ?? ? ? (2.55) with ()()()()()[] 2 1 2 2 22 cossin2sincos ssobsobssobs RRrrxrr +?????+??= ???? (2.56) ( ) ( ) r rM M obss r ?? sinsin ??? = (2.57) 0 a R M s s ? = ? (2.58) ??? ?= (2.59) () () () () (){}???? ?? sinsinsincos sin 2 3 ????+??? ??? = ? ? sobs obssr Rrr r rM t M (2.60) where M s is the Mach number of a point on the quarter chord in the spanwise direction, r is the distance between the source and observer, r obs is the distance of the observer to the origin, ? is the angular velocity, ? is the retarded time, R s is the distance between the source and the origin, ? is the angle of the observer with the x axis, and ? is the angle of the source with respect to the y axis. A geometric illustration is shown in Figure 2.33. 56 The volume displacement of the propeller is simulated by a Rankine [5] body with a single source sink pair on each spanwise collocation position. The source strengths are determined by solving for the Rankine body in the freestream with the major axis length equal to the local propeller chord. Thus the maximum diameter, h, of the Rankine body is determined from: 2 ht c t S ?=?? ? ? ? ? ? ?? ? (2.61) with ?S the spanwise width of the panels and (t/c) the local thickness ratio. 57 Figure 2.33 Acoustic signature analysis geometry Since the aerodynamic propeller analysis of this work does not consider airfoil thickness, a generic thickness profile based on thin bladed propellers is applied, which is shown in Figure 2.34. . ? ? x y z x s x y observer, z=0 source . F T F R F ? ? yxr ?= () () ()?? sin,cos, 0,, ??= = sss RRxy yxx R s r obs 58 The point mass flux, q m , is obtained by solving the source/sink position from ( ) 01 22 2 2 =+??? ???? (2.62) where h c =? 2 ? (2.63) and chordmidfrompositionk sourcec ____sin2 =?? (2.64) Once ? is known, q m is obtained from 22 2 2 inf0 2 ?? ? ? ?? +? ? ? ? ? ? ? ? ? ?? ? ? ? ? ? ???= c Vq m (2.65) To determine the overall pressure perturbation the time of each spanwise pressure perturbation to reach the observer must be known as the propeller turns. This is done by initially calculating the distance each blade pressure perturbation has to the observer and dividing it by the speed of sound. The propeller blade is then rotated backwards in time in 1? increments and the observer time for each spanwise pressure perturbation is now determined by computing the times the pressure perturbations take to reach the observer minus the time the blade took to rotate 1?. Then, by matching the times of each pressure perturbation, the total pressure perturbation is the sum of all the local pressure perturbations with equal observer times. 59 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 012345678 Propeller Span [inch] Th i c k n e s s to Ch ord Ra t i o t/ c [-] Figure 2.34 Linearized airfoil thickness distribution of a APC 14x7 propeller The test case for verifying the acoustic noise signature is the data obtained by Succi [67]- [70] . The simulated propeller is a ? scale model of a Cessna 172 Skyhawk propeller with the type description 1C160. The geometry of the 1C160 quarter scale propeller is taken from Succi who describes the chord and propeller pitch in the spanwise direction. In Appendix G the propeller geometric description is given. Based on this information, a flat plate propeller was generated with 18 panels spanwise and 3 panels chordwise which is illustrated in Figure 2.35. In the spanwise direction are placed sectional forces which are located at the quarter chord position. The point mass sources and sinks are located with equal distance from the middle chord of each spanwise panel center. 60 Figure 2.35 1C160 quarter scale Cessna 172 propeller For the test case verification, the observer was located 0.5 meters from the propeller. The wave form is for an observer in the plane of rotation behind the propeller ? = 130?. The propeller speed was 10000 rpm and the free stream velocity 29 [m/s]. In Figure 2.36 the experimental and theoretical results of Succi as compared to the results of Miller. The implemented acoustic signature model of his work is shown in Figure 2.37 and agrees well when compared to the results of Succi and Miller. 61 Figure 2.36 Pressure perturbation comparison experiment and theory from Miller [5] 62 Time [ms] P r es s u r e p e r t u b a t i o n [ P a] -1.5 -1 -0.5 0 0.5 1 -40 -30 -20 -10 0 10 20 Figure 2.37 Pressure perturbation of 1C160 quarter scale propeller 63 3 PROPELLER AERODYNAMIC PERFORMANCE PREDICTION PROGRAM The propeller performance program is the principal portion of the objective function for the GA. The program consists of nine major routines. ? Airfoil and propeller geometry generation ? Influence coefficient determination ? Gauss Seidel solver subroutine to solve for the unknown circulations ? Main program to connect the subroutines and to compute the lift, drag and power of the propeller ? Compressibility correction subroutine ? Viscous model analysis ? Pressure perturbation determination ? Flow separation model The geometry of the propeller is defined by 40 variables: ? 22 variables describe the airfoil geometries ? 5 variables describe the spanwise angle of attack function ? 5 variables describe the spanwise chord length variation function ? 5 describe the propeller sweep function ? 1 for the number of propeller blades ? 1 for the propeller radius ? 1 is related to the position about which the pitch of the propeller is changed 64 The 40 propeller geometry design variables are allowed to change during the optimization process and other fixed parameters such as the free stream velocity, rotational speed, number of panels and propeller hub diameter are constants. Figure 3.1 illustrates the program structure. A description of the variables can be found in Section 2 and Appendix B. The open program structure allows for optional activation of certain individual subroutines to account for compressibility effects, viscous and acoustic signature determination. In addition to single (design and operating) point analysis the program is capable of analyzing a propeller over a range of free stream velocities and propeller rotational speeds. Additionally variable pitch propellers can be designed by adding variables which offset the angle of attack function of the propeller by a constant value. The program starts by loading the propeller geometry input parameters from the file ?Snglerun.dat.? The subroutine ?coordgen? generates the propeller shape based on the input parameters. The influence coefficients are determined and the unknown circulations are returned by a Gauss-Seidel solver. The accuracy of this basic inviscid solution can be increased by activating a viscous, compressibility model. Additionally an acoustic signature model can be used to determine propeller pressure perturbations. Based on the optimization requirements either thrust, power input, efficiency or acoustic signatures can be returned and used as an optimization objective. 65 Figure 3.1 Flow chart propeller aerodynamic performance program Input Parameters: pr = (x(1)) nxbl = (x(2)) ua_1 = (x(3)) ua_2 = (x(4)) ua_3 = (x(5)) ua_4 = (x(6)) ua_5 = (x(7)) ua_6 = (x(8)) ua_7 = (x(9)) ua_8 = (x(10)) ua_9 = (x(11)) ua_10 = (x(12)) ua_11 = (x(13)) ua_12 = (x(14)) ua_13 = (x(15)) ua_14 = (x(16)) ua_15 = (x(17)) ua_16 = (x(18)) ua_17 = (x(19)) ua_18 = (x(20)) ua_19 = (x(21)) ua_20 = (x(22)) fau_1 = (x(23)) fau_2 = (x(24)) ? Compressibility Model ? Viscous Model ? Flow Separation Model ? Noise Model ?Propeller airfoil and Geometry ?Grid Discretization ?Compressibility Domain stretching Influence Coefficient Determination Gauss Seidel Solver Write Output 66 4 CODE VERIFICATION The aerodynamic propeller performance code was developed without significant use of existing programs. Thus the numerical scheme was verified using experimental published and existing numerical data. In an effort to reduce optimization computation times the Gauss Seidel subroutine which solves for the unknown circulations was investigated to determine the minimum requirement for the convergence of solution. The minimum propeller wake extension and the minimum number of chordwise and spanwise panels which still produce accurate results were also investigated to further reduce computation times. 4.1 Gauss Seidel Solver convergence After the determination of the influence coefficients in the x, y and z-direction and the application of the BC that the flow is parallel to the propeller surface 0? =?nV r , a system of linear dependent equations is obtained. ( )() ()() ()() ? ? ? ? ? ? ? ? ? ? ? ? ?+?? ?+?? ?+?? = ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? mmmmmm m Vr Vr Vr b bb bbb ??? ??? ??? ? cossin ... cossin cossin ...... ... 4 1 222 111 2 1 2221 11211 (4.1) 67 with ),(),(inf_),(),(inf_),(),(inf_ jinjicoeffjinjicoeffjinjicoeffb zzyyxxij ?+?+?= (4.2) To solve for the unknown circulations ?, an implicit or explicit solver can be used. Due to the large number of panels on the propeller (up to 5000) an explicit Gauss Seidel solver [71] is used to determine the individual circulations of each quadrilateral and horseshoe element. The advantage of the explicit solver is that it is much faster when solving a system with a large number of equations. [71] The Gauss Seidel solver begins with the first linear equation, solves it for the first unknown ? as shown in Equation 4.3. The second linear equation is then solved after the circulation for the second panel and is updated with the vortex strength of the first panel. This procedure is repeated until all equations are solved once. The process is then reversed and all linear equations are solved again starting at the bottom of the matrix. [ ] [] [] mm m n mm n m n mn m n mm nn n n mm nn n b RHSbbb b RHSbbb b RHSbbb +??+??+?? =? +??+??+?? =? +??+??+?? =? + ?? ++ + + + + 1 11 1 22 1 111 22 22323 1 1211 2 11 113132121 1 ... .......... .......... ... ... (4.3) The convergence of the solution is measured by the summation of the errors between the left hand sides (LHS) and the right hand sides (RHS) of the linear equation after each matrix sweep. Equation 4.4 defines the determination of the related errors. 68 () ? = +++ +??+??+??= m i i n mim n i n i RHSbbberror 1 11 22 1 11 ... (4.4) A two bladed propeller with panel numbers ranging from 32 to 1568 has been investigated to determine the minimum number of required matrix sweeps which still produce accurate results. Figure 4.1 shows the convergence of the solution over the number of matrix sweeps. Noticeable is that after about 250 matrix sweeps the summation of the errors for all panel densities reach the numerical accuracy of the computer. Number of Matrix Sweeps [-] S u mma t i o n o f E r r o r [ m/ s ] 200 400 600 10 -16 10 -14 10 -12 10 -10 10 -8 10 -6 10 -4 4x4 8x8 12 x 12 16 x 16 20 x 20 24 x 24 28 x 28 Figure 4.1 Summation of LHS and RHS difference after matrix sweeps 69 To evaluate the convergence of the solution the thrust of a two bladed propeller using 512 and 1565 panels respectively is plotted over the number of matrix sweeps. These results of this study are shown in Figure 4.2 below. 1 10 100 1000 10000 11010 Number of Sweeps [-] Thr ust [ N ] 0.00001 0.0001 0.001 0.01 0.1 1 10 100 % o f th r u s t di ffer e nce w i th r esp . to 1000 sw eeps Thrust [N] 16 x16 Thrust [N] 28 x 28 % difference to 1000 16 x 16 % difference to 1000 28 x 28 Figure 4.2 Accuracy of the solution after matrix sweeps by Gauss Seidel solver A single matrix sweep is defined by solving all equations in the matrix once downstream and once upstream while updating the circulations continuously. The results shown in Figure 4.2 are referenced against the solution obtained after 1000 matrix sweeps. It is observed that with increasing panel numbers the number of matrix sweeps must be increased to maintain a prescribed numerical accuracy. Figure 4.2 shows that at about 70 sweeps, the solutions are within 1%. Thus for all further propeller performance 70 computations the number of matrix sweeps is set to 100, which will provide an accuracy of 0.01 in this case while maintaining computational efficiency. 4.2 Wake investigation To minimize the computation times during the optimization process, a wake analysis was done to determine the convergence behavior of the aerodynamic performance parameters. The goal is to reduce the number of downstream rotations of the propellers rigid helical vortex sheet, while maintaining accuracy of the solution. The wake analysis was done on optimized two and four bladed propellers with a diameter of 0.32 meters, a freestream velocity of 21 [m/s] and a propeller speed of 5800 RPM. The wake vortex sheet was extended in steps from 1/8 to 5 rotations downstream and propeller performance parameters were noted. An illustration of the wakes is shown in Figures 4.2-4.4: 71 X Z Y X Z Y Figure 4.2-4 Wake vortex sheet extensions for ?, 2 and 5 rotations downstream 72 The investigated performance parameters are thrust and power since in most propeller optimizations these are of special interest. Figure 4.5 shows the thrust generated form two and four bladed propellers, depending on the wake extension in numbers of 360 o rotations. Also shown is the percent difference of the propeller thrust with respect to the max vortex sheet (5 rotations downstream). More details of the wake convergence investigation are given in Appendix C. 2.5 3.5 4.5 5.5 6.5 7.5 0123456 Number of wake rotations [-] Thr us t T [ N ] 0 1 2 3 4 5 6 7 8 9 10 Di f f e r e n c e i n Thr us t [ % ] Thrust 2 bladed propeller Thrust 4 bladed propeller % difference in thrust 2 bladed propeller % difference in thrust 4 bladed propeller Figure 4.5 Propeller thrust variation due to wake vortex sheet extension In both, the two and four bladed propeller cases, the thrust converges within 0.3% after only 2 full downstream rotations of the wake vortex sheet. This result is confirmed by Cheung [23] who showed that, even with a limited wake, good agreements with 73 experimental data is obtained. Thus in all further optimizations the propeller wake is extended only to two 360? rotations downstream. 4.3 Number of Chordwise and Spanwise Panels To minimize the optimization computation times, a grid convergence study has been done. The goal is to reduce the numbers of chordwise and spanwise panels while maintaining solution accuracy. Two different propeller grid discretizations have been investigated as shown in Figure 4.6 below. The chosen schemes are based on the work done by S. Werner [24] , Fiddes & Gaydon [28] and Labrujere & Zandbergen [72] . The investigated schemes are 1. Full cosine distribution in the chordwise direction and even spacing in the spanwise direction 2. Full cosine distribution in the chordwise direction and half cosine in spanwise direction Figure 4.6 Propeller Discretization schemes Propeller root Propeller tip 1 2 74 Figure 4.7 shows that the discretization scheme with a full cosine distribution in the chordwise direction and the half cosine in the spanwise direction converge the fastest. The accuracy when compared to the converged thrust value is within 0.5% with a 7x14 (7 panels in chordwise and 14 panels in spanwise direction) panel arrangement. This agrees with the findings form Cheung [23] , who found good agreement with experimental data with a blade discretization of 7 chordwise and 10 spanwise panels. 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0 25 50 75 100 125 150 175 200 Number of Panels [-] Thr u st T [ N ] Full cosine chordwise 3x20 Full cosine chordwise 4x20 Full cosine chordwise 5x20 Full cosine chordwise 6x20 Full cosine chordwise 7x20 Full cosine chordwise 7x5 Full cosine chordwise 7x10 Full cosine chordwise 7x25 Full cosine chordwise 7x30 Full cosine chordwise 7x n half cosine spanwise Figure 4.7 Lift force convergence by discretization scheme and panel number 75 4.4 Validation of the numerical method 4.4.1 Comparison to lifting line and lifting surface method The VLM by lifting surface was applied to a simple test case to validate the numerical method. The validation includes the accuracy of lift coefficient and drag coefficient for a flat plate wing with zero thickness and an aspect ratio of AR=6 at different angles off attack. The wing is simulated with a single layer of vortex ring elements and is discretized in four even spaced panels spanwise and two even spaced panels chordwise. Figure 4.8 illustrates the wing geometry and the horseshoe trailers. Even though the lifting line method agrees well with results obtained by the lifting surface solution when compared to flat plates, the lifting line method is not capable of predicting more complex geometries accurately. 76 X Y Z Figure 4.8 Wing and wake discretization The results are compared to the solutions obtained from lifting line method by Kuethe and Chow [18] , which is illustrated in Figure 4.9. 77 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 2 4 6 8 10 12 14 Angle of Attack [deg] L i ft Coeffi c i en t cl [-] cl from Kuethe & Chow (lifting line method) cl from Propeller Program (lifting surface method) 0 0.01 0.02 0.03 0.04 0.05 0.06 2 4 6 8 101214 Angle of Attack [deg] Dr a g Coef fic i e n t c d [-] cd from Propeller program (lifting surface method) cd from Kuethe & Chow (lifting line method) Figure 4.9 Lifting line and lifting surface method cl and cd of flat plate wing AR=6 78 Good agreement between the lifting surface solution and the lifting line method is obtained for both, the lift and the drag coefficient for the investigated flat plate wing with AR=6. 4.4.2 NACA 109622 Propeller In a second validation case the lifting surface method was used to predict the performance of a NACA109622 straight blade propeller with three blades. The blades are modeled by ten spanwise and two chordwise panels with a blade angle ?=45.4? at the 75% propeller radius. The geometry of the propeller is given in Appendix E. Figure 4.10 shows the propeller geometry with the horseshoe trailers extending two full propeller rotations downstream. 79 X Z Y Figure 4.10 NACA109622 propeller geometry with wake The results are plotted in Figure 4.11-13 and incorporate experimental data as well as data obtained from lifting line method by Cheung [23] with and without 2D drag data. The lifting line method with no drag data and the lifting surface method agree very well. This is primarily due to the fact that the propeller blade is not cambered, thus the lifting surface method does not improve the accuracy of the solution. 80 0 0.2 0.4 0.6 0.8 1 1.2 1.5 1.6 1.7 1.8 1.9 2 2.1 2.2 2.3 2.4 2.5 Advance Ratio J [-] E f ficie n c y ? [- ] Experiment Lifing Line Method no Thickness and cd=0 Lifting Line Method no Thickness with Drag Data Lifting Surface Method no Thickness and cd=0 Figure 4.11 Propeller efficiency comparison at different advance ratios Figures 4.12 and 4.13 show the thrust and power coefficients obtained from experimental data from a lifting line solution with and without drag data included and from the lifting surface solution developed in this effort. The lifting line and the lifting surface solution agree very well for both the thrust and the power coefficient. Discrepancies between the experimental data are mainly due to the simplified model which does not include thickness or Reynolds number effects. 81 0 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 1.5 1.6 1.7 1.8 1.9 2 2.1 2.2 2.3 2.4 2.5 Advance Ratio J [-] Th rust Coeffici ent c T [- ] Experiment Lifting Line Method no Thickness and cd=0 Lifting Line Method no Thickness with Drag Data Lifting Surface Method no Thickness and cd=0 Figure 4.12 Propeller thrust coefficient at various advance ratios 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 1.7 1.8 1.9 2 2.1 2.2 2.3 2.4 2.5 Advance Ratio J [-] P o wer Coefficien t c P [-] Experiment Lifting Line Method no Thickness and cd=0 Lifting Line Method no Thickness with Drag Data Lifting Surface Method no Thickness and cd=0 Figure 4.13 Propeller power coefficients at various advance ratios 82 5 GENETIC ALGORITHM The optimization of designs in early stages of product development has become more popular in recent years. The goal is to find an optimal solution for any given design problem. This effort is driven by cost savings in the development stage and product competitiveness. Real world applications are found today in every sector of life, which is mainly due the introduction of very powerful and cheap computers. This lead to the development of various mathematical optimization theories which range from simple methods using first derivatives to advanced evolutionary algorithms. 5.1 Optimization methods Even though optimizations existed before the 17 th century, with the development of Newton?s calculus, initially differentiation was used with simple functions to find local minimums and maximums. Problems emerged when multiple functions were combined for a system description and nonlinear algebraic equations were produced which were not easy to solve. In the 19 th century gradient based methods emerged, which used the steepest ascent or descent (steps are proportional to the gradient) to find the functions maximum or minimum. The disadvantages of this method are that the function must be differentiable in the design space with respect to all independent variables and only local minimums or maximums can be found, which requires a good initial guess of the solution 83 before starting the optimization. Linear programming as stated by Gabasov and Kirillova [73] was first formulated by a Soviet mathematician L.V. Kantorovich in the 1930?s in which an objective function with linear constraints was used. The disadvantage of this method is that it cannot handle complex problems with multiple variables since the objective function and the constraints must be linear. Another method developed around 1960 by Hooke and Jeeves is pattern search optimizations which looks at the behavior of a function to be analyzed and based on this information generates new points which are closer to the optimum. The disadvantages of pattern methods are that they may find a local minimum or maximum depending on the initially guessed starting solution. The concept of ?Design of Experiments? was initially developed in the 1920s with the requirement to analyze vast amounts of experimental data. Weber and Skillings [74] defined steps for the experimenter to define the structure of the experiment and to find the proper testing setup. In this method, goals are defined, experiment variables and ranges are set, the experiment is run and the results are analyzed. A linear statistical model is then used to describe the structure of the data resulting from the experiment. Monte Carlo methods are also popular methods since they are able to solve complex physical and mathematical problems with many variables. They are computationally expensive since they randomly ?walk? through the design space by changing variable numbers to find a minimum or maximum. An incorporated gradient based method can facilitate the determination of the optimum. A relatively new method is the Particle Swarm Optimization method which was developed by Kennedy and Eberhart [75] in the mid 1990s. It is a population based stochastic optimization technique inspired by social behavior of bird flocking or fish 84 schooling. PSO shares many similarities with evolutionary computational techniques such as Genetic Algorithms (GA?s). The system is initialized with a population of random solutions and searches for optima by updating generations. However, unlike GA?s, PSO?s have no evolutionary operators such as crossover and mutation. In PSO, the potential solutions, called particles, fly through the problem space by following the current optimum particles. Each particle keeps track of its coordinates in the problem space which are associated with the best solution (fitness) it has achieved so far and stores the fitness value. Another "best" value that is tracked by the particle swarm optimizer is the best value, obtained ?to date? by any particle in the neighborhood of the particle. If a particle takes all the population as its topological neighbors, the best value is a global best. The particle swarm optimization concept consists of an evaluation of the performance of each particle at each time step and based on the particles and the global best solution a change in the velocity of (accelerating) each particle is determined. Acceleration is weighted by a random term, with separate random numbers being generated for acceleration toward best (individual particle and global) locations. In the past several years, PSO has been successfully applied in many research and application areas Ref [76]- [81]. It has been demonstrated that PSO gets results in a faster, cheaper way as compared with Genetic Algorithms [82]- [84] . Another reason that PSO is attractive is that there are few parameters to adjust. One version, with slight variations, works well in a wide variety of applications. The drawbacks of PSO?s are that they may prematurely converge and are not able to improve on the solution. The second drawback is that stochastic approaches have problem dependent performance. This dependency normally results from the parameter settings of each algorithm. 85 5.2 Genetic Algorithm Mitchell [85] states that ?Genetic Algorithms were invented by John Holland in the 1960?s and were developed by Holland and his colleagues at the University of Michigan in the 1960?s and 1970?s?. Holland wanted to take the evolutionary processes that are hypothesized to occur in nature (adaptation, survival of the fittest, etc) and incorporate them into a computer system. Holland?s classic book ?Adaptation in Nature and Artificial Systems? formulated evolutionary/population based algorithms that can be used to optimize a variety of real world systems. According to Coley [86] , Genetic Algorithms are numerical optimization algorithms built around natural selection and natural genetics. Sometimes GA?s are referred to as evolutionary optimizers because they mimic some of the process in Darwinian evolution theory. Darwinian evolution theorizes that totally new species can be produced via mutation and random chance. However, unlike evolution, GAs operate more on the principal of improvement of an initial population of solutions to a design problem rather than pure optimization. Coley [86] describes the make-up of a typical GA in a number of ways. First, a GA can be described as a number, or population, of guesses of the solution to the problem. Second, The GA employs a method of calculating how good or bad the individual solutions within the population are. This method is called determining the fitness of the solution. Additionally, a method for mixing fragments of the better solutions to form new and on average, even better solutions can be used. This method is called crossover. Finally, the GA is able to use a mutation operator to avoid permanent loss of diversity within the solutions. 86 One of the main benefits of using a GA is the fact that the algorithm can start without a single point, or initial guess, to get the optimization running. The previously described direct methods all, except the design of experiments and PSO, need an initial guessed solution to ?march? toward the desired optimized result, which can be a maximum or a minimum value. The GA uses a population of guesses that are random and spread throughout the search space. Powerful operators such as selection, crossover and mutation help direct members of each population toward the desired goal(s) of the problem. A binary encoding system allows for a host of variables to be manipulated by the GA and then used in a suite of performance codes. These codes analyze the performance of each member of the population and the GA ranks each one according to how well the member of the population meets the desired goal(s). In GA terminology, the objective function is the function that determines the performance of a particular chromosome (i.e. member of the population). [87] Thus, in this study, the suite of performance codes, grouped together as a whole, represent the objective function. At the same time the GA has some disadvantages. In using the GA there is a greater likelihood that a global optimum solution will be found. However, finding this global optimum is not guaranteed. Even if the GA is in the neighborhood of the global optimum, there is a possibility, that through crossover and mutation the global optimum may not be selected. Also the GA does not address the robustness of the individual design solutions it creates. The GA simply attempts to meet the desired goals and will adjust the design parameters accordingly. Thus, it is up to the user to ensure the proper operation of the GA and to verify that the results are genuine. Finally, the satisfactory operation of the GA 87 relies on the accuracy of the system models that make up the objective function. The accuracy of the system and overall objective function is described in Chapter 2 and 3. In this effort a binary encode tournament based GA which was developed by Anderson [87] was used to drive the design optimization. The population size of each generation was 400 and the optimizations were run for about 700 generations. Results from initial GA propeller optimizations showed premature convergence behavior. A second real coded GA, which was recently developed by John [88] , was implemented and showed better convergence to an optimum solution when compared to the binary encoded GA. 88 6 IMPLEMENTATION OF THE PROPELLER PERFORMANCE PROGRAM INTO THE GENETIC ALGORITHM The propeller aerodynamic performance program is integrated in the GA as a subroutine which, when required, is called by the GA. The GA design variables are listed in a separate file with min and max value limitations as well as step resolution. This input file also contains the GA parameters, population size, number of generations and advanced features such as elitism, niching and pareto. The GA starts by creating a population of different propellers based on a random selection of the blade geometry variables. All members are then sent to the propeller program where the performance parameters lift, drag, power and efficiency are determined. In the final step the value of the objective parameter or function is send back to the GA for evaluation. The best performing members of the first generation are selected and, based on the values of these variables, a new generation of propellers is created. Figure 6.1 shows the optimization structure. 89 Figure 6.1.Flow chart propeller performance program implementation into GA The objective function returns a parameter which can be minimized or maximized by the GA. Based on the requirements of the user, different propeller performance parameters can be integrated into the objection function. When a combination of parameters is included in one single objective function, close attention must be given to how each parameter is weighted in the function. To avoid multiple design variables Genetic Algorithm (binary) Generation with 400-propeller members Best performers are chosen and new population is generated based on Gannl.dat Optimization parameters Population size Number of generations Optimization characteristic (niche, pareto, elitist?) Propeller Parameters (max, min, resolution) Number of propeller blades Propeller span Propeller 3D geometry function (22) AOA function (5) Sweep function (5) Chord function (5) Propeller hinge point Propeller Performance 3D propeller geometry (with Chord, AOA and Sweep function) Grid generation Results: Vortex strength Normal force Viscous effects Compressibility Flow separation model Propeller acoustic signature Aerodynamic coefficients, power, efficiency Objective function 90 within a function, the GA allows for simultaneous optimization of multiple objective functions. This simplifies the goal function and leads to more diverse goal options. 91 7 OPTIMIZATION RESULTS The motivation behind this effort was to investigate the feasibility of designing propellers using a binary coded GA and a real encoded GA. The goal has been to speed up the design process while optimizing the propeller for a specific objective composed of a single or multiple goals. For the purpose of finding an optimal propeller the design space of available geometries has been developed with a wide range of design parameters. Throughout all optimizations, the number of chordwise and spanwise panel elements (7 and 14) has been fixed based on the findings of Section 4.3. In addition, the number of propeller blades was limited in the program to a maximum of four blades. The cases are broken down into single and multiple point optimizations as well as inclusions of compressibility, viscous effects and propeller acoustic signature. The objective was to design a propeller with better performance than existing propellers. The propeller operating conditions are defined by the cruise power, the flight speed and the maximum propeller diameter. These propeller operating conditions were set to match a small generic UAS like the Raven which is currently used by the US Army and which was developed by AeroVironment. Therefore the power was set to 200-300 Watts in cruise, the maximum propeller diameter was 0.25 meters and the flight speed 36 [m/s]. The test cases considered are described in the following. 92 7.1 Comparison of propeller optimizations with and without viscous and compressibility effects taken into account For comparison of viscous and compressibility effects, three different propeller optimization cases at different operational points have been investigated. In the first case a fixed pitch propeller was optimized for a cruise airspeed of 36 [m/s]; in the second case a fixed pitch propeller is optimized for a launch condition of 9 [m/s]; and in the third case a fixed pitch propeller is optimized using both cruise and launch conditions as the objective. The parameters of interest are the propeller thrust (the first goal to be matched) and power (second goal to be minimized to maximize propeller efficiency). 7.1.1 Single point Optimization for Cruise Condition: Binary encoded GA Three different single point optimization runs have been conducted to determine the influence of viscous and compressibility effects. The first case did not include compressibility or viscous effects, the second one considered only compressibility, and the last took only viscous effects into account. All of the cases had the same fixed input parameters as shown in Table 7.1. Table 7.1 Fixed input parameters cruise condition Free stream velocity 35 [m/s] Propeller speed 7500 [rpm] Propeller max diameter 0.25 [m] Propeller hub diameter 0.06 [m] 93 The binary encoded GA was set up to drive the design using two objective functions which must be minimized. Within the GA the two objective functions are combined to one single function which is stored in the ?GAout1? file and used to determine convergence behavior. The first objective function contains the required thrust in cruise, which was set in this design effort to 7.0 [N]. The cruise thrust selected was a guess, based on available data for small, UAS and does not represent a particular system. Subtracted from the required thrust is the generated thrust determined by the aerodynamic prediction program. Equation 7.1 shows the function which is minimized when the propeller thrust reaches the required thrust. Thrustobj ?= 0.71_ (7.1) In the second objective function the power was set as a goal to be minimized which drives the propeller efficiency. The power was divided by 300 to equalize the weighting of the objective functions. 300 2_ Power obj = (7.2) Table 7.2 shows the results of the three optimization cases which have been run for about 200 generations with each generation consisting of 400 members. In all of the cases the thrust required was matched. The propeller efficiencies range between 78 and 87% depending if viscosity or compressibility are taken into account. A variable pitch propeller with similar diameter was tested in the wind tunnel at 18 [m/s] free stream 94 velocity and was rated with a max efficiency of 62 %. Wind tunnel data are given in Appendix H. Table 7.2 Optimized propellers performance parameter results cruise condition Case 1 Inviscid Incompressible Case 2 Inviscid Compressible Case 3 Viscous Incompressible Thrust [N] 7.0 7.0 7.0 Power [W] 306.9 290.2 321.3 Efficiency [%] 82.0 86.8 78.4 In Figures 7.1, 7.2 and 7.3 the optimized propellers for inviscid, compressible and viscous solutions are illustrated. All blades are two bladed with the inviscid and the viscous solutions showing some form of a proplet. Noticeable in Figures 7.1-7.3 are also the leading edge panels of the propeller blades are bent downwards in the tip region which would explain why the efficiencies in some cases are lower. 95 Figure 7.1 Optimized cruise propeller geometry from binary encoded GA; inviscid, incompressible 96 Figure 7.2 Optimized cruise propeller geometry from binary encoded GA; inviscid, compressible 97 Figure 7.3 Optimized cruise propeller geometry from binary encoded GA; viscous, incompressible 98 Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 inviscid, incompressible inviscid, compressible viscous, incompressible Figure 7.4 Optimized propeller spanwise circulation distributions for cruise The spanwise circulations for the propellers considered in this section are presented in Figure 7.4. All the cases have discontinuities in the spanwise circulation, which reduces propeller efficiencies and is a sign that only a local optimum has likely been found. The viscous incompressible case in Figure 7.3 developed a noticeable large proplet at the tip, which allows for a higher tip loading which is shown in Figure 7.4 by an increased circulation. 99 Number of Generations [-] Number of Evaluations [-] O b je c t i v e F u n c t io n [ - ] 0 50 100 150 200 250 0 50000 100000 1 1.5 2 inviscid, incompressible inviscid, compressible viscous, incompressible Figure 7.5 Objective function convergences binary encoded propeller Since the three optimizations have been run for only 200 generations a separate run has been done to investigate the long term convergence behavior. Therefore the inviscid incompressible case has bee run for 1300 generations which relates to 520,000 propeller evaluations and took about 6 days of computation. The convergence result is illustrated in Figure 7.6 and demonstrates that after about 100 generations no significant improvements are achieved even though a gain of about 3-4 % in propeller efficiency is possible. 100 Number of Generations O b j e c t i v e F unc t i on 0 500 1000 1 1.2 1.4 1.6 1.8 2 2.2 inviscid, incompressible Figure 7.6 Objective function convergence inviscid, incompressible, for binary encoded GA 7.1.2 Single point Optimization for Cruise Condition: Real encoded GA The three different cases presented in section 7.1, the inviscid, incompressible, the inviscid, compressible and the viscous incompressible were run again using a real encoded GA to investigate convergence behavior to determine whether on not the objective functions are minimized to an lower value. All input parameters and objective functions are identical to the cases shown in Section 7.1.1. In table 7.3 are illustrated the results obtained from the optimizations of the real coded GA. The thrust objective in all cases is satisfied exactly with the inviscid incompressible 101 case having the highest efficiencies. The compressible case lags the efficiency by about 7 % when compared to the results obtained be the binary encoded GA, which indicates a premature convergence thus a global optimum is not reached. Table 7.3 Optimized propellers performance parameter results real encoded GA Case 1 Inviscid Incompressible Case 2 Inviscid Compressible Case 3 Viscous Incompressible Thrust [N] 7.0 7.0 7.03 Power [W] 288.9 320.5 320.4 Efficiency [%] 87.2 79.2 79.0 The propeller geometries generated by the real coded GA are shown in the Figures 7.7 ? 7.9. It is again characteristic that two propeller geometries have proplets which would indicate a loading shift to the tip. Though the circulation distributions in Figure 7.10 shows that the loading is actually more spread out when compared with the results obtained by the binary encoded GA. The lag in efficiency of the viscous incompressible case is also seen in the discontinuous circulation distribution illustrated in Figure 7.10. 102 Figure 7.7 Optimized cruise propeller geometry from real encoded GA; inviscid, incompressible 103 Figure 7.8 Optimized cruise propeller geometry from real encoded GA; inviscid, compressible 104 Figure 7.9 Optimized cruise propeller geometry from real encoded GA; viscous, incompressible 105 Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.2 0.4 0.6 0.8 1 1.2 inviscid, incompressible inviscid, compressible viscous, incompressible Figure 7.10 Optimized circulation distributions for cruise propeller, real encoded GA When comparing Figure 7.5 and 7.10 it is noticeable that even though both GA versions seem to level off in the convergence behavior around 100 generations, the real encoded GA shows a higher gradient in the objective functions after the first 100 generations. Thus for all the following optimizations the real coded GA is used. The case where viscous effects are taken into account is the only one used in the remaining optimizations since the subroutine for compressibility effects should only be used in determining the performance of a propeller when compared to an inviscid or viscous solution. 106 Number of Generations [-] Number of Evaluations [-] O b j e c t i v e F unc t i on [ - ] 0 50 100 150 200 250 0 50000 100000 0.8 1 1.2 1.4 1.6 1.8 2 inviscid, incompressible inviscid, compressible viscous, incompressible Figure 7.11 Objective function convergences real encoded GA, cruise condition To verify that the optimized propeller has better performance a CAM 13x7 propeller, set up with a variable pitch mechanism, was tested in a wind tunnel at about 19 [m/s] free stream velocity. The power input was kept constant while the propeller pitch was changed gradually. Figure 7.12 shows that the maximum propeller efficiency of the CAM 13x7, at 19 [m/s], is 63% while the optimized propeller operates at about 80% efficiency. 107 0 10 20 30 40 50 60 70 0.06 0.08 0.1 0.12 0.14 0.16 0.18 0.2 Propeller Torque [Nm] Pro p e l l e r Ef f i ci en cy [ % ] CAM 13x7 Propeller Figure 7.12 CAM 13x7 propeller efficiency at 19 [m/s] free stream velocity 7.2 Single point Optimization Launch Condition A single point optimization was conducted to design a propeller for a launch condition of 9 [m/s]. The investigation included the case where viscous effects are taken into account as described in section 7.1.1. The only difference in the fixed input parameters is the launch airspeed which simulates a hand launch of a UAS by a human and the thrust requirement of 14 [N]. The input parameters for this case are shown in Table 7.3. The thrust required for the launch was set to a high value since the operational mode of this 108 system does not include a mechanical launch by catapult or bungee cord. Thus the optimization objective functions are Thrustobj ?= 0.141_ (7.3) 300 2_ Power obj = (7.4) Table 7.3 Fixed input parameters for launch condition Free stream velocity 9 [m/s] Propeller speed 7500 [rpm] Propeller max diameter 0.25 [m] Propeller hub diameter 0.06 [m] The optimized propeller satisfies the thrust requirement by generating exactly 14 [N] of thrust. The efficiency of 65 percent seems initially low, but this is due to the characteristic of the efficiency equation which only produces high efficiencies at higher free stream velocities when propeller size is limited. Results from wind tunnel tests of a CAM 13x7 propeller reached an efficiency of 58 %. The data of the wind tunnel test are shown in Appendix H. Table 7.4 Optimized propellers performance parameter results for launch condition Case: Viscous, Compressible Thrust [N] 14.0 Power [W] 193.6 Efficiency [%] 65.1 109 Figure 7.13 Optimized launch propeller geometry: viscous, incompressible case The optimized propeller has an efficiency of 65 % which is 7% better than the propeller tested during wind tunnel experiments. The propeller which has four blades due to the thrust requirement of 14 [N] is shown in Figure 7.13. 110 Nondimensioal Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.2 0.4 0.6 0.8 1 1.2 viscous, incompressible Figure 7.14 Optimized propeller spanwise circulation distribution viscous, incompressible case The spanwise circulation of the launch propeller shown in Figure 7.14 is smooth in the tip region, but shows breaks toward the propeller root. The reason is due to either premature convergence or convergence that has not yet been reached. 111 Number of Evaluations [-] O b je c t iv e F u n c t io n [ - ] 0 10000 20000 30000 40000 0.5 1 1.5 2 2.5 viscous, incompressible Figure 7.15 Objective function convergence from real encoded GA: viscous, incompressible case The convergence of the objective function in Figure 7.15 is obtained after about 32,000 evaluations. 7.3 Dual point Optimization Cruise and Launch Condition In this case a propeller is optimized for both, the launch and the cruise conditions. The only parameter which changes is the free stream velocity which, for the launch condition is equal to the single point optimization of section 7.1.2 and, for the cruise condition is 36 [m/s] as described section 7.1.1. The first objective function (Equation 7.5) includes the 112 absolute value of required thrust for launch minus the generated propeller thrust and the second objective function (Equation 7.6) contains the cruise thrust requirement as well as the actual cruise power. The weighting of the propeller thrust and power in the second objective function (Equation 7.6) was normalized by dividing the power by 300. The first objective function did not contain any power term since it is assumed that the launch power setting is only a very short time period thus it is negligible in the overall flight endurance. launch Thrustobj ?= 0.141_ (7.5) 300 0.72_ cruise cruise Power Thrustobj +?= (7.6) The fixed input parameters shown in Table 7.5 are the same as from section 7.1.1 and 7.1.2. Table 7.5 Fixed input parameters for dual point optimization Free stream velocity 9 / 36 [m/s] Propeller speed 7500 [rpm] Propeller max diameter 0.25 [m] Propeller hub diameter 0.06 [m] The thrust required in both the launch and cruise condition is well matched. When compared to the single point optimizations the efficiencies in the dual point optimizations are lower which is expected since the propeller has to satisfy both operating conditions. Noticeable is that the geometry of the optimized two bladed propeller appears to be very 113 smooth even in the tip region where the single point optimized propellers exhibits small ripples in the geometry. Table 7.6 Optimized propeller performance parameters: cruise and launch Case: Viscous, Incompressible Launch Cruise Thrust [N] 14.01 6.97 Power [W] 373.9 396.6 Efficiency [%] 33.7 63.7 114 Figure 7.16 Optimized cruise and launch propeller; viscous, incompressible case The circulation distributions for both the launch and cruise conditions shown in Figure 7.17 have smooth distributions, with the launch case showing a higher circulation due to the larger thrust requirement. 115 Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 launch condition v=9 [m/s] cruise condition v=36 [m/s] Figure 7.17 Objective function convergences real encoded GA, cruise condition The convergence behavior of the goal for the two advance ratio optimization is shown in Figure 7.18. The objective function convergence appears to occur at about 60,000 case evaluations which are about 20,000 evaluations higher than the single point optimizations done before. 116 Number of Evaluations [-] O b j e c t i v e F unc ti on [ - ] 0 50000 100000 0.8 1 1.2 1.4 1.6 1.8 2 2.2 2.4 2.6 2.8 3 Figure 7.18 Objective function convergences real encoded GA, cruise condition 7.4 Single Point Optimization in Cruise considering Propeller Noise The first case was set to optimize the cruise propeller as described in 7.1.1 with the same fixed input parameters. The objective functions contained the total pressure perturbation (divided by 5 to obtain about equal weighting) and the thrust required. Thus the two objective functions are Thrustobj ?= 0.71_ (7.7) 117 ( ) 0.5 2_ 0 total pp obj ? = (7.8) The propeller power was not part of the optimization goal since the propeller efficiency was not of interest in this optimization. The optimized propeller satisfies the thrust requirement by producing exactly the required thrust as shown in Table 7.7. Figure 7.19 illustrates the loading shift of the circulation to the inboard section which relates to a noise reduction. Also the span of the propeller is reduced, which results in lower noise due to lower tip Mach numbers. Finally the optimization lead to an increase from two to four blades, shown in Figure 7.20, which is known to lower the overall propeller noise. All of these findings follow the conclusions made by Miller [5] who investigated propeller noise optimization. The optimized propeller efficiency, shown in Table 7.7, is low since the power minimization was not an optimization objective. Table 7.7 Optimized propellers performance parameters cruise condition with noise considered Case: Viscous, Incompressible Thrust [N] 7.0 Power [W] 509.4 Efficiency [%] 49.5 118 Nondimensioal Propeller Span [-] Ci r c u l at i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 1.2 1.3 Noise included in objective function Noise not included in the objective function Figure 7.19 Spanwise circulation distribution of the optimized propeller in cruise condition 119 Figure 7.20 Optimized cruise propeller from real coded GA; viscous, incompressible with noise considered 7.5 Single/Dual Point Optimization in Cruise/Launch Condition considering Propeller Noise and Performance The optimizations of three different propeller operational points are considered in this section. A cruise propeller at 36 [m/s] free stream velocity as described in section 7.1.1, a launch propeller at 9 [m/s] free stream velocity from section 7.1.2 and the combined case with cruise and launch condition described in section 7.1.3. For this analysis the propeller acoustic signature is included in the objective functions to reduce the overall propeller noise. The acoustic computation is done with the model described in section 2.8 without 120 the thickness noise taken into account. This is done because the thickness distribution cannot be changed during the optimization and acts only as a fixed source of noise. Additionally the thickness noise model only applies to propellers with a high aspect ratio and does not predict the thickness noise well for propellers with large chords. Thus the objective functions for the three different propeller operational points are as follows: For cruise condition: Thrustobj ?= 0.71_ (7.9) ( ) 0.3000.10 2_ _0 Power pp obj cruisetotal + ? = (7.10) For launch condition: Thrustobj ?= 0.141_ (7.11) ( ) 0.3000.10 2_ _0 Power pp obj launchtotal + ? = (7.12) For launch / cruise condition: launchcruise ThrustThrustobj ?+?= 0.140.71_ (7.13) ( ) ( ) 0.100.3000.10 2_ _0_0 cruisetotal cruise launchtotal pp Power pp obj ? ++ ? = (7.14) 121 The results from the three different optimizations are included in Table 7.8. The thrust requirements are all matched exactly and the propeller efficiencies seemed close to the ones where noise was not considered. Table 7.8 Propellers performance parameter results with noise optimization Case 1 Cruise 36 [m/s] Case 2 Launch 9 [m/s] Case 3 Launch / Cruise 9 / 36 [m/s] Thrust [N] 7.0 14.0 13.7 /6.62 Power [W] 322.9 267.8 366.8 / 350.9 Efficiency [%] 78.0 47.0 30.0 / 70.2 The geometries of the optimized propeller blades, illustrated in Figures 7.21-7.23 share similar shapes. The root chord is small and increases to the propeller mid section with a decrease in chord length towards the tip. The propeller blade geometries have a smooth shape with no geometric perturbation visible. 122 Figure 7.21 Optimized cruise propeller; viscous, incompressible case with noise and performance considered 123 Figure 7.22 Optimized launch propeller; viscous, incompressible case with noise and performance considered 124 Figure 7.23 Optimized launch / cruise propeller; viscous, incompressible case with noise and performance considered 125 Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.2 0.4 0.6 0.8 1 1.2 1.4 36 [m/s] Free Stream Velocity, Noise considered 36 [m/s] Free Stream Velocity, No noise considered Figure 7.24 Spanwise circulation distribution of the optimized propeller in cruise The optimization of the cruise propeller produced a significant reduction in the circulation strength and a shift towards the propeller mid section, which is illustrated in Figure 7.24. The noise optimized propeller shown in Figure 7.24 is normalized with respect to the propeller span of the optimized propeller where no noise was considered. Besides the reduction in circulation the noise optimized propeller shows also a smaller span which lowers the tip Mach number and as a result reduces the acoustic signature. [5] Even though the circulation distribution in the case of the noise optimized propeller is more uniform, the propeller efficiency is a little lower. The reasons for that are the 126 smaller propeller diameter which lowers the efficiency and the wider blade chord which generates higher viscous drag. Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.2 0.4 0.6 0.8 1 9 [m/s] Free stream Velocity, Noise considered 9 [m/s] Free Stream Velocity, No noise consiedered Figure 7.25 Spanwise circulation distribution of the optimized propeller at launch The circulation distribution of the launch condition propeller does not show any significant shift or reduction in the circulation. This may be explained by the significant objective function improvement toward the end of the optimization, where the design was driven by improving the efficiency and not by the acoustic propeller signature. 127 Nondimensional Propeller Span [-] Ci r c u l a t i o n [ m 2 /s ] 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1.1 0.5 1 1.5 2 36 [m/s] Free Stream Velocity, Noise considered 9 [m/s] Free Stream Velocity, Noise considered 36 [m/s] Free Stream Velocity, No noise considered 9 [m/s] Free Stream Velocity, No noise considered Figure 7.26 Spanwise circulation distribution of the optimized propeller in cruise and launch condition The last case considered both acoustic signature at launch and the acoustic signature at cruise condition. The noise optimized propeller shows a reduction in the circulation strength at both operating points as shown in Figure 7.26. The loading seemed to shift towards the tip region which is an indication that power and thrust required in the objective function are dominant. 128 Number of Generations [-] Number of Evaluations [-] O b j e c t i v e F unc ti on [ - ] 0 50 100 150 200 250 0 50000 100000 1 1.5 2 2.5 Free Stream Velocity 36 [m/s] Free Stream Velocity 9 [m/s] Free Stream Velocity 9 / 36 [m/s] Figure 7.27 Objective function convergences, cruise, launch and combined case The convergence of the objective functions of the three investigated propeller operating points is illustrated in Figure 7.27. The single point optimizations for launch and the cruise conditions show a convergence after about 30,000 evaluations, but the case considering the launch and the cruise condition together shows after an initial convergence tendency a steep gradient in the objective function around 90,000 evaluations which indicates that an optimum has not yet been achieved. 129 8 CONCLUSION & RECOMMENDATIONS The goal of this research was to develop a tool for performance and acoustic noise signature optimization applicable to small thin-airfoil propellers with complex geometries. The optimization scheme consists of a propeller performance program, an acoustic signature subroutine and a GA for the optimization. A major thrust of this work was the development of a general geometric description of propeller blades which allows for a wide design space while generating smooth, viable blade geometries. The VLM using a lifting surface scheme with vortex ring elements was used to determine the aerodynamic propeller performance parameters. This allows for the computation of swept propellers including dihedral and camber variation. All propeller blades were considered to be thin and modeled using a single layer of vortex elements to describe the geometry. The research includes an investigation to reduce computational times by finding the minimum number of required panel elements and the minimal wake extension which produces accurate propeller performance results. The performance analysis scheme was verified by comparing it with flat plate wing results as well as a straight-blade NACA109622 propeller. The implemented acoustic signature model to determine propeller noise is based on the work of Miller, [5] and during the optimization, only loading noise was considered. 130 The optimizations were performed to generate a propeller for a small electric powered UAS. Fixed input parameters such as the maximum chord, the maximum propeller span and the propeller speed were kept constant to allow for comparisons of the optimizations. Three different propeller operational conditions were considered. One launch condition with a free stream velocity of 9 [m/s], one cruise condition with 36 [m/s] free stream velocity and one case were launch and cruise condition were considered as the same time to demonstrate multiobjective optimization. The analysis included inviscid incompressible, compressible and viscous flows. Additionally, a binary and a real encoded GA have been used to investigate the method that converges quicker and which method produces better performing propellers. Results showed that even though both GA?s converged within the same number of case evaluations, the real encoded GA was generally more successful over longer optimization runs. Some optimizations show premature convergence which can be improved upon by changing ranges and limitations on GA variables. The method produced trade-offs between best performing single run and dual-point optimizations. Some of the propellers showed the development of a proplet with a loading shift to the blade tip section. Since the optimization did not include a propeller structure model which would eliminate propellers with unrealistic high bending moments, propellers with proplets are competitive candidates. Noticeable is that the circulation distributions for the single point optimizations showed unevenness whereas the case with multiple objective goals generated a smoother circulation distribution. The inclusion of noise in the objective function resulted in a reduction of the propeller radius for some cases and a loading shift of the circulation towards the propeller mid section. When multiple parameters were used in a single objective function the 131 weighting of the individual goals within the objective functions must be carried out judiciously to achieve a propeller design which satisfies all goals. Further research on the effects of weighting of individual propeller parameters within an objective function should be investigated. Recommendations ? The recommendations for improvement of this design tool are concerned mainly with premature convergence during the optimizations. Further investigations must be completed to fully understand why in some cases, the GA converges to a local and not a global minimum. ? The implemented empirical flow separation and viscous model could be replaced by a viscous/inviscid interaction method to describe flow separation and viscous effects better. ? To improve the noise prediction for blades with large chords, multiple chordwise noise sources are required. ? To allow for the simulation of thick airfoil propellers, a second vortex layer must be added to describe the bottom airfoil geometry. ? The implementation of a structure model based on the two vortex layers to describe the propeller geometry would generate a realistic thickness distribution. ? The extension of the aerodynamic propeller performance program to include a duct around the propeller would allow for the optimization of ducted propeller/fan combinations. ? Finally, for highly loaded and low advance ratio conditions, a free wake model can improve the aerodynamic prediction. 132 REFERENCES [1] Lynam, F. J. H. and Webb, H. A. ?the Emission of Sound by Airscrews.? Technical Report of the Advisory Committee if Aeronautics for the Year 1918-1919 Vol. 2 (London) Technical Reports and Memoranda No. 624, March 1919, pp. 792-801. [2] J. Grasmeyer, ?Application of a Genetic Algorithm with Adaptive Penalty Functions to Airfoil Design?, AIAA 1997-7, Aerospace Science Meeting Exhibit, 35 th , Reno, NV, Jan. 6-9, 1997 [3] Fanjoy, D. W., Crossley W. A., ?Aerodynamic Shape Design for Rotor Airfoils via Genetic Algorithm?. American Helicopter Society 53 rd Annual Forum, Virginia Beach, VA, April 19-may 1, 1998. [4] Gardner A. B., Selig M. S., ?Airfoil Design using a Genetic Algorithm and an Inverse Method?. AIAA 2003-43, 41 st Aerospace Science Meeting and Exhibit, 6-9 January, 2003. [5] Miller C. J., ?Optimally Designed Propellers Constrained by Noise?, Purdue University, Dissertation, December 1984. [6] Powell, M. J. D., ?An Efficient Method for Finding the Minimum of a Function of Several Variables without Calculating Derivatives,? Computer Journal, Volume 7, 1964, pp. 155-162. [7] Anders Smaerup Olsen, ?Optimization of Propellers using the Vortex Lattice Method?, Technical University of Denmark, Dec 2001. 133 [8] Hampsey M., ?Multiobjective Evolutionary Optimization of Small Wind Turbine Blades.? Dissertation, Department of mechanical Engineering, University of Newcastle, 2002. [9] Rankine, W. J. M., ?On the Mechanical Principles of the Action of propellers?, Transactions of the Institute of Naval Architects, Volume 6, 1865. [10] Froude R. E., ?On the Part Played in Propulsion by Difference in Fluid Pressure?, Transactions of the Institute of Naval Architects, Volume 30, 1889. [11] Vogeley A. W., ?Axial Momentum Theory for Propellers in Compressible Flow?, NACA TN No. 2164, August, 1950. [12] Froude W., ?On the elementary Relation between Pitch, Slip and Propulsive Efficiency?, Transaction of the Institute of Naval Architects, Volume 19, 1878. [13] Drzewiecki S., ?Theorie Generale de l?H?eleice?, Paris, 1920. [14] Anderson J.D., Jr. ?Fundamentals of Aerodynamics?, 3rd Edition, McGraw-Hill, New York 2001. [15] Glauert H., ?An Aerodynamic Theory of the Airscrew?, Grear Britain Aeronautical Research Committee, Reports & Memoranda No. 786, 1922. [16] Katz J., Plotkin A., ?Low Speed Aerodynamics? 2nd Edition, Cambridge University Press, New York, 2005. [17] Prandtl L., ?Tragfluegel Theorie? 1 und 2 Mitteilung. Nachrichten von der Koeniglichen Gesells. D. Wissens. Math.-Phys. Klasse 451 (1918) und 107 (1919). Reprinted in Abhandlungen zur Hydrodynamik und Aerodynamik. Goettingen 9(1927). [18] Kuethe A. M., Chow C. Y., ? Foundations of Aerodynamics, Bases of Aerodynamic Design? 5th Edition, John Wiley & Sons, Inc., New York, 1998 134 [19] Purvis J. W. and Burkhalter J. E., ?Simplified Solution of the Compressible Lifting Surface Problem?, AIAA Journal, Volume 20, No. 5, Page 589, May 1982. [20] Burkhalter J. E., Harris M., ?Non-linear aerodynamic analysis of grid fin configurations?, AIAA Applied Aerodynamics Conference, AIAA-1995- 1894, 13th, San Diego, CA, June 19-22, 1995. [21] Brooks R. A., Burkhalter J. E., ?Experimental and analytical analysis of grid fin configurations?, Journal of Aircraft, 0021-8669, vol.26, No.9 (885-887), 1989. [22] Burkhalter J. E. (U.S. Air Force Academy, Colorado Springs, CO), ?Downwash measurements on a pitching canard-wing configuration?, Journal of Aircraft, 0021-8669, vol.30, No.6 (1005-1008), 1993. [23] Shiu Wing Cheung, R. ?Numerical Prediction of Propeller Performance by Vortex Lattice Method?, UTIAS Technical Note No.265 Nov 1987. [24] Sofia Werner, ?Application of the Vortex lattice Method to Yacht Sails?, University of Aukland, Master Thesis, July 2001. [25] Goldstein, S., ?On the Vortex Theory of Screw Propellers?, Proceedings of the Royal Society of London, Series A, Vol. 123, 1929. [26] Christopher J. Fisichella, ?A case for improving prescribed wake analysis?, AIAA ASME Wind Energy Symposium, 19th, Aerospace Science Meeting and Exhibit, 38th NV, Jan 10-13 2000. [27] A. J. Chorin and P.S. Bernard, ?Discretization of a vortex sheet, with an example of roll-up, J. Comp. phys, (13):423, 1973. [28] S.P. Fiddes, and J. Gaydon, ?A vortex lattice method for calculating the flow past yacht sails?, Journal of Wind Engineering and Industrial Aerodynamics, 63:35-59, 1996. 135 [29] Karamcheti K. ?Principles of Ideal-Fluid Aerodynamics?, Rep. Edition, John Wiley & Sons, Inc., New York, 1980. [30] D. F. Thrasher, D. T. Mook, and A. H. Nayfeh, ?A computer based method for analyzing the flow over sails?, Chesapeake Sailing Yacht Symposium, pages119-127, SNAME, 1983. [31] J. J. Bertin and M. L. Smith, ?Aerodynamics for Engineers?, Prentice-Hall, Inc., Englewood Cliffs, 2nd Edition, 1989. [32] A. H. Day, ?Steps toward an optimal yacht sail?, Technical Report, The Royal Institution of Naval Architecture, 1993. [33] J. L. Guermond, ?Collocation methods and lifting surfaces?, European Journal of mechanics, B/Fluids, 8(4):283-305, 1989. [34] Sobieczky, H., ?Parametric Airfoils and Wings?, Notes on Numerical Fluid Mechanics, Vol. 68, pp.71-88, Vieweg Verlag, 1988. [35] Samareh, J.A., ?Survey of Shape Parameterization Techniques for High-Fidelity Multidisciplinary Shape Optimization?, AIAA Journal Vol. 39, No.5, May 2001. [36] Robinson, G. M., Keane, A. J., ?Concise Orthogonal Representation of Supercritical Airfoils?, Journal of Aircraft, Vol. 38, No.3. [37] Song, W., Keane, A. J., ?A Study of Shape Parameterization Airfoil Optimization?, AIAA-2004-448 10th AIAA, ISSMO Multidisciplinary Analysis and Optimization Conference, Albany, New York, Aug. 30-01, 2004. [38] Padula, S., Li, W., ?Options for Robust Airfoil Optimization Under Uncertainty?, 9th AIAA Multidisciplinary Analysis and Optimization Symposium, 4-6 September2002. 136 [39] Hicks, R. M., Henne, P.A.,? Wing design by numerical optimization?, Journal of Aircraft, Vol.15, pp. 407-412, 1978. [40] Kulfan B., M., Bessoletti J., E., ?Fundamental Parametric Geometry Representation for Aircraft Component Shapes?, 11th AIAA/ISSMO Multidisciplinary Analysis and Optimization Conference, 6-6 September 2006, Portsmouth, Virginia. [41] Kulfan B., M., ?A Universal Parametric Geometry Representation Method ? CST?, 45th AIAA Aerospace Sciences Meeting and Exhibit, 8-11 January, Reno, Nevada, 2007. [42] Gennaretti M. and Morino L., ?A Boundary Element Method for the Potential, Compressible Aerodynamics of Bodies in Arbitrary Motion?, Aeronautical Journal, Vol. 96, Jan. 1992, pp15-19. [43] Curri I. G., ?Fundamental Mechanics of Fluids?, McGraw-Hill, 1993. [44] Fung Y. C., ?An Introduction to the Theory of Aeroelasticity?, John Wiley & Sons, Inc., New York, NY, 1955. [45] C. J. Szymendera, ?Computational Free Wake Analysis of a Helicopter Rotor?, Master Thesis, Pennsylvania State University, 2002. [46] Tromp, Jeffrey C., ?Numerical Study of Three Viscous/Inviscid Interaction Methods?, ADA202575, Defense Technical Infromation Center, Sep. 1988. [47] Cebeci, T. and Besnard, E. ?An efficient and accurate approach for analysis and design of high lift configurations.? Canadian Aeronautics and Space Journal, 44(4):1-17 1998. [48] John Dreese, http://www.dreesecode.com/primer/airfoil5.html, 2004 137 [49] Borst, H. V. and Associates, ?Summary of Propeller Design.? Summary of Propeller Design Procedures and Data, Volume 1 Aerodynamic Design and Installation, AD-774 831, NTIS, 1973. [50] Leishman, J. G., ?Challenges in Modeling the Unsteady Aerodynamics of Wind Turbines,? Wind Energy, Vol. 5, 2002, pp. 86?132. [51] Snel, H., ?Review of the Wind Turbine Aerodynamics,? Wind Energy, Vol. 6, No. 3, July/September, 2003, pp. 203?211. [52] Hansen A. C. and Butterfield C. P., ?Aerodynamics of Horizontal-Axis Wind Turbines,? Annual Review of Fluid Mechanics Vol. 25, 1993, pp. 115?149. [53] Fingersh, L. J., Simms, D., Hand, M., Jager, D., Cotrell, J., Robinson, M., Schreck, S. and Larwood, S. ?Wind Tunnel Testing of NREL?s Unsteady Aerodynamics Experiment,? In proceedings of 39th Aerospace Sciences Meeting, AIAA Paper 2001?0035, Reno, NV, Jan. 8?11, 2001. [54] Corrigan, J. J., and Schilings, J. J., ?Empirical Model for Stall Delay Due to Rotation,? American Helicopter Society Aeromechanics Specialists Conference, San Francisco, CA, 1994. [55] Du, Z. H. and Selig, M. S., ?A 3-D Stall Delay Model for Horizontal Axis Wind Turbine Performance Prediction, ?ASME/AIAA Joint Wind Energy Symposium, AIAA Paper 1998?0021, Reno, NV, Jan. 1998, pp. 9?19. [56] Pierce, K. and Hansen, A. C., ?Prediction of Wind Turbine Rotor Loads Using the Beddoes?Leishman Model for Dynamic Stall,? Journal of Solar Energy Engineering, Vol. 117, 1995, pp. 200?204. [57] Beddoes, T. S., ?Representation of Aerofoil Behaviour,? Vertica, Vol. 7, 1983, pp. 183?197. [58] Leishman, J. G., and Beddoes, T. S., ?A Semi-Empirical Model for Dynamic Stall,? Journal of American Helicopter Society, Vol. 34, No. 3, 1989, pp. 3? 17. 138 [59] Schreck, S. J., and Robinson, M. C., ?Dynamic Stall and Rotational Augmentation in Recent Wind Turbine Aerodynamic Experiments,? In Proceedings of 32nd AIAA Fluid Dynamics Conference and Exhibit, St. Louis, MS, June 24?26, 2002. [60] James L. Tangler, J. David Kocurek, ?Wind Turbine Post-Stall Airfoil Performance Characteristics Guidelines for Blade-Element Momentum Methods?, 43rd AIAA Aerospace Sciences Meeting and Exhibit 10-13 January 2005, Reno, Nevada. [61] Viterna, L .A., and Corrigan, R. D., ?Fixed Pitch Rotor Performance of Large Horizontal Axis Wind Turbines, ?DOE/NASA Workshop on Large Horizontal Axis Wind Turbines, Cleveland, Ohio, July 1981. [62] Viterna, L .A., and Janetzke, D. C., ?Theoretical and Experimental Power from Large Horizontal-Axis Wind Turbines,? NASA TM-82944, Sept 1982. [63] Robert E. Sheidahl R., E., Klimas P., C., ?Aerodynamic Characteristics of Seven Symmetrical Airfoil Sections through 180-Degree Angle of Attack for use in Aerodynamic Analysis of Vertical Axis Wind Turbines?. Sandia National Laboratories, SAND80-2114, 1981. [64] F. B. Metzger, ?A Review of Propeller Noise Prediction Methodology 1919- 1994?, NASA Contractor Report 198156, National Aeronautics and Space Administration, Langley Research Center, June 1995. [65] Kenneth R. Fidler, ?Subsystem acoustic testing of a Vertical Takeoff and Landing Ducted Propeller Unmanned Aerial Vehicle? Technical Report, AMR-SS-04- 05. [66] Lowson, M. V., ?The sound field for singularities in motion.? Proceedings of the Royal Society (London), Volume A286, 1965, pp. 559-572. [67] Succi G. P., ?Noise and Performance of General Aviation Aircraft: A Review of the MIT Study,? SAE Technical Paper No. 810586, April 1981. 139 [68] Succi G. P., Munro, D. H., Zimmer, J. A., Dunbeck, P. D., Lasabee E. E., K. U. and Kerrenbrock J. L., ?Noise and Performance of Propellers for Light Aircraft,? NASA CR-165732, 1981. [69] Succi G. P., Munro, D. H. and Zimmer J. A., ?Experimental verification of Propeller Noise Prediction,? AIAA Journal, Volume 2, No. 11, pp.1483- 1491., November 1982. [70] Sullivan J. P., ?The Effect of Blade Sweep on Propeller Performance,? AIAA paper No. 77-716, 1977. [71] Anderson, Jr. J. D., ?Computational Fluid Dynamics, The basics with Applications?. McGraw-Hill, Inc., New York, 1995. [72] Th. E. Labrujere and P.J. Zandenbergen. ?On the application of a new version of lifting surface theory to non-slender and kinked wings?. Journal of Engineering mathematics, 7:85-96, 1973. [73] Gabasov, R.F. and Kirillova, F.M., Methods of Optimization, Optimization Software Inc., New York, NY, 1988. [74] Weber, D.C. and Skillings, J.H., A first Course in the Design of Experiments: A linear Models Approach, CRC Press, Boca Raton, FL, 2000. [75] Kennedy, J. Eberhart, R. ?Particle Swarm Optimization?, Proceedings of the IEEE International Conference on Neural Networks, Perth, Australia 1995, pp 1942-1945. [76] Fukuyama, Y., Takayama, S., Nakanishi, Y., and Yoshida, H. ?A particle swarm optimization for reactive power and voltage control in electric power systems. Proceedings of the Genetic and Evolutionary Computation?, Conference 1999 (GECCO 1999), Orlando, Florida, USA. pp. 1523-1528, 1999. 140 [77] Carlisle, A. and Dozier, G. ?Adapting particle swarm optimization to dynamic environments?. Proceedings of International Conference on Artificial Intelligence (ICAI 2000), Las Vegas, Nevada, USA. pp. 429-434, 2000. [78] Venter, G. and Sobieszczanski-Sobieski, J., "Multidisciplinary optimization of a transport aircraft wing using particle swarm optimization," Structural and Multidisciplinary Optimization, pp. preprint, 2004. [79] Onwubolu, G. C. and Clerc, M., "Optimal path for automated drilling operations by a new heuristic approach using particle swarm optimization," International Journal of Production Research, vol. 4 pp. 473-491, 2004. [80] Zhou, Y., Zeng, G., and Yu, F., "Particle swarm optimization-based approach for optical finite impulse response filter design," Applied Optics, vol. 42, no. 8, pp. 1503-1507, Mar.2003. [81] Zheng, Y., Ma, L., Zhang, L., and Qian, J. Robust PID controller design using particle swarm optimizer. Proceedings of IEEE International Symposium on Intelligence Control 2003, pp. 974-979, 2003. [82] Jenigiri, S. ?Comparative study of efficiency of genetic algorithms and particle swarm optimization technique to solve permutation problems?. 1996. [83] Angeline, P. J. Evolutionary optimization versus particle swarm optimization: philosophy and performance differences. Evolutionary Programming VII: Proceedings of the Seventh Annual Conference on Evolutionary Programming, 1998. [84] Eberhart, R. C. and Shi, Y. Comparison between genetic algorithms and particle swarm optimization. Evolutionary Programming VII: Proceedings of the Seventh Annual Conference on Evolutionary Programming, San Diego, CA. 1998. [85] Mitchell, M., An Introduction to Genetic Algorithms, MIT Press, Cambridge, MA, 1998. 141 [86] Coley, D.A., An Introduction to Genetic Algorithms for Scientists and Engineers, World Scientific Publishing Co., Singapore, 1999. [87] Anderson, M.B., ?Design of a Missile Interceptor Using a Genetic Algorithm?, Doctoral Dissertation, Auburn University, 1998. [88] Dyer John, personal communication, Auburn University, 2007 142 APPENDIX A: BIOT-SAVART LAW The following describes the derivation of the Biot-Savart law. Biot-Savart Law: The continuity equation in non-conservation form is 0 q =???+? ? Dt D (A.1) with ? being the density and q the velocity. For incompressible fluids (?=const.) the continuity equation reduces to 0q 0q = ? ? + ? ? + ? ? =?? =?? z w y v x u (A.2) Angular velocity and vorticity The motion of a fluid element consists of translation, rotation, and deformation. The angular velocity of a fluid element is due to the velocity variation within the element. As a result the fluid element may deform and rotate. The derivation of the angular velocity is described in Ref. [16]. Thus the vector form of the angular velocity is 143 V???= 2 1 ? (A.3) Vorticity ? is defined as twice the angular velocity. curlVV =??=?? ?2? (A.4) To obtain the velocity field as a function of the vorticity distribution, Equation (A.4) must be inverted. As described in Ref. [16] the velocity field may be expressed as the curl of a vector field B, such that Bq ??= (A.5) Since the curl of a gradient vector field is zero, B is indeterminant to within the gradient of a scalar function of position and time, and B can be selected such that 0=?? B (A.6) The vorticity then becomes ( ) ( ) BBB 2 q? ?????=????=??= (A.7) The application of Equation A.6 reduces this to Poisson?s equation for the vector potential B: B 2 ? ??= (A.8) The solution of this equation, using Green?s theorem (see Ref. [29]) is 144 dVB V ? ? = 10 rr ? 4 1 ? (A.9) Here B is evaluated at point P (which is a distance r 0 from the origin, shown in Figure A.1) and is a result of integrating the vorticity ? (at point r 1 ) within the volume V. The velocity field is then the curl of B: dV V ? ? ??= 10 rr ? 4 1 q ? (A.10) Figure A.1 Velocity at point P due to a vortex distribution Before proceeding with the integration, an infinitesimal piece of the vorticity filament ?, as shown in Figure A.2 is considered. The cross-sectional area dS is selected such that it P . r 0 -r 1 ? dV V r 1 r 0 145 is normal to ?, and the direction dl on the filament is dld ? ? l = (A.11) Also, the circulation ? is dS?=? (A.12) and dldSdV = (A.13) Figure A.2 Induced velocity at point P by a vortex segment so that dl ds ? P r 0 r 1 Origin r 0 ? r 1 ? 146 1010 rr l rr ? ? ???= ? ?? d dV (A.14) Carrying out the curl operation while keeping r 1 and dl fixed we get ( ) 3 10 10 10 rr rrl rr l ? ?? ?= ? ??? d d (A.15) Substitution of this result back into Equation A.6 results in the Biot-Savart law, which states ( ) ? ? ?? ? = 3 10 10 rr rrl 4 q d ? (A.16) or in differential form ( ) 3 10 10 rr rrl 4 q ? ?? ? =? d ? (A.17) 147 APPENDIX B: PROPELLER GA VARIABLES Shown below are the propeller variables which define the propeller geometry, the number of propeller blades as well as the max propeller span. If a variable pitch propeller is considered, the optional variables, which offset the angle of attack by a constant value are activated. Propeller GA Variables pr xray(1) Propeller span nxbl xray(2) Number of propeller blades ua_1 xray(3) BP Shape Function Parameter airfoil upper side ua_2 xray(4) BP Shape Function Parameter airfoil upper side ua_3 xray(5) BP Shape Function Parameter airfoil upper side ua_4 xray(6) BP Shape Function Parameter airfoil upper side ua_5 xray(7) BP Shape Function Parameter airfoil upper side ua_6 xray(8) BP Shape Function Parameter airfoil upper side ua_7 xray(9) BP Shape Function Parameter airfoil upper side ua_8 xray(10) BP Shape Function Parameter airfoil upper side ua_9 xray(11) BP Shape Function Parameter airfoil upper side ua_10 xray(12) BP Shape Function Parameter airfoil upper side ua_11 xray(13) BP Shape Function Parameter airfoil upper sid ua_12 xray(14) BP Shape Function Parameter airfoil upper side ua_13 xray(15) BP Shape Function Parameter airfoil upper side ua_14 xray(16) BP Shape Function Parameter airfoil upper side ua_15 xray(17) BP Shape Function Parameter airfoil upper side ua_16 xray(18) BP Shape Function Parameter airfoil upper side ua_17 xray(19) BP Shape Function Parameter airfoil upper side ua_18 xray(20) BP Shape Function Parameter airfoil upper side ua_19 xray(21) BP Shape Function Parameter airfoil upper side ua_20 xray(22) BP Shape Function Parameter airfoil upper side fau_1 xray(23) Class Function upper Airfoil Parameter fau_2 xray(24) Class Function upper Airfoil Parameter aoa_1 xray(25) BP Angle of Attack Function Parameter aoa_2 xray(26) BP Angle of Attack Function Parameter 148 aoa_3 xray(27) BP Angle of Attack Function Parameter aoa_4 xray(28) BP Angle of Attack Function Parameter aoa_5 xray(29) Angle of Attack max min value chord_1 xray(30) BP Chord Length Function Parameter chord_2 xray(31) BP Chord Length Function Parameter chord_3 xray(32) BP Chord Length Function Parameter chord_4 xray(33) BP Chord Length Function Parameter chord_5 xray(34) Chord Length max min sweep_1 xray(35) BP Propeller Sweep Function Parameter sweep_2 xray(36) BP Propeller Sweep Function Parameter sweep_3 xray(37) BP Propeller Sweep Function Parameter sweep_4 xray(38) BP Propeller Sweep Function Parameter sweep_5 xray(39) Propeller Sweep in multiples of max chord length cprp xray(40) Propeller rotational point 149 APPENDIX C: PROPELLER WAKE INVESTIGATION Propeller Wake Investigation of a two and four Bladed Propeller Below are shown the results of the thrust and power obtained from the propeller performance program based on the wake extension. The propeller speed was 5800 rpm, the free stream velocity V inf = 21.0 [m/s] and propeller diameter was d=0.32 [m]. Wake extenstion 2 bladed propeller 4 bladed propeller in revolutions Thrust [N] Power [W] Thrust [N] Power [W] 0.125 3.97 104.07 7.70 203.43 0.250 3.93 103.22 7.45 198.42 0.500 3.87 102.01 7.22 193.83 1.000 3.83 101.28 7.10 191.36 2.000 3.82 100.99 7.05 190.48 3.000 3.81 100.94 7.04 190.23 4.000 3.81 100.92 7.03 190.16 5.000 3.81 100.91 7.03 190.12 150 APPENDIX D: PROPELLER PANEL DISCRETIZATION Propeller Panel Discretization and Panel Density Below are the results of the panel discretization and panel density investigation. This was done to determine the least number of panels which still produce accurate aerodynamic performance prediction in an effort to reduce computation times during the optimization process. The investigation contains two panel arrangement types. The first one uses a full cosine distribution of the panels in the chordwise direction. This assures that the panels are denser at the leading and training edge. The second panel arrangement has besides the full cosine distribution of the panels in the chordwise direction a half cosine distribution in the spanwise direction. In addition to that the Table below shows also the effects on the thrust when the number of spanwise and chordwise panels are changed. Full cosine panel distribution in the chordwise direction Full cosine panel distribution in the chordwise and half cosine panel distribution in the spanwise direction Panel Numbers Chordwise x Spanwise Thrust [N] Panel Numbers Chordwise x Spanwise Thrust [N] 3x20 2.8144 7x7 3.9418 4x20 2.9766 7x10 3.8095 5x20 3.3385 7x12 3.7544 6x20 3.6052 7x14 3.7281 7x5 3.936 7x16 3.7195 7x10 3.8262 7x18 3.7183 7x20 3.8162 7x20 3.7203 7x25 3.8152 7x24 3.7269 7x30 3.8145 7x26 3.7316 151 APPENDIX E: NACA 109622 PROPELLER GEOMETRY NACA 109622 Propeller Geometry Data The propeller performance program was validated using a NACA 109622 straight blade propeller. The data below are from Ref. [23], where a three bladed propeller was simulated with a blade angle of ? .75 =45.4?.at the 75% spanwise position. The x, y and z coordinates are given of the spanwise geometric description with C/r being the section chord to span ratio. x y z c/r ? 0.0 0.275 0.0 0.2415 1.15581 0.0 0.3475 0.0 0.2415 1.07155 0.0 0.42 0.0 0.2415 1.0034 0.0 0.4925 0.0 0.2415 0.94246 0.0 0.565 0.0 0.2415 0.88481 0.0 0.6375 0.0 0.2415 0.83769 0.0 0.71 0.0 0.2415 0.77661 0.0 0.7825 0.0 0.2415 0.74869 0.0 0.855 0.0 0.2415 0.70855 0.0 0.9275 0.0 0.2415 0.67713 0.0 1 0.0 0 0 152 The propeller geometry is generated based on the data from Table above. The propeller has 10 spanwise and two chordwise panels. Thus the individual nodal and collocation points of the first propeller blade are Nodal Points Leading Edge Center Trailing Edge x y z x y z x y z 0.0000 0.2750 0.0000 0.1113 0.2750 0.0449 0.2225 0.2750 0.0899 0.0000 0.3475 0.0000 0.1071 0.3475 0.0542 0.2142 0.3475 0.1083 0.0000 0.4200 0.0000 0.1032 0.4200 0.0613 0.2063 0.4200 0.1226 0.0000 0.4925 0.0000 0.0992 0.4925 0.0675 0.1984 0.4925 0.1350 0.0000 0.5650 0.0000 0.0952 0.5650 0.0731 0.1903 0.5650 0.1462 0.0000 0.6375 0.0000 0.0916 0.6375 0.0775 0.1832 0.6375 0.1550 0.0000 0.7100 0.0000 0.0867 0.7100 0.0829 0.1734 0.7100 0.1659 0.0000 0.7825 0.0000 0.0844 0.7825 0.0853 0.1687 0.7825 0.1707 0.0000 0.8550 0.0000 0.0809 0.8550 0.0886 0.1618 0.8550 0.1773 0.0000 0.9275 0.0000 0.0781 0.9275 0.0911 0.1561 0.9275 0.1823 0.0710 1.0000 0.0700 0.0779 1.0000 0.0913 0.0860 1.0000 0.1000 153 Collocation Points Leading Edge Panels Trailing Edge Panels xc yc zc xc yc zc 0.0556 0.3113 0.0225 0.1669 0.3113 0.0674 0.0535 0.3838 0.0271 0.1606 0.3838 0.0812 0.0516 0.4563 0.0307 0.1547 0.4563 0.0920 0.0496 0.5288 0.0337 0.1488 0.5288 0.1012 0.0476 0.6013 0.0365 0.1428 0.6013 0.1096 0.0458 0.6738 0.0387 0.1374 0.6738 0.1162 0.0434 0.7463 0.0415 0.1301 0.7463 0.1244 0.0422 0.8188 0.0427 0.1266 0.8188 0.1280 0.0404 0.8913 0.0443 0.1213 0.8913 0.1330 0.0390 0.9450 0.0456 0.1171 0.9450 0.1367 154 APPENDIX F: NACA 0009 AIRFOIL DATA The lift and drag coefficient data of a NACA 0009 airfoil section. [63] are taken to generate the empirical stall model for the propeller. 155 156 APPENDIX G: PROPELLER GEOMETRY OF SKYHAWK 172 The geometric description of the ? scale Cessna 172 Skyhawk propeller as described by Succi [69] and the generated geometry for the numerical verification of the test case is given in the following table. Geometric propeller description by Succi: Radius r/R Chord C/R Thickness t/C Pitch ? 0.133 0.149 0.43 31.5 0.244 0.153 0.204 30.9 0.32 0.153 0.156 29.5 0.4 0.154 0.127 27.1 0.48 0.153 0.109 24.6 0.64 0.141 0.092 20.3 0.8 0.119 0.092 17.3 0.88 0.103 0.08 16.1 0.96 0.079 0.08 15.1 1 0 0 14.3 157 The collocation pints for the for the Cessna 172 1C160 blade model are x-Control Point y-Control Point z-Control Point 0.00471 0.03879 0.00773 0.00473 0.05218 0.00786 0.00469 0.06346 0.00795 0.00459 0.07263 0.00801 0.00447 0.08204 0.00809 0.00431 0.09169 0.00821 0.00413 0.10135 0.00830 0.00394 0.11100 0.00836 0.00361 0.12548 0.00829 0.00317 0.14478 0.00809 0.00274 0.16408 0.00770 0.00233 0.18339 0.00714 0.00203 0.19787 0.00663 0.00182 0.20752 0.00619 0.00160 0.21717 0.00563 0.00136 0.22682 0.00494 0.00104 0.23406 0.00389 0.00064 0.23889 0.00247 0.01412 0.03879 0.02318 0.01419 0.05218 0.02357 0.01407 0.06346 0.02385 0.01378 0.07263 0.02402 0.01340 0.08204 0.02428 0.01293 0.09169 0.02463 0.01240 0.10135 0.02491 0.01182 0.11100 0.02509 0.01084 0.12548 0.02488 0.00951 0.14478 0.02426 0.00822 0.16408 0.02310 0.00699 0.18339 0.02142 0.00609 0.19787 0.01990 0.00547 0.20752 0.01858 0.00480 0.21717 0.01689 0.00408 0.22682 0.01483 0.00311 0.23406 0.01167 0.00192 0.23889 0.00740 158 The nodal pints for the for the Cessna 172 1C160 blade model are x-Nodal Point y-Nodal Point z-Nodal Point 0.00000 0.03209 0.00000 0.00000 0.04549 0.00000 0.00000 0.05888 0.00000 0.00000 0.06805 0.00000 0.00000 0.07722 0.00000 0.00000 0.08687 0.00000 0.00000 0.09652 0.00000 0.00000 0.10617 0.00000 0.00000 0.11582 0.00000 0.00000 0.13513 0.00000 0.00000 0.15443 0.00000 0.00000 0.17374 0.00000 0.00000 0.19304 0.00000 0.00000 0.20269 0.00000 0.00000 0.21234 0.00000 0.00000 0.22200 0.00000 0.00000 0.23165 0.00000 0.00000 0.23647 0.00000 0.00000 0.24130 0.00000 0.00939 0.03209 0.01533 0.00944 0.04549 0.01558 0.00948 0.05888 0.01584 0.00929 0.06805 0.01595 0.00909 0.07722 0.01607 0.00878 0.08687 0.01631 0.00846 0.09652 0.01654 0.00807 0.10617 0.01667 0.00768 0.11582 0.01678 0.00677 0.13513 0.01639 0.00590 0.15443 0.01596 0.00505 0.17374 0.01485 0.00427 0.19304 0.01371 0.00385 0.20269 0.01283 0.00345 0.21234 0.01194 0.00295 0.22200 0.01057 0.00248 0.23165 0.00920 0.00167 0.23647 0.00636 0.00089 0.24130 0.00351 0.01879 0.03209 0.03066 0.01887 0.04549 0.03117 0.01896 0.05888 0.03168 0.01857 0.06805 0.03191 0.01818 0.07722 0.03213 0.01756 0.08687 0.03261 159 0.01693 0.09652 0.03308 0.01615 0.10617 0.03333 0.01537 0.11582 0.03357 0.01355 0.13513 0.03278 0.01180 0.15443 0.03191 0.01011 0.17374 0.02970 0.00854 0.19304 0.02742 0.00770 0.20269 0.02565 0.00689 0.21234 0.02388 0.00591 0.22200 0.02115 0.00497 0.23165 0.01840 0.00334 0.23647 0.01272 0.00179 0.24130 0.00701 160 APPENDIX H: WIND TUNNEL DATA CAM 13x7 PROPELLER Propeller wind tunnel data for a CAM 13x7 propeller. The propeller hub allowed for blade pitch changes to maximize the thrust. Motor efficiencies are about 72% for a power input of 67 [W]. The first Table below contains thrust and torque data of the CAM 13x7 propeller at 60 feet free stream velocity and a constant power input of 67 Watts, while changing blade pitch. Pitch change 60 feet/sec free stream velocity Thrust [lbf] Torque [lb inch] 0.032 0.685 0.282 0.882 0.335 0.987 0.345 1.035 0.348 1.072 0.37 1.24 0.362 1.308 0.345 1.369 0.339 1.534 0.27 1.565 The Table below contains the thrust data of the 13x7 propeller in a fixed pitch mode as it changes with increased free stream velocity and with a constant power input of 67 Watts. Free Stream Velocity [ft/sec] Fixed pitch propeller thrust [lbf] 0 1.07 20 0.92 30 0.7 40 0.524 50 0.26 60 0.22