OPTIMIZATION OF FUEL-AIR MIXING FOR A SCRAMJET COMBUSTOR GEOMETRY USING CFD AND A GENETIC ALGORITHM Except where reference is made to the work of others, the work described in this thesis is my own or was done in collaboration with my advisory committee. This thesis does not include proprietary or classified information Vivek Ahuja Certificate of Approval: John Burkhalter Roy Hartfield, Chair Professor Emeritus Professor Aerospace Engineering Aerospace Engineering Rhonald Jenkins Andrew Shelton Professor Emeritus Assistant Professor Aerospace Engineering Aerospace Engineering George T. Flowers Dean Graduate School Auburn Alabama December 19, 2008 OPTIMIZATION OF FUEL-AIR MIXING FOR A SCRAMJET COMBUSTOR GEOMETRY USING CFD AND A GENETIC ALGORITHM Vivek Ahuja A Thesis Submitted to the Graduate Faculty of Auburn University in Partial Fulfillment of the Requirements for the Degree of Master of Science iii OPTIMIZATION OF FUEL-AIR MIXING FOR A SCRAMJET COMBUSTOR GEOMETRY USING CFD AND A GENETIC ALGORITHM Vivek Ahuja Permission is granted to Auburn University to make copies of this thesis at its discretion, upon request of individuals or institutions and at their expense. The author reserves all publication rights. Signature of Author Date of Graduation iv THESIS ABSTRACT OPTIMIZATION OF FUEL-AIR MIXING FOR A SCRAMJET COMBUSTOR GEOMETRY USING CFD AND A GENETIC ALGORITHM Vivek Ahuja Master of Science, December 19, 2008 (B.S., University of Pune, 2007) 124 Typed pages Directed by Roy Hartfield A new methodology for the optimization of fuel-air mixing in a scramjet combustor using integrated Genetic Algorithms and Computational Fluid Dynamics is presented. A typical combustor design involving Mach 2 crossflow over a rearward facing step with staged normal injection is considered for study and is optimized using this method. The CFD results are validated against experimental results prior to optimization to allow for grid refinement and high accuracy of results. Quantification of typical v combustor performance and design parameters is discussed and adaptation for use with CFD grids is presented. An integrated system of computers and software designed for fast computation times has been created. Correlations between variations in physical geometry and optimization of fuel-air mixing are presented. vi ACKNOWLEDGMENTS The author would like to thank Dr. Roy Hartfield for providing the encouragement, support and a free-thinking environment for conducting this research. Thanks are also due to Dr. Rhonald Jenkins for helping improve the mathematical models used in the study and to Dr. John Burkhalter, Dr. Andrew Shelton and Scott Thomas for helping with the development of the integrated software used in the research. Thanks are also due to my family members for their back stage support during a wonderful year of work. vii Style manual or journal used: American Institute of Aeronautics and Astronautics Journal Computer software used: FLUENT 6.2.17, GAMBIT, TECPLOT 10/360, FORTRAN, AUTOCAD 2008 viii TABLE OF CONTENTS LIST OF FIGURES .................................................................................................. ix NOMENCLATURE ................................................................................................xii 1. INTRODUCTION ............................................................................................. 01 2. ANALYSIS AND RESEARCH ......................................................................... 04 2.1 Combustor Section Geometry and Design .................................................... 04 2.2 Grid Development and Refinement .............................................................. 07 2.3 Computational Methods and Solvers ............................................................ 14 2.4 Boundary Conditions .................................................................................... 21 2.5 Experimental Validation ............................................................................... 25 2.6 Flow without injection .................................................................................. 26 2.7 Flow with injection ...................................................................................... 33 2.8 Mass Fraction analysis ................................................................................. 42 2.9 Genetic Algorithm and CFD Integration ....................................................... 44 2.10 Design Parameters and Design Space ......................................................... 48 2.11 Design Goals .............................................................................................. 51 2.12 The Networked Environment ...................................................................... 56 3. THE RESULTS ................................................................................................. 61 3.1 Case Run-1................................................................................................... 61 3.2 Case Run-2................................................................................................... 75 3.3 Case Run-3................................................................................................... 82 3.4 Case Run-4................................................................................................... 89 4. CONCLUSIONS AND RECOMMENDATIONS ............................................. 96 REFERENCES ........................................................................................................ 98 APPENDIX ............................................................................................................ 101 ix LIST OF FIGURES Figure 2.1.1. Combustor section geometry in 3-D ...........................................................5 Figure 2.2.1. Flow field features for the non-injection case .............................................7 Figure 2.2.2. Two dimensional grid refinements .............................................................8 Figure 2.2.3. Grid refinement in greater detail and Zone interactions ..............................9 Figure 2.2.4. Combustor section grid for the non-injection case .................................... 11 Figure 2.2.5. Grid structure around the Injection surfaces ............................................. 13 Figure 2.2.6. 3D visualization of combustor section for the injection case ..................... 13 Figure 2.3.1. Residual plots for the Hybrid solver model............................................... 19 Figure 2.6.1. Pressure field comparison for the non-injection case ................................ 27 Figure 2.6.2. Static pressure profile comparison at the primary injector ........................ 28 Figure 2.6.3. Simulation results for the pressure contours in a 3D view ........................ 29 Figure 2.6.4. Comparisons for the temperature field for the non-injection case ............. 30 Figure 2.6.5. 3D visualization of the temperature field ................................................. 30 Figure 2.6.6. Static temperature profile comparison at the primary injector ................... 31 Figure 2.6.7. Comparisons for the velocity field for the non-injection case ................... 32 Figure 2.7.1. Flow features at the centerline section for the injection case ..................... 34 Figure 2.7.2. 3D visualization for the pressure field in the injection case....................... 34 x Figure 2.7.3. Nature of the injectant flow through the combustor .................................. 35 Figure 2.7.4. Comparisons for the pressure field for the injection case .......................... 36 Figure 2.7.5. Comparisons for the static pressure profile at the primary injector ............ 37 Figure 2.7.6. Comparisons for the static pressure profile at the secondary injector ........ 37 Figure 2.7.7. Comparisons for the temperature field for the injection case..................... 38 Figure 2.7.8. Comparisons for the static temperature profile at the primary injector ...... 39 Figure 2.7.9. Comparisons for the static temperature profile at the secondary injector ... 40 Figure 2.7.10. 3D visualization of the temperature field...................................................39 Figure 2.8.1. Comparisons for the injectant mole fraction distribution.......................... 43 Figure 2.8.2. Mass Fraction distributions along the combustor length. ......................... 43 Figure 2.10.1. The Design Space for the optimization of the Scramjet Combustor ......... 50 Figure 2.12.1. An overview of the network structure at Auburn University .................... 57 Figure 2.12.2. A concept of operations view for a GA run ............................................. 58 Figure 2.12.3. Computational time versus number of processors .................................... 59 Figure 3.1.1. Total pressure ratio versus the generational best and worst performers .... 63 Figure 3.1.2. Mixing efficiency versus the generational best and worst performers ...... 64 Figure 3.1.3. Primary injector axial variation for the best and worst performers ........... 66 Figure 3.1.4. Secondary injector axial variation for the best and worst performers ....... 67 Figure 3.1.5. The primary and secondary injector variation for the best performers ...... 67 Figure 3.1.6. Optimizing the injector Axial (X) location .............................................. 69 Figure 3.1.7. Primary injector transverse variation for best and worst performers ......... 70 Figure 3.1.8. Secondary injector transverse variation for best and worst performers ..... 70 Figure 3.1.9. Effect of bringing the injectors closer to each other. .................................. 73 xi Figure 3.2.1. Single Dimension Longitudinal Design Space for Case Run 2 ................... 76 Figure 3.2.2. Total Pressure Ratio for best and worst performers versus generations ...... 78 Figure 3.2.3. Mixing Efficiency for best and worst performers versus generations ......... 78 Figure 3.2.4. Movement of the primary and secondary injectors versus generations ....... 80 Figure 3.2.5. Relative distance between the injectors versus GA Generations ................ 81 Figure 3.3.1. Single Dimension Longitudinal Design Space for Case Run 3 ................... 84 Figure 3.3.2. Mixing Efficiency for best and worst performers versus generations ......... 85 Figure 3.3.3. Total Pressure Ratio for best and worst performers versus generations ...... 86 Figure 3.3.4. Primary and secondary injectors versus generation numbers ..................... 88 Figure 3.4.1. Mixing Efficiency for best and worst performers versus generations ......... 90 Figure 3.4.2. Total Pressure Ratio for best and worst performers versus generations ...... 91 Figure 3.4.3. Movement of the primary and secondary injectors versus generations ....... 92 Figure 3.4.4. Movement of the primary and secondary injectors versus generations ....... 94 Figure 3.4.5. Relative distances between the two injectors for Case Run 4 ..................... 94 Figure 3.4.6. Final results for the optimized injector locations ....................................... 95 xii NOMENCLATURE u Local velocity at any point in a grid ?g2206 Distance between two grid nodes in a grid space ?g2202 Time step g2247 Courant number N Total number of grid points in a grid surface g2251s Average static density for a grid surface g2251? Average total density for a grid surface g2251?REF Reference total density g2761g2201g4666g2191g4667 Static density at any nodal point ?i? in the grid surface Ps Average static pressure for a grid surface P? Average total pressure for a grid surface P?REF Reference total pressure g1790g2201g4666g2191g4667 Static pressure at any nodal point ?i? in the grid surface Ts Average static temperature for a grid surface T? Average total temperature for a grid surface T?REF Reference total temperature g1794g2201g4666g2191g4667 Static temperature at any nodal point ?i? in the grid surface xiii M Mach number g2169g2206g4666g2191g4667 Mach number at any nodal point ?i? in the grid surface calculated using the ?X? axis component of velocity g2204g2206 Velocity component of the fluid along the ?X? axis direction g2204g2206g4666g2191g4667 Velocity component of the fluid at any grid point ?i' along the ?X? axis direction g2183g2174g2161g2162g4593 Reference sonic velocity of the fluid g2237 Specific heat ratio for the fluid R Universal gas constant 1 1. INTRODUCTION For more than three decades extensive effort has been directed towards deploying an operational Scramjet. However, an operational Scramjet is viable only if every subsystem works at near optimum levels. Optimization of the subsystems is now possible using Genetic Algorithms (GA) in conjunction with CFD solvers. This effort presents the methodology for such optimization processes 1-6. The study of fuel-air mixing in a supersonic cross-flow has been looked at as potential scramjet combustor geometry. With the development of computing technology, it is possible to develop optimized preliminary designs for scramjet combustors using a CFD solver and a GA 7. Experimental results from research conducted in the 1990s have been used for the validation of CFD solutions 8. The experiments were highly focused on developing accurate data sets for a single case flow situation 10-11. This effort builds on this single validated case by considering geometric variations of the combustor design, solving for the flow and arriving at a geometry which is optimized for mixing efficiency with minimum total pressure loss. 2 This thesis includes validation efforts and the results obtained from the optimization effort. The fuel-air mixing phenomena inside a scramjet combustor section are modeled as a function of injector locations and geometric aspect ratio changes. The GA is a tournament based, binary encoded FORTRAN code and the CFD solver is a Reynold?s Averaged Navier Stokes code developed by FLUENT. The combustor section is taken to be a rectangular cross-sectioned type with a rearward facing step 8, 10-13. The inlet velocity at the entry plane of the combustor section upstream of the rearward facing step is taken to be Mach 2 throughout the course of the optimization effort, in accordance with experimental results 8. Two injectors downstream of the rearward facing step which are moved in both the axial (X) direction and the transverse (Z) directions independently during the course of the optimization. In effect, this characterizes the fuel injection into the combustor as being a ?Staged Normal Injection? type. For the purpose of demonstrating the viability of the GA driven CFD approach, this is the only geometric variation considered during the present effort. The initial flow simulation results obtained from this effort were compared with the results obtained from the experiments9-11 to validate the accuracy of the simulation result. Upon completion of this validation effort, the next step involved the automation and integration of the above process with a Genetic Algorithm. A centralized program structure was created to distribute the computation effort over a large computer network available at Auburn University. Further studies were then conducted to reduce the 3 computational time at both the CFD level as well as the external hardware level to make this optimization effort feasible. To summarize, the objective of this effort was to develop the methodology for optimizing a scramjet combustor fuel-air mixing performance through the combined use of a GA and CFD as enumerated below: A) Development of the flow-field grid and the validation of the accuracy of the CFD results using experimental data. B) Integration of the CFD and GA software codes into a networked system of computers and improving the robustness, reliability and computational efficiency. C) Conducting case runs on the above system with standard scramjet combustor geometry and presenting the optimization results. D) Providing the interpretation of the GA results in terms of actual physical changes in the combustor geometry to determine trends and correlations and also creating the basis for future work. 4 2. ANALYSIS AND RESEARCH 2.1 COMBUSTOR SECTION GEOMETRY AND DESIGN The combustor section in this optimization effort is based on the generic rectangular flow cross-section design with a rearward facing step as shown in Figure 2.1.1(a) and Figure 2.1.1(b). The values of the various geometries of the combustor design are presented in Table 2.1.1. The locations of the injection points are also shown in Figure 2.1.1(a) and 2.1.1(b). Known as the primary and secondary injectors, they are placed at three and seven step heights downstream of this step respectively. The availability of experimental results dictated the initial values of these dimensions to enable a comparative study of the simulation and experimental results 8. The nature of the optimization objectives meant that these locations would thereafter vary. These initial values represent a validated starting point for the optimization effort. 5 Figure 2.1.1(a) Combustor section geometry in 3-D Table 2.1.1 Combustor Geometry Injection surface diameter (D) 1.93 mm Combustor section height 21.29 mm Combustor section width 30.48 mm Step height 3.18 mm First injector location (Initial) 9.54 mm Second injector location (Initial) 22.26 mm Step location along axial direction 11.02 mm 6 Figure 2.1.1 (b) Combustor section geometry in the plan view 7 2.2 GRID DEVELOPMENT AND REFINEMENT The first grid refinement objective was the refinement of a grid for the two dimensional flow field using the axial (X) and vertical (Y) directions. To facilitate this process, the first flow simulation was conducted using no fuel injection, and was used to refine the grid in the X-Y plane in the region downstream of the rearward step where the shear layer was expected to develop as a result of the rearward facing step. In the absence of injection of fuel, this shear layer is seen to be symmetrical in all X-Y planes along the transverse (Z) axis except in the near wall regions. Figure 2.2.1 presents the details of the flow-field in the case of no injection. Figure 2.2.1. Flow field features for the non-injection case 8 8 The grid was divided into two regions: above the step and behind the step. The region above the step is where the expansion fan and the single oblique shock downstream of this fan were expected to be seen. It was acknowledged that this is a relatively simple flow field to be simulated in FLUENT, and that the grid need not be highly refined. Hence, a relatively coarse grid was developed to maintain computational efficiency. A more refined mesh zone surrounding the expected region of the shear layer for the given flow conditions was developed. The main effort at this stage was to achieve good fidelity between the experimental data set and the simulations with a minimal number of grid points. The grid for the X-Y plane is shown in Figure 2.2.2. The refined region and the zone segregations are shown in Figure 2.2.3(a) and Figure 2.2.3(b). This decentralized grid structure for the X-Y plane is held constant for the remainder of the analysis. Figure 2.2.2. Two dimensional grid refinements of the Combustor X-Y plane and various nodal zones 9 Figure 2.2.3(a) (left) and Figure 2.2.3(b) (Right). Grid refinement in greater detail and Zone interactions (Combustor X-Y plane) The exact values of the nodes varied for each grid zone. Table-2.2.1 presents the number of nodes for the 2-dimensional flow case for the various grid zones. As can be seen, using double the number of zones as normally required helped to make the grid very simple to construct and visualize. This concept was then carried into the 3- dimensional cases as well. One of the focal points during this effort was fixing the exact values of the node numbers within the structured grid. Accuracy and computational time are always directly proportional to each other. The higher the demands for accuracy, the higher were the computational time for getting the CFD solutions. And the higher the computational times, the higher was the run times. At the same time accuracy could not be sacrificed to save on computational time. Unique methods were developed to reduce computational time from the networked structure and this is discussed later. Suffice to say that high 10 accuracy could be maintained at this early stage of the effort. An incremental approach was initiated to determine the value of the nodal densities beyond which no significant change in accuracy was noticed. This value of the nodal density for the various zones was then fixed and maintained throughout the analysis. Table 2.2.1 Node points per zone Zone Number of Node points (X*Y) 1 275 (25x11) 2 2500 (25x100) 3 3000 (30x100) 4 2750 (55x50) Given the rectangular nature of the flow field geometry, and the initial case of no injection surfaces in the X-Z plane, the 3-dimensional grid was relatively simple to develop and merely extended the existing grid network into the Z coordinates. This gave the required 3-dimensional mesh as shown in Figure 2.2.4. All existing zones were extended into their respective volumes and correspondingly meshed. Table-2.2.2 presents the number of nodes for each of the volumes within this 3-dimensional environment. The total number of nodes for the 2-dimensional case amounted to 8,525 for a single 11 combustor X-Y plane and those for the 3-dimensional case amounted to 170,500 for the entire combustor volume. Figure 2.2.4. Three dimensional visualization of combustor section geometry for the non-injection case At no point has the axial-length of the combustor section been stated. One of the reasons for this was that this value needed to be determined from external data and is very much a part of the optimization process and thus not arbitrary. However, not having a fixed combustor axial length presented a unique problem during the construction of the mesh. This problem was a result of the fact that an output value of the flow parameters was required in terms of pressure and mass flow rate. However, this aspect is more clearly discussed in later sections of this thesis. 12 Table 2.2.2 Nodes per grid volume Volume Number of Nodes No. (X*Y*Z) 1 5500 (25x11x20) 2 50000 (25x100x20) 3 60000 (30x100x20) 4 55000 (55x50x20) The second step in the grid development involved refinement around the injection locations on the floor of the combustor in the X-Z plane and the results were again compared with the experimental data. For this flow case, the simulations were done 3- dimensions with the X-Y grid being maintained the same as the one achieved during the first step of refinement discussed above. The X-Z plane mesh was extensively altered once the treatment of the injection surfaces was considered. The circular nature of these surfaces complicates the mesh and is shown in Figure 2.2.5. The grid zone structures remained roughly the same, though the minor differences are highlighted in Figure 2.2.6. 13 Figure 2.2.5. Grid structure around the Injection surfaces on the floor of the Combustor (X-Z plane) Figure 2.2.6. 3-Dimensional visualization of combustor section geometry for the injection case with mesh distortions 14 2.3 COMPUTATIONAL METHODS AND SOLVERS When considering the combined effects of supersonic flow over a rearward facing step coupled with the shock interactions of a staged normal injection process, it quickly became obvious from experimental data that the required flow-field was highly complex. It was also clear that since the scope of this analysis moves far beyond the single case runs with no geometrical variation, the flow field could not be predicted in great detail in advance. This necessitated the need for a robust grid structure and a similarly robust CFD solver model. High speed was also necessary for the large number of runs required by the Genetic Algorithm. Consideration of the available solvers showed that the two available choices with the FLUENT software were as follows: a) Segregated Solver Model b) Coupled Solver Model The segregated solver represented a fast solution model that fit well with regard to the time constraint issues as discussed above, giving significantly faster solutions to a given flow field as opposed to the coupled model (both explicit and implicit versions). 15 However, a significant disadvantage was the segregated solver?s inability to handle highly complex flow fields as needed for the scramjet combustor analysis. Indeed, initial trials even with the non-injection cases using this model resulted in diverging solutions within a handful of initial calculations, necessitating a shutdown of the analysis. The problem only got worse when injection cases were attempted. On the other hand, the coupled solver proved to be a far more robust model at the fundamental level as opposed to the segregated model. This model was better suited for difficult and complex flow situations. The difference between segregated solvers and coupled solvers is in their method of solving the continuity, momentum and energy equations. While the segregated model solves these equations sequentially, the coupled solver solves them simultaneously, thus maintaining superior control over solution divergence. For these initial validation studies, the use of the implicit coupled model (which solves for all variables in all cells at the same time) was employed. The disadvantage of this model was efficiency due to a mathematical choke point built into the model that provides both the robustness of the model and its slow run times. To understand this more clearly, an understanding of the fundamental differences between the two solution models is in order. Essentially, the iteration scheme consists of the following steps: a) Fluid properties are initialized for each calculation based on the previous solution or the initial condition, whichever is the case. 16 b) The continuity, momentum and energy equations are solved simultaneously. Here, the system of governing equations is cast in an integral, Cartesian form for an arbitrary control volume ?V? with differential surface area ?dA? as follows: (2.3.1) Where the vectors ?W?, ?F? and ?G? are defined as: (2.3.2) and the vector ?H? contains source terms such as energy sources. All variables in the above matrices are defined on a ?per unit mass? basis. c) Equations for turbulence are solved based on the selected model. d) A convergence check is conducted and the cycle is repeated. The coupled set of governing equations as given above is discretized in time for steady and unsteady calculations. In this case, with the flow situations being steady, the time marching scheme proceeds until a steady state solution is reached. Temporal discretization of the coupled equations is accomplished by an implicit time-marching scheme. 17 The time-step is calculated from the CFL (Courant-Friedrichs-Lewy) condition (a condition for algorithms solving partial differential equations to be convergent; named after Richard Courant, Kurt Friedrichs and Hans Lewy who described it in their 1928 paper) 15, defined by the parameter: (2.3.3) This parameter is known as the Courant Number and is the critical parameter to be determined when using the coupled solver. Since the coupled solver is a time marching iteration technique, this number connects the refinement of the mesh with the time iteration values (and hence solution divergence issues). Essentially this works around the idea that for a time marching iteration technique, the time step must be less than the time taken by an event to cross between two successive grid points. This discussion involves the flow of fluids from one grid point to another grid point. This is an extremely important issue for achieving convergence during the iterations since the value of the Courant number varies for each grid. For all grids the Courant number was determined after an incremental advancement scheme to be 0.1, the maximum value to avoid solution divergence 12. One note to be made here is that the grid developed previously did not have a uniform grid structure. Indeed, there were several grid zones with each zone having its own nodal densities. Some zones, such as the one handling the shear layer, possessed varying nodal densities even within the zone thanks to the use of successive ratio schemes for deliberately increasing nodal densities near the walls and reducing density further 18 downstream to allow the grid to interact with the other zones. With such a structure, the Courant number must be set for the smallest grid distance for the entire grid. To improve computational efficiency, a unique method was developed to combine the segregated and coupled solvers into a single system designed to take advantage of the robustness of the coupled solver with the speed of the segregated solver during various stages of the calculation process. In essence, the idea was to use the coupled solver model initially, and have it run through a number of time steps to allow the flow field to develop around the rearward facing step; however, the coupled solver was not allowed to run through the length of the combustor section. The flow section covered by the coupled solver was the section of the combustor that the segregated solver had proven to be incapable of solving. Outside of this region, the flow field is relatively simple, and thus within the capabilities of the segregated solver (at much higher speeds). It is here that the model is shifted from the coupled solver to a segregated solver. In other words, the coupled solver sets up the flow field to a point where the segregated solver can pick up and converge at much faster speeds. This hybrid CFD solver model was found to effectively combine the unique advantages of both solvers into a single system to allow converged solutions for all geometries and boundary conditions in small time periods. The above process is effectively visualized when reviewing the convergence data plots as shown in Figure 2.3.1. It can be seen that the coupled solver initially engaged to run through a set number of iterations (each of which represents a time step for the flow 19 field) provides very slow solution convergence owing to the manner in which the building blocks have been set up. Once the model is shifted over to the segregated solver, the combination of its building blocks and the ?initial guess? given to it by the coupled solver allows it to converge almost instantly to very low values, and allows the completion of the case run for post processing in a much smaller number of iterations than would have been required had the coupled solver been used exclusively. Typical numbers for single case runs for the scramjet combustor geometry are as included in Table-2.3.1. Figure 2.3.1. Residual plots for a single case run using the Hybrid solver model 20 Table 2.3.1 CFD Solver computational times Solver Time to solution convergence (<10-6) for all flow parameters (single case run, single processor) Coupled solver only Convergence achieved in 20+ Hours Segregated solver only Solution Divergence Combined Coupled (Implicit) and Segregated Solvers Convergence achieved in 6+ Hours 21 2.4 BOUNDARY CONDITIONS All fluid entry surfaces were designated as pressure-far-fields with all in-flow boundary values fixed. All fluid exit planes were designated as pressure-outlets while all remaining faces were designated as wall boundaries. The no-slip condition was applied for all walls along with zero normal pressure and temperature gradients. The K-epsilon Turbulence model was utilized with enhanced wall treatment and viscous heating effects enabled. The density model for the air was defined as the ideal-gas model. All operating pressures were referenced against the standard sea-level atmospheric values for the center point in the entry far-field plane for the air. Finally, at the injector, all dependant variables are specified to be physical values of the temperature, pressure and velocity and the injector surface itself was treated as planar. No wall radiation effects were simulated. The mass flow rate of air was 0.20 kg/sec at nominal stagnation conditions. The nominal Mach-2 cross-flow entry stagnation conditions were 300K and 274kPa, and maintained at the entry plane (Y-Z). The injector stagnation conditions were 263kPa and 300K at Mach 1 with a uniform entry. 22 The injectant mass flow rate was maintained at 1.64 g/sec. These values conform to those used for the experiments against which the initial simulations were validated 7, 8. Mass flow rate tracking was enabled. FLUENT uses the wall functions as defined by Launder and Spalding which are quite commonly used throughout the industry. These are the default options with this software with the additional option of Enhanced Wall Treatment options. This latter option uses what are known as Enhanced Wall Functions. These functions extend the applicability of throughout the near wall region including the laminar sublayer, buffer region and fully turbulent outer region by formulating the law-of-the-wall as a single wall law for the entire wall region. As such it is possible to monitor the accuracy of the grid refinement process in FLUENT using one of two ways: manual or automatic adaptations. In both cases fixing a limit on the y+ values is extremely important. The Log-Law model is found to be effective beyond values of y+ = 30 and less than y+ = 60. Below this range of values the Enhanced Wall Treatment of the concerned region is necessary to accurately model the flow. The manual adaptation process involves plotting out the y+ values for a given grid during a grid refinement process and adjusting the refinement process until the values falls within the right range. This process is advantageous to use when single case runs are being conducted. However, with the proposed GA driven CFD case runs with changing geometric features, this process runs into limitations in that over-refining the grid might 23 be necessary to ensure that the wall region is being adequately modeled despite the changing geometries. The automatic grid adaptation process bypasses this problem by conducting automatic, localized grid adaptations in the wall regions using the user-defined input values for y+. This may involve coarsening the grid in some regions and refining the grid in others. Overall, the process ensures good wall accuracy for a given geometry at the cost of reduced grid uniformity. Given the networked structure of this optimization effort, it was clearly advantageous to use the self adaptation features of the CFD solver. However, it was noticed that the reduced grid uniformity could lead to unexpected solution divergence issues. This is a consequence of the hybrid solver setup as defined in the previous topic. Since the coupled solver was only defined to resolve a certain region of the flow with the given value of the courant number, which in turn is based on the grid refinement by its definition, any changes initiated by the self adaptation scheme led to changes in the requirement of the value of the courant number. Although the changes were usually small, it often led to a situation where the solution diverged for a given case rather unexpectedly because of the changes in the grid adaptation scheme. This posed a problem that was overcome by using the grid zones structure defined in Section 2.1. The zones were created with the issue of y+ values in mind so that only the zones affected by the movement of the injector locations were tagged for the self adaptation process and the value of the courant number was set to compensate for the 24 uneven changes in the grid. The Enhanced Wall Treatment option was not used. Given the nature of the grid with different zones of varying nodal densities, the y+ values also vary, but have been maintained to values of ~30. 25 2.5 EXPERIMENTAL VALIDATION Hartfield, Hollo, Fletcher and McDaniel 10 extensively investigated non-reacting staged transverse injection behind a rearward facing step. The simulation data was compared against a spatially complete data set of the various flow parameters obtained experimentally using optical techniques based on Laser-Induced-Iodine fluorescence and supersonic wind tunnel facilities 8. The validation efforts were divided into two parts: a) Flow without injection b) Flow with injection The above was necessary since the small diameter injection holes on the floor of the combustor meant that the shock interactions and flow perturbations as a result of the injection process were unlikely to affect the nature of the flow along the entire width of the combustor except in the immediate region near the injection surfaces and downstream of the rearward facing step. As a result, regions remaining unaffected by the injection process would continue to exhibit the nature of flow as if there were no injection whereas, areas near the injection would exhibit complex and highly variable three dimensional flow features. 26 2.6 FLOW WITHOUT INJECTION The supersonic cross-flow over a rearward facing step meant that the focus in this section was on developing the accuracy of the expansion fan, the shear layer and the oblique shock downstream of the step. The experimental flow field can be characterized by using the flow diagram shown in Fig.2.2.1. The typical flow from upstream of the rearward step is roughly maintained except for the boundary conditions and the flow expands over the step leading to the formation of the classic expansion waves. The presence of the step causes the flow to behave as if it is a mixing case of two parallel flow fields having different velocities. As a result, with the flow above the step moving at Mach 2, the flow just behind the step is characterized as a low speed recirculation. The formation of a shear layer having an axial length of multiple step heights takes place. Figure 2.6.1 compares the experimental results for the static pressure with those from the simulations. As can be seen, a close compatibility is maintained between the theoretical calculations, the experimental results and the simulation results for all flow features. 27 Figure 2.6.1. Experimental 7 and simulation comparisons for the pressure field for the non-injection case The above plots are for the X-Y plane at the centerline axis of the combustor. Nevertheless, for the non-injection cases, the symmetry of the flow is maintained in all subsequent X-Y planes except for those near the wall edges where the boundary wall effects cause some perturbations and hence an a slight asymmetry in the X-Y flow fields. In order to be able to compare the results more clearly, a comparison of the experimental and simulation results for the pressure profile has been shown in Figure 2.6.2. A close compatibility of the simulation and experimental results is clearly seen. 28 Figure 2.6.2. Comparison of the experimental 7 and simulated static pressure profile along the line from the combustor floor to the top of the combustor at the center of the primary injector One prominent effect noted in the comparison of the experimental and simulated results was the inability of the CFD solvers to accurately depict the expansion fan region in the vicinity of the upper edge of the rearward step. Experimentally, the very low pressure region just downstream of the rearward step is seen to provide a suction effect that causes the expansion effects to occur upstream of the physical step location 8. However, repeated CFD efforts failed to reproduce this effect to the extent as that seen in the experimental data even after the use of localized ultra-high resolution grid domains near the step location. This is found to concur with other independent efforts using the SPARK three dimensional Navier-Stokes codes wherein the same problem was encountered and remained unresolved 13. 0 0.2 0.4 0.6 0.8 1 1.2 Sta tic Pr es su re (P/ Pin f) Experimental Simulation 29 Figure 2.6.3. Simulation results for the pressure contours in a 3-Dimensional view in the FLUENT display screen Similarly, the static temperature results are shown in Figure 2.6.4. Here the results compare equally well as far as the flow field was concerned. However it was noticed that the temperatures obtained experimentally were somewhat higher than those obtained via calculations in the recirculation region behind the shear layer. Figure 2.6.5 displays the wall temperature distribution for the entire combustor section as a result of the flow and effectively shows the high temperature region on the top of the combustor section where the stagnation conditions are nearly reached. 30 Figure 2.6.4 Experimental 7 and simulation comparisons for the temperature field for the non-injection case Figure 2.6.5. Three dimensional visualization of the temperature field in the FLUENT display mode for the test section 31 Figure 2.6.6. Comparison of the experimental 7 and simulated static temperature profile along the line from the combustor floor to the top of the combustor at the center of the primary injector One of the main points of focus for the validation was the accuracy with which the shear layer was modeled immediately downstream of the rearward facing step. This was the driving factor behind the grid development. The details of this region are shown in Figure 2.6.7. The close compatibility between the results obtained through the simulation and the experimental results is clearly visible in these two images. 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 0 5 10 15 20 Sta tic Te mp era tu re (T/ Tin f) Experimental Simulation 32 Figure 2.6.7. Experimental 7 and simulation comparisons for the velocity field for the non-injection case The recirculation region is reproduced accurately from a geometrical standpoint. The length of the shear layer along the axial (X) direction is seen to conform to experimental data as does the near planar nature of the shear plane itself. 33 2.7 FLOW WITH INJECTION The flow with injection must be evaluated at several sections for different theoretical explanations for the flow features. Nevertheless, some distinct features remain constant throughout. The main features of this flow are the formation of detached bow shocks upstream of the injection points. These are a result of the injection stream acting as a blunt body being moved through a supersonic cross-flow as shown in Figure 2.7.1 and the under expanded jet cores of the two injectors. The shock is detached in all three dimensions and this can be visualized in Figure 2.7.2. The shock essentially moves around the disturbance created by the injection. The flow upstream of the first injection point is uniform causing the detached shock to be very clearly defined. The same is not true of the second injection point, which faces severe disturbances in the upstream flow as a result of the first injection point. This can be seen in Figure 2.7.3. 34 Figure 2.7.1 Flow features at the centerline section for the injection case 7 Figure 2.7.2. Visualization for the pressure field in the injection case to highlight the 3-Dimensional nature of the shock patterns 35 The second injector is seen to penetrate further into the free-stream than the first due to the disturbed nature of the free-stream downstream of the first injector. This entails lower Mach numbers and higher pressure fluctuations. More shear layers are seen to develop in the region downstream of the second injector. Streamwise vortex formations are seen in the Y-Z planes. The nature of the flow of the injectant along the length of the combustor downstream of the injection locations is shown in Figure 2.7.3. The injectant plane is seen to be in the shape of long cylindrical sections with low dispersion into the cross flow even at substantial distances downstream of the injection locations. Notice the long, narrow, cylindrical shape of the mass fraction distribution of the fuel. The divergence is achieved mainly in the immediate downstream region of the injection locations. Figure 2.7.3 Nature of the injectant flow through the combustor. 36 Further, the expansion fan faces severe disturbances as a result of this injection as well. The detached shock is seen to effectively destroy the lower expansion shocks (near the bottom of the fan) upstream of the injection location while the overall angle of the expansion fan is reduced as a result of reduced pressure gradients upstream and downstream of it. The experimental and simulated results compare favorably in all these aspects, as shown in Figure 2.7.4 for the pressure contours. Figure 2.7.4 Experimental 7 and simulation comparisons for the pressure field for the injection case 37 Figure 2.7.5 Comparison of the experimental 7 and simulated static pressure profile along the line from the combustor floor to the top of the combustor at the center of the primary injector Figure 2.7.6 Comparison of the experimental 7 and simulated static pressure profile along the line from the combustor floor to the top of the combustor at the center of the secondary injector 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 5 10 15 20 Sta tic Pr es su re (P/ Pin f) Experimental Simulation 0 0.2 0.4 0.6 0.8 1 1.2 0 5 10 15 20 Sta tic Pr es su re (P/ Pin f) Experimental Simulation 38 Figure 2.7.4 also illustrates the formation of the barrel shocks just above the injection surfaces. This shock pattern is seen to be deflecting in the downstream direction. A better realization of this is achieved by considering the temperature plots shown in Figure 2.7.7. Figure 2.7.7 Experimental 7 and simulation comparisons for the temperature field for the injection case 39 Figure 2.7.8 and Figure 2.7.9 shows the comparison between the experimental and simulated results in the form of line profiles for the primary and secondary injector centers. Note that these lines are perpendicular to the floor of the combustor and move vertically upwards towards the roof of the combustor. Along this line various data points have been extracted for comparison. Figure 2.7.10 shows the complex nature of the three dimensional temperature distributions in the combustor cross section and the much higher temperature gradients obtained in the injection case as compared to the non-injection case. Figure 2.7.8 Comparison of the experimental 7 and simulated static temperature profile along the line from the combustor floor to the top of the combustor at the center of the primary injector 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 0 5 10 15 20 Sta tic Te mp era tu re (T/ inf ) Experimental Simulation 40 Figure 2.7.9 Comparison of the experimental 7 and simulated static temperature profile along the line from the combustor floor to the top of the combustor at the center of the secondary injector Figure 2.7.10 Three dimensional visualization of the temperature field in the FLUENT display screen The velocity vector plots shown in Figure 2.7.11 also present a favorable comparison between the experimental and simulation results. Some unique flow features are revealed in the velocity vector plots. As the streamlines show, the flow immediately 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 0 5 10 15 20 Sta tic Te mp era tu re (T/ Tin f) Experimental Simulation 41 above the injectors is seen to be compressing in the area upstream of the injection point while in the downstream area, the flow is seen to be expanding. This is also consistent with the deflection of the barrel shock structures discussed previously. Note that the streamlines shown in the figures are shown in a two- dimensional plane, and does not fully represent the three dimensional velocity fields. Figure 2.7.11 Experimental 7 and simulation comparisons for the velocity field for the injection case 42 2.8 MASS FRACTION ANALYSIS Figure 2.8.1 compares the experimental and simulation results for the mole fraction distribution of the injectant in the X-Y plane (single diameter transverse offset). It should be noted that the simulation studies were compared with experimental data that dealt with injecting air into air. In the case of the CFD simulations, the injectant air was ?tagged? and tracked in the flow-field as if it were a different fluid altogether. As will be described in later sections of this thesis, the fuel-air mixing phenomenon was suitably modified as well. Note also that since air is being injected into air, the mole and mass fraction distributions are the same. Good accuracy is seen to be maintained for the mole fraction distributions. Figure 2.8.2 presents a three dimensional plot of the mass fraction distribution for the injectant along the length of the combustor and the evolution of the distribution field at the exit of the combustor can be clearly visualized. 43 Figure 2.8.1 Comparison of the experimental 7 and the simulation results for the mole fraction distribution of the injectant along the centerline axis Figure 2.8.2. Simulation results for the Mass Fraction distribution for various Y-Z planes along the combustor length. 44 2.9 THE GENETIC ALGORITHM AND CFD INTEGRATION A Genetic Algorithm is an evolutionary biologically inspired optimization technique which uses various types of evolutionary mechanisms such as mutation, selection, inheritance etc, to move toward an optimized solution for a given function or mathematical process 1-6. For this effort, the IMPROVE V2.8 binary GA based on the Fortran 77 software base was used. A more detailed description of the operational procedure is discussed in the results section given the direct connection between the results obtained and the GA setup. In this section a brief overview of a typical run is presented for creating a fundamental base for the more detailed discussions that follow. IMPROVE is the acronym for Implicit Multi-Objective PaRameter Optimization Via Evolution 16. Here the term ?evolution? is a misnomer given that the GA is more adaptive than evolutionary in the absolute biological sense. The GA adapts existing genes to give better performance but does not create new genes and is thereby restricted to a given gene pool as defined by the design space. For complex problems with non- definable objective functions, GA?s have been found to be much better at estimating the global maximum or minimum than gradient based methods that are prone to finding local maxima or minima. 45 A binary GA operates on the binary encoded input parameters. In this sense the input must include a range of values for each parameter that defines the design space as well as a field resolution value for the said parameter. Since the length of a chromosome in the binary GA is defined by the number of bits for the given design space (sum of the number of bits for each design parameter), the resolution value, the parameter design space and the number of parameters define the size of the optimization problem as shown in equation 2.9.1. g1840g1873g1865g1854g1857g1870 g1867g1858 g1854g1861g1872g1871 = uni2211 [g2922g2924g4672 g3171g3159g3182g4666g3284g4667g3127 g3171g3167g3172g4666g3284g4667 g3293g3280g3294g3290g3287g3296g3295g3284g3290g3289 g4673 g2922g2924g4666g2870g4667 ] g3015g3048g3040g3029g3032g3045 g3042g3033 g3017g3028g3045g3028g3040g3032g3047g3032g3045g3046 g3036g2880g2869 (2.9.1) A typical run with this version of the GA therefore involves conversion of the input parameters into binary forms and the generation of an initial population of members whose size is fixed by the user before the activation of the run. This initial ?generation? is composed of a random selection of parameters from the design space. The design values are then passed to a code to evaluate the performance of a given engineering system (in our case the initiation and evaluation of CFD data) and allows the GA to stack the members of the population in order of their individual performances. After the performance evaluation, the adaptation process begins with a tournament selection process in which the members of the population are chosen for reproduction. This is followed by processes such as probability-based crossover and mutation to develop the next generation of members for evaluation. Since these members have been constructed based on the results of the previous generation best performers, the 46 best performers of the new generation are statistically found to be an improvement over those of the previous generation. This effort was a multi-goal analysis which required a critical review of the various GA fitness stacking selection schemes such as strict, relaxed and apportioned Pareto methods. The use of these schemes evolved with the analysis of each run and their description fits better within the results section. However, there are other crucial GA schemes like micro GA and steady state GA versus conventional GA selection, creep selection, uniform versus non-uniform crossover etc. that defined the final results but unlike the goal evaluation schemes, were common throughout the various runs and can be described here in detail. The steady-state GA selection is referred to as a population elitist strategy in that a steady state GA keeps only the best performers from the current generation and the previous generation. As such, it can lead to a loss of genetic diversity. It also means that this process can be used with a higher mutation rate. The mutation rate defines the rate at which mutation occurs within the chromosome while crossover probability is the probability of crossover when two parents come together after the end of the tournament selection process. Higher values of these schemes along with the steady-state GA can lead to quicker solution convergence. Elitism was another scheme that required usage for cases when the total pressure ratio was being optimized (discussed later). It is used to preserve the best performer for each goal. Creep mutation is used to ?fine tune? the GA once significant hyperplanes have 47 been found by randomly varying the least significant bit in a random parameter. Uniform Crossover was kept inactive to allow for random matchup of parent genes within the children members of the next population. As a result of the above procedure, the following requirements must be met when attempting to use the GA in coordination with CFD solvers for the optimization of the Scramjet Combustor Geometry: a) Selecting the variable design parameters b) Defining the design space c) Selecting the optimization goals d) Automating the grid generation process e) Automating the setup of the grid and the case parameters within the CFD software f) Integrating the CFD software with the GA code g) Modification of the GA to allow parallel processing of numerous members of a generation over multiple computer processors simultaneously to reduce calculation time h) Modification of the CFD solver setup to allow parallel processing of each given case over multiple computer processors for each member of the generation. i) Transfer of the above structure to a networked system of computers 48 2.10 THE DESIGN PARAMETERS AND DESIGN SPACE The overall geometry and boundary conditions of the scramjet combustor for this analysis remains fixed throughout the course of the optimization study. The physical aspect of the combustor varied in the optimization is the location of the two injectors on the combustor floor (X-Z plane). This can then be subdivided into four variables as discussed below: a) The axial (X) location of the first injector b) The transverse (Z) location of the first injector c) The axial (X) location of the second injector d) The transverse (Z) location of the second injector Having fixed the parameters which are varied during the course of the optimization, the range and resolution within which the GA can select the values must be specified based on the concept mentioned in the previous topic. It should be noted here that this design space is not fixed as further runs were conducted based on the evaluation of trends and analysis of results from the initial results. 49 As such, the design space evolved or was restricted depending on the requirement of further runs as will be discussed in the results section. Therefore, the values mentioned in this section correspond only to the initial runs and are designed to give an idea of the typical design space and values that will be used throughout the rest of this analysis. The range of values for each of the above parameters was fixed to be within one step height of its original experiment values and to have a resolution of selection of 0.1. The step height is 3.18 mm, and therefore the value for each of the parameters can change from +3.18 to -3.18 mm beyond the original position. The overall design space then is the area swept by these limiting values and in this case takes the form of a square for each of the two injector locations as shown in Figure 2.10.1 50 Figure 2.10.1. The Design Space for the optimization of the Scramjet Combustor 51 2.11 DESIGN GOALS The primary performance parameter of interest in this effort are the total pressure loss and the mass fraction of the injectant, with both being defined for a cross sectional (Y-Z) plane of the combustor. The total pressure loss represents overall propulsion performance degradation while the mass fraction plots are considered as a means of quantifying the mixing phenomena 17. Generally a tradeoff exists between the total pressure decrease and the degree to which injectant is mixed into the free stream. The optimization includes the total pressure loss versus what is referred to as the mixing efficiency 18. The term ?mixing efficiency?, as defined by Anderson 18, is an empirical, one dimensional measure of the degree of mixing completeness which takes into account both the near-field mixing (initial macro-mixing phase) and far-field mixing (molecular diffusion/micro-mixing). For an overall fuel-lean mixing location, it is ?the amount of fuel that would react if complete reaction occurred without further mixing divided by the amount of fuel that would react if the mixing were uniform? 18. 52 In order to be able to define the mixing efficiency at any given Y-Z plane downstream of the second injector location, the mixing phenomena must be quantified against a reference standard. In this context, the mixing efficiency is defined as the ratio of the number of cells on the required Y-Z grid surface which possess a fuel-air ratio within a given range of values where the combustion process can occur to the total number of cells on the entire surface. In other words, beyond the standard stoichiometric fuel-air ratio, there exists a fuel-rich value and a fuel-lean value which represent the practical range of values for proper combustion. Any value of the fuel-air ratio above or below this range will have insufficient reactants to combust properly and thus is disqualified from the selection procedure. Any value in between these two is selected as having desirable combustion conditions as far as mixing of the reactants is concerned. The total grid nodes are fixed to be the same for a given Y-Z plane at the end of the combustor where all calculations for the exit conditions are conducted. This Y-Z plane exists at an axial distance of ten step heights beyond the second injector location in the original setup and thereafter remains fixed at that distance regardless of the change in the second injector location during the optimization process. This plane is therefore located at a distance of 65.08 mm downstream of the inlet of the combustor for all calculations. The number of cells for this exit plane of the combustor is also fixed at 1890. 53 The flammability limits for the fuel-air ratio above and below the stoichiometric ratio must be defined in order to define the selection criteria for any given cell in the exit plane. However, so far all validation studies discussed here have involved injection of air into air. Therefore there is no fuel involved. To move beyond this condition would involve changing the injectant from air to some fuel. In order to provide some use for the validation studies involving gaseous injectant discussed in previous sections, hydrogen was the obvious choice as the fuel. It was further determined that the flammability limits in terms of the fuel-air ratio of hydrogen to air per mass basis varies from 0.0520 (fuel rich) to 0.002774 (fuel lean). The stoichiometric value is 0.0294. A validation study of hydrogen injection into air was initially bypassed and the simulations were conducted using the previous limits for hydrogen while keeping the injectant as air. In this way the flow field retains the nature discussed in the validation study in previous sections and yet provides some realistic values for the calculation of the mixing efficiency. The mixing efficiency at this point is therefore defined as the ratio of the number of cells in the exit plane having a fuel-air ratio between the flammability limits of hydrogen mentioned above to the total number of cells for the entire exit plane. The total pressure ratio in three dimensional flow fields must also be quantified. The concept of single dimension analysis employed here treats the total pressure ratio as the ratio of the total pressure at the exit of a component to that at the entry of the component. It logically follows that in order to do the same with regard to entry and exit planes in a three dimensional flow case, one must determine the averaged total pressure for both the entry and the exit planes of the given scramjet combustor. At the same time, 54 the validity of the governing equations of the flow must be maintained. These requirements are met with the development of the averaging equations using the mass- momentum averaging technique which is a modification of the Stewart mixing analysis23. The equations used in this analysis have been customized to handle rectangular geometries and meshed surfaces. This derivation for the mass-momentum averaging technique for CFD grids is explained in detail in Appendix A. The resulting equations needed for quantifying the total pressure ratio for the scramjet combustor are given in equations 2.11.1 to 2.11.4. The average total pressure across a grid plane in a flow field is found to be: g1842? = g4678g3330g3267g3269?g3276 g3267g3254g3255g4594 g4679uni2211 g3017g3294g4666g3036g4667g3436g3262g3299g4666g3284g4667g3118 g3297g3299g4666g3284g4667 g3440g3002g4666g3036g4667g3263g3284g3128g3117 g3344g3294 g3344? g3297g3299 g3276g3267g3254g3255g4594 g3002 (2.11.1) Where all parameters are defined in the nomenclature. As can be seen, except for the denominator terms, the entire equation is composed of terms that can be determined from the CFD grid data. As for the denominator terms, three additional equations are required to quantify them in terms of the grid data as well. These equations are: g1855 = g3437 g3211RTg4594 g3276g3267g3254g3255? g3118 g3441uni2211 g3427g3017g3294g4666g3036g4667g3014g3299g3118g4666g3036g4667g3002g4666g3036g4667g3431g3263g3284g3128g3117 g2878 g4672g3211g3126g3117g3118g3211 g4673uni2211 [g3017g3294g4666g3036g4667g3002g4666g3036g4667]g3263g3284g3128g3117 g4678g3330g3267g3269?g3276 g3267g3254g3255g4594 g4679uni2211 g3017g3294g4666g3036g4667g3436g3262g3299g4666g3284g4667g3118 g3297g3299g4666g3284g4667 g3440g3002g4666g3036g4667g3263g3284g3128g3117 (2.11.2) 55 g3096g3294 g3096? = g3428g883 g3398 g4672 g3082g2879g2869 g3082g2878g2869g4673g484 g3049g3299g3118 g3028g3267g3254g3255? g3118 g3432 g4672 g3117g3330g3127g3117g4673 (2.11.3) g3049g3299 g3028g3267g3254g3255g4594 = g4672 g3082g3030 g3082g2878g2869g4673g3398 g3495g4666 g3082g3030 g3082g2878g2869g4667 g2870 g3398 g883 (2.11.4) These equations then allow for the accurate calculation of the total pressure ratio over the inlet and outlet planes (Y-Z) of the scramjet combustor and then take the ratio of the total pressure at the outlet to that at the inlet as a goal for the optimization. With the mixing efficiency and the total pressure ratio now defined based on the flow parameters at the inlet and the outlet, the optimization goals for the GA can be evaluated. 56 2.12 THE NETWORKED ENVIRONMENT Automation is an important requirement for optimization studies. Another requirement was that the entire program code be initiated on as many processors as possible to allow for rapid case runs. This meant that the CFD solver and the GA must both run in parallel. The computer network available at Auburn University consisted of a total of sixty available processors for CFD runs on thirty computer nodes. All processes are initiated on a central node referred to as the head node. In addition, a total of four CFD solvers can be run simultaneously on the above cluster. This means that there were a total of ten processors available for each CFD run for use as parallel processors. All CFD simulations are run through the FLUENT software, whereas all grid generation requirements are met by the use of GAMBIT software with either software being remotely initiated and controlled in a non GUI environment by the GA. Figure 2.12.1 shows the overview of how GAMBIT as a grid generator is integrated into the code to allow for automated grid generations and also its interaction with FLUENT. GAMBIT can be initiated in a batch mode with the creation of a journal file from a FORTRAN subroutine. This is done only on the head-node of the cluster to prevent increase in computational time on the nodes for the CFD runs. Once initiated in a non-GUI mode on the head-node, GAMBIT creates a mesh file using the set of 57 commands written in the journal file and this is combined with a similar data file for FLUENT. This allows a mesh file and a command data file to be passed to the nodes where the CFD cases are being run. Figure 2.12.1. An overview of the network structure involving GAMBIT and FLUENT on the Computer cluster at Auburn University Once the CFD run for a given case is completed, case and data files are automatically written based on the commands in the data file (.dat) and passed back to the head-node. Here they are reopened in a non-parallel mode by the Fortran program and an ASCII file is written containing the required flow field parameters. This is necessarily done on the head-node since FLUENT cannot execute the formation of ASCII files (which are exported) while in parallel mode on the cluster nodes. 58 Once the formation of the ASCII file is complete, it is opened by the FORTRAN subroutine designed to process the flow field parameters for the quantification of the required goal parameters and used by the Genetic Algorithm to proceed forward. Figure 2.12.2 shows the nature of the interaction of the head-node and the cluster nodes. In effect, the mesh files generated for each member of a generation are passed on from the head node to each node on the cluster as per the requirement. Since the amount of time required for the case runs is usually large, the time required for the generation and transfer of the mesh files sequentially and the transfer of the calculated data back to the head-node does not significantly affect the overall computational time adversely. Figure 2.12.2. A generalized concept of operations view for a given generation in a GA run 59 An important issue is the computational time versus number of processors being used. The maximum number of processors available for each CFD run on the network is ten. With a single processor, the time required for each case run is several hours, but as the number of processors is increased, the time for overall calculation decreases dramatically from around six hours for a single processor to around two hours for four processors or about one and a half hours for six processors. Figure 2.12.3 shows the nature of the curve for calculation time versus the number of processors. It was noted that after a maximum of around eight processors, the curve levels off, and further decrease in calculation time by adding even more processors is offset by the increased I/O time between the large numbers of processors. Figure 2.12.3. Computational time for a single FLUENT case run versus number of parallel computer processors used 0 1 2 3 4 5 6 7 0 5 10 15 Ho ur s Number of Processors Calculation time 60 Based on the above results, the structure for this optimization analysis conducted four simultaneous CFD case runs (i.e. four members of a generation in parallel mode) with a parallel processor structure involving a total of eight processors leading to the overall use of 32 processors at any given time. The calculation time for each case run for such a structure was found to be around one hour on average. This meant that an entire generation of 20 members of a GA run with full CFD results could be completed in as little as five to six hours allowing for the completion of around four generations per day and the completion of a single optimization run to be completed in less than a week. 61 3 THE RESULTS 3.1 CASE RUN - 1 The first optimization run was conducted using the following parameters: a) Optimization goals: maximize total pressure as primary goal and maximize mixing efficiency as secondary goal using the strict pareto methodology. Both parameters are calculated from the CFD data as discussed above. b) Design space: movement of the two injection surfaces in the X-Z direction on the floor of the combustor within one step height of their original experimental location with a resolution of 0.1mm. c) Genetic Algorithm parameters: a) A 20 generation run b) Each generation has 20 members c) Two goal optimization run d) Total pressure ratio has priority over mixing efficiency as optimization goals 62 A brief discussion of the pareto methodology is presented here to explain the use of the strict pareto scheme used for this run. The Pareto scheme, in general, differs from the conventional GA setup when multiple goals (>1) are to be simultaneously optimized. The goals are optimized individually based on the input parameters common to all goals but the fitness selection differs in the way goals are classed above or below one another. This is known as ?goal domination? and defines the selection process in that one goal set, which is a collection of individual goal performances for the said input, must dominate the others if the said parameters are to survive. In this way all goal sets are compared to determine which parameter set is chosen over the other sets. A ?strict pareto? is one that conducts this domination statistics analysis but with a relative measure of significance between the goals themselves in that the total pressure ratio was stacked as higher priority than the mixing efficiency. By default, the GA stacks the results for all members of a given generation in the order of increasing performance (looking top to bottom), with the performance of the priority goal taking precedence over the secondary goal. As a result, for each generation there exists a best performer and a worst performer in its list of members. All results presented in this paper include plots of the best and worst performers of a given generation versus the optimization goals or design parameters over the entire number of generations available. This allows one to follow the trend that the optimization efforts tend to take and provide valuable insights into the inner workings of the optimization process. 63 Figure 3.1.1 shows the results obtained for the total pressure ratio over the course of the GA run. It is seen that the best performer of each generation continually improves over that of the previous generation until a period of stabilization is achieved after which the best performer does not practically improve further. Another feature visible from Figure 3.1.1 is that the worst performer undergoes a lot of fluctuations in its values over the course of the optimization process as each worst performer is replaced with a new result based on the previous generation as the GA tries to improve performance with each generation. Occasionally a result of such manipulations is that a member having better performance than the best performer of the previous generation, in which case the latter gets replaced by the new best performer and the GA then builds up on it. Figure 3.1.1 The total pressure ratio versus the generational best and worst performers 0.775 0.78 0.785 0.79 0.795 0.8 0.805 0.81 0 5 10 15 20 To tal Pr es su re Ra tio Generations Best Performers Worst Performers 64 Figure 3.1.2 shows the results obtained for the mixing efficiency over the course of the GA run. It is clearly visible that the mixing efficiency has been treated as a secondary goal by the GA over the pursuit of the total pressure ratio optimization. In other words, when referring to the results for the mixing efficiency, it is clearly seen that the absolute best performer in terms of just the mixing efficiency has been classed as the worst performer for a given generation because of the total pressure loss. Clearly, the best performer in terms of the overall analysis objective is not necessarily the best performer in terms of the mixing efficiency alone, and vice versa. Figure 3.1.2. The mixing efficiency versus the generational best and worst performers 0.175 0.18 0.185 0.19 0.195 0.2 0.205 0.21 0.215 0 5 10 15 20 Mi xin g E ffic ien cy Generations Worst performers Best performers 65 It is seen that the total pressure ratio (as well as mixing efficiency) can be compared to a number of reference points within the system. The reference points are the best and worst performers of the first generation. That is, the best performer of the last generation can be compared to the best and worst performers of the first generation to compare the overall optimization effort. However, in this analysis, we have an experimental reference point for the two injection locations and they will serve as the ?global? reference points for the entire analysis. The value of the total pressure ratio obtained for the experiment based setup was found to be 0.8043 (or 80.43%) based on the equations set up previously and the mixing efficiency was found to be 0.1852 (or 18.52%). As such then one can calculate the improvement in both the total pressure ratio and the mixing efficiency against the above mentioned values. For the total pressure ratio, this value comes out to be 0.11% (80.52% from 80.43%) and for the mixing efficiency there is an associated decrease in value from 18.52% to 18.50% for a total loss of 0.10%. the improvement of total pressure ratio can be justified against this loss in fuel- air mixing. However, this is the first case run, and as we shall see, the results obtained with improvements in our understanding of what the system trends are will lead to improved results. The plots shown in Figure 3.1.3 and Figure 3.1.4 are similar in nature to that presented previously in that all results are shown versus the generation number on the X- axis and design range for each parameter in the Y-axis of the plot. Both the best and the worst performers for each generation are presented. 66 Figure 3.1.3 and Figure 3.1.4 present the results for the combustor axial (X) variation of the two injector locations. It is immediately seen that the GA moves towards the inner limits of the design space for these two parameters beyond generation number six. Figure 3.1.3 The primary injector axial (X) variation for the generational best and worst performers 2 2.2 2.4 2.6 2.8 3 3.2 3.4 3.6 3.8 4 0 5 10 15 20Prim ary In jec to r X Va ria tio n ( Ste p H eig ht s) Generations Best Performers Worst Performers 67 Figure 3.1.4 The secondary injector axial (X) variation for the generational best and worst performers Figure 3.1.5 The primary and secondary injector variation for the best performers 6 6.2 6.4 6.6 6.8 7 7.2 7.4 7.6 7.8 8 0 5 10 15 20Sec on da ry Inj ec to r X va ria tio n ( Ste p h eig ht s) Generations Best Performers Worst Performers 68 Figure 3.1.5 shows the physical presence of both of the injectors and their axial (X) locations within the design space. It explains more clearly the idea that the design space is seen to be limiting the results obtained as far as the axial variation of the injectors is concerned. This is related to the type of definition used for the mixing efficiency and the physical nature of the injectant flow. From Figure 2.5.10 it can be seen that the nature of the injectant flow after entering the combustion chamber is in the form of a thin long frustum cone. Any section in the Y-Z plane of the combustor taken from this conical flow structure has a different fuel-air dispersal contour set. In order for the maximum dispersal to take place, the conical structure must be pulled as far upstream within the combustor as possible as shown in Figure 3.1.6. This is exactly what the GA tries to do. This suggests that either the design space should increase along the axial direction of the combustor or a different definition of the mixing efficiency which is independent of the length of the combustor is needed. 69 Figure 3.1.6 optimizing the injector Axial (X) location Figure 3.1.7 and Figure 3.1.8 display the results corresponding to the best and worst performers for the variation in the combustor transverse (Z) direction versus the generation numbers. Here the results are more promising in that optimized results are obtained well within the current design space. One very early conclusion suggested is that the design space for transverse variation of the injector locations does not need to be increased beyond the current range for any future runs with the current setup. However, using different fuels or inlet conditions will render this assumption invalid. 70 Figure 3.1.7. The primary injector transverse (Z) variation for the generational best and worst performers Figure 3.1.8. The secondary injector transverse (Z) variation for the generational best and worst performers -3 -2 -1 0 1 2 3 0 5 10 15 20 Pri ma ry Inj ec to r Z va ria tio n ( Ab ou t Ce nt erl ine , m m) Generations Best Performers Worst Performers -3 -2 -1 0 1 2 3 0 5 10 15 20 Se co nd ary In jec to r Z va ria tio n ( Ab ou t ce nt erl ine , m m) Generations Best Performers Worst Performers 71 For both plots of injector transverse variations, it is seen that the location for the injectors in the combustor ?Z? direction is almost constant after the fifth generation. However, the worst performers show more disturbances even after the fifth generation and it is these disturbances that provide some very useful insights into the optimization process. From Figure 3.1.8 it can be seen that the secondary injector has its transverse position nearly fixed for both the best and worst performers after the fifth generation with very few perturbations between the best and worst performers. At the same time, however, the primary injector displays large fluctuations for the worst performers even after the fifth generation. Here an interesting trend that is noticed is that all fluctuations for the first injector seem to take place closer towards the centerline of the combustor floor and almost entirely fail to cross the best performer position. That is, once the GA has started to optimize the mixing efficiency after having finished the optimization of the total pressure ratio to the best available value, there is a trend towards moving the two injectors closer towards each other in the transverse direction. After the fifth generation, even the axial variations of the two injectors were limited to their best possible locations within the design space as discussed previously. As a result, with the second injector remaining relatively fixed in the transverse direction as well, the only parameter varying noticeably after the fifth generation is the ?Z? variation of the first injector. 72 From Figure 2.5.10 and in light of the above discussions, some interesting results are observed. There appears to be a direct relationship between the mixing efficiency for the worst performer and its corresponding transverse location around the centerline of the combustor. Even more interesting is that the best result for mixing efficiency is obtained by bringing the two injectors closer to each other along the transverse direction rather than spreading them out. This has a lot to do with the type of fuel used and its flammability limits given the current definition of the mixing efficiency. Figure 3.1.9(a), Figure 3.1.9(b) and Figure 3.1.9(c) represent this more visually with a simple example of concentric mass fraction contours. As is seen from these figures, the nature of the movement clearly affects the mixing efficiency directly, and is therefore dependent on the type of fuel and the values used. While it is clear from the above discussion that, given the current type of fuel and the corresponding flammability limits, the mixing efficiency is seen to improve as the two injectors are brought closer. The overlap of the mass fraction contours clearly presents an advantage initially to try and increase the mixing efficiency. It is also clear that bringing them too close to each other creates a reverse effect as seen in Figure 3.1.9(c). Once the two injectors are brought nearly in line with each other, with little or no transverse distance between them, the initial effect of increasing mixing efficiency is lost. 73 Figure 3.1.9(a) (left) and 3.1.9(b) (right). Effect of bringing the injectors closer to each other. Figure 3.1.9(c). Effect of bringing the injectors closer to each other. 74 This result is also seen from Figure 3.1.7 and Figure 3.1.8 at generation 18. At this point, the two injectors are closer to each other than at any other point. It shows both a drop in the total pressure ratio and a drop in the mixing efficiency. This suggests that there exists a certain optimum location between these two transverse location values depending on the values put forward for a given fuel?s flammability limits. 75 3.2 CASE RUN - 2 One of the interesting results obtained from the first case run was the movement of the injector locations in the longitudinal direction. It was clearly seen that the GA had been restricted in its movement of the two injector locations in its effort to improve the mixing efficiency without losing out on total pressure ratio. Further, by considering the worst performers for that run, it can also be seen that the movement of the two injector locations closer to each other might yield a better total pressure ratio and mixing efficiency. However, it should be noted that moving the injectors too close to each other would prove detrimental to total pressure ratio but prove highly beneficial to mixing efficiency thanks to the turbulence created within the combustor due to close Mach 1 injections. In this sense there was clearly a need to further explore this longitudinal variation in more detail and a more aggressive design space. This therefore formed the basis for Case Run 2. 76 Figure 3.2.1 The Single Dimension Longitudinal Design Space for Case Run 2 Several changes were made to design space based on the trends observed from Case-1. The transverse variations were completely removed and gave way to a single dimension design space along the centerline of the combustor floor as shown in Figure 3.2.1. Further, the design space was extended for the primary injector to bring them closer to the rearward facing step in accordance with previously observed data. The new design space therefore allowed a variation in step heights from 1.0 to 4.0 for the primary injector and 6.0 to 8.0 for the secondary injector. The strict pareto scheme was removed and the number of goals reduced to one by combining the two goals (both of same 77 magnitude) with suitable domination factors. For this case run, the total pressure ratio was again kept dominant over the mixing efficiency. The results obtained are shown in Figure 3.2.2 and Figure 3.2.3 for the total pressure ratio and the mixing efficiency respectively. The GA performs as expected and the plots are typical of GA optimization results. The total pressure ratio has clearly been improved over the course of the optimization and so has mixing efficiency, though the latter is not close to its full potential independently as observed from the worst performer data in Figure 3.2.3. Given that the transverse variation in the injector movements have been eliminated for this run, it is easy to compare the results obtained from the CFD data with those obtained from the experimental data. Namely, the total pressure ratio has been seen to improve by 0.11% (80.52% from 80.43%) while the mixing efficiency has improved by 10.69% (20.5% from 18.52%). These numbers are encouraging but not especially impressive. The trends are clear in that the mixing efficiency has a potential to improve dramatically given the chance to be dominant in the optimization even without the need for the transverse variation improvement discussed in Case-1 results. That is encouraging enough to be the core of the setup for the next case run, i.e. Case-3. 78 Figure 3.2.2 The Total Pressure Ratio for best and worst performers for Case Run 2 versus GA Generations Figure 3.2.3 The Fuel-Air Mixing Efficiency for best and worst performers for Case Run 2 versus GA Generations 0.8025 0.803 0.8035 0.804 0.8045 0.805 0.8055 0 5 10 15 20 To tal Pr es su re Ra tio Generations Best Performer Worst Performer 0.17 0.18 0.19 0.2 0.21 0.22 0.23 0.24 0.25 0 5 10 15 20 Mi xin g E ffic ien cy Generations Best Performer Worst Performer 79 Further analysis of the results of Case-2 helps validate the suspected trends as noticed during Case-1. Figure 3.2.4 shows the details of injector movement for the best performer data along the combustor floor. Without any transverse movement, this analysis is relatively simple to visualize. It can be clearly seen that the GA attempts to find an optimum variation between the two injectors as also their global location on the combustor floor. As such, the analysis can be broken up into two parts: global injector locations and relative separation between the two. In the global sense, the two injectors are seen to move modestly closer towards the base of the rearward facing step. It is especially interesting to note that despite the increased design space that moved the possible primary injector location within one step height of the base of the rearward step and the trend noticed in Case-1 regarding the global movement of the two injectors and their movement upstream, the GA has found a result that does not move beyond our initial design space for Case-1. This result makes sense if one notices the pressure contours as presented in the CFD grid validation studies at the beginning of the thesis. Given a set of inlet conditions for the combustor and a fixed geometry of the combustor section, the formation of the expansion fan at the edge of the rearward facing step is fixed despite the movement of the injectors, within certain limits. It is this ?limit? that gets broken if the primary injector were to move further upstream. In other words, if the primary injector were to move upstream, the associated bow shock would move into the expansion fan leading to increased total pressure losses. 80 Case-2 was set up to preserve total pressure ratio, and the GA successfully achieves that objective by limiting the location of the primary injector location. In retrospect the Case-1 setup was not necessarily restrictive as earlier thought for the goal of preserving total pressure ratio. Figure 3.2.4 The movement of the primary and secondary injectors versus Generation numbers for Case Run 2 Figure 3.2.5 displays the results for the relative distances between the two injectors for best performers of Case-2. Again the trends observed in Case-1 are validated in that it can be clearly seen that the general trend in preserving total pressure ratio by moving the two injector locations away from each other but closer to the forward edge of the design space (just like in Case-1) provides increased total pressure ratios with lower mixing efficiencies. 0 1 2 3 4 5 6 7 8 9 0 5 10 15 20 Inj ec to r L on git ud ina l M ov em en t Generations 81 Figure 3.2.5 The relative distance between the two injectors for Case Run 2 versus GA Generations 0 1 2 3 4 5 6 0 5 10 15 20 Ste p He igh ts Generations 82 3.3 CASE RUN - 3 The results of Case-1 and Case-2 have suitably established some trends regarding the optimization of total pressure ratios for the given combustor geometry. Based on the discussion highlighted in the previous topic it was noticed that extending the design space in the longitudinal direction provided no visible improvement (and therefore the need for the extension) as far as the total pressure ratio was concerned. It was also noticed that despite efforts to create a suitable enviornment for improving the total pressure ratio (including the high amount of damage to mixing efficiency), the amount of improvement in the value of the total pressure ratio above the experimental values was hardly justifiable for the losses incurred in other goals. It was also noticed that the mixing efficiency was showing great potential for improvement, given the chance. Three separate trends had been noticed that individually provided improvement in the mixing efficiency almost independently of each other. It was therefore possible that if combined, the overall improvement in mixing efficiency might completely sweep away the accompanying losses in the total pressure ratio. 83 But before moving on the case run with a full two dimensional design space with mixing efficiency as the dominant goal, it was necessary to validate the claims seperately. In conjunction with the results of Case-1 where the trends in the transverse direction regarding the mixing efficiency optimization were distinctively shown, it was necessary to analyze the longitudinal optimization potential for mixing eficiency. This then formed the basis for Case-3. The changes in the GA setup for Case-3 relative to Case-2 were minor. The same longitudinal single dimension design space was maintained as that for Case-2. However, the domination factors were changed to make the mixing efficiency as the dominant goal over the total pressure ratio. Figure 3.3.1 shows the details of the design space for Case-3 optimization run and is in fact indentical to that for Case-2. 84 Figure 3.3.1 The Single Dimension Longitudinal Design Space for Case Run 3 The results obtained are shown in Figure 3.3.2 and Figure 3.3.3 for the mixing efficiency and the total pressure ratio for the best and worst performers of the run respectively. The mixing efficiency has shown a drastic improvement over that of Case-1 and Case-2 results and the experimental results. At a current best value of 24.82% it has improved by 34.01% over the experimental setup. Clearly this was the dominant result beyond Generation-3 as seen in Figure 3.3.2. The reason why this value remains the same beyond Generation-3 all the way to the final Generations is explained by the setup of the GA schemes. 85 It was mentioned previously in this thesis that the steady state scheme was used. This scheme is basically an elitist population scheme and therefore preserves the best performers to the detriment of reducing the gene diversity. In this case the steady state scheme has preserved this parameter set throughout and indeed the diversity of the parameter sets was seen to reduce dramatically as the optimization took its course. The majority of the members for the final generations showed very little deviation from this best value towards the end of the run. Figure 3.3.2 The Fuel-Air Mixing Efficiency for best and worst performers for Case Run 2 versus GA Generations The total pressure ratio obviously took a hit as a result of the improvement in the mixing efficiency. As can be seen from Figure 3.3.3, the value of the total pressure ratios corresponding to the best performers for Case-3 clearly show reduced values as compared 0.15 0.17 0.19 0.21 0.23 0.25 0.27 1 6 11 16 Mi xin g E ffic ien cy Generations Best Performer Worst Performer 86 to the total pressure ratio of the worst performers in a trend similar to that observed for Case-1 and Case-2 with respect to the Mixing efficiency. The reduction in the total pressure ratio for Case-3 is 0.11% (80.34% from 80.43%). Given the enormous increase in mixing efficiency, these losses can be best considered negligible. Figure 3.3.3 Total Pressure Ratio for best and worst performers for Case Run 2 versus GA Generations However, the physical trend behind this drastic increase in fuel-air mixing efficiency still needs further investigation. Figure 3.3.4 is a plot of the injector locations in the global setup on the floor of the combustor as a function of step heights away from the rearward step through the course of the Case-3 run (best performers only). 0.8028 0.803 0.8032 0.8034 0.8036 0.8038 0.804 0.8042 0.8044 0.8046 0 5 10 15 20 To tal Pr es su re Ra tio Generations 87 It can be clearly seen that the dominant trend is to push the primary injector as far forward as possible and to push the secondary injector as far back as possible. This is consistent with the discussion for the Case-2 injector longitudinal movements in that moving the primary injector as far forward as possible leads to the destruction of the expansion fan region and leads to increase in turbulence in the flow at the cost of total pressure ratio. While this led to the GA not using it for Case-2 in order to preserve total pressure ratio, the GA did in fact use it in Case-3 to improve the mixing efficiency drastically. Further, it should be noted that from the analysis of the experimental validation studies, that the recirculation region just downstream of the rearward facing step was an interesting possibility in terms of spreading the injectant further along the transverse direction. This has been noticed by other researchers as well 12. Indeed, it should be mentioned that even with the experimental validation studies using CFD, it was noticed that a small quantity of the fuel from the primary injector had a tendency to get dragged into this recirculation region and get spread out along the floor of the combustor. In Case-3 runs, this phenomenon has been exploited by the GA to further improve mixing efficiency. As the primary injector was moved forward with the aim of destroying the expansion fan near the rearward facing step, it also moved the injector well within the recirculation region which in conjunction with the increased turbulence and erratic fluid motions caused a significant portion of the fuel from the primary injector to be dispersed far more effectively than would have been possible otherwise. 88 Figure 3.3.4 The movement of the primary (blue) and secondary (red) injectors versus Generation numbers for Case Run 3 It can now be appreciated that the trends observed in the Case-3 run regarding improvement of mixing efficiency are relatively independent of the transverse movement trend of Case-1 to the level that if combined, they provide a possibility for further increase in the mixing efficiency of the combustor. This forms the basis of Case-4, which is then our final combination study for the trends of Case-1 and Case-3. This is discussed in detail in the next topic. 1 2 3 4 5 6 7 8 0 5 10 15 20 Inj ec to r X va ria tio n i n s tep he igh ts fro m rea rw ard fa cin g s tep Generations 89 3.4 CASE RUN ? 4 As mentioned previously, Case-1 and Case-3 provided the basis for the setup for Case-4. The modification for this case run was the change from the one dimensional design space back to a two dimensional design space allowing both transverse and longitudinal variation in the injector locations. The mixing efficiency was now set as the dominant goal just like for Case-3. The remaining setup is the same as that for Case-3. Figure 3.4.1 shows the results obtained for the mixing efficiency versus the generation numbers. The mixing efficiency is found to increase above the experimental mixing efficiency (which was found to be 18.52% as previously explained) by 37.70% to a value of 25.55%. More interestingly, the value of the Case-4 best performer has increased beyond that of the Case-3 Best performer by 2.74%. 90 Figure 3.4.1 The Fuel-Air Mixing Efficiency for best and worst performers for Case Run 4 versus GA Generations Figure 3.4.2 shows the results for total pressure ratio versus the generation numbers for the Case-4 best and worst performers. There are two observations that are made for Figure 3.4.2. Firstly, it is found that the total pressure ratio for the best performers has been maintained at very nearly the same value as that of the Case-3 best performers (a reduction of 0.09% has been noted) while improving the mixing efficiency by 2.74% as mentioned previously. 0.15 0.17 0.19 0.21 0.23 0.25 0.27 0 5 10 15 20 Mi xin g E ffic ien cy Generations 91 Figure 3.4.2 The Total Pressure Ratio for best and worst performers for Case Run 4 versus GA Generations Secondly, it can be seen from Figure 3.4.2 that as the Case-4 run has progressed, not only has the mixing efficiency been improved beyond that of the generational worst performers as expected, it has also improved the total pressure ratio beyond that of the worst performers. This is significant in that in all previous runs, the total pressure ratio and the mixing efficiency had always been found to be inversely proportional to each other so that if one of these two goals was made dominant, the other goal almost always performed poorly to the level that the generational worst performers for that goal outperformed the best performers. 0.8026 0.8028 0.803 0.8032 0.8034 0.8036 0.8038 0.804 0.8042 0.8044 0.8046 0.8048 0 5 10 15 20 To tal Pr es su re Ra tio Generations 92 Figure 3.4.3 The movement of the primary (blue) and secondary (red) injectors with respect to the rearward facing step versus Generation numbers for Case Run 4 Figure 3.4.3 displays the longitudinal movement of the two injectors on the combustor floor for the best performers and the result is almost exactly the same as that obtained in Case-3. This validates the idea of the independency of the trends of Case-1 and Case-3 as originally expected. Figure 3.4.4 displays the relative transverse movements of the two injectors on the combustor floor for the best performers of each generation. Figure 3.4.4 displays the movement of the two injectors versus the run progress expressed in generation numbers. Clearly the trend is not entirely clear from this plot and 0 1 2 3 4 5 6 7 8 0 5 10 15 20 Dis tan ce in st ep he igh ts fro m th e R ea rw ard ste p Generations 93 the data seems to be behaving erratically. In order to find the underlying trend behind this data the relative distances between the two injectors are presented in Figure 3.4.5. It can be seen in Figure 3.4.5 that the GA attempts to find the optimum transverse difference between the two injectors as the run proceeds. Here it is noticed that the data seems to be oscillating about the expected optimum position when the data is compared to the minute differences in the mixing efficiency values from Figure 3.4.1. As the run proceeds, however, the oscillations are seen to be reducing and the GA is seen to be converging on the optimum location. Figure 3.4.6 shows the final locations of the optimized injector locations at the end of Case-4. This is in accordance with the trend from Case-1 where the GA finds the optimum transverse distance between the two injectors based on the value of the flammability limits of the fuel being injected into the combustor. The results from Case-1 to Case-4 therefore allow one to evaluate the optimization trends observed for the staged normal injection combustor very effectively. 94 Figure 3.4.4 The movement of the primary (blue) and secondary (red) injectors with respect to the Combustor Centerline versus Generation numbers for Case Run 4 Figure 3.4.5 The relative distance between the two injectors for Case Run 4 versus GA Generations -4 -3 -2 -1 0 1 2 3 0 5 10 15 20 Dis tan ce in m m ab ou t c en ter lin e Generations 0 0.5 1 1.5 2 2.5 0 5 10 15 20 Tra ns ve rse di sta nc e b etw ee n i nje cto rs in ste p h eig ht s Axis Title 95 Figure 3.4.6 The final results for the optimized injector locations 96 4 CONCLUSIONS AND RECOMMENDATIONS This effort involved the development of a unique way of using Genetic Algorithms and Computational Fluid Dynamics in an automated and integrated optimization structure designed to optimize a scramjet combustor for improved total pressure ratio and mixing efficiency by varying the injector locations on the floor of the combustor. Validation studies were designed to utilize the available experimental data for single case runs to develop the geometry to conduct the optimization analysis. Fuel injection was simulated so that actual values of hydrogen flammability limits were chosen for the calculations but the analysis remained close to its experimental origins by simulating air-on-air injection as fuel-on-air injection by tagging the injectant flow in the simulations. Also presented was the description of the use of the binary GA with the tournament selection method and its use in a multi-goal optimization process. A design space was created and the design parameters were specified. Also specified was the determination process for the design goals from the available three-dimensional flow simulation data. 97 The overall structure behind the automation of the entire process was presented. Computational times involved and the numerous methods used over the course of various sections of the research were discussed. Successful results were obtained in the above efforts and the current hardware and software setup available at Auburn University is now a fast, proven, robust and effective system and one that was used in the optimization process successfully. Finally, discussions of the results obtained from the optimization runs were presented. Successful results were obtained for the total pressure ratio and the mixing efficiency and the methodology was proven as being effective. Improvements of mixing efficiency in excess of 30% above those obtained in experiments with minimum losses in total pressure ratio effectively established the potential for the use of GA in propulsion system designs. Future work should include moving towards more practical situations by simulating combustion between the fuel and air within the same combustor, variations in fuel input to the system and variations in the input flow parameters. In addition, CFD methods could be made more robust with the use of Full-Multi-Grid (FMG) options now becoming available in the industry along with the use of more advanced and robust solvers based on the hybrid solver developed for this effort. Further work should also be directed towards other combustor designs including ramp based injector modeling of full scramjet engines using the system developed for this thesis for optimization efforts. 98 REFERENCES 1. Gage, P., and Kroo, I., ?A Role of Genetic Algorithms in a Preliminary Design Environment,? AIAA Paper No. 93-3933. 2. 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. 3. Burger, C. and Hartfield, R.J., ?Propeller Performance Optimization using Vortex Lattice Theory and a Genetic Algorithm?, AIAA-2006-1067, presented at the Forty-Fourth Aerospace Sciences Meeting and Exhibit, Reno, NV, Jan 9-12, 2006. 4. Doyle, J., Hartfield, R.J., and Roy, C. ?Aerodynamic Optimization for Freight Trucks using a Genetic Algorithm and CFD?, AIAA 2008-0323, presented at the 46th Aerospace Sciences Meeting and Exhibit, Reno, NV, January 2008. 5. 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 6. Hartfield, R.J., Burkhalter, J. E., and Jenkins, R. M., ?Scramjet Missile Design Using Genetic Algorithms?, Applied Mathematics and Computation, in press. 7. Hartfield, R. J., Jr., ?Planar Measurement of Flow Field Parameters in Nonreacting Supersonic Flows with Laser-Induced Iodine Fluorescence?? PhD. 99 Dissertation, Department of Mechanical and Aerospace Engineering, University of Virginia, May 1991. 8. David M. Peterson,_ Pramod K. Subbareddy_ and Graham V. Candler, ?Assessment of Synthetic Inflow Generation for Simulating Injection Into a Supersonic Cross-flow?, University of Minnesota, Minneapolis, MN, 55455 9. J. McDaniel, D. Fletcher, R. Hartfield Jr and S. Hollo, ?Staged Transverse Injection Into Mach-2 Flow Behind a Rearward-Facing step: a 3-D Compressible Test Case for Hypersonic Combustor Code Validation?, AIAA Third International Aerospace Planes Conference, 3-5 December 1991/ Orlando, FL 10. Eklund, D. R., Fletcher, D. G., Hartfield, R. J., Jr., Northam, G. B., and Dancey, C. L., ?A Comparative Computational/Experimental Compressible Flow Field Investigation. Mach 2 Flow Over a Rearward-Facing Step,? Computers in Fluids, Vol. 24, No. 5, pp. 593-608. 11. Eklund, D. R., Fletcher, D. G., Hartfield, R. J., Jr., McDaniel, J. C., Northam, G. B., Dancey, C. L., and Wang, J. A., ?A Comparative Computational/Experimental Compressible Flow Field Investigation. Staged Normal Injection into a Mach 2 Flow Behind a Rearward-Facing Step,? AIAA Journal, Vol. 32, No. 5, May 1994, pp. 907-916. 12. Dean R. Eklund, Douglas G. Fletcher, G. Burton Northam, J. C. McDaniel and Roy J. Hartfield, Jr., ?Mach-2 Flow over a Rearward-Facing step: a comparison between calculation and measurement? 13. K. Uenishi and R. C. Rogers, ?Three-Dimensional Computation of Mixing of Transverse injector in a ducted supersonic airstream?, AIAA/SAE/ASME 22nd Joint Propulsion Conference 14. Abbitt, J. D., Hartfield, R. J. and McDaniel, J. C., ?Mole Fraction Imaging of Transverse Injection in a Ducted Supersonic Flow?, AIAA Journal, Vol. 29, No 3, March 1991, pp. 431-435. 100 15. R. Courant, K. Friedrichs and H. Lewy, ?On the partial difference equations of mathematical physics?, IBM Journal, March 1967, pp. 215-234, English translation of the 1928 German original 16. Murray B. Anderson, ?An Optimization Software Package based on Genetic Algorithms?, Sverdrup Technology Inc./TEAS Group. 17. Anderson, G. Y., ?An examination of Injector/Combustor Design effects on Scramjet Performance?, Proceedings of the 2nd International Symposium on Air- Breathing Engines, Sheffield, England, March 1974. 18. William H. Heiser and David T. Pratt, with Daniel H. Daley and Unmeel B. Mehta, ?Hypersonic Airbreathing Propulsion?, AIAA Education series publication 19. Pai, S. I., ?The Fluid Dynamics of Jet?, D. Van Nostrand, New York, 1954 20. A. Abdelhafez, A. Gupta, R. Balar and K. Yu, ?Evaluation of Oblique and Traverse Fuel Injection in a Supersonic Combustor?, AIAA Paper 5026, 43rd AIAA/ASME/SAE/ASEE Joint Propulsion Conference & Exhibit, 2007 21. C. Kim, E. Jeong, J. Kim and I. Jeung, Seoul National University, Seoul, South Korea , ?Mixing and Penetration Studies of Transverse Jet into a Supersonic Crossflow?, AIAA Paper 5420, 43rd AIAA/ASME/SAE/ASEE Joint Propulsion Conference & Exhibit, 2007 22. Rogers, R. C., ?Mixing of Hydrogen injected from multiple injectors Normal to a supersonic airstream?, NASA TN D-6476, September 1971. 23. Amjad Ali Sheikh, ?A study of overall mixing loss downstream of a blunt trailing edge in a subsonic stream?, Tuskegee Institute, Alabama 36088, May 1982. 101 APPENDIX MATHEMATICAL DERIVATION OF THE AVERAGE TOTAL PRESSURE ACROSS A FLOW PLANE The equations used in this analysis have been customized to handle rectangular geometries and meshed surfaces. The mass-momentum averaging technique is explained by the further development of the following fundamental equations as detailed in the following sections. All notations remain as detailed in the main body of the thesis. A. THE CONTINUITY EQUATION: The continuity equation for a given grid plane in the three dimensional flow field can be represented as shown in equation (1): uni2211 g2025g3046g4666g1861g4667g1874g3051g4666g1861g4667g1827g4666g1861g4667 = g2025g3046g3015g3036g2880g2869 g1874g3051g1827 (1) Equation (1) can be modified by dividing both terms with the termg4666g2025g3019g3006g3007g4593 g1853g3019g3006g3007g4667g4593 ) to get the following result: uni2211 g3096g3294g4666g3036g4667g3096? Rg3137g3138g3028g3267g3254g3255g4594 g1874g3051g4666g1861g4667g1827g4666g1861g4667g3015g3036g2880g2869 = g3096g3294g3096 g3267g3254g3255? g3049g3299 g3028g3267g3254g3255g4594 g1827 (2) 102 Now, before going any further, it should be noted that for the defined reference and total conditions, using the ideal gas equation, assuming adiabatic conditions and using isentropic relationships, the following results are obtained: a) Adiabatic flow conditions: g1846g4593 = g1846g3019g3006g3007g4593 (3) b) Ideal gas law: g3096g4594 g3096g3267g3254g3255g4594 = g3017g4594 g3017g3267g3254g3255g4594 (4) c) Definition for sonic velocity: g1853g3019g3006g3007g4593 = g3493g2011g1844g1846g3019g3006g3007 = g3495g4672g2870g3082g3019g3021g3267g3254g3255g4594g3082g2878g2869 g4673 (5) d) Isentropic Relationships: g4672g3021g3294g3021?g4673 = g4678g883g3398g4672g3082g2879g2869g3082g2878g2869g4673g3436 g3049g3299g3118g3028 g3267g3254g3255g4594 g3118g3440g4679 (6) These relationships will prove highly useful in the mathematical manipulations that follow. Returning to equation (2), the right hand side of the equation can be modified in the following manner: 103 uni2211 g3096g3294g4666g3036g4667g3096? Rg3137g3138g3028g3267g3254g3255 g4594 g1874g3051g4666g1861g4667g1827g4666g1861g4667g3015g3036g2880g2869 = g3096g3294g3096? g4672 g3096g4594g3096 g3267g3254g3255? g4673 g3049g3299g3028 g3267g3254g3255g4594 g1827 (7) Substituting equation (4) in the right hand side of equation (7), yields: uni2211 g3096g3294g4666g3036g4667g3096? Rg3137g3138g3028g3267g3254g3255 g4594 g1874g3051g4666g1861g4667g1827g4666g1861g4667g3015g3036g2880g2869 = g3096g3294g3096? g4672 g3017g4594g3017 g3267g3254g3255? g4673 g3049g3299g3028 g3267g3254g3255g4594 g1827 (8) Let either side of equation (2) be represented by a function ?I?. B. MOMENTUM EQUATION: Similar to the continuity equation, the momentum equation can be represented for a given grid surface as follows: uni2211 [g2025g3046g4666g1861g4667g1874g3051g2870g4666g1861g4667g1827g4666g1861g4667g3397 g1842g3046g4666g1861g4667g1827g4666g1861g4667]g3015g3036g2880g2869 = g2025g3046g1874g3051g2870g1827g3397g1842g3046g1827 (9) Let equation (9) be represented as: J = J1 + J2 = J3 + J4 (10) Where ?J? signifies either of the two sides of the momentum equation whereas the other terms are defined as: g13g883 = uni2211 [g2025g3046g4666g1861g4667g1874g3051g2870g4666g1861g4667g1827g4666g1861g4667]g3015g3036g2880g2869 (11) 104 g13g884 = uni2211 [g1842g3046g4666g1861g4667g1827g4666g1861g4667]g3015g3036g2880g2869 (12) g13g885 = g2025g3046g1874g3051g2870g1827 (13) g13g886 = g1842g3046g1827 (14) Dividing the entire expression by g2025g4593g2902g2889g2890g1853g3019g3006g3007? g2870 and simplifying the terms individually, the following sets of results are obtained. From equation (11) g13g883 = uni2211 [g2025g3046g4666g1861g4667g1874g3051 g2870g4666g1861g4667g1827g4666g1861g4667]g3015g3036g2880g2869 g2025g4593g2902g2889g2890g1853g3019g3006g3007? g2870 Applying the ideal gas law to equation (11): g13g883 = g4678 g2902g2904g4594 g3017g4594Rg3137g3138g3028g3267g3254g3255? g3118 g4679uni2211 g4674g4672 g3017g3294g4666g3036g4667g3019g3021 g3268g4666g3036g4667 g4673g1874g3051g2870g4666g1861g4667g1827g4666g1861g4667g4675g3015g3036g2880g2869 (15) But g1839g3051g4666g1861g4667 = g3049g3299g4666g3036g4667g3493g3082g3019g3021 g3268g4666g3010g4667 , so that a modification of equation (16) is: g13g883 = g4678 g2963g2902g2904g4594 g3017g4594Rg3137g3138g3028g3267g3254g3255? g3118 g4679uni2211 g3427g1842g3046g4666g1861g4667g1839g3051g2870g4666g1861g4667g1827g4666g1861g4667g3431g3015g3036g2880g2869 (16) 105 From equation (12) g13g884 = uni2211 [g1842g3046g4666g1861g4667g1827g4666g1861g4667] g3015g3036g2880g2869 g2025g4593g2902g2889g2890g1853g3019g3006g3007? g2870 Substituting equation (5) in above equation: g13g884 = uni2211 [g1842g3046g4666g1861g4667g1827g4666g1861g4667] g3015g3036g2880g2869 g2025g4593g2902g2889g2890 g3436g884g2011g1844g1846g3019g3006g3007 g4593 g2011 g3397g883 g3440 This on further substitution with equation (4) yields: g13g884 = g4672g2963g2878g2869g2870g2963 g4673uni2211 [g3017g3294g4666g3036g4667g3002g4666g3036g4667]g3263g3284g3128g3117g3017g4594 Rg3137g3138 (17) From equation (13), g13g885 = g3096g3294g3049g3299g3118g3002 g3096g4594Rg3137g3138g3028g3267g3254g3255? g3118 (18) Applying some simple modifications to equation (18): g13g885 = g3096g3294g3049g3299g3118g3002 g3096g4594g3028g3267g3254g3255? g3118 g3436 g3096g4594g3096g4594 Rg3137g3138 g3440 (19) 106 Substituting equation (4) in equation (19), g13g885 = g3096g3294g3049g3299g3118g3002 g3096g4594g3028g3267g3254g3255? g3118 g4672 g3017g4594g3017g4594 Rg3137g3138 g4673 (20) g13g886 = g3017g3294g3002 g3096g4594Rg3137g3138g3028g3267g3254g3255? g3118 (21) Substituting equation (5) in equation (21): g13g886 = g3017g3294g3002 g3096g4594Rg3137g3138g3118g3330g3267g3269g3267g3254g3255 g4594 g3330g3126g3117 (22) Now substituting equation (4) in equation (22), we get: g13g886 = g4672g2963g2878g2869g2870g2963 g4673g4672g3017g3294g3017?g4673g4672 g3017g4594g3017 g3267g3254g3255g4594 g4673g1827 (23) Substituting equations (15), (17), (20) and (23) in equation (10), the following expression for the momentum equation is obtained: g13 = g4678 g2963g2902g2904g4594 g3017g4594Rg3137g3138g3028g3267g3254g3255? g3118 g4679uni2211 g3427g1842g3046g4666g1861g4667g1839g3051g2870g4666g1861g4667g1827g4666g1861g4667g3431g3015g3036g2880g2869 g3397 g4672g2963g2878g2869g2870g2963 g4673uni2211 [g3017g3294g4666g3036g4667g3002g4666g3036g4667] g3263g3284g3128g3117 g3017g4594Rg3137g3138 = g3096g3294g3049g3299g3118g3002 g3096g4594g3028g3267g3254g3255? g3118 g4672 g3017g4594g3017g4594 Rg3137g3138 g4673g3397 g4672g2963g2878g2869g2870g2963 g4673g4672g3017g3294g3017?g4673g4672 g3017g4594g3017 g3267g3254g3255g4594 g4673g1827 (24) 107 Now, from the right hand side terms of equations (8) and (24): g3011 g3010 = g3344g3294g3297g3299g3118g3250 g3344g4594g3276g3267g3254g3255? g3118 g3436 g3265g4594g3265g4594 Rg3137g3138 g3440g2878 g4672g3211g3126g3117g3118g3211 g4673g4672g3265g3294g3265?g4673g4678 g3265g4594g3265 g3267g3254g3255g4594 g4679 g3344g3294 g3344?g4678 g3265g4594 g3265g3267g3254g3255? g4679 g3297g3299 g3276g3267g3254g3255g4594 (25) Equation (25) on further simplification gives us the following expression: g3011 g3010 = g3344g3294g3297g3299g3118g3250 g3344g4594g3276g3267g3254g3255? g3118 g2878 g4672g3211g3126g3117g3118g3211 g4673g4672g3265g3294g3265?g4673 g3344g3294 g3344? g3297g3299 g3276g3267g3254g3255g4594 (26) Now defining a variable ?c? so that: g1855 = g3011g3010 (27) Substituting equation (27) in equation (26) and modifying the equation: g3436 g3049g3299g3118g3028 g3267g3254g3255g4594 g3118 g3440 g514 g1855g4672 g3049g3299 g3028g3267g3254g3255g4594 g4673g3397 g4672 g2963g2878g2869 g2870g2963 g4673g4672 g3017g3294 g3017?g4673g4672 g3096g4594 g3096g3294g4673 = g882 (28) Substituting equation (4) in equation (28): 108 g3436 g3049g3299g3118g3028 g3267g3254g3255g4594 g3118 g3440 g514 g1855g4672 g3049g3299 g3028g3267g3254g3255g4594 g4673g3397 g4672 g2963g2878g2869 g2870g2963 g4673g4672 g3021g3294 g3021?g4673 = g882 (28) Substituting equation (6) in equation (28): g3436 g3049g3299g3118g3028 g3267g3254g3255g4594 g3118 g3440 g514g4672 g2870g3030g3082 g3082g2878g2869g4673g4672 g3049g3299 g3028g3267g3254g3255g4594 g4673g3397 g883 = g882 (29) Equation (29) is a quadratic equation for the variable g4672 g3049g3299g3028 g3267g3254g3255g4594 g4673 and can be solved as: g3049g3299 g3028g3267g3254g3255g4594 = g4672 g3082g3030 g3082g2878g2869g4673g3398 g3495g4666 g3082g3030 g3082g2878g2869g4667 g2870 g3398 g883 (30) Also, from isentropic relationships, the following relationship for the density ratio is obtained: g3096g3294 g3096? = [g883 g3398 g4672 g3082g2879g2869 g3082g2878g2869g4673g484 g3049g3299g3118 g3028g3267g3254g3255? g3118 ]g4672 g3117 g3330g3127g3117g4673 (31) Now, from equation (8): g1835 = g3533 g2025g3046g4666g1861g4667g2025? g2902g2889g2890g1853g3019g3006g3007 g4593 g1874g3051g4666g1861g4667g1827g4666g1861g4667 g3015 g3036g2880g2869 = g2025g3046g2025? g4678 g1842 g4593 g1842g3019g3006g3007? g4679 g1874g3051 g1853g3019g3006g3007g4593 g1827 109 Substituting equation (4) in the left hand side term of equation (8): g1835 = uni2211 g4672g3096g3294g4666g3036g4667g3019g3021?g3017 g3267g3254g3255 ? g4673g4672g3049g3299g4666g3036g4667g3028 g3267g3254g3255g4594 g4673g1827g4666g1861g4667g3015g3036g2880g2869 (32) On further manipulations of equation (35): g1835 = uni2211 g4672g3017g3294g4666g3036g4667g3082g3019g3021?g3017 g3267g3254g3255 ? g4673g4672 g3014g3299g4666g3036g4667g3118g3028 g3267g3254g3255g4594 g3049g3299g4666g3036g4667 g4673g1827g4666g1861g4667g3015g3036g2880g2869 (33) Or on simplifying the equation: g1835 = g4672 g3082g3019g3021?g3017 g3267g3254g3255 ? g3028g3267g3254g3255g4594 g4673uni2211 g1842g3046g4666g1861g4667g4672g3014g3299g4666g3036g4667g3118 g3049 g3299g4666g3036g4667 g4673g1827g4666g1861g4667g3015g3036g2880g2869 (34) Substituting equations (24) and (34) into equation (27): g1855 = g3437 g3211RTg4594 g3265g4594Rg3137g3138g3276g3267g3254g3255? g3118 g3441uni2211 g3427g3017g3294g4666g3036g4667g3014g3299g3118g4666g3036g4667g3002g4666g3036g4667g3431g3263g3284g3128g3117 g2878 g4672g3211g3126g3117g3118g3211 g4673uni2211 [g3265g3294g4666g3284g4667g3250g4666g3284g4667] g3263g3284g3128g3117 g3265g4594Rg3137g3138 g4678 g3330g3267g3269?g3265 g3267g3254g3255 ? g3276g3267g3254g3255g4594 g4679uni2211 g3017g3294g4666g3036g4667g3436g3262g3299g4666g3284g4667g3118 g3297g3299g4666g3284g4667 g3440g3002g4666g3036g4667g3263g3284g3128g3117 (35) 110 As can be seen, the g1842g4593g2902g2889g2890 term cancels from both the numerator and the denominator of equation (35), and on further simplification: g1855 = g3437 g3211RTg4594 g3276g3267g3254g3255? g3118 g3441uni2211 g3427g3017g3294g4666g3036g4667g3014g3299g3118g4666g3036g4667g3002g4666g3036g4667g3431g3263g3284g3128g3117 g2878 g4672g3211g3126g3117g3118g3211 g4673uni2211 [g3017g3294g4666g3036g4667g3002g4666g3036g4667]g3263g3284g3128g3117 g4678g3330g3267g3269?g3276 g3267g3254g3255g4594 g4679uni2211 g3017g3294g4666g3036g4667g3436g3262g3299g4666g3284g4667g3118 g3297g3299g4666g3284g4667 g3440g3002g4666g3036g4667g3263g3284g3128g3117 (36) Equation (36) provides the value for ?c? in terms of the grid data. Finally, equating the right hand term of equation (7) and (34): g3096g3294 g3096? g4672 g3096g4594 g3096g3267g3254g3255? g4673 g3049g3299 g3028g3267g3254g3255g4594 g1827 = g4672 g3082g3019g3021? g3017g3267g3254g3255 ? g3028g3267g3254g3255g4594 g4673uni2211 g1842g3046g4666g1861g4667g4672 g3014g3299g4666g3036g4667g3118 g3049g3299g4666g3036g4667 g4673g1827g4666g1861g4667 g3015g3036g2880g2869 (37) Equation (40) when further simplified, gives: g1842? = g4678g3330g3267g3269?g3276 g3267g3254g3255g4594 g4679uni2211 g3017g3294g4666g3036g4667g3436g3262g3299g4666g3284g4667g3118 g3297g3299g4666g3284g4667 g3440g3002g4666g3036g4667g3263g3284g3128g3117 g3344g3294 g3344? g3297g3299 g3276g3267g3254g3255g4594 g3002 (38) Equation (38) provides the final form of the equation which can be solved for the average total pressure g1842? across a given grid plane in the three dimensional flow field. The terms in the denominator for equation (38) can be found from equations (30), (31) and (36) 111 which are themselves a result of flow field grid data. Note that g1842g1314 and g1874g3051, while average values, simultaneously satisfy the continuity and momentum equations. They thus have physical as well as statistical meaning.