Experimental and Computational Investigation of a Turbofan Inlet Duct by Zachary Mitchel Hall A thesis submitted to the Graduate Faculty of Auburn University in partial fulfillment of the requirements for the Degree of Master of Science Auburn, Alabama December 18, 2009 Approved by Anwar Ahmed, Chair, Professor of Aerospace Engineering Roy Hartfield, Professor of Aerospace Engineering Andrew Shelton, Assistant Professor of Aerospace Engineering ii Abstract The flow field of an axisymmetric turbofan inlet was investigated experimentally and computationally. Flow visualization results were obtained using LIF, hydrogen bubbles, and PIV. Incidence angle was determined to be a primary factor in the variation of pressure coefficients near the leading edge of the inlet geometry. Formation of an inlet vortex was observed when the inlet was in the proximity of the ground. An optimization study for a non-axisymmetric turbofan inlet duct using Genetic Algorithms and Computational Fluid dynamics was also conducted. The CFD model included flow conditions of the CFM56-5B turbofan inlet duct geometry to provide baseline results for comparison with the optimized geometry. A grid refinement study was conducted to ensure accuracy of results. This grid generation and computational model was used in the optimization study. iii Acknowledgements The author would like to thank Dr. Anwar Ahmed for instilling an interest in aerodynamics and for guidance and support throughout the course of this effort. The author would also like to thank his committee: Dr. Roy Hartfield for his guidance in aerospace system optimization methods and Dr. Andrew Shelton for his support in computational aerodynamics. Thanks are also due to Vivek Ahuja for his help with the development of the optimization code used in this effort and Hamza Ahmed for his training in experimental methods. The author would especially like to thank his fianc? for her enduring support and patience during this effort. iv Table of Contents List of Figures ................................................................................................................... vii List of Tables ..................................................................................................................... xi Nomenclature .................................................................................................................... xii Introduction ......................................................................................................................... 1 Methodology ..................................................................................................................... 19 1 Experimental ...................................................................................................... 19 1.1 Model Geometry ............................................................................... 19 2 Experimental Data Acquisition .......................................................................... 22 2.1 Test Facilities .................................................................................... 22 2.2 Experimental Setup ........................................................................... 22 2.3 Tests Performed ................................................................................ 24 2.3.1 Flow Visualization ........................................................................... 24 2.3.2 Particle Image Velocimetry (PIV) ................................................... 25 3 Computational Fluid Dynamics ......................................................................... 28 3.1 Grid Generation ................................................................................ 28 3.2 CFD Solver Model ............................................................................ 32 3.3 Boundary Layer Resolution: Near ? Wall Treatment ...................... 35 3.4 Flight Conditions .............................................................................. 37 v 3.5 Grid Refinement................................................................................ 38 3.5.1 y + Sensitivity Study ......................................................................... 38 3.5.2 Far-Field Boundary Refinement ...................................................... 47 3.5.3 Grid Spacing Refinement ................................................................ 51 4 Inlet Optimization .............................................................................................. 56 4.1 Genetic Algorithms ........................................................................... 56 4.1.1 Tournament Selection ...................................................................... 56 4.2 Integrated CFD and GA Networks ................................................... 58 4.3 Optimization Goals ........................................................................... 60 4.4 Design Space ..................................................................................... 63 Results & Discussion ........................................................................................................ 66 5 Experimental Results.......................................................................................... 66 5.1 Effect of Reynolds Number .............................................................. 70 5.2 Effect of Angle of Incidence ............................................................. 74 5.3 Formation of an Inlet Vortex ............................................................ 78 6 Computational Fluid Dynamic Simulations ....................................................... 82 6.1 Simulation of CFM56-5B ................................................................. 82 6.2 Optimized Geometry ......................................................................... 85 Conclusions ....................................................................................................................... 89 References ......................................................................................................................... 90 Appendices ........................................................................................................................ 93 1 CST Transformation & Bernstein Polynomials ................................................. 93 vi 2 Example GA Curve Fitting Design Space.......................................................... 99 3 Example Curve Fitting of Bernstein Polynomials to CFM56-5B Cowl .......... 100 4 PIV Uncertainty Analysis................................................................................. 103 vii List of Figures Figure 1. Inlet Locations ? Integrated engines [1] .............................................................. 2 Figure 2. Inlet Locations ?Podded engines [1] ................................................................... 3 Figure 3. Thin lip inlet auxiliary devices at low speed [5] ................................................. 4 Figure 4. Subsonic Inlet Nomenclature [1] ......................................................................... 5 Figure 5. Typical Streamline Patterns for Subsonic Inlets [4] ............................................ 7 Figure 6. Possible locations of boundary layer separation [4] ............................................ 8 Figure 7. Example Total-pressure distortion contour plot [5] .......................................... 10 Figure 8. Sources of distortion [5] .................................................................................... 12 Figure 9. Ground vortex visualization [16]....................................................................... 13 Figure 10. Inlet vortex flow model (inlet vortex source) .................................................. 14 Figure 11. Potential flow streamlines a) quiescent b) low speed c) high speed [13] ........ 16 Figure 12. Axisymmetric Model Geometry (inches) ........................................................ 20 Figure 13. CFM56-5B Duct and Engine Layout .............................................................. 21 Figure 14. Water tunnel experimental flow visualization diagram ................................... 23 Figure 15. Water tunnel flow visualization setup ............................................................. 24 Figure 16: CFM56-5B mounted on wing.......................................................................... 29 Figure 17: 3-D geometry generation using circular revolution of 2-D CFM56-5B ......... 30 Figure 18. CFM-56-5B CFD grid ..................................................................................... 34 viii Figure 19: Universal wall law plot for turbulent boundary layers on smooth, solid surfaces [23] .................................................................................................... 36 Figure 20: Grid comparison of area-weighted y+ solutions on inner duct surface ........... 40 Figure 21: Grid comparison of area-weighted y+ solutions on outer duct surface ........... 40 Figure 22: Comparison of coarse and fine mesh for y+ sensitivity study ........................ 42 Figure 23: y+ grid comparison of non-dimensionalized total drag on inner duct surface 44 Figure 24: y+ grid comparison of non-dimensionalized total drag on nose cone surface 45 Figure 25: y+ grid comparison of non-dimensionalized fan-face total pressure .............. 46 Figure 26: y+ grid comparison of non-dimensionalized fan-face mass flow rate ............ 46 Figure 27: Comparison of coarse and coarser mesh for far-field boundary study............ 48 Figure 28: Far-field grid comparison of non-dimensionalized total drag on inner duct surface ............................................................................................................. 48 Figure 29: Far-field grid comparison of non-dimensionalized total drag on nose cone surface ............................................................................................................. 49 Figure 30: Far-field grid comparison of non-dimensionalized fan-face total pressure .... 49 Figure 31: Far-field grid comparison of non-dimensionalized fan-face mass flow rate ... 50 Figure 32: Comparison of coarse and fine mesh for grid spacing study .......................... 51 Figure 33: Grid spacing comparison of percent total drag on inner duct surface ............. 54 Figure 34: Grid spacing comparison of percent total drag on inner duct surface ............. 54 Figure 35: Grid spacing comparison of non-dimensionalized fan-face total pressure ..... 55 Figure 36: Grid spacing comparison of non-dimensionalized fan-face mass flow rate ... 55 Figure 37: Example tournament selection and creation of new generations .................... 58 ix Figure 39. Generalized concept of operations for a given generation GA run [34] ......... 59 Figure 40. Network structure of GAMBIT and FLUENT [34] ........................................ 60 Figure 41: Design space for optimization study ............................................................... 63 Figure 42: Maximum and minimum input geometry variables for GA design space ...... 65 Figure 43. Example LIF Data ........................................................................................... 67 Figure 44. Example PIV Data ........................................................................................... 68 Figure 45. Re? = 3700 4.5 hz a) LIF b) PIV ..................................................................... 69 Figure 46. Sub-critical operating condition Re? = 3700 4.5 hz a) LIF b) PIV ................ 70 Figure 47. Critical operating condition Re? = 11000 13.4 hz a) LIF b) PIV ................... 71 Figure 48. Super-critical operating condition Re? = 15000 17.9 hz a) LIF b) PIV .......... 72 Figure 49. Super-critical operating condition Re? = 18000 a) LIF b) PIV ....................... 73 Figure 50. Super-critical operating condition Re? = 22000 a) LIF b) PIV ....................... 73 Figure 51: PIV Results Re? = 3700 Incidence Angle = a. -10? b. 0? c. 10? ..................... 77 Figure 52 PIV Results Re? = 15000 Incidence Angle = a. -10? b. 0? c. 10? .................... 77 Figure 53. Quiescent flow PIV measurements Re? = 0 .................................................... 79 Figure 54. Inlet Vortex Formation Re? = 0, dt = 20/32 seconds ..................................... 80 Figure 55. Inlet Vortex Stagnation Region, Re? = 0 ....................................................... 81 Figure 56. Vector Flow Field of CFM56-5B at an Angle of Attack of two degrees ........ 82 Figure 57. Flow field pressure coefficient contour plots for the CFM-56-5B .................. 83 Figure 58. Flow field Mach number contour plots for the CFM-56-5B ........................... 84 Figure 59. Boundary Layer Growth of CFM56-5B AoA: 2? ........................................... 84 Figure 60. Flow uniformity: best & worst performers vs. GA Generations ..................... 86 x Figure 61. Total pressure ratio: best & worst performers vs. GA generations ................. 87 Figure 62. Comparison of the original and the optimized CFM-56-5B Geometry .......... 88 Figure 63. General equation representation for a round LE, sharp TE airfoil [32] ......... 93 Figure 64. Class Function 2D Design Space .................................................................... 94 Figure 65. CFM56-5B Duct Geometry to Equation Fit .................................................... 97 Figure 66. GA Design Space ............................................................................................ 99 Figure 67. GA Equation Fitting after 50 Generations..................................................... 100 Figure 68. GA Equation Fitting after 200,000 Generations ........................................... 101 Figure 69. Equation Fitness with Increasing Generations of GA Run ........................... 102 Figure 70. Example standard deviation, Re? = 22000 .................................................... 103 Figure 71: PIV Standard deviation: Variation in angle of incidence and Re.................. 104 xi List of Tables Table 1: Initial boundary layer aspect ratio comparison for y+ sensitivity study ............. 39 Table 2: Mesh sizes used for y+ sensitivity study ............................................................ 42 Table 3: Solution comparison for y+ sensitivity study ..................................................... 43 Table 4: Solution comparison for Far-Field sensitivity study .......................................... 47 Table 5: Grid comparison of initial and maximum grid spacing ...................................... 51 Table 6: Grid comparison of y+ solutions for grid refinement study ............................... 52 Table 7: Solution comparison for grid refinement study .................................................. 52 Table 8. Fitness Ratio Data of Increased Generations .................................................... 102 xii Nomenclature A shaping function variable ( )iA Cell area at ?ith? location on grid surface AoA Angle of attack REFa? Reference speed of sound of the fluid BP Bernstein Polynomials C Class function c chord length CFD Computational Fluid Dynamics pC Pressure coefficient CST Class-Shape-Transformation d Fan-blade diameter dt Time step GA Genetic Algorithm i Variable counter k Total number of velocity values lbs Pounds LE leading edge of an airfoil M Mach number xiii ( )iM x x-axis component of Mach number at ?ith? location on grid surface N Class function variable N? Total number of grid points on a grid surface n Order of Bernstein polynomial p Local pressure PIV Particle image velocimetry ( )iPS Static pressure at ?ith? location on grid surface ?p Freestream pressure R Universal gas constant RANS Reynold?s Averaged Navier Stokes r Variable counter R Radius Re Reynolds number Re? Reynolds number S Wing area Sr,n(x) Shape function T? Average total temperature for a grid surface TE Trailing edge of an airfoil U Freestream velocity u Velocity ? Averaged mean velocity ui Velocity at the ith location xiv u* Friction velocity V Local velocity ?V Freestream velocity v Kinematic viscosity ?x Velocity component in the x-axis direction ( )ivx Velocity component at grid point ?i? in the x-axis direction x Variable location on x-axis y Variable location on y-axis y+ Dimensionless distance from wall z Variable location on z-axis ? Boatail angle ? Change in ? Index of flow uniformity ? Specific heat ratio ? Symbolic variable ? Courant number ? Ratio of x-axis location to chord length ? Ratio of y-axis location to chord length S? Static density for a grid surface ?? Average total density for a grid surface ? Degrees ? Freestream 1 Introduction Efficient inlets are an integral part of a propulsion system because the performance of an engine depends greatly on the characteristics of flow provided by the inlet. The purpose of an engine inlet or diffuser is to match the mass flow requirements of the compressor for efficient operation at a given flight condition. In doing so, a diffuser slows down the flow to a lower Mach number with the best pressure recovery so as to provide the highest stagnation pressure upstream of the compressor face for improved thrust generation. Diffusion of flow is accomplished by decelerating the flow such that the steep stream wise pressure gradient does not result in boundary layer separation rapidly and flow unsteadiness and consequently the distortions in the velocity profiles at the compressor inlet that can lead to compressor stall and surge. Inlets are categorized as either integrated or podded based upon where they are located on the aircraft. Integrated inlets are a part of the body of the aircraft, commonly on the fuselage and are also known as buried inlets. Common locations include on the nose of the fuselage (chin), above or below the fuselage (chin or over-fuselage), underneath the wings (armpit), on the fuselage before the tail (tail root), or on the wing leading edge, Figure 1. 2 Figure 1. Inlet Locations ? Integrated engines [1] Integrated inlets are commonly used in military applications for their improved aerodynamic characteristics in transonic and supersonic flight. However, efficient and stable supersonic diffusion for a wide range of supersonic Mach numbers is difficult to achieve. Supersonic integrated inlets must capture the engines required mass flow rate while managing shocks in the inlet. To ensure this, variable inlet geometry may be required to minimize inlet loss and drag [2]. Integrated inlets are more constrained by internal aerodynamics due to engine placement and therefore require longer ducts that bend and turn the flow. Therefore, internal flow separation and secondary flow are the primary design concerns [3]. Unlike integrated inlet, podded inlets are placed at a distance from the aircraft and are commonly used in subsonic commercial applications due to two important 3 considerations: 1) the engines are farther from the fuselage, therefore, the engines radiate less noise to the cabin, and 2) maintenance is generally easier. Common locations include under and below the wings, on the tail, and even on the wingtips, Figure 2 . Figure 2. Inlet Locations ?Podded engines [1] Since podded inlets are placed away from the aircraft, the engine ingests undisturbed freestream flow. They also have short inlet ducts and therefore have the most direct route possible from freestream flow conditions to the engine [3][4] and a near- isentropic internal diffusion can be achieved [2]. Variable inlet geometries are used to adjust the inlet for optimal flow conditions at multiple operating conditions. Variable inlet lip geometry is not required but is commonly used to improve performance at low speeds, such as during takeoff [5]. A few examples of variable inlet geometry devices are shown for reference in Figure 3. 4 Figure 3. Thin lip inlet auxiliary devices at low speed [5] The primary aerodynamic considerations for a podded inlet design are the distance air travels from freestream conditions to the fan-face, the shape of the inlet, and the type of attachment to the wing. A subsonic, podded inlet was chosen for this investigation; therefore further literature on podded inlet aerodynamics is presented here. Podded inlets consist of a nacelle or duct that surrounds the engine. The cross section of a typical subsonic inlet nacelle and its major geometry parameters is presented in Figure 4. The inlet area, Al, sometimes referred to as the inlet highlight, is the area of the forward most cross section of the duct. The maximum internal area of the duct is at the fan-face. The main purpose of the inlet, to increase the pressure by decelerating the flow through subsonic diffusion, is achieved by an increase in area from the inlet to the fan-face, Af. 5 Figure 4. Subsonic Inlet Nomenclature [1] The area of flow that is ingested into the engine, at a given condition, is known as the capture area. This area is normally non-dimensionalized by the inlet area, and is referred to as the capture area ratio. The magnitude of the capture area ratio is dependent on the freestream flow conditions and the mass flow requirements of the engine. To understand how capture area ratios may be greater or less than one, Seddon?s [3] conclusions analysis of incompressible flow through an inlet is presented here. ( )( ) peefe CqqAAk ?=+ ? 11 2qe (6) Using the definition of dynamic pressure: ( ) 21 21 22 ? ? ? ? =??? ? ??? ?== ? qqqqVVAA e e e ee ?? (7) Substituting Equation 6 in Equation 7 obtains: 6 ( )21 1 fe pe e AAk CAA + ?= ? (8) For a real world engine, the area at the fan-face is the exit area: k CAA pf f + ?= ? 1 1 (9) The final result, Equation 9, shows that, assuming incompressible flow, that the capture area is determined primarily by the area at the fan face. The inlet area was simplified out of the equation; this proves that the capture area is independent of the inlet area. Therefore, the flow at the entry to the inlet adapts to the required mass flow value set by the fan-face area whether or not the inlet area is large or small. Equation 9 was obtained by assuming incompressible, subsonic flow, but Seddon asserted that this conclusion may be applied to compressible subsonic flow as well. Qualitative results for compressible and incompressible subsonic flow are the same, though quantitative results differ. Therefore a qualitative understanding of compressible flow through an inlet can be gained from incompressible flow. During low-speed high-thrust operations (take-off and climb) the inlet may be at a sub-critical operating condition, and not initially ingest enough air to meet the mass flow requirement. If this occurs, the engine will attempt to ingest the air required [1]. A streamline pattern of additional suction is presented in Figure 5a. The inlet area is less than the freestream area and flow accelerates into the inlet. During level-flight at cruise velocities, the inlet may be at a super-critical operating condition, and ingest more air than is required. During this condition inlet must 7 spill the excess flow out the front lip of the inlet; referred to as spillage [1] [4]. A streamline pattern of spillage is presented in Figure 5b. When spillage occurs, the inlet area is greater than the freestream area, and flow decelerates into the inlet. Figure 5. Typical Streamline Patterns for Subsonic Inlets [4] The freestream capture area is commonly non-dimensionalized by the area of the inlet and referred to as the capture area ratio. The capture area ratio is greater than one during sub-critical operating conditions, less than one at super-critical conditions. The critical condition is when the capture area ratio is equal to one, and is commonly referred to as full flow. For all capture area ratio values, the flow bifurcates at the duct leading edge of the inlet and either is ingested into the engine or flows externally around the duct. Freestream flow that is not ingested accelerates externally over the surface of the duct and causes a region of lower pressure near the external surface of the duct lip. At high-speed conditions, the flow may accelerate and become supersonic, resulting in the formation of shock waves. The ingested flow decelerates in the diffuser and, depending on the rate of the pressure rise, can be a source of local flow separation. Flow separation also occurs on 8 the external and internal surfaces of the inlet: (1) on the exterior surface, on the (2) interior lip, (3) or nose cone, Figure 6. Figure 6. Possible locations of boundary layer separation [4] During high alpha operations, the likelihood of local flow separation and formation of separation bubbles increase [6][7]. This possibility is manifested in cases of nacelles with sharp lips. Sharp lips are designed to have optimal flow for a specific flight condition, at the cost of performance at off-design conditions. Rounded lips are commonly used to ensure well behaved internal flow characteristics at most flight conditions. Internal flow separation may occur if the area of the inlet increases too quickly. Limits exist on the radii of curvature of the leading edge and internal diffuser section of the duct which maximize the pressure recovery with minimal to no flow separation [9]. Internal flow separation has large adverse effects on engine performance. Many research efforts have focused on inlet design to predict and ultimately prevent flow separation through the use of wind tunnel testing, empirical models, and numerical analysis [8][10]. The purpose of the inlet is to increase the pressure of the flow by lowering the Mach number. Because of this, during compressible flow, an inlet is a form of a compressor. Like a compressor, the inlet ingests free stream flow, at a high Mach number and low pressure, and increases the local pressure by lowering the Mach number of the 9 flow as it passes through the inlet. (The inlet?s ability to increase the pressure is known as the pressure recovery and is considered the best gauge of inlet efficiency.) The pressure recovery is the ratio of the work done in compression to the total available kinetic energy of the flow. In high freestream flows, this ratio is equivalent to the total pressure ratio, or the ratio of the total pressure at the exit of the diffuser to the freestream pressure. As previously stated, diffusion in podded inlets is near isentropic, and therefore are quite efficient, commonly having pressure recovery values over 90%. The pressure recovery is critical to engine performance because approximately a 1% reduction in inlet pressure recovery will reduce thrust by about 1.3% [1]. Reductions in pressure recovery can be caused by turbulence due to flow separation, losses due to boundary layers on the diffuser walls, and external and internal shock waves. The total effect of these disturbances on pressure recovery is known as distortion. Distortion is used to describe the magnitude of how non-uniform and non- ideal a flow is across the face of the engine: at the aerodynamic interface plane (AIP) between the inlet and the engine [5]. Distortion refers to an adverse distribution of flow at the API and is described in terms of the variation of total pressure across the engine face [3]. Distortion levels are traditionally visualized by plotting total-pressure contours in polar coordinate form at the AIP [5]. The data is plotted as if looking upstream, to view how the engine sees the oncoming flow. A typical total-pressure distortion plot comparing experimental and analytical data is presented in Figure 7. 10 Figure 7. Example Total-pressure distortion contour plot [5] A certain level of radial non-uniformity will always exist, even without flow separation, due to the boundary layer growth on the inlet walls, that can cause major engine aerodynamic failures [11],[12]. The most severe malfunctions caused by distortion are compressor stall and engine surge. Compressor stall is when airflow separates from the airfoils in the compressor, decreasing engine performance. Flow separation, if severe enough, may cause engine surge. Flow separation can prevent the compressor from increasing the local pressure inside the engine, causing the pressure inside of the engine to be higher than the pressure in the compressor and inlet. If this occurs, the high- pressure flow from the combustor surges back through the front of the inlet and can cause severe structural damage. Malfunctions due to distortion can not be completely prevented through careful design, but inlets can be designed to limit their occurrences to a minimum that can be considered negligible. This is accomplished by maintaining a tolerable amount of 11 distortion. This tolerance limit is the amount of distortion an engine can experience without surge occurring. 12 Distortion can be caused by a number of sources, depending on the type of inlet and location on the aircraft. Presented here is an overview of common aerodynamic effects that increase distortion, depicted in Figure 8 [5]. Distortion can be caused by: (a) Flow separation due to the interaction between the inlet shock and the compression surface boundary layer; (b) Spillage of the boundary layer from the front of the fuselage that can be ingested into the inlet; (c) Vortices formed upstream of the aircraft can be ingested: vortices formed by the aircraft or vortices generated by other aircraft if the aircraft crosses the path of another aircraft before the vortices have dissipated; (d) Separated flow in a curved or S-duct diffuser; (e) Separation at the duct lip due to sideslip or high incidence angle. Figure 8. Sources of distortion [5] 13 Other types of ingestion induced distortion include ingestion of weapon exhaust, re-ingestion of hot-exhaust during vertical and short take-off and landing (VSTOL) flight and use of thrust reversers due to, and the ingestion of inlet vortices near the ground. Inlet vortices are vortices that form between the inlet and the ground; and are also known as ground vortices [13]. They form when the engine is operating at high power near the ground while at low speeds or stopped [14]. Normally the vortex core cannot be seen with the naked eye, but condensation can cause the vortex core to become visible [15]. Figure 9. Ground vortex visualization [16] In the absence of the ground, inlet vortices do not occur. In freestream flight, streamlines of ingested flow are axisymmetric .At low speeds, the capture area is near infinite and freestream flow may be ingested from a considerable distance away from the engine centerline. In the presence of the ground, the ground impedes these streamlines below the inlet and the flow lines become asymmetric; causing a stagnation point to form on the ground [17]. 14 However, inlet vortices do not occur solely from an inlet being near the ground. Inlet vortices require a non-ideal flow field to form; therefore only in the presence of flow disturbances will inlet vortices form. A common misconception is that the vortex is formed by the gas turbine rotor, but this is incorrect [13]. Disturbances such as wake vortices from other aircraft and cross winds cause horizontal flow to be introduced to the vertical up-flow of air from the ground to the inlet. A vortex immediately forms when this occurs [14]. High axial velocities in the vortex core occur and cause a large drop in pressure [18]. A cartoon of the vortex core flow field is presented in Figure 10. Figure 10. Inlet vortex flow model (inlet vortex source) Inlet vortices are a problem to the aerospace industry because of the drop in pressure in the vortex core. This pressure drop causes the vortex to become like a ?vacuum? and to pick up rocks, sticks, dust, etc. and create foreign object debris [14][19]. 15 The vortices are not strong enough to carry debris the entire distance into the engine, however, the vortices can be strong enough to toss objects upward, and the inlet flow field ingests the object. In the absence of an inlet vortex, the ingested flow field is not strong enough to suck objects on the ground into the engine [17]. Therefore if the formation of inlet vortices can be prevented, damage due to FOD will be prevented. Inlets can be designed to minimize the formation of inlet vortices but they can never be wholly prevented. Some novel concepts have been developed to periodically break up the inlet vortex. However these have been found to stir up as much debris as the inlet vortices that they prevent [14][19]. A common method implemented by the aerospace industry to minimize ingestion of objects into the engine is to keep the runways clean by sweeping and or vacuuming the runways on a regular basis [13]. Inlet vortex suppression through design requires a strong understanding of the characteristics of their formation. The main characteristics of the vortex are its stagnation location on the ground relative to the engine and its strength. The strength of the vortex and location of the stagnation point depend on the cross wind conditions and are not fixed [17]. To gain a better understanding of this, Colehour?s analytical flow model is presented here [13]. Colehour divided the flow into inviscid potential flow and viscous flow. He concluded that the potential flow model allows the deduction of a fairly complete flow model in the region of the ground plane. He asserted that the complete three-dimensional flow field can be studied using the method of Rubbert by assuming the flow is inviscid and incompressible [20]. With these assumptions, his potential flow field results are shown here. 16 A simple sink flow near a ground plane was used to model flow with quiescent ambient conditions, Figure 11a. With the addition of tangential flow to the sink flow, a model of flow with some ambient wind condition is obtained, Figure 11b. It was found that the addition of a head-wind caused the stagnation point to move in the ambient flow direction. A large enough addition of a tangential flow was found to remove the stagnation point from the ground, Figure 11c. His results support the assertion that the stagnation point location is variable and dependent on the ambient flow conditions. Figure 11. Potential flow streamlines a) quiescent b) low speed c) high speed [13] a b c 17 Genetic Algorithms (GA) and Motion Adjoint methods have been used in the past to optimize aerodynamic shapes for freight trucks [25], gas turbine engines [26], missiles [27][28], and aircraft wings and airfoils [29], [30].An efficient tool used for optimization of aerodynamic shapes is the Class-Shape-Transformation (CST) method with Bernstein Polynomials (BP) [31],[32]. This method has been successfully used with the GA to create equations for and to optimize airfoil and propeller [33] geometries. Optimization of an inlet duct can be used to increase pressure recovery and ensure flow uniformity at the inlet to an engine, resulting in better performing and efficient engine inlet geometries. However, little research has been done on optimization of engine inlet geometries using Bernstein Polynomials and GA. The CFM56-5B turbofan inlet duct was used as the baseline geometry for this effort. It was chosen because it is widely used around the world. It is the primary engine for the A320 family of aircraft, and any implemented improvements to the inlet geometry would have a significant economic impact on the airline industry. The CFM56 engine is a single stage, high-bypass ratio turbofan engine designed for commercial and military use. The CFM56-5B has the highest fan pressure ratio of all CFM56 engines, providing between 22,000 and 33,000 lbs of thrust, and is the first commercial engine to use an ultra-low emission combustor. A binary encoded GA was used to drive Computational Fluid Dynamics (CFD) that ran three dimensional flow simulations with the aim of optimizing engine geometry. The inlet geometry was defined as multiple Bernstein Polynomials whose coefficients were varied using a binary encoded GA over a wide design space to optimize the inlet 18 geometry for cruise conditions. The system created for this optimization effort was thus both versatile and modular and allowed for designs that, in future, could enable further modular additions including engine pylons and wing interfaces. 19 Methodology 1 Experimental 1.1 Model Geometry The axisymmetric model?s geometry was determined from available two- dimensional CFM56-5B specifications using a genetic algorithm and the Class Shape Transformation method (CST) [31]. Exact data points were obtained for each of the duct surfaces (two internal and two external duct surfaces) using AUTOCAD. This was done to obtain equations to represent each surface by curve fitting the data point. Bernstein polynomials have been proven accurately to represent airfoils and nacelles [31]. Therefore, a real-coded genetic algorithm was used to randomize these Bernstein polynomials to match the twenty points obtained from the figure for each curve. This method is presented in detail in the appendices. Using this method four equations were found, with a maximum percent difference from the data points of less than 0.1 %. This was assumed to be a negligible error. A 3-Dimensional model was created from these 2-D equations with circular swept revolution from the upper surfaces to the lower. Ideally, a non-axisymmetric model would have been used for experimental tests. However an axisymmetric model was chosen over the real world non-axisymmetric CFM56-5B geometry due to machining limitations using a CNC lathe. The axisymmetric geometry was obtained by revolving the 20 equations representing the real world duct upper surfaces 360? about the engine centerlines. The lower portion of the duct is blunt, which was assumed to minimize the formation of inlet vortices when near the ground. The upper airfoil has a larger camber and a shape that matches more conventional duct geometry. The upper duct surfaces were chosen to create the model because it was assumed that they were more critical to cruise condition design. A comparison of the model geometry used and the CFM56-5B is presented in Figure 12 and Figure 13. Figure 12. Axisymmetric Model Geometry (inches) 21 Figure 13. CFM56-5B Duct and Engine Layout Photos courtesy of CFM International, a 50/50 joint company between Snecma and General Electric The model was cut from a solid block of acrylic using a CNC lathe from a STEP file created in SOLIDEDGE. Acrylic was chosen to be able to permit light to shine through the model to illuminate clearly the flow field inside the model. This was necessary to be able to conduct flow visualization inside the model cavity. The nose cone was painted black to reduce reflections that would adversely affect the captured images. Initially after machining, the model had a turned finish and was opaque. The model was polished using NOVUS polish remover level 2 and then level 1 to obtain a clear finish. 22 2 Experimental Data Acquisition Flow visualization was conducted on scale model geometry of a high-bypass turbofan inlet duct. 2.1 Test Facilities All tests were conducted in the Auburn University Aerospace Engineering 45 cm x 45 cm cross-section test section water tunnel. The test section was 2 m long and transparent, which allowed for flow visualization and quantitative measurements. The water tunnel was capable of maximum velocity of 1.1 m/s and had a turbulence intensity of 1% at maximum velocity. 2.2 Experimental Setup Tests were conducted in the water tunnel and the model was sting mounted on a specially designed support system. The model was connected to a constant volumetric flow rate water pump to achieve a favorable pressure gradient in the cavity of the model. The extracted flow was deposited back into the water tunnel through a hose, at approximately thirty times the body diameter downstream of the model support system, which was assumed to have a negligible effect on the flow field. The vertical support was attached to the sting mount at approximately four body diameters behind the trailing edge of the model. The entire support assembly rested on the top of the tunnel walls. This distance was chosen to limit the moment arm caused by the model on the support system, while also minimizing disturbance caused by the cylindrical vertical support strut. The mass flow rate was determined from the measured volumetric flow rate using a Fill-RITE digital flow meter attached directly after the rate water pump. A laser light sheet was 23 created using an argon laser and reflected using a mirror to illuminate regions of interest, as can be seen in Figure 14. Figure 14. Water tunnel experimental flow visualization diagram The laser sheet was set at an obtuse angle to prevent shadows caused by refraction in the leading edge of the lower surface from being present in the PIV images. A dark colored material was placed on the opposite side of the water tunnel from the camera to absorb reflected light that would interfere with flow visualization results, as can be seen in Figure 15. Note that the inlet model is not attached to the mounting system in the image shown. 24 Figure 15. Water tunnel flow visualization setup 2.3 Tests Performed The tunnel was run at Reynolds numbers of 3700, 7000, 11000, 15000, 18000, and 22000 based on the diameter of the inlet near the fan-face (2 inches). 2.3.1 Flow Visualization Two flow visualization methods were used in this effort: hydrogen bubbles and laser induced fluorescence (LIF); simultaneously and alone. Hydrogen bubbles were produced through electrolysis by applying voltage across a platinum wire and anode placed in the water tunnel. The electrolysis caused the hydrogen and oxygen in the water to separate and form visible hydrogen bubbles that followed accurately the path of the flow. A constant voltage of 34 V was applied through the probe for all tests to create optimum size and amount of bubbles for flow visualization. A bubble wire probe that created bubbles vertically on the X-Y plane was used to obtain data on the vertical 25 symmetry plane of the model. The probe consisted of three equally spaced 32 Swg platinum wires. The probe was mounted onto a traversing system that rested on the tunnel side walls. A 2 mm thick laser light sheet was placed in the plane of symmetry to illuminate the hydrogen bubbles. For the LIF tests, a solution of sodium flourescein salt and water was injected into the freestream flow at approximately two model body diameters away from the leading edge of the model. An airfoil shroud was placed onto the dye injection tube to minimize interference caused by the introduction of the dye injection system to the flow. Flow visualization results were recorded using a high speed CCD camera. 2.3.2 Particle Image Velocimetry (PIV) Particle Image Velocimetry is a non-intrusive technique that tracks the motion of seeding particles in the flow. The particles are illuminated by a dual-pulsed laser light sheet in the plane of interest. The laser light is reflected off of the particles during the duration of the laser pulse. A CCD camera is positioned perpendicular to the plane of interest to capture a series of still images from the light reflected from seeding particles. The camera is synchronized with the dual-pulse laser to obtain pairs of images with a known time difference between the images of each pair. Each seeding particle displacement distance is determined by calculating the change in location of individual seeding particles between the first and second image of each pair. The time between each dual pulse laser is set according to the freestream velocity of the flow field for each experiment and the cross-correlation camera is synchronized with the dual-pulse laser. Therefore the time and displacement of the seeding particles between each image is 26 known. From each pair, the velocity vector magnitude and direction of the seeding particles can be obtained from the known displacements of the seeding particles between each image. Measurements were statistically averaged from 200 pairs of images to remove any outlying disturbances in the flow as well as obtain the major trends of the flow. Each pair of images may have not caught every intricacy of the flow with the same fidelity. Averaging multiple pairs of images ensures that vector data of high fidelity is obtained for the entire flow field. Increasing the total number of images increases quality of the averaged data, however a limit exists after which increasing total pairs of images gives diminishing returns. The number of pairs of images was chosen to obtain a well defined flow field, while not requiring an unreasonable amount of time to compute the statistical average of the pairs. Rifki determined that ?the dual-pulsed laser and cross-correlation camera allow the system to significantly reduce images acquisition rate between a pair of consecutive images in the illuminated observation plane [34].? An increased acquisition rate decreases the time required for each series of image capturing, thus allowing for a larger number of total image pairs to be taken and statistically averaged in the same amount of time. A Dantec Dynamics PIV system, consisting of a Highsense 1k x 1k cross- correlation CCD camera and a PIV 2100 processor, and New Wave Research 50 mJ dual- pulse ND: YAG laser were used to acquire PIV images that were processed using Dantec Flow Manager Software. The interrogation window used for cross-corexhaurelation computation was 32x32 pixels with 50% overlap. PIV data was post-processed using 27 Tecplot 360. The flow was seeded with silver coated hollow glass spheres with an average diameter of 20 microns. These were chosen because they are highly reflective, with specific gravity approximately equal to one, and follow the streamlines of the flow closely. 28 3 Computational Fluid Dynamics The CFM56-5B geometry used in the experimental investigation was also used for the CFD solution. 3.1 Grid Generation One of the major requirements of the automated optimization discussed later in this effort was the need for fully autonomous and computerized grid generations driven by a GA. Such an automated grid generation program was created for the CFM56-5B engine geometry. The program created a journal file based on the geometry of the inner and outer surfaces of the duct. The journal file was run by the grid generator code. The commercial code GAMBIT was used to develop an unstructured mesh from the journal file defined in one zone. The input to the grid generation code included values for the constants in the Bernstein Polynomials as well as the chord lengths of the 2D duct geometry on the vertical X-Y plane A 3D grid was also generated based on the CFM56-5B engine duct. Assumptions were made about the engine geometry to obtain a 3D grid from available 2D data of the duct. It was assumed that the upper surface of the duct CFM56-5B is revolved over a circular path to the lower surface. As can be seen in Figure 16, the CFM56-5B duct is nearly circular. Therefore this assumption was deemed to be reasonable. 29 Figure 16. CFM56-5B mounted on wing Photo courtesy of Gary Brossett via the Aircraft Engine Historical Society A geometric representation of the circular revolution used to model the engine is shown in Figure 17. The image on the left is the circular revolution of the interior surface points and the image on the right is the circular revolution of the exterior surface points. The concentrated blue regions near the leading edge are representative of the clustering of points near the leading edge as previously mentioned. 30 Figure 17. 3-D geometry generation using circular revolution of 2-D CFM56-5B In the figure, the vertical lines that connect the upper and lower surfaces are the circular revolutions, while the vertices on the symmetry plane, representative of the 2-D data for the CFM56-5B, are shown as vertices. Center points for each pair of upper and lower surface points were obtained and the distances from the center to the surface were determined. This distance was set as the radius of the circular revolution. The vertices where the radii are extended on the Y-Z plane are shown as black vertices. The virtual line that connects these points does not appear ?flat? because y-axis location of the center points between the upper and the lower surface points was not constant as the upper and lower surfaces varied differently from each other along the chord. The assumption of a circular revolution introduced a known discrepancy in design between the optimized geometry and the real world engine because the CFM56-5B duct is not a perfect circle. However, this assumption was required for computational efficiency.. 31 The fuselage, wings, and pylon were not modeled in this effort. The engine exhaust was not an area of interest for this work and CFD modeling of the exhaust would have greatly increased the complexity of the investigation. Therefore the exhaust volume was replaced by a sting of five fan-face diameters in length. The engine was assumed to be symmetric about the X-Y plane to reduce computational time and allow for a finer mesh in areas of high gradients. A cylindrical shaped grid was used instead of a rectangular grid, with the engine centerline placed along the centerline of the cylinder to further reduce the computational time. The grid was generated with a sufficient number of nodal points to resolve the flow field around the CFM56-5B geometry accurately. A successive ratio was applied to the growth of the mesh relative to defined edges and surfaces. The ratio was chosen to place the smallest x-axis vertex spacing closest to the leading edge and increase this spacing as the vertices were placed further along the chord length. This successive ratio was applied to the inner and outer surfaces of the duct as well as to the tip of the nose cone. The primary region of interest was near the cowling and inside the duct inlet before the fan face. The mesh was refined to model accurately the pressure gradients in this region and is discussed later. The nodal points were carefully developed to create a well-defined mesh that adequately modeled the flow while not requiring excessive computational times. 32 3.2 CFD Solver Model The CFD solver used for this effort was the commercially available Reynold?s Averaged Navier Stokes (RANS) code, FLUENT. An implicit pressure-based solver was used for this effort. The pressure-based solver was chosen because it maintained accuracy and robustness of the solution model given a large number of geometries prescribed by the GA. Convergence was defined to be when all of the residuals dropped below a 10-6 in magnitude. A Mach number ramp-up was used to ensure the stability of the solution during initial solution convergence. A first-order accuracy solver was used until the solution converged. To increase the accuracy of the solution, the converged first-order solution was used in a second-order accuracy solver. The Spalart-Allmartas (SA) turbulence model was used because it was designed specifically for aerospace applications involving wall-bounded flows [21]. The SA model is a one-equation turbulence model and was chosen over more complex two-equation turbulence models (ie. k-epsilon, k-omega) or Reynolds stress models due to reduced computational times. Turbulence models are accurate only for attached flows. The flow field around the CFM56-5B was assumed to remain attached during cruise. Detached flow would be avoided during the original design process of the engine. However, for the optimization study, if the GA generated a geometry that caused flow separation to occur, the inaccurate separation modeling would not affect the optimization study. The turbulence model would inaccurately model the flow separation, however the incorrect solution 33 would accurately solve for increased losses in total pressure, compared to geometries with no detached flow. The separated flow solution would not be chosen as a ?best geometry? by the GA using the tournament selection method, explained later, and therefore would be ignored. Therefore an attached flow solver was deemed to be adequate. The boundary conditions for the solution were set as far-field pressure, pressure outflows, or no-slip walls. The surfaces of the inner and outer duct, sting and nose cone were set as no-slip walls. The entry plane upstream of the engine and curved plane of the cylindrical grid were defined as pressure far-field boundaries to model the flow coming into the control volume. A pressure outflow boundary was used for the fan-face to model air ingestion. The exit plane of the grid downstream of the engine was set as a pressure outflow. For the pressure outflow boundaries, the mass flow rate was not enforced, and therefore was allowed to change during convergence. The grid used is presented in Figure 18 below. 34 Figure 18. CFM-56-5B CFD grid 35 3.3 Boundary Layer Resolution: Near ? Wall Treatment Boundary layer growth on the interior surfaces of the duct and nose cone causes losses in total pressure. To better understand this growth, and to ensure that the grids used in the y+ sensitivity study accurately resolved the boundary layer, a brief overview of boundary layer profiles by Frank White is presented below [22]. The velocity profile of a boundary layer consists of three regions: the viscous sub- layer (y+ < 5), the outer region (350 ? y+ > 30), and the overlap region (5 ? y+ ? 30). In the viscous sub-layer, viscous shear is dominant and the velocity profile is linear. However, in the outer region, turbulent shear is dominant and the velocity profile is logarithmic. In the blending region, both viscosity and turbulence are important in resolving the boundary layer, and the velocity profile is represented by a blending of the linear and logarithmic profiles. Presented in Figure 19 is the experimental turbulent modeled in wall coordinates. As can be seen in the figure, the law of the wall fits the experimental data. The law of the wall however for y+ > 1000, the data departs the wall law due to an adverse pressure gradient. 36 Figure 19. Universal wall law plot for turbulent boundary layers on smooth, solid surfaces [23] To accurately model the boundary layer and near-wall region, multiple approaches can be used. Given here is a description of the two near-wall treatments predominantly used: near-wall models or wall functions [21]. For near-wall modeling, a mesh is developed to resolve the viscosity-dominated region completely to the wall, which includes the viscous sub-layer and overlap region. However, to properly resolve this region a considerably fine mesh is required. Wall-functions do not resolve the viscous region near the wall. Instead, semi- empirical ?wall functions? are used to model, as opposed to fully resolve, the near-wall region. Because wall functions do not require a fine mesh all the way to the wall, turbulence models that use wall functions do not require grids to be as fine, compare to near-wall models. Therefore, solutions using wall functions are faster and therefore more 37 practical. This is true only for high-Reynolds-number flows because of the thin boundary layer. As previously stated, the turbulence model used was the SA model. The SA model uses an enhanced wall function. Instead of modeling the near-wall region as three separate wall laws (viscous sub-layer, outer layer, and overlap layer), enhanced wall functions formulate the law-of-the-wall as a single wall law for the entire near-wall region. Because of the lack of universally accepted models in the overlap region, grid solutions with y+ values in the blending region were avoided for the optimization study. 3.4 Flight Conditions The goal of the optimization study was to demonstrate a methodology for optimizing a turbofan inlet for cruise conditions. Exact flow parameters for the CFM56- 5B engine duct during cruise are proprietary information of GE and Snecma, and had to be approximated for this effort. The engine is the primary engine for the A320 family of aircraft. Therefore the far-field conditions were based upon the A320 cruise conditions at a Mach number of 0.8 at 38,000 ft altitude. A static pressure was estimated at the fan-face by assuming a face Mach number of 0.4 and assuming a stagnation pressure at the constant freestream value. This is the proper boundary condition for inlet modeling. The static pressure was set, while the total pressure at the fan-face was free to vary due to loses in the boundary layer. 38 3.5 Grid Refinement A grid refinement study was conducted on the CFM56-5B geometry to ensure that variation in solutions of different sized meshes based on y+ solutions and grid spacing were negligible. 3.5.1 y + Sensitivity Study Before conducting the grid mesh size refinement study, a grid sensitivity study of the aspect ratio of the first boundary layer mesh was conducted to verify that variations in y+ did not affect the solution. The goal of the optimization was primarily concerned with the interior of the duct. Therefore conditions at the fan-face and on the interior walls of the duct and nose cone versus number of grid cells are presented to compare results. Initial y+ estimations were based on a boundary layer on a flat plate for freestream and fan-face Reynolds numbers to determine approximate heights of the first layer of the boundary layer mesh. The initial height of a boundary layer mesh is commonly non-dimensionalized by the tangential grid spacing and given as a percentage; this is referred to as the boundary layer aspect ratio. The first boundary layer mesh aspect ratio was varied to generate coarse and fine grids to determine the aspect ratios ranges that corresponded to each boundary layer region. As can be seen in Table 1, aspect ratios of 200% and higher corresponded to y+ values in the outer region (y+ > 30). Aspect ratios between 25% and 200% corresponded to y+ values in the blended region (5 ? y+ ? 30). Aspect ratios of 12.5% and below corresponded to y+ values in the viscous sub-layer (y+ < 5). 39 Table 1. Initial boundary layer aspect ratio comparison for y+ sensitivity study Area-Weighted y+ Initial BL Aspect Ratio Grid Points Inner Duct Surface Outer Duct Surface 300% 911097 61.70 122.03 200% 996064 40.44 79.46 100% 1164641 19.16 37.98 50% 1328332 9.03 17.94 25% 1450240 4.81 8.70 12.5% 1616464 2.46 4.47 The y+ solutions of the inner duct and outer duct surface for each grid are presented in Figure 20 and Figure 21. Horizontal lines are overlain onto the figures to show the boundaries between the three regions of the boundary layer. As can be seen in the figures, the majority solutions of the grids generated were in the blending region, which is known to be the least accurate. 40 1 10 100 0.9 1.1 1.3 1.5 1.7A rea -W eig hte d y + I nn er Du ct Su rfa ce # of Cells x 106 Figure 20. Grid comparison of area-weighted y+ solutions on inner duct surface 1 10 100 0.9 1.1 1.3 1.5 1.7 Ar ea- We igh t y + O ute rr Su rfa ce # of Cells x 106 Figure 21. Grid comparison of area-weighted y+ solutions on outer duct surface (Outer Region) (Blending Region) (Viscous Sub-layer) (Outer Region) (Blending Region) (Viscous Sub-layer) 41 For comparison of solutions, a course mesh and fine mesh were chosen from the grids generated based on their solutions of the area-weighted y+ values on the interior and exterior surface of the duct. The aspect ratio for the course mesh was chosen such that the solution y+ values were in the outer region of the boundary layer while the aspect ratio for the fine mesh was chosen such that the solution y+ values were in the laminar sub-layer of the boundary layer. The mesh sizes and y+ solutions of both grids are shown in Table 2. Shown in Figure 22 below is a comparison of the symmetry plane mesh near the upper surface duct leading edge also for both grids. 42 Table 2. Mesh sizes used for y+ sensitivity study Initial BL Aspect Ratio Area-Weighted Interior Duct Surface y+ Area-Weighted Exterior Duct Surface y+ Number of Cells Course Mesh 200% 40.44 79.46 996064 Fine Mesh 12.5% 2.46 4.47 1616464 Figure 22. Comparison of coarse and fine mesh for y+ sensitivity study Table 3, the total drag values for the coarsest grid, for an aspect ratio of 300 %, were significantly less than the solutions using the less coarse grids. This was assumed to be because the boundary layer mesh was too coarse to properly resolve the boundary layer on the duct surface. As previously stated, Sparta-Allmartas solutions in the blending region are less accurate. However, the solutions presented in the table below, which includes solutions in all three regions of the boundary layer, can be seen to be nearly identical, excluding the coarsest grid. Therefore it was concluded that that the solutions were insensitive to changes in y+. 43 Table 3. Solution comparison for y+ sensitivity study Total Drag Initial BL Aspect Ratio Grid Points Inner Duct Surface (N) Nose Surface (N) Fan-Face Area- Weighted Total Pressure (N/m2) Mass Flow Rate (kg/s) 300% 911097 162.87 314.41 12407.02 51.2436 200% 996064 2205.30 1370.14 12426.50 51.3746 100% 1164641 2208.58 1369.00 12405.34 51.2235 50% 1328332 2199.96 1369.43 12424.21 51.3647 25% 1450240 2191.46 1371.70 12402.54 51.1904 12.5% 1616464 2204.07 1370.36 12400.2 51.1800 Plots of the interior solution values are presented below to see more easily the negligible differences between the solutions. The plotted data is non-dimensionalized and presented as percentages of the coarse mesh solution. Solution data of the coarsest grid, aspect ratio of 300%, is not presented. Figure 20 and Figure 21 present the total drag of the inner duct surface and nose cone. The solutions for drag on these surfaces are presented because drag prediction is sensitive to the boundary layer resolution. Therefore, if the differences in drag between each grid were small, it could be concluded that the boundary layer for each grid was properly resolved. 44 Figure 23. y+ grid comparison of non-dimensionalized total drag on inner duct surface 45 Figure 24. y+ grid comparison of non-dimensionalized total drag on nose cone surface The goals of the optimization study were to maximize flow uniformity and total pressure ratio at the fan-face. The area-weighted total pressure and the mass flow rate at the fan-face are presented in Figure 25 and Figure 26 to prove that the solutions of the grid chosen for the optimization study were insensitive to changes in y+ and valid. 46 Figure 25. y+ grid comparison of non-dimensionalized fan-face total pressure Figure 26. y+ grid comparison of non-dimensionalized fan-face mass flow rate 47 As can be seen in the figures, the differences between solutions were minimal, and therefore were considered to be negligible. It was concluded that the solutions were insensitive to the values of y+. 3.5.2 Far-Field Boundary Refinement A study was conducted on the sensitivity of the solutions to the minimum distance from the edge of the finite grid to the engine model, in terms of fan-face diameters. The coarse mesh chosen for the y+ sensitivity study was chosen as the baseline for comparison for this study. The minimum distances of the grids generated are shown in Table 4. Shown in Figure 27 below is a comparison of the mesh on the symmetry plane seen from afar. Table 4. Solution comparison for Far-Field sensitivity study Total Drag Distance (# Fan Diameters) Grid Points Inner Duct Surface (N) Nose Surface (N) Fan-Face Area- Weighted Total Pressure (N/m2) Mass Flow Rate (kg/s) 40 752591 2201.01 1369.57 12394.30 51.131 20 745101 2198.62 1369.72 12379.59 51.0293 10 737694 2198.92 1370.41 12425.34 51.3697 48 Figure 27. Comparison of coarse and coarser mesh for far-field boundary study Figure 28. Far-field grid comparison of non-dimensionalized total drag on inner duct surface 49 Figure 29. Far-field grid comparison of non-dimensionalized total drag on nose cone surface Figure 30. Far-field grid comparison of non-dimensionalized fan-face total pressure 50 Figure 31. Far-field grid comparison of non-dimensionalized fan-face mass flow rate As can be seen in the figures, the differences between solutions were minimal, and therefore were considered to be negligible. It was concluded that the solutions were insensitive to the far-field boundary, on the range of distances investigated in the far-field sensitivity study. The far-field distance of forty times the fan-face diameter was chosen for the optimization study due to increased solution accuracy at little cost due to the small increase in total grid size. This slight increase in total grid cells was due to the cells furthest from the geometry were considerably large. 51 3.5.3 Grid Spacing Refinement The coarse mesh chosen for the y+ sensitivity study was chosen as the baseline for comparison for the grid spacing refinement study. The initial and maximum mesh sizes were varied to obtain a more coarse mesh and fine mesh. The initial and maximum grid spacing, dictated by the size function starting at the leading edge of the duct and affecting the inner and outer duct surfaces, of the grids generated are shown in Table 5. Shown in Figure 32 below is a comparison of the symmetry plane mesh near the upper surface duct leading edge also for both grids. Table 5. Grid comparison of initial and maximum grid spacing Interior Duct Surface Initial % Fan- Diameter Maximum % Fan-Diameter Number of Cells 200.00% 0.290 14.496 510045 Coarse Mesh 100.00% 0.145 7.248 996064 66.66% 0.096 4.832 1511458 Fine Mesh 50.00% 0.072 3.624 2048880 Figure 32. Comparison of coarse and fine mesh for grid spacing study 52 The solutions were found to be insensitive to variations in y+, however y+ solutions are presented in Table 6 to determine which region of the boundary layer each mesh solved. The data shows that the y+ solutions of all grids were above 30, and therefore were in the outer region of the boundary layer. Table 6. Grid comparison of y+ solutions for grid refinement study Area-Weighted y+ Grid Spacing Grid Points Inner Duct Surface Outer Duct Surface 200% 510045 47.96 172.67 Coarse Mesh 100% 996064 40.43 79.46 66.66% 1511458 36.08 56.39 Fine Mesh 50% 2048880 31.90 46.46 The solutions of the y+ sensitivity study were chosen again compared but between solutions of the grids generated in the grid refinement study. Presented in Table 7 are these solutions. As can be seen in the table, the solutions of all of the grids were close to each other. The coarsest mesh, which had nearly four times less number of cells than the fine grid, had almost the same values. Table 7. Solution comparison for grid refinement study Total Drag Grid Spacing Grid Points Inner Duct Surface (N) Nose Surface (N) Fan-Face Area-Weighted Total Pressure (N/m2) Mass Flow Rate (kg/s) 200% 510045 2210.35 1363.07 12412.92 51.2248 100% 996064 2205.30 1370.14 12426.50 51.3746 66.66% 1511458 2203.17 1370.63 12408.47 51.2698 50% 2048880 2199.88 1371.38 12399.50 51.2117 53 The data in Table 7 is plotted in the following figures. The plotted data is again non-dimensionalized and presented as percentages of the coarse mesh solution, which had a grid spacing of 100%. As can be seen in the figures, the differences between the solutions are negligible. It was therefore concluded that the solutions were insensitive to variations in grid spacing. 54 Figure 33. Grid spacing comparison of percent total drag on inner duct surface Figure 34. Grid spacing comparison of percent total drag on inner duct surface 55 Figure 35. Grid spacing comparison of non-dimensionalized fan-face total pressure Figure 36. Grid spacing comparison of non-dimensionalized fan-face mass flow rate 56 4 Inlet Optimization 4.1 Genetic Algorithms A genetic algorithm is an optimization technique that is based upon the biological process natural selection through successive reproduction generations of an evolving population. A member of a population is a sample design for the system of interest. Members with beneficial genetic traits survive and pass on these traits to their children. Members without these positive traits generally do not survive as often to reproduce and pass on their genes, and genes that are not beneficial are gradually removed, or naturally selected out, from the population. After multiple generations, the result is a population that is better suited to their design goals than the original population. 4.1.1 Tournament Selection The genetic algorithm employed in this investigation is modeled after natural selection through the use of a tournament selection. Selection is based on a designated fitness to the desired qualities, representative of a natural environment, developed by Murray Anderson. Initially a population is created from a range of traits, representative of genetic chromosomes. New populations are created from parents chosen randomly from the population. As previously mentioned, the main goals of the optimization study were to maximize flow uniformity and total pressure ratio at the fan-face. To do this, the GA used a best performer selection method known as tournament selection. The tournament selection method picks, at random, two pairs of members in the total population of each generation. The ?better? member of each pair is chosen based on the quality value 57 calculated from the weight single representative values for each optimization goal. The next generation is created from the two ?better? members by ?mating? the two geometries. Each member of the new generation is created by randomly splicing the codes of each of the ?better? members. A flow chart of an example GA run including the tournament selection process for three generations of an eight member population is shown in Figure 37. The roman numerals above the flow chart represent the steps of the process. For the first generation, generation a, the GA randomly creates eight members based upon input parameters given to the GA. Two pairs of members are randomly chosen from the total population for the tournament selection process, [I]. The first step of the tournament selection process is the fitness calculation and determination of the ?best? fitness member for each pair, [II]. A member?s fitness is a measure of a desired characteristic of the member. The ?best? member of each pair is then mated by randomly splicing their characteristics to a member of the next generation b, [III]. This step is repeated until the entire population of the new generation is created. The entire process, steps [I] [II] and [II] are repeated to create the third generation, generation c. 58 Figure 37. Example tournament selection and creation of new generations The process is not limited to only eight members and three generations; these values were chosen only as an example. This process can be used for any combination of total members in each generation and number of generations. 4.2 Integrated CFD and GA Networks The network of integrated computers set up to allow multi-node parallel CFD case ran under the direction of the GA had the capability to drive up to four simultaneous CFD case runs with one smaller network of computers to allow fully automated, rapid three dimensional grid development under the direction of the GA. This computer network consisted of a maximum of sixty available processors for CFD runs on thirty computer 59 nodes. All processes were initiated on a central node referred to as the head node. In addition, a total of four CFD solvers were run simultaneously on the above cluster. All CFD simulations were run with the FLUENT software, whereas grid generation requirements were met using the GAMBIT based grid generator as discussed previously with both being remotely initiated and controlled in a non-GUI environment by the GA. The current network uses the FLUENT V6.2 version software as the CFD solver and was upgraded to use the latest FLUENT V6.3 version software that served as the basis for the CFD runs conducted for this optimization effort. A brief overview of the networked environment is presented in Figure 38. The relationship between single nodes in the computer cluster with the head node that drove the entire networking operation along with the GA is presented in Figure 39. The file structure shown here was used for each of the processors running in parallel for each of the CFD cases being run. Figure 38. Generalized concept of operations for a given generation GA run [34] 60 Figure 39. Network structure of GAMBIT and FLUENT [34] 4.3 Optimization Goals The main optimization goals for the engine cowling were identified as maximizing the total pressure ratio and flow uniformity at the engine fan face. The goal emphasis for the first run was on maintaining the total pressure ratio while improving the flow uniformity. This was achieved by adding weighting factors on the two goals to ensure that one goal acted as the primary goal (flow uniformity) while the other acted as the secondary goal (total pressure ratio). Improvements in the primary goal had precedence over the secondary goal but still were not set so dominant as to improve at the cost of the secondary goal. The method for quantification of the two goals is important. For both the total pressure ratio and flow uniformity, the calculations must be averaged for the entire fan 61 inlet grid surface so that net goal values are sent back to the GA for evaluation. This averaging process, however, must take into account the validity of the governing equations of the flow. This constraint was met by using a modified version of the Stewart Mixing Analysis designed for gridded surfaces [34], which is applicable to CFD grids. The flow uniformity was quantified as the standard distribution of the velocities at all the cells of the grid about the mean value of the flow velocity as calculated by the modified Stewart Mixing Analysis. Equation 10 presents the resulting form of the mean velocity of the flow along the axial direction of the engine at the fan inlet surface. 111 2 ??? ? ? ??? ? + ??? ??? ? ??? ? + ??= ? ? ? ? ? REF x a v (10) Where, ( ) ( ) ( )[ ] ( ) ( )[ ] ( ) ( )( ) ( )? ? ? = = = ? ? ? ? ? ? ? ? ??? ? ??? ? ??? ? ??? ? ?? ? ???? ? ??? ? ++ ??? ? ??? ? ?? ? =? N i x x S REFREF REF N i SN i xS REFREF iAiv iMiPaP TR P iAiP iAiMiPaP TR ? 1 2 ? 1 ? 1 2 2 2 1 ? ? ?? (11) The flow uniformity is therefore defined as: ( )( )? = ?= N i xx N viv? 1 2 ?flow ofDeviation Standard (12) The average total pressure is similarly defined using the mixing analysis as: 62 ( ) ( )( ) ( ) Aav iAiv iMiPaRT P REF xS N i x x S REF REF ?? ? ? ? ? ? ? ? ? ??? ? ??? ? ??? ? ??? ? ? =? ? = ? ? ? ? 1 2 (13) The total pressure ratio is then evaluated by dividing the pressure obtained from Equation 13 by the free stream total pressure value. 63 4.4 Design Space A constrained design space for a limited scope optimization study was used to demonstrate the optimization abilities of the genetic algorithm. Variations of two- dimensional geometric parameters of the three-dimensional duct were used to create the limited design space. The geometry of the duct was allowed to vary by changing the coefficients of the Bernstein Polynomials, representing the interior and exterior duct surfaces, about the baseline CFM56-5B geometry. For the optimization study, the representative curves of the CFM56-5B that were varied can be seen below in Figure 40 as dashed lines, while the geometries that were held constant can be seen as solid lines. Figure 40. Design space for optimization study 64 The upper and lower leading edges were varied by varying the exterior chord lengths of the inlet. The exterior chord lengths were set as independent input, designed by the GA, with a variation of 10% about the CFM-56-5B values. The range of the possible leading edges can be seen in brackets in Figure 40. The interior chord lengths were set to be dependent on the exterior chord lengths. This was accomplished by the trailing edge locations being held constant; as the CFM56- 5B trailing edge locations, seen as black circles in Figure 40. The varied quantities created a design space comprised of a total of twelve variables. The two optimization goals plus the twelve variables constituted the overall input and output quantifications for the GA runs. A matrix of possible geometries at the maxima and minima of the design space is presented in Figure 41. The CFM56-5B geometry is represented by thin lines, while the geometries created by the GA are represented by thick lines. 65 Figure 41. Maximum and minimum input geometry variables for GA design space In the figure, the upper plots have larger Bernstein Polynomial coefficients, and therefore have thicker chords, while the lower plots have smaller Bernstein Polynomial coefficients. The plots on the right have longer chord lengths, while the plots on the left have shorter chord lengths. 66 Results & Discussion 5 Experimental Results PIV was found to be a more efficient for external flow field investigation compared to dye injection and bubble wire because of entrainment. Streamlines for the entire flow field were obtained simultaneously using PIV. However, only one streamline was obtained at a time using LIF. Bubble wire was able to obtain streamlines for the entire flow field simultaneously as well, but an image of a single frame did not convey the flow field as clearly. However, it was found that LIF and bubble wire were better flow visualization techniques for obtaining the internal flow field due to limitations in PIV data acquisitions. PIV data is not presented above the upper geometry surface or inside of the geometry due to these limitations. Data in regions where the laser passed through the model was unreliable due to the light having been internally reflected and refracted off of the acrylic model geometry. In the flow visualization images, a second streamline was observed that had less intensity than the main streamline, Figure 42. This second streamline was found to be a reflection of the laser used to illuminate the flow field and was ignored. 67 Figure 42. Example LIF Data PIV measurements are presented using streamlines of the flow field superimposed onto pressure coefficient contour plots, as can be seen in Figure 43 below. Pressure coefficients were calculated from PIV velocity vector field measurements using the incompressible pressure coefficient equation: 2 1 ???????= ?V VCp (14) The freestream pressure coefficient, where the pressure coefficient is equal to zero, is seen as orange or green in the figures. Regions of lower pressure coefficients, where the local flow velocities are greater than the freestream velocity, are seen as colder colors in the figures and conversely, the regions of higher-pressure coefficients, where the local flow velocities are less than the freestream velocity, are seen as warmer colors in the figures. Reflection 68 Figure 43. Example PIV Data Contour plots reveal that the pressure coefficient contours in the vicinity of the inlet were circular in trend. Multiple contour levels can be seen near the model inlet area. The contour level spacing closest to the inlet area is relatively small, and moving further upstream and away from the inlet, the spacing becomes larger, until the contour level values are equal to freestream. The small contour spacing is indicative of a high pressure gradient, representative of the simulated ingestion of flow into the inlet. Higher pressure coefficients were measured near the upper and lower leading edges of the duct due to stagnation. However, the contour level spacing was greater than the contour level spacing of ingestion. Therefore the pressure gradient due to stagnation was concluded to be less than the pressure gradient due to suction. At a distance further from the duct surface, and further from the engine center line, the contour level values return to freestream values. A close investigation of the external flow field upstream of the model centerline revealed that the freestream flow was at a slight positive angle from the horizontal. This 69 was due to the upper boundary of the water being a free surface, instead of bounded by a wall. Therefore the water level height was free to rise and fall due to the water displaced by the model obstructing the flow. Dye injected upstream and outside of the external surface of the model was ingested into the inlet at a freestream Reynolds number of 3700, as can be seen in Figure 44a. The dye streamline can be seen continuing to the porous screen, representative of the engine fan face. Because the capture area was greater than one, it was concluded that the test was analogous to a sub-critical operating condition. PIV results corroborated the flow visualization findings. Streamlines obtained by PIV closely match streamlines of flourescein dye used for flow visualization, Figure 44b. The streamlines allow one to visualize the capture stream tube decreasing in area, normal to the plane of view. Figure 44. Re? = 3700 4.5 hz a) LIF b) PIV a b 70 5.1 Effect of Reynolds Number As an aircraft goes through its flight envelope, the engine inlet encounters multiple flow fields, with varying freestream flow Reynolds numbers and angles of attack. The most common of these flight conditions are take-off, climb-out, cruise, and approach. The effect of varying the freestream Reynolds number, while the incidence angle was held constant, on the external flow field around the model geometry is presented in this section. This was done to determine the effect of Reynolds number on the flow field separately from changes in incidence angle. It was determined that at a Reynolds number of 3700, the capture area was larger than the inlet area, therefore the capture area ratio was greater than one, Figure 45. The capture area can be seen as the two opposite streamlines furthest from the engine centerline that are ingested into the inlet. The distance between these streamlines converge as they are ingested indicative of a sub-critical operating condition. Figure 45. Sub-critical operating condition Re? = 3700 4.5 hz a) LIF b) PIV a b 71 Regions of pressure gradients were observed to be approximately elliptical. The ellipses of pressure gradients all had the same semi-major axes: fixed as the distance between the upper and lower inlet leading edges. The size of regions of constant pressure gradients increased further from the inlet leading edges. The elliptical regions also increased in eccentricity further from the leading edge. As the freestream Reynolds number was increased, the streamline of the ingested dye moved closer to the internal surface of the duct geometry. At a critical freestream Reynolds number of 11000, the streamline of the injected dye stagnated at the leading edge of the duct, as can be seen in Figure 46a. After stagnation, the streamline of the bifurcated and a fraction of the streamline was ingested into the engine while the rest of the streamline spilled over the duct surface and flowed externally around the model, Figure 46a. The capture area was found to be approximately equal to one, indicative of a nearly critical operating condition. Figure 46. Critical operating condition Re? = 11000 13.4 hz a) LIF b) PIV a b 72 Pressure coefficient contours revealed that external flow near the surface of the duct had lower pressure coefficients. This was due to the flow accelerating over the surface of the duct. Higher pressure coefficients were measured near the duct leading edges, due to stagnation. The pressure coefficient contours near the model inlet were only slightly higher than the freestream. As the freestream Reynolds number was increased past the critical number, spillage occurred. At a Reynolds number of 15000, dye injected at a height inside of the duct, closer to the model centerline, was found to flow externally around the duct, as can be seen in Figure 47a, indicative of a super-critical operating condition. The external flow near the duct was found to have even greater decreased pressure coefficient contours, Figure 47b. The pressure coefficient contours near the model inlet had a significant increase compared to the freestream. Figure 47. Super-critical operating condition Re? = 15000 17.9 hz a) LIF b) PIV a b 73 Further increasing the freestream Reynolds number further increased the differences in local and freestream pressure coefficients. This can easily be seen by comparing Figure 47b, Figure 48b, and Figure 49b below. Figure 48. Super-critical operating condition Re? = 18000 a) LIF b) PIV Figure 49. Super-critical operating condition Re? = 22000 a) LIF b) PIV External flow was found to accelerate and cause lower pressure coefficients. Ingested flow decelerated and caused increased pressure coefficients. Increasing the Reynolds number was found to amplify these trends. a b a b 74 It was therefore concluded that the capture area of the inlet is sensitive to the freestream Reynolds number. This corroborated Seddon?s assumptions of flow through a duct, as described earlier. 5.2 Effect of Angle of Incidence A quantitative analysis, using PIV, of the effect of varying the incidence angle of the model on the external flow field around the model geometry is presented in this section. The incidence angle was varied while the Reynolds number was held constant. The geometry was set at positive ten degrees, and zero degrees incidence angle. These angles were chosen because they represent approximately angles of attack experienced during take-off and cruise conditions. PIV measurements were not obtained for the flow field above the model due to refraction of the laser sheet inside of the geometry. The duct was also set at negative ten degrees. The model geometry was axisymmetric; therefore the flow field below the model, when set at a negative incidence angle, was equivalent to the flow above the model, when set at a positive incidence angle. It was determined earlier that there were differences between the flow fields, due to the presence of the vertical support in the experimental setup, however the flow fields were observed to be comparable, and Figure 51 below. Pressure gradients were observed to be closest together near the inlet duct leading edge, indicative of a strong potential. The regions of pressure gradient variations for each angle of incidence were observed to be similar. Elliptical regions of pressure gradients were observed for all angles of incidence. As previously seen, the ellipses semi-major axes were fixed at the upper and lower inlet 75 leading edges. Therefore the angle of the external regions of pressure coefficient was found to be equal to the incidence angle of the model. However, increasing the incidence angle was found to not increase the eccentricity of the elliptical distributions. Compared to an incidence angle of zero, the region near the forward most duct leading edge, when at a non-zero incidence angle, experienced increased pressure coefficient values. Simultaneously, the region near the rearward most duct leading edge experienced decreased pressure coefficient values. For the sub-critical operating condition, an increase in incidence angle was found to increase the pressure coefficient near the stagnation region, Figure 50. However, for the super-critical operating condition, an increase in the incidence angle further decreased the pressure coefficient of the flow field near the external surface of the inlet. This decrease in pressure coefficient is indicative of higher velocities, assumed to be caused by the flow accelerating over the external inlet surface. The region of variation to the external flow field varied with changes to angle of incidence. Streamlines near the lower lip, at negative angles of incidence, would be ingested into the inlet. However, the same streamlines, when the model was placed at positive angles of incidence, flowed externally around the geometry. The streamlines near the external surface of the duct were affected by the duct impeding their flow when the model was at an incidence angle. However, external streamlines followed closely the curvature of the external surface for all angles of incidence. 76 The freestream streamlines were found to be not significantly affected by variations in the incidence angle, and therefore capture area did not change as well. The ingestion or spillage of flow streamlines of the external flow field was found to be dependent on angle of incidence. However, the major trends were concluded to be independent of the incidence angle of the model. It was therefore concluded that the capture area of the inlet is independent of the incidence angle. 77 Figure 50. PIV Results Re? = 3700 Incidence Angle = a. -10? b. 0? c. 10? Figure 51. PIV Results Re? = 15000 Incidence Angle = a. -10? b. 0? c. 10? a b c a b c 78 5.3 Formation of an Inlet Vortex An interesting phenomenon known as an inlet vortex was captured during the course of the experiment. The vortex formed when the model was placed two diameters from the bottom wall the water tunnel, which served as the ground plane. The camera was placed 45 degrees to the tunnel wall to obtain a cross-planar view of the inlet. Previous studies used a similar setup to capture the formation of an inlet vortex [13]. The mass flow rate was held constant for all tests. To gain an understanding of quiescent flow, PIV measurements were taken with the model placed within two body diameters from the bottom surface of the water tunnel. As can be seen in Figure 52, the magnitude of the velocity increased as the radial distance from the inlet decreased. The streamlines can be seen to be ingested from both tangential directions from the inlet leading edge. Away from ground plane, the streamlines would be axisymmetric. However, the ground plane impeded the streamlines below the model and caused the flow to be asymmetric. This caused a stagnation region to form below the inlet; the approximate stagnation streamline can be seen in Figure 52 directly below the model leading edge. Pressure coefficients could not be calculated because the freestream velocity was zero. This would cause a division by zero in the pressure coefficient equation. Instead, to gain a quantitative understanding of the flow field, the velocity magnitude was calculated and superimposed onto streamlines in Figure 52 using the following equation: 22 vuv ?= (15) 79 Figure 52. Quiescent flow PIV measurements Re? = 0 Inlet vortices were known to require a disturbance to the flow field, therefore a disturbance in the form of a cross-face motion was applied to the flow field near the inlet leading edge. The formation of the vortex in quiescent flow is presented in Figure 53 as still images each with a 20/32 second time difference between them. Flourescein dye was injected upstream of the inlet near the wall of the tunnel. The dye traveled along the bottom surface until it was ingested into the inlet. The vortex formation can be seen in frame a of Figure 53. The vortex formation was found to be a transient, unstable process, as can be seen from frames a ? d. After frame d, the vortex became relatively stable, as can be seen in as frames e ? h being almost identical. 80 Figure 53. Inlet Vortex Formation Re? = 0, dt = 20/32 seconds Injection of dye into the flow field introduced a small, but non-negligible amount of perturbations to the otherwise laminar freestream flow in the x-axis direction. This addition of momentum was assumed to cause slight adverse effects on the flow field. As an alternative method to capture the inlet vortex, a flourescein salt crystal was sprinkled on the ground plane in the vicinity of the region where the inlet vortex was observed to form. As the crystal dissolved, a solution of fluorescent dye was produced and contained in the termination region of the vortex on the ground plane. The fluorescent dye contained in this region can be seen as an intense bright region in Figure 54. a b c d e f g h 81 Figure 54. Inlet Vortex Stagnation Region, Re? = 0 The large axial velocities produced by the vortex caused a drop in pressure in the vortex core. This drop in pressure entrained the fluorescent dye into the vortex core. The entrained dye allowed the tight vortex to be seen. The inlet vortex originated on the ground plane and terminated at the porous screen, representative of the fan-face, as can be seen in Figure 54. 82 6 Computational Fluid Dynamic Simulations 6.1 Simulation of CFM56-5B A vector plot on the symmetry plane is presented in Figure 55. The flow bifurcates at the stagnation points of the duct leading edges, which is represented by dark blue, and either flows externally around the duct, or is ingested into the engine. The flow that is ingested decelerates, seen as colder coloration, while the flow that goes externally around the duct accelerates, seen as warmer coloration. Figure 55. Vector Flow Field of CFM56-5B at an Angle of Attack of two degrees Deceleration causes the pressure coefficient to increase, while acceleration causes the pressure coefficient to decrease. This can be seen in pressure coefficient contour plots with the grid mesh superimposed onto the plots, Figure 56a. and Figure 56b. The pressure 83 coefficient is referenced from the free-stream pressure and is defined by the equation for pressure coefficient: ? ??= q ppC p (16) Therefore, a pressure coefficient of zero represents regions where the pressure is equal to the free-stream pressure, as can be seen in Figure 56a. as areas of bright green. Figure 56. Flow field pressure coefficient contour plots for the CFM-56-5B Mach number contours are plotted in Figure 57a and Figure 57b to show the manner in which the geometry of the duct affects the mach number of the flow field (Note that the two figures do not have the same scale). The fan inlet for the engine in both figures can be seen to have a Mach number very near 0.4 as specified. Since the pressure value set at the pressure outflow boundary condition was simply an estimate assuming the freestream stagnation pressure, the solution confirms that this was in fact a reasonable assumption for the fan face. a b 84 Figure 57. Flow field Mach number contour plots for the CFM-56-5B As previously stated, the boundary layer on the surface of the duct was resolved using an aspect ratio boundary layer mesh. A closer examination of the vector flow field of Figure 58 reveals that the boundary layer growth along the lengths of the surfaces was indeed resolved and the growth of the boundary layer is evident. Figure 58. Boundary Layer Growth of CFM56-5B AoA: 2? a b 85 6.2 Optimized Geometry The initial GA run was initiated for a total run of 10 generations with each generation having 20 members. The run was initiated over the computer cluster as defined previously. The GA found and selected the optimum goal values at the end of this run and the results are shown in Figure 59 and Figure 60. Shown in Figure 60 is the generational improvement of the total pressure ratio at the fan inlet of the engine, with the baseline geometry results presented as a solid black line. It should be noted that for this run, the total pressure ratio had a lower weighting factor than the flow uniformity, which was the primary goal. As such, only marginal improvement could be expected for this parameter as the run progressed. This is evident from results. The flow uniformity showed a marked improvement over the course of the GA run. Figure 59 shows the results for this parameter. Being the primary goal for this run, it was found to improve by 25% over the baseline CFM-56-5B case. The improvement over the best performer for the initial generation was 5% and the improvement over the worst performer of the initial generation was found to be 9.5%. 86 Figure 59. Flow uniformity: best & worst performers vs. GA Generations The total pressure improvement over the initial generation best performer was found to be equal to 0.08%, shown in Figure 60. The improvement over the worst performer of the initial generation was found to be equal to 0.2%. The GA generated best performer however had a 1% decrease in total pressure compared to the baseline CFM56- 5B geometry. This was assumed to be due to the larger weighting given to the flow uniformity. 87 Figure 60. Total pressure ratio: best & worst performers vs. GA generations The method used to obtain an average of each design goal at the fan face was determined to correctly model the fitness improvement trends. However, it could not be assumed that that the values were quantitatively accurate. This was determined to be adequate to prove the abilities of the optimization method and setup. A comparison of the optimized and baseline geometries is presented in Figure 61. As can be seen in the figure, the optimization decreased both the external chord lengths from the baseline geometry. 88 Figure 61. Comparison of the original and the optimized CFM-56-5B Geometry 89 Conclusions The capture area was found to be directly proportional to the freestream Reynolds number, but independent of the incidence angle of the model. An inlet vortex was determined to be a steady well behaved vortex, caused by cross-flow disturbances introduced to the flow field. Computational modeling of a full scale turbofan geometry results suggest that the assumptions of grid spacing and y+ values for the mesh used were reasonable. The computational model of the CFM56-5B engine was used as the baseline for a geometry optimization study to improve total pressure and flow uniformity at the fan-face. The improved turbofan inlet geometry, compared to the CFM56-5B model, had a 25% increase in flow uniformity while at the cost of a 1% decrease in total-pressure at the fan- face. The hardware and software setup was used for the optimization study was found to be a fast, proven, robust, and effective system. 90 References [1] Raymer, D. P. (1989). Aircraft Design: A Conceptual Approach, AIAA Education Series, American Institute of Aeronautics and Astronautics, Inc, Washington, D.C. [2] Mattingly, J.D., Heiser, W.H., Daley, D.H., Aircraft Engine Design, AIAA Education Series, AIAA, New York, 1987, pp. 354-393 [3] Seddon, J., and Goldsmith, E. L., 1999, "Inlet aerodynamics," 2nd edition, AIAA education series, AIAA [4] Hill, P., Peterson, C., Mechanics and Thermodynamics of Propulsion, Addison- Wesley Publishing, Massachusetts, 1992. [5] Goldsmith. E.L., Seddon, J. 1993, ?Practical Inlet Aerodynamic Design?, AIAA Educational Series [6] Jakubowski, A.K., Luidens, R.W., ?Internal cowl-separation at high incidence angles?, AIAA 13th Aerospace Science Meeting, Pasadena, California, January 20-22, 1975 [7] Hurd, R., ?Subsonic pitot inlets ? high-speed high-incidence performance?, Rolls- Royce (Bristol) Report PD 2029, 1976 [8] Albers, J.A. and Miller, B.A. (1073) ?Effect of subsonic inlet lip geometry on predicted surface and flow Mach number distributions?. NASA TN D7446 [9] Fox, R.W. and Kline, S.J. (1962) ?Flow regimes in curved subsonic diffusers?. J. Basic Eng. Trans ASME. [10] Carter, E.C (1972) ?Experimental determination of inlet characteristics and inlet and airframe interference?. AGARD, LS 53 [11] Aulehla, F. (1982) ?Inlet swirl, a major disturbance parameter in engine/inlet compatibility. 13th Congress of ICAS/AIAA, Seattle, August 1982. [12] Hercock, R.G. and Williams, D.D. (1974) ?Distortion-induced engine instability aerodynamic response?. AGARD,LS72-Paper No. 3. [13] Colehour, J.L. and Farquhart B.W. Inlet Vortex. J. Aircraft, The Boeing Company, Seattle, Wash. Vol. 8. no. 1. Jan. 1971 91 [14] Smith, David, Dorris III John, inventors; McDonnel Douglas Corporation. Aircraft Engine Apparatus with Reduced Inlet Vortex. United States patent 6,129,309. 2000 Oct 10. [15] Yadlin, Y. ?Simulation of Vortex Flows for Airplanes in Ground operations?, AIAA Paper 2006-56,2006 [16] Laurent, R. "Jet Engine Ground Vortex Studies", Master?s Thesis, Cranfield University Sept. 2007. [17] Rodert, L. & Garrett, F. ?Ingestion of Foreign Objects into Turbine Engines by Vortices?, NACA TN 3330, 1955 [18] Kendall, J. M. Jr., ?Experimental Study of a Compressible Viscous Vortex,? TR No. 32-290, June 1962, Jet Propulsion Laboratory, C.I.T., Pasadena, Calif. [19] Johns, C. ?The Aircraft Engine Inlet Vortex Problem?, AIAA?s Aircraft Technology, Integrations and Operations (ATIO) 2002 Technical, 1-3rd 2002, Los Angles, California, AIAA 2002-5894 [20] Rubbert, P. E., et al., ?A General Method for Determining the Aerodynamics Characteristics of Fan-In Wing Configurations,? D-6-13476-1, Nov. 1967, The Boeing Company, Seattle, Wash. [21] ANSYS, Inc., Fluent 6.3 User?s Guide, 2008. http://fluent.com. [22] White, F.M., Viscous Fluid Flow, 3rd ed. McGraw Hill, New York, 1991 [23] Clauser, F. H., ?The Turbulent Boundary Layer,? in Advances in Applied Mechanics, Vol. IV, Academic Press, New York, 1956 [24] Dyer, J.D., ?Aerospace Design Optimization Using a Real coded Genetic Algorithm?, Master?s Thesis, Auburn University, Alabama, May 2008. [25] Doyle, J., Hartfield, R.J., and Roy, C. ?Aerodynamic Optimization for Freight Trucks using a Genetic Algorithm and CFM?, AIAA 2008-0323, presented at the 46th Aerospace Sciences Meeting and Exhibit, Reno, NV, January 2008. [26] Torella, G., Blasi, L., ?The Optimization of Gas Turbine Engine Design by Genetic Algorithms?, AIAA Paper 2000-3710, 36th AIAA/ASME/SAE/ASEE Joint Propulsion Conference and Exhibit, July 2000. [27] J.E. Burkhalter, R.M. Jenkins, and R.J. Hartfield, M.B. Anderson, G.A. Sanders, ?Missile Systems Design Optimization Using Genetic Algorithms?, AIAA Paper 2002-5173, Classified Missile Systems Conference, Monterey, CA, November, 2002 92 [28] Hartfield, Roy J., Jenkins, Rhonald M., Burkhalter, John E., ?Ramjet Powered missile Design Using a Genetic Algorithm?, AIAA 2004-0451, presented at the forty- second AIAA Aerospace Sciences Meeting, Reno NV, January 5-8, 2004. [29] Anderson, M.B., ?Using Pareto Genetic Algorithms for Preliminary Subsonic Wing Design?, AIAA Paper 96-4023, presented at the 6th AIAA/NASA/USAF Multidisciplinary Analysis and Optimization Symposium, Bellevue, WA, September 1996. [30] Perez, R.E., Chung, J., Behdinan, K., ?Aircraft Conceptual Design Using Genetic Algorithms?, AIAA Paper 2000-4938, presented at the 8th AIAA/USAF/NASA/ISSMO Symposium on Multidisciplinary Analysis and Optimization, Bellevue, WA, September 2000. [31] 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. [32] Kulfan B., M., ?A Universal Parametric Geometry Representation Method ? CST?, 45th AIAA Aerospace Sciences Meeting and Exhibit, 8-11 January, Reno, Nevada, 2007. [33] Burger, C. and Hartfield, R.J., ?Propeller Performance Optimization using Vortex Lattic Theory and a Genetic Algorithm?, AIAA-2006-1067, presented at the Forty- Fourth Aerospace Sciences Meeting and Exhibit, Reno, NV, Jan 9-12, 2006. [34] Rifki, R., ?Flow Around Axisymmetric and Two-Dimensional Forward-Facing Cavities?, Master?s Thesis, Auburn University, Alabama, May 2006 [35] Ahuja, V., ?Optimization of Fuel-Air Mixing for a Scramjet Combustor Geometry using CFD and a Genetic Algorithm?, Master?s Thesis, Auburn University, Alabama, Dec. 2008. [36] Ahuja, V. and Hartfield, R.J., ?Optimization of Fuel-Air Mixing for a Scramjet Combustor Geometry using CFD and a Genetic Algorithm?, Multidisciplinary Analysis and Optimization Conference (MDO), Victoria City, Vancouver, Canada, September 2008. [37] Anderson, M. B., ?Using Pareto Genetic Algorithms for Preliminary Subsonic Wing Design,? AIAA Paper 96-4023, presented at the 6th AIAA/NASA/USAF Multidisciplinary Analysis and Optimization Symposium, Bellevue, WA, September 1996. 93 Appendices 1 CST Transformation & Bernstein Polynomials There are a number of possible methods to approximate equations of airfoil geometries. One of the most efficient methods is the class function, shape function transformation (CST). Presented here is the theory developed by Kulfan [31][32]. The general equation of an airfoil can be represented as ( ) ( ) ( ) TESn ??????? ?+???= 1 (17) The first two terms of Eq. (14) form the class function, Sn is the shape function, and the final term is the trailing edge thickness. In the physical space Eq. (14) becomes: ( ) ( ) ( ) ( )TEczcxcxSncxcxcxcz ?+???= 1 (18) A graphical representation of Eq. (14) is shown in Figure 62 . Figure 62. General equation representation for a round LE, sharp TE airfoil [32] 2?zTE ? = x/c ?= z/c ?TE = ?zTE/c 94 The class function controls the LE radius, airfoil thickness, and boatail angle while the shape function determines the geometry between the LE and TE. Depending on the shape function chosen, a number of possible shapes in the design space can be modeled. The class function in the design space is defined as [31][32]. ( ) ( ) [ ] 211 2 1 NNNNC ??? ??= (19) with ? being the fraction of the chord. In the physical space the unit class function is: ( ) ( ) [ ] 21 1 NN cxcxcxC ??= (20) where c is the chord length and x is the x position along the chord. The first term of Eq. (17) affects the shape of the leading edge and the last term affects the shape of the trailing term, as defined by N1 and N2. An example of possible shapes applicable to airfoils can be seen in Figure 63. N1 = 1.0, N2 = 1.0 N1 = 0.5, N2 = 1.0 N1 = 0.5, N2 = 0.5 Figure 63. Class Function 2D Design Space The entire airfoil was represented by the combination of two unit class functions, for the upper and the lower surfaces, multiplied by a unit shape function. 95 The unit shape function is defined by a BP of the order n with the variable x/c ranging from 0-1.0. The BP?s were chosen due to the mathematical property of ?Partition of Unity? which states that all of the terms of a BP sum to equal 1.0 [32]. This is a desired quality because all powers of BP are equal. So, any order of n for BP can be used. If BP of different orders of n were not all equal, then BP could not be used for this method. The shape function is defined as [32] ( ) ( ) ( )? = ??? ?= n r rn nr cxrnr nxS 0 , 1!! ! (21) The first term of the equation defines the binomial coefficients with increasing order of n for the BP. For each order of n, there exists n + 1 number of total terms in the BP. The first term of Equation (18) defines the leading edge radius by the equation: First Term cRLE2= (22) The second term of Eq. (18) defines the boatail angle by the equation: Last Term czTan ?+= ? (23) The remaining terms, those between the first and last, are known as shaping terms, which affect the geometry of the airfoil between the leading and trailing edge. 96 The components of the peaks of the shape functions are equally spaced along the chord, for i ranging from 0 to n, which is determined by the equation: ( ) nicx S =max (24) The components of the peaks of the component airfoils are also equally spaced along the chord, for i ranging from 0 to n, which is determined by the equation: ( ) ( ) ( )nNNiNcx Z +++= 211max (25) Increasing the order of n increases the ability of the shape function to equate airfoil geometries. This is because as the order of n increases, the number of peaks to be able to match the desired geometry increases [32]. A negative effect of increased order of n is increased computational times. Eventually increasing the order of n past a critical number will have diminishing returns; improvements in accuracy between orders of n will be small enough to be negligible and will not be worth the increased computational times. Values of n = 6-9 for BP accurately modeled airfoil geometries [32]. A Bernstein Polynomial of order n = 6 was chosen because it is a high enough order of n to accurately model an airfoil9 and has decreased computational times compared to higher orders of n. The entire 2D duct geometry that is to be curve fitted, with the exception of the turbine blade geometry and nose cone, can be seen in Figure 64. 97 Figure 64. CFM56-5B Duct Geometry to Equation Fit The entire airfoil can be created with two unit class functions, one for the lower surface and one for the upper surface. The entire 2D duct, consisting of 2 airfoils, can be created with 4 class functions. Equations 23 & 24 define the entire duct geometry. The first two terms of each equation are the class function. The last terms represent the shape function, BP of order n = 6. The variables N1 and N2 were set equal to 1.0 and 0.5 to create a class-function with a rounded leading edge and a sharp trailing edge8. The variables of A were modified [33] to apply to the CFM56-5B coordinates and preliminary bounds were set at -100 to 100 to allow for a larger design space. This 2D class/shape function for each curve consists of 9 variables. Each airfoil consists of an upper and lower surface, which are defined by two different sets of BP equations. The entire 2D duct consists of two airfoils, with two surfaces each, for a total of four surfaces and 36 variables for all of the surfaces. 98 Airfoil, Outer Surface definition: ( ) ( )[ ] ( ) ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( )[ ] zcx cxAcxcxA cxcxA cxcxA cxcxA cxcxAcxA cxcxz NNOuter *)/( 1*1***6 1***15 1***20 1***15 1***6* *1* 6 7 51 6 42 5 33 4 24 3 15 2 6 1 21 + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?+?+ ?+ ?+ ?+ ?+ ?= (26) Airfoil Inner Surface definition: ( ) ( )[ ] ( ) ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( ) ( )[ ] ( )[ ] zcx cxAcxcxA cxcxA cxcxA cxcxA cxcxAcxA cxcxz NNInner *)/( 1*1***6 1***15 1***20 1***15 1***6* *1* 6 14 51 13 42 12 33 11 24 10 15 9 6 8 43 + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?+?+ ?+ ?+ ?+ ?+ ?= (27) 99 2 Example GA Curve Fitting Design Space Having fixed the input parameters, the GA design space could finally be defined by fixing the range of values above and below the reference values. The resulting design space is shown in Figure 65. In Equation 23, the variables A1-7 were given bounds of -200 and 200. One can see in Figure 65 that the boundary conditions of the Genetic Algorithm used to fit the data. The comparison of the maximum and minimum lines to the known data points of the duct geometry proved that a mathematical solution existed in these boundaries. -50 -40 -30 -20 -10 0 10 20 30 40 50 0 20 40 60 80 100 120 140 x-axis location, in z-a xis lo ca tio n, in Data to Fit Max Bound Min Bound Figure 65. GA Design Space 100 3 Example Curve Fitting of Bernstein Polynomials to CFM56-5B Cowl The number of generations that the GA was run for was increased from 50 generations to 200,000 generations to determine the ability of the GA to match the data. It can be seen that the more generations that the GA was run for curve fitting the upper surface of the top airfoil, the closer the equations came to matching the data. The curve fitting of the upper surface of the top airfoil is represented in Figure 66 and Figure 67. All four curves that create the duct geometry were curve fit, but only the fitting of the upper surface are presented here. 0 2 4 6 8 10 12 0 20 40 60 80 100 120 140 x-axis location, in z-a xis lo cat ion , in Data to Fit GA Output Figure 66. GA Equation Fitting after 50 Generations 101 0 2 4 6 8 10 12 0 20 40 60 80 100 120 140 x-axis location, in z-a xis lo cat ion , in Data to Fit GA Output Figure 67. GA Equation Fitting after 200,000 Generations At approximately 10,000 generations, the returns of increasing generations were greatly diminished (Table 1). In Figure 68 a hyperbolic relationship between generations and fitness can be seen, with an asymptote at approximately fitness equal to 0.2. For a perfect curve fit, the fitness would be equal to 0. This discrepancy can be attributed to errors in recording geometry data from the available CFM56-5B duct geometry as well as limits of BP?s ability to fit the curve. However, the maximum distance between the CFM56-5B duct geometry and the BP representation of the curve is .0944 in, less than a tenth of an inch. 0 0.5 1 1.5 2 2.5 3 3.5 0 20,000 40,000 60,000 80,000 100,000 120,000 140,000 160,000 180,000 200,000 # of Generations Fit ne ss 102 Figure 68. Equation Fitness with Increasing Generations of GA Run Table 8. Fitness Ratio Data of Increased Generations Generations Fitness ?Fitness % Difference 500 4.800183773 0 0.0000% 1,000 2.863285065 1.936898708 40.3505% 2,000 1.610103011 1.253182054 43.7673% 5,000 0.860965073 0.749137938 46.5273% 10,000 0.39534381 0.465621263 54.0813% 20,000 0.386857092 0.008486718 2.1467% 50,000 0.258321494 0.128535599 33.2256% 100,000 0.224858895 0.033462599 12.9539% 150,000 0.224622801 0.000236094 0.1050% 200,000 0.224622801 0 0.0000% 103 4 PIV Uncertainty Analysis To ensure confidence in the data obtained using PIV, an uncertainty analysis was conducted. Standard deviation is a gauge of the uncertainty of the data. Therefore the standard deviation for each dataset was obtained to conduct the analysis. Shown below is an example standard deviation contour plot. Figure 69. Example standard deviation, Re? = 22000 The standard deviation was obtained for all test cases. Shown in Figure 70 is the standard deviation for each test case. As can be seen in the figure, the uncertainty increased as Reynolds number increased for all angles of incidence. Interestingly, the PIV measurements taken when the model was placed at an angle of incidence of zero had the largest uncertainty of all measurements. PIV measurements at a negative angle of attack had the least uncertainty. 104 Figure 70: PIV Standard deviation: Variation in angle of incidence and Re The average of the standard deviations for all datasets was found to be % 2.42. This uncertainty value was assumed to negligible. Therefore the measured data was concluded to be precise.