NUMERICAL ANALYSIS OF HIGH PRESSURE INJECTION OF NATURAL GAS INTO DIESEL ENGINE COMBUSTION CHAMBERS By Paul Walsh B.A.Sc., The University of British Columbia, 1989 A THESIS SUBMITTED IN PARTIAL FULFILMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE in THE FACULTY OF GRADUATE STUDIES (MECHANICAL ENGINEERING DEPARTMENT) We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA December 1991 ©PaulWsh,19 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. (Signatu Department of if&C /1.0171/cg ( The University of British Columbia Vancouver, Canada Date DE-6 (2/88) vtAevty^/ 92_ ABSTRACT The behaviour of a transient turbulent jet of natural gas as it is injected into a simulated combustion chamber of a diesel engine was investigated using numerical techniques. The TEACH code developed by A.D. Gosman of Imperial College, London, was used to investigate the influence of parameters such as injection angle, engine speed, and reservoir tank-to-chamber pressure ratio on the development of the jet. It has been shown that the TEACH code is fully capable of predicting details of jet behaviour such as radial and axial velocity and concentration profiles when compared to known data. The code has been modified to use a compressing and expanding grid to simulate the effects of piston motion. A model of a fixed geometry combustion chamber revealed that the most influential parameter on jet behaviour is the injection angle. The jet had a tendency to adhere to either the top wall of the chamber near the injector tip or to the bottom wall directly opposite depending on the injector angle. The compressing grid simulation showed that the presence of piston motion combined with other parameters such as injection angle and pressure ratio, produced jet characteristics that were dissimilar compared to the fixed boundary model. In general it was shown that the jet was less sensitive to injection angle and strongly influenced by increased pressure ratio as a result of the moving boundary. ii TABLE OF CONTENTS ABSTRACT^ ii TABLE OF CONTENTS^ iii LIST OF FIGURES^ vi ACKNOWLEDGEMENTS^ xi NOMENCLATURE^ xii CHAPTER 1: INTRODUCTION 1.1 General^ 1 1.2 Gas Jets^ 3 1.3 Numerical Analysis^ 8 1.4 Thesis Objectives and Structure^ 13 CHAPTER 2: MATHEMATICAL FORMULATION 2.1 General^ 16 2.2 Governing Equations^ 16 2.3 k-e Turbulence model^ 19 2.4 Discretization of Governing Equations ^ 22 2.5 The SIMPLE Algorithm ^ 26 CHAPTER 3: JET SIMULATION 3.1 General^ 30 3.2 Configuration and Assumptions^ 30 3.3 Moving Boundary Modification ^ 35 3.4 Density Change Due to Piston Motion and Gas Injection ^ 40 3.5 Boundary Conditions^ iii 42 CHAPTER 4: ROUND FREE JET 4.1 General^ 49 4.2 Steady-State Jet^ 50 4.3 Grid Size and Time Step Dependence ^ 51 4.4 Comparison with Analytical and Experimental Results ^ 53 4.5 Radial Profiles and Jet Penetration Data ^ 57 CHAPTER 5: FIXED PISTON RESULTS 5.1 General^ 64 5.2 Comparison to Experiment ^ 65 5.3 Injection Angle^ 74 5.4 Pressure Ratio^ 81 5.5 Combustible Mixture Region^ 87 5.6 Clearance Gap^ 91 CHAPTER 6: MOVING PISTON MODEL 6.1 General^ 98 6.2 Effect of Cylinder Pressure and Temperature Level ^ 99 6.3 Pressure Ratio Effects^ 100 6.4 Injection Angle Effects^ 105 6.5 Engine Speed^ 110 6.6 Combustible Mixture Region ^ 115 CHAPTER 7: CONCLUSIONS AND RECOMMENDATIONS ^118 REFERENCES^ 120 iv APPENDIX A: GAS JET ANALYSIS ^ v 123 LIST OF FIGURES Figure 2.4.1 Standard control volume configuration ^ 21 Figure 3.2.1 Conical gas jet in relation to piston and cylinder ^ 30 Figure 3.2.2 Details of the conical gas jet^ 31 Figure 3.2.3 Computation grid inside combustion chamber^ 32 Figure 3.2.4 Detailed computation grid^ 33 Figure 3.3.1 Details of compressing and non-compressing grid regions. ^35 Figure 3.5.1 Injection speed^ 42 Figure 3.5.2 Underexpanded jet details^ 43 Figure 3.5.3 Boundary control volumes ^ 44 Figure 4.2.1 Axial velocity along the jet axis at various times; 10 m/s initial jet speed 48 Figure 4.2.2 Radial concentration profiles at various times; 0.46 m from jet orifice; 10 m/s initial jet speed^ 49 Figure 4.2.3 Radial velocity profiles^ at various times; 0.46 m from jet orifice; 10 m/s initial jet speed 50 Figure 4.3.1 Radial velocity profiles for various time step sizes; 0.435 m from jet ^ orifice; 10 m/s initial jet speed 51 Figure 4.3.2 Radial concentration profiles at various time steps; 0.46 m from jet ^ 52 orifice; 10 m/s initial jet speed; mass fraction in kg/kg Figure 4.3.3 Radial concentration profiles at various grid sizes; 0.46 m from jet ^ 53 orifice; 10 m/s initial jet speed; C m is centerline concentration Figure 4.4.1 Comparison of numerical and analytical axial velocity profiles; U. ^ 54 is 10 m/s initial jet speed; U m is the centerline velocity Figure 4.4.2 Corrected prediction of the axial velocity profile; U. ^ 55 is 10 rri/s initial jet speed; Um is the centerline velocity Figure 4.4.3 Comparison of numerical and experimental axial concentration profiles for a methane jet; 10 m/s initial jet speed; C o is the jet ^ outlet concentration vi 56 Figure 4.5.1 Comparison of velocity profiles of Tolhnien, Schlichting and the numerical analysis; 10 m/s initial jet speed; U. is the jet centerline velocity 57 Figure 4.5.2 Comparison of numerical radial concentration profile with Tollmien's data ; 0.46 m from jet orifice; 10 m/s initial jet speed; C. is centerline concentration 58 Figure 4.5.3 Comparison of a numerical model and Birch's data of a radial concentration profile of a methane jet; 0.46 m from jet orifice ; 10 m/s initial jet speed; C. is centerline concentration 59 Figure 4.5.4 Comparison of modelled jet penetration times with that of Kuo and Bracco; 10 m/s initial jet speed.^ 60 Figure 5.1.1 Concentration field for case 1 showing jet development with time, clearance gap 0.3125 inches, pressure ratio 2:1, jet angle 10 deg. contour interval 0.04 (kg/kg), injector tip in top left corners. ^63 Figure 5.1.2 Concentration field for case 8 showing jet development with time, clearance gap 0.6250 inches, pressure ratio 5:1, jet angle 30 deg. contour interval 0.04 (kg/kg), injector tip in top left corners. ^64 Figure 5.1.3 Velocity field for case 8 showing jet development with time, clearance gap 0.6250 inches, pressure ratio 5:1, jet angle 30 deg. , injector tip in top left corners. ^ 65 Figure 5.2.1 Schlieren image compared to case 8; jet angle 30 deg; pressure ratio 5:1, clearance gap 0.6250 inches, concentration contours 0.08 (kg/kg) , injector tip in top left corners. ^ 66 Figure 5.2.2 Schlieren image compared to case 7; jet angle 10 deg; pressure ratio 5:1, clearance gap 0.6250 inches, concentration contours 0.08 (kg/kg) , injector tip in top left corners. ^ 67 Figure 5.2.3 Schlieren image compared to case 1; jet angle 10 deg; pressure ratio 2:1, clearance gap 0.3125 inches, concentration contours 0.08 (kg/kg) ^ , injector tip in top left corners. 68 Figure 5.2.4 Schlieren image compared to case 6; jet angle 30 deg; pressure ratio 5:1, clearance gap 0.3125 inches, concentration contours 0.08 (kg/kg) , injector tip in top left corners. ^ 69 Figure 5.3.1 Concentration fields for case 1 and 2 showing the jet angle effect, pressure ratio 2:1, clearance gap 0.3125 inches, contour interval ^ 0.04 (kg/kg), injector tip in top left corners. vii 73 Figure 5.3.2 Velocity fields for case 1 and 2 showing the jet angle effect, pressure ratio 2:1, clearance gap 0.3125 inches, injector tip in top left corners. Figure 5.3.3 Concentration fields for case 7 and 8 showing the jet angle effect, pressure ratio 5:1, time 3.59 ms, contour interval, 0.04 (kg/kg), injector tip in top left corners. Figure 5.3.4 Velocity fields for case 7 and 8 showing the jet angle effect, pressure ratio 5:1, time 3.59 ms, clearance gap 0.6250 inches, injector tip in top left corners. Figure 5.3.5 Concentration fields showing jet angle sensitivity, pressure ratio 5:1, time 2.13 ms, contour interval, 0.08 (kg/kg), injector tip in top left corners. Figure 5.4.1 Concentration fields for case 1 and 5 showing pressure ratio effects, jet angle 10 deg., clearance gap 0.3125 inches, contour interval, 0.04 (kg/kg), injector tip in top left corners. Figure 5.4.2 Velocity fields for case 1 and 5 showing pressure ratio effects, jet angle 10 deg., clearance gap 0.3125 inches, injector tip in top left corners. Figure 5.4.3 Concentration fields for case 8 and 3 showing pressure ratio effects, jet angle 30 deg.,time 2.13 ms, clearance gap 0.6250 inches, contour interval 0.08 (kg/kg), injector tip in top left corners. Figure 5.4.4 Velocity fields for case 8 and 3 showing pressure ratio effects, jet angle 30 deg., time 2.13 ms, clearance gap 0.6250 inches, injector tip in top left corners. Figure 5.4.5 Pressure fields for case 8 and 3 showing pressure ratio effects, jet angle 30 deg., time 2.13 ms, clearance gap 0.6250 inches, injector tip in top left corners. Figure 5.5.1 Combustible region development for case 7 ,pressure ratio 5:1, jet angle 10 deg., clearance gap 0.6250 inches, contours given in mass fraction (kg/kg), injector tip in top left corners. Figure 5.6.1 Concentration fields for case 6 and 8 showing the effect of clearance gap, pressure ratio 5:1, jet angle 30 deg.,time 2.13 ms, contour interval 0.08 (kg/kg), injector tip in top left corners. viii Figure 5.6.2 Concentration fields for case 6 and 8 showing the effect of clearance gap, pressure ratio 5:1, jet angle 30 deg.,time 3.59 ms, contour interval 0.08 (kg/kg), injector tip in top left corners. ^ 90 Figure 5.6.3 Concentration fields for case 1 and 4 showing the effect of clearance gap, pressure ratio 2:1, jet angle 10 deg.,time 2.13 ms, contour interval 0.08 (kg/kg), injector tip in top left corners.^ 92 Figure 5.6.4 Concentration fields for case 1 and 4 showing the effect of clearance gap, pressure ratio 2:1, jet angle 10 deg.,time 3.59 ms, contour interval 0.08 (kg/kg), injector tip in top left corners.^ 93 Figure 5.6.5 Velocity fields for case 6 and 8 showing the effect of clearance gap, pressure ratio 5:1, jet angle 30 deg.,time 3.59 ms, injector tip in top left corners.^ 94 Figure 6.2.1 Concentration fields for case 7 and 1 comparing actual engine and STP initial conditions, pressure ratio 2:1, jet angle 30 deg., engine speed 1200 RPM, contour interval 0.08 (kg/kg), injector tip in top left corners. ^ 98 Figure 6.3.1 Concentration fields for case 1 and 2 showing pressure ratio effects, jet angle 30 deg., engine speed 1200 RPM, contour interval 0.08 (kg/kg), 99 injector tip in top left corners.^ Figure 6.3.2 Velocity fields for case 1 and 2 showing pressure ratio effects, jet angle 30 deg., engine speed 1200 RPM, injector tip in top left corners.^ 100 Figure 6.3.3 Concentration fields for case 8 and 10 showing pressure ratio effects, jet angle 10 deg., engine speed 1200 RPM, contour interval 0.08 (kg/kg), 102 injector tip in top left corners.^ Figure 6.4.1 Concentration fields for cases 8,1 and 3 showing injection angle effects, pressure ratio 2:1, engine speed 1200 RPM, crank angle 6.8 deg. BTDC, contour interval 0.08 (kg/kg), injector tip in top left corners. 104 Figure 6.4.2 Velocity fields for cases 8,1 and 3 showing injection angle effects, pressure ratio 2:1, engine speed 1200 RPM, crank angle 6.8 deg. BTDC, contour interval 0.08 (kg/kg), injector tip in top left corners. ^105 Figure 6.5.1 Concentration fields for cases 2 (fixed boundary),6 and 1 showing effects of engine speed, pressure ratio 2:1, jet angle 30 deg., time 5.3 ms after injection start, contour interval 0.08 (kg/kg), injector tip in top 107 left corners.^ ix Figure 6.5.2 Velocity fields for cases 2 (fixed boundary),6 and 1 showing effects of engine speed, pressure ratio 2:1, jet angle 30 deg., time 5.3 ms after injection start, injector tip in top left corners.^ 107 Figure 6.5.3 Concentration fields for cases 1 (fixed boundary),9 and 8 showing effects of engine speed, pressure ratio 2:1, jet angle 10 deg., time 5.3 ms after injection start, contour interval 0.04 (kg/kg), injector tip in top 109 left corners.^ Figure 6.5.4 Velocity fields for cases 1 (fixed boundary),9 and 8 showing effects of engine speed, pressure ratio 2:1, jet angle 10 deg., time 5.3 ms after injection start, injector tip in top left corners. ^ 110 Figure 6.6.1 Combustible region development for case 1 ,pressure ratio 2:1, jet angle 30 deg., engine speed 1200 RPM, contours given in mass fraction 113 (kg/kg), injector tip in top left corners. ^ x ACKNOWLEDGEMENT Th , author would like to express his sincere appreciation to his supervisor Prof. P.G. hill ,, for all his expert advice and encouragement which made the completion of this thesis possible in a timely and orderly fashion. I would also like to extent my gratitude to my co-supervisor Prof. M. Salcudean for her technical expertise and generous advice. I would also like to thank Mr. Patric Ouellette for the use of his experimental results related to this work. The author is also indebted to his family for their advice and support, and would like to express his appreciation. I would finally like to express my appreciation to Ms. Dorothy Miyake for her assistance on matters related to computer software operation. xi NOMENCLATURE Ae^- mi.,. of control volume face A,ap.e„ - neighbouring coefficients of linearized equations C^- concentration d,D^- diameter de^- effective/equivalent diameter (high pressure ratio) Fr^- body force term F,D^- convective,diffusive components G^- rate of strain term k^- turbulent kinetic energy Pa,Po^- ambient,tank pressure Pe^- Peclet number ( ratio of F over D) Red^- Reynolds number ( diametral) r,y^- radial coordinate S^- source term t,t*^- time, nondimensional time U,u^- axial velocity ( mean,fluctuating) V,v^- radial velocity ( mean,fluctuating) y+^- nondimensional distance from wall z,x,x^* - axial coordinate ( * - nondimensional ) 4^- piston face distance from injector tip K,a,d,ce^- constants xii Greek - turbulent energy dissipation • nondimensional axial coordinate ( z/zh) r^- diffusivity (p/a) - dependant variable ( U,V,k,e,C) p ^- density Peff - dynamic viscosity ( turbulent,laminar,effective) a^- Prandtl number - kinematic viscosity ti W^- wall shear stress Superscripts - future point in time o - previous point in time - guessed value - averaged value Subscripts o - jet inlet value m^- jet centerline value p,n,e,s,w^- neighbouring nodes 1 CHAPTER 1 : INTRODUCTION 1.1 GENERAL Natural gas has been used as an automotive fuel for many years. Only recently has it come under close study as an alternative to conventional automotive fuels. As a new environmental awareness occupies public attention, the quest for automotive fuels that offer the performance of standard fuels but with far fewer harmful emissions has expanded. The public has already demanded that their commercial diesel carriers meet strict new guidelines through legislation to be introduced in the next three years. The new laws will limit the amounts of nitrogen oxides (NOx) and particulate matter to levels considered to be restrictive by the diesel engine manufacturing establishment. Particulate levels must be reduced to 1/6 their 1989 levels while NOx must be reduced by half. It may be that the only feasible manner in which these new standards can be met is with the use of alternative fuels such as natural gas. Unlike gasoline and standard diesel fuel, natural gas is a simple hydrocarbon resulting in fewer combustion byproducts, giving it the potential to be a clean air fuel [34]. With careful mixture and timing control, emissions of carbon monoxide, particulates, NOx and reactive hydrocarbons can be reduced to levels below those of standard spark ignition and diesel engines. Additionally, natural gas has a low photochemical reactivity meaning that unburned fuel will not contribute significantly to smog problems plaguing urban centres. Clearly, natural gas has the potential to run cleaner than standard diesel fuel. Although natural gas has been used in spark ignition engines for some time, its use as a diesel fuel raises certain problems. The flame speed of natural gas is low, especially at high 2 compression ratios. This will lead to poor ignition and burning characteristics causing lower power output. Adding a small amount of 'pilot' diesel fuel to the natural gas to promote ignition is one possible remedy to this problem. Several diesel engine strategies are available related to the timing and to hardware requirements of dual fuel systems. Early cycle injection places the natural gas charge into the cylinder before the pilot diesel fuel, early in the compression stroke when cylinder pressures are low, eliminating the need for a high pressure injection system. However, if the fuel is well mixed, a detonation limit may be reached where the fuel pre-ignites once certain conditions are achieved. This will restrict compression ratios to values below which preignition will not occur. This limit can hamper the achievement of the full operating potential of the engine. The fumigation process injects the natural gas fuel into the intake air stream before the compression stroke commences. The natural gas will displace approximately 10% of the intake air flow, reducing the maximum amount of intake air and performance compared to late cycle injection systems. A detonation limit will restrict the operating range of an engine configured in such a manner. Such systems do not perform as well as late cycle direct injection engines. High pressure, late cycle, direct injection with electronic control of injection duration and timing offers technical advantages. Since injection occurs late in the compression stroke, none of the intake air is displaced by the natural gas. No fuel of any kind is present in the cylinder during the compression stroke; removing any possibility of pre-ignition. Therefore, a late cycle injection can preserve the diesel cycle efficiency. Even though a late cycle high 3 pressure system appears to offer the best chance for a clean efficient system, the proper timing and mixing is still required to achieve peak performance. The Alternative Fuels Laboratory at the University of British Columbia is developing a dual fuel, late cycle injection system with an innovative injector design. In contrast to traditional dual fuel systems which have one injector for each fuel, the U.B.C. project combines the two fuels in one injector. The aim is to introduce the pilot diesel directly into the stream of natural gas as it is injected into the cylinder. It is thought that this process will result in a uniform distribution of fine diesel droplets which should promote quicker ignition and faster burning. In spite of the fact that late cycle direct injection systems can offer superior performance compared to other systems, very few diesel engine researchers have attempted to implement it since it is technically the most challenging. Late cycle injection requires both a high pressure gas supply and fast acting injectors, which represents a considerable development effort. Once constructed, the injector jet characteristics will need to be investigated to allow optimization for peak engine performance. To meet such goals, both experimental and numerical modelling techniques are called upon to investigate the jet performance. 1.2 GAS JETS The injector tip is designed to produce a continuous sheet of gas instead of several small round jets as in the traditional liquid diesel fuel injector. Initially the jet is very thin, approximately one third of a millimetre in thickness. It emerges from the injector tip as a high speed axially symmetric conical sheet. The injector tip is centrally located in the engine cylinder giving axial symmetry to the jet-cylinder system. This provides the opportunity to 4 assume two dimensional flow for any model of this system. The pressure ratio between the fuel tanks and the engine cylinders may be as high as 5 to 1. The jet may therefore be under-expanded, meaning that the pressure at its exit plane will be higher than the cylinder pressure into which it emerges. The jet will expand to cylinder pressure as it leaves the injector, forming a region where compressibility effects will be significant. Ewan and Moodie [11] have made a study of such jets. Their results show that for a fuel tank to cylinder pressure ratio of 5 to 1 the length of this expansion region will be approximately 1.3 times as large as the exit diameter. Considering that the jet velocity drops to one third its original velocity within approximately 18 exit diameters of this initial expansion, any compressibility effects will therefore be limited to a very small region ( one tenth the cylinder diameter). It is not unreasonable then to assume that the flow field inside the cylinder is mainly incompressible. The gas jet behaviour depends on the injection angle, the injection pressure, the chamber geometry and the engine speed. The effects of these parameters can best be assessed through the determination of the jet velocity and concentration fields. These fields can be determined at least qualitatively by either analytical, experimental or numerical means. Simplified models can be solved analytically to predict flow fields and jet behaviour in an attempt to understand mixing and penetrating character. Experimental analysis uses physical representations either full size or scaled models in order to simulate the actual in cylinder processes . Much analytical and experimental work has been conducted on the subject of free incompressible steady-state turbulent jets. A review of some of these results will allow a comparison to be made to a steady-state version of the code used in this thesis. If the 5 numerical code can successfully predict turbulent jet behaviour, then this will be a positive assessment of its analysis capabilities and will give confidence for its use in in-cylinder transient jet investigations. Schlichting [29] gives an analytical derivation created by Tollmien for the velocity profile of a circular turbulent jet emptying into free space. The governing equations are approximated using momentum integral techniques assuming that the centerline velocity is inversely proportional to the axial distance. The solution gives an estimate of the axial velocity as a function of the radial and axial coordinates. u - 3 K^1 87c c ox 1 = 1 \ 3 ,,rfe y I -ie.x 11 -. (1+ _,,, 1 2) 2 4" where K= 2n f u 2 ydy 0 Here, y and x are the radial and axial jet coordinates respectively, and K is the specific momentum. The virtual kinematic viscosity for round free jets has been determined experimentally by Riechard [29] to be e ° = 0 . 0161 . 1,8? This value is a constant for round free jets as a result of the centreline velocity , u m , being inversely proportional to the axial distance x. The results computed using this analysis agree well with experimental results presented by Schlichting. The universality of the velocity profile for a turbulent jet is demonstrated by the reduction of the experimental data at various axial locations in the jet and of the predicted values onto one standard curve of iliu m vs. T. ^ 6 Abramovich [2] also quotes Tollmien and gives tabulated results of velocity, concentration and axial velocity decay for round free turbulent jets. Tollmien's axial velocity decay formula for round jets begins with u m = C/x. Using momentum conservation it can be shown that urn ^u _ 0 . 96 0^ax 1?„, u. and R. are the initial jet axial velocity and radius respectively. The constant 0.96 is the dimensionless distance (ax/R o = 0.96) from the jet origin where x=0, to the transition cross section where u m/u. = 1. This theory agrees well with data taken by Trupel and Gottingen [2], who found the constant a to be 0.066. A study of the concentration fields of steady-state turbulent round free jets of methane was conducted by Birch et al [7]. They found that methane jets can be analyzed using incompressible jet theory so long as a correction is made to the starting diameter to account for the density difference between the gas and the ambient air. The adjusted diameter was referred to as the effective diameter, and was computed from similarity conditions by Thring and Newby (1953). A free jet can be characterized by its initial diameter d o , its total momentum M, and its total initial mass flow rate m. For an initially uniform velocity these are: 1 1 n do2 p o u 0. M = --4 ndo2 p o u o2^m = — 4 From these Thring and Newby deduced that an incompressible jet of initial density p o , differing from the ambient density p a , will behave as a jet having p o equal to p a but having 7 a different starting diameter. For the same mass and momentum fluxes the effective diameter (dE ) will be: dE = d (p o /p a ) 1 / 2 The initial jet diameter is d, and the gas and air densities at ambient pressure and temperature are respectively p o and p a . Their experiments led them to conclude that the centerline concentration of a methane jet decayed with distance according to the following relation; C _ kide Co^(z +a) where C is concentration, C o is inlet concentration, k 1 and a are constants, and z is the axial coordinate. Concentration profiles in the radial direction as a function of axial position were also investigated and were found to obey the relation C = exP( D(r/ z) 2 ) . - m Here, D is a constant, found experimentally to be 73.6. Cm is the centreline concentration, and r is the radial coordinate. Kuo and Bracco [21] investigated the transient nature of incompressible jets numerically. They compared their results produced by the implicit finite volume TEACH code, which employs the k-E turbulence model, to experiments of air injected into air to test the performance of the code. Care was taken in the experiments to ensure a uniform velocity profile at the nozzle outlet in order to justify such a velocity profile in the numerical 8 code. Time and axial position were non-dimensionalized for comparing cases. t* - t Llin^ x* — ^X D (41 °53 )^D (4' °53 ) Here, t indicates time, x axial distance from the orifice, D the nozzle diameter, and D the kinematic viscosity. They observed the time required for a position on the centreline downstream to reach 70.0% of its steady state velocity. The centreline velocity decay with axial position for the steady jet was tested as well. The results of the numerical analysis appeared to agree fairly well with the experimental data, validating the use of the TEACH code for transient turbulent jet studies. The penetration time of their transient jet investigation is given in the following relation. t*=0 .235x* 2^x*k7 With the experimental information on jet behaviour just described, the validity of using a numerical code for its analysis can then be assessed. 1.3 NUMERICAL ANALYSIS The purpose of the research described in this thesis is to obtain insight into jet-cylinder interactions through numerical means. Numerical analysis attempts to investigate jet behaviour through an iterative solution of the governing equations of mass and momentum conservation by computer algorithm given a set of inlet and boundary conditions. 9 Numerical models of fluid dynamic phenomena, have become an increasingly popular research tool for a number of reasons. First, a more comprehensive insight of the phenomena under study can be obtained . The effects of numerous parameters such as injection pressure and spray angle can be investigated quickly and with minimal effort. Answers to important questions such as how piston motion affects the jet behaviour and the subsequent in-cylinder flow field may be obtained. Secondly, complete representations of concentration, velocity and pressure fields at any point in time are possible. Finally, new design ideas can be tested thoroughly over the complete range of engine operation before costly engine experiments are performed. Numerical fluid flow models approximate the governing equations through various means. Numerous methods have been devised to reduce the governing equations to a series of linear algebraic equations in a process known as discretization. The control volume method is one discretization scheme used frequently in engine modelling. The region of interest is subdivided into a grid of nodes each surrounded by a control volume . The governing equations are integrated over each control volume to produce one linear equation per dependent variable per control volume. The values of dependent variables are calculated at the node point within each control volume from the corresponding linear equation produced by the discretization process. Values at adjacent nodes are related through a number of possible weighting schemes, which will form the basis of the numerical code. The finite difference method which is also common in engine modelling derives the linearized governing equations in a different manner. The terms appearing in the governing equations are approximated through a Taylor series expansion of the dependent variable. 10 The work described in this thesis is the development and application of a two-dimensional finite volume model for the study of transient injections of natural gas into a diesel engine cylinder. Codes such as this have only appeared in the last two decades with the introduction of high speed computation devices. The first true two-dimensional code to analyze fluid motion in engines was written by Watkins in 1973 [32]. His code was a finite difference model that solved only mass and momentum conservation without combustion or compressibility effects. Geometry was restricted to a rectangular enclosure without the complexities of inlet ports or piston bowls. Piston motion was simulated by a transformation of the axial coordinate to a non-dimensional time independent form, which allowed the grid to expand and contract. The flow field was further simplified by neglecting turbulence and considering only laminar viscous flow. The collapsing grid to simulate piston motion was the most significant contribution of his analysis. Soon afterwards Gosman et al. [13] developed a similar code called Reciprocating Piston Motion ( RPM) based on the same geometry as Watkins. Gosman's code included valve geometry and used the k-E turbulence model. When compared to available engine data, their results were only qualitatively accurate. The most notable success was the code's ability to predict recirculating flow patterns seen in experiments when air was drawn through an open valve. He also was able to compare his results to data taken in a motored engine. In some cases the discrepancy between the predicted and measured values of air motion was considerable. The researchers felt that their biggest source of errors was the lack of knowledge of the boundary conditions and the possible inability of the turbulence model to account for larger turbulent structures. 11 Later, Gosman developed the TEACH code which is the base program upon which the analysis in this study is built. TEACH is a two dimensional finite volume code originally written to analyze steady-state fluid flow inside circular ducts. It uses the k-e turbulence model and assumes the fluid flow is incompressible. Later, Watkins [32] built on his previous work by incorporating into his piston motion code the k-E turbulence model and a rudimentary compressibility model. His results showed that flow behaviour was strongly dependent on initial conditions, such as at the valve inlet. A considerable portion of his work concerned the development of adequate boundary conditions at valve inlets. In the past decade a considerable number of engine simulation codes have appeared. , Spalding [23] created the PHOENICS code as a fluid flow simulator, and it was one of the first codes to appear with the finite volume discretization technique. Its usefulness as an engine modelling code arises from its ability to analyze multiphase flow, three-dimensional phenomena, and chemical kinetics. Cloutman [23] introduced the CONCHAS program in 1984. It was a finite difference code that possessed a simple combustion model based on chemical kinetics and a simple Arrhenius relation for flame speed. He simulated the combustion chamber of a spark ignition engine. His results showed that the numerical calculations could qualitatively predict the heat transfer and nitrogen monoxide characteristics, but prediction of other variables was less accurate. This was only accomplished after the combustion model was adjusted by the variation of several parameters. Cloutman found that the optimal parameter settings varied with equivalence ratio, volumetric efficiency, spark timing and even the choice of engine. He 12 showed a simple combustion model was incapable of the generality that was desired in a numerical combustion model. Some of the most advanced codes produced to analyze internal engine phenomena are the KIVA codes created by Amsden et al. [23]. The KIVA I code appeared in 1985 and was capable of computing three dimensional flow fields with chemical reaction. It also possessed a liquid fuel spray break up and evaporation model for use with diesel engine configurations. The KIVA II code replaces the explicit nature of the original KIVA with an implicit finite volume method. Liquid fuel jets are computed by a discrete particle technique, where a number of droplets are represented by one computation particle. Fuel droplet size and coalescence are determined statistically. The KIVA codes have gained popularity for modelling liquid diesel sprays in conventional diesel engines. Little work has been done on the numerical modelling of gaseous hydrocarbon injections into engine combustion chambers. However, one study by Gaillard [12] is worth noting. Gaillard simulated the injection of four methane jets directed radially into a two dimensional round combustion chamber using a semi-implicit fmite difference scheme. Initially, Gaillard compared computed values of transient jet penetration with known results of Kuo and Bracco. He found that the jet penetration and concentration contours agreed well with experiment. Gaillard's intention was to study the effects that swirling air motion inside the cylinder would have on jet development. His model showed the effects of swirl motion at various speeds in a two- dimensional plane. However, Gaillard did not study the effects of jet angle and penetration in an axial direction. Many numerical codes exist that are capable of modelling internal combustion engines. 13 Virtually all codes that have been written model diesel engines have spent considerable effort analyzing the break up of the liquid diesel fuel jet. Codes, such as KIVA II, have been used in numerous studies on liquid fuels in diesel engines. Very little has been done on the study of gaseous fuels injected in diesel engines. The work of Gaillard is the only study conducted in this area. Although his work is insightful on the effects of swirl on the transient jet, it does not consider the effects of jet angle, piston motion and pressure ratio. To answer these questions a two dimensional numerical model in radial and axial coordinates has been constructed and is presented in this thesis. 1.4 THESIS OBJECTIVES AND STRUCTURE Numerical analysis was used to investigate the injection of natural gas fuel into diesel engines for several reasons. Little experimental or numerical data exists on the behaviour of such jets, making optimization of the dual fuel diesel injection system difficult. Apart from being less expensive and time consuming than experiment, numerical analysis may be able to provide insight into areas where experiment cannot. Piston motion may have a strong influence on the jet behaviour, which is easily investigated numerically. The original TEACH I code was written to calculate steady-state incompressible pipe flow. It was selected for its simplicity and flexibility, and its reputation as a proven universal flow analysis program. The code must be modified to perform transient computations through the introduction of time-dependent source terms. Piston motion is simulated through a transformation of the axial coordinate to accommodate a compressing and expanding grid. The accumulation of mass inside the closed cylinder as gas is injected is accounted for 14 through a correction to the density field. Two models are used in this thesis: a fixed piston model will have a stationary lower boundary that can be set at various positions; and a moving piston model, with a moving lower boundary that can be made to move at various speeds. In each model the jet angle and pressure ratio may be altered. The objectives of this thesis are stated now in the form of questions to be addressed. 1). Can the TEACH code, with adjustment of the initial conditions, simulate a two dimensional transient jet ? 2). Does the jet angle influence the jet penetration and development in both moving and fixed piston models, and if so in what way ? 3). Is there a change in jet penetration when a higher injector tip pressure ratio is used ? 4). Does the chamber size influence the jet development ? 5). How does engine speed influence the jet development ? 6). Is there a difference in the jet character between the fixed and moving piston models ? The chapters in this thesis are outlined. Chapter two goes over the TEACH code, which was modified to perform this analysis. The governing equations and their discretization are all described. Chapter three describes the assumptions used in developing the models and their configuration. A derivation of the general and local compressibility model as well as an introduction to the moving boundary transformation is presented. Chapter four compares the free jet calculations to known analytical and experimental 15 results. The comparison will validate the use of the TEACH code for analyzing the turbulent jet under consideration. Chapter five gives the results for the analysis of the transient jet erupting into a fixed, closed environment. This study will identify the parameters that appear to affect the jet behaviour the most. Jet angle, pressure ratio, and clearance will be considered. Chapter six gives the results of the moving boundary analysis. The effects of jet angle, pressure ratio and engine speed will be investigated. Chapter seven draws conclusions from the results, and recommends a course for further study. 16 CHAPTER 2 : MATHEMATICAL FORMULATION 2.1 GENERAL The motion of fluid inside the combustion chambers of diesel engines can best be described with the aid of the mass conservation and Navier-Stokes equations. The NavierStokes equations express the conservation throughout the computational domain of the fluid momentum. The TEACH fluid flow analysis code , which was used in this investigation, solves the mass and momentum conservation relations. The second section ( 2.2 ) of this chapter reviews the basic governing equations used in the TEACH code. The effect of turbulent fluctuations is introduced through a redefinition of the velocity field into mean and fluctuating components in a manner proposed by Reynolds ( 1895 ). The development of the k-c turbulence model using the Reynolds form of the Navier-Stokes equations is outlined in the next section of the chapter ( 2.3 ). Reduction of the generalized governing equation through discretization to a linear algebraic form is reviewed in section 2.4, along with a description of the hybrid differencing scheme used in the process. The SIMPLE algorithm used to solve the momentum equations is briefly reviewed in section 2.5. 2.2 GOVERNING EQUATIONS The basic equations governing fluid motion are the mass conservation and the NavierStokes equations of momentum conservation, derived in detail in Schlichting [291. The conservation of momentum provides one equation each for the radial and axial directions. The mass conservation equation provides a way to determine the pressure field, which will be shown in a later section. The mass and momentum conservation equations in these ^ 17 coordinates are shown below. The variables u and v are the axial and radial velocity components respectively, z and r are the axial and radial coordinates respectively, t is time, p is pressure , p is density, and F r and Fz are body force terms, in the radial and axial directions per unit volume respectively. continuity: ^ 4-Ea (Pry) + az (Pu) =° (2.1) momentum conservation: (2.2) radial direction:^ ay^,,av‘^ap., 1 a ' r ay\ y .321 ar r^ar • 7r ark tlat • var • az) ' r 2 az2 ^n iav . 4. - - - (2.3) axial direction:^ pl a2 u) pi au^au^u^F^ap^( 1 a (..au) _ r ar ar^a z 2 ^at^ar^az) -^- az •^ L -r These equations imply constant viscosity p. For variable p and allowing the symbol 4) to stand for either the u or v velocity components, eqs. (2.2) and (2.3) can be shown to be special cases of the following general equation. 18 General Momentum Equation: ^ (2.4) a as ^s vti, _^(^) a (4) p ) — aaz ( Pu. ) +^ (rP^az 11—E^r ar^ar at The general momentum equation is valid for laminar and turbulent flows. The dependent variables of the general momentum equation are instantaneous values, meaning they are defined at every point in time and can fluctuate with turbulent motion. The dependent variables can be separated into mean ( indicated by an overline) and fluctuating 4)' components, (2.5) = The mean values are determined through phase averaging, described by Reynolds (1980) [27]. The cycle-to-cycle variations in the cylinder of the moving piston model are removed by averaging over a great number of samples of a dependent variable taken at a specific point inside the engine and in its cycle, shown in equation (2.6). The engine rotates with a period of To . 1°+1) r, t)^Lim I 1 \ 10 (z r t 1-nTo ) (2.6) If the fluctuating and mean components of the dependent variables are inserted into the general momentum equation (2.4), the following equation is obtained by averaging. The overlines having been omitted on the mean values. 19 at (4) a^a + az (P14) —r— ar (rPvC — 2.7 = 1( 11 1 -^4-1-(r11:4- rP7V) S• The averaged product of the fluctuating components of the velocities (for example pu'v'), is the Reynolds stress term, which represents the diffusion of momentum by turbulent fluctuations. 2.3 k - e TURBULENCE MODEL The presence of the Reynolds stresses complicates the solution of the general equation (2.7) by introducing further variables into the equations. With the Reynolds stress, the number of dependent variables will exceed the number of equations that can be used to solve them. The Reynolds stress terms must them be modelled in terms of the other dependent variables in order to invoke closure. Boussinesq, in 1877, developed the concept of an eddy viscosity to approximate the Reynolds stress terms. He suggested that the turbulent shear stress ( Reynolds stress ) is equivalent to the product of the mean strain rate and the eddy viscosity ( PT ). -puivl = µT ay^ (2.8) The turbulence model used in the TEACH code is the two equation k-E model. The eddy viscosity is derived from the turbulence kinetic energy k, and the turbulent energy dissipation rate e. Both turbulent kinetic energy and turbulent energy dissipation rate are determined through their own transport equation. The governing equations for k and e and the derivation 20 of the eddy viscosity from these parameters are given by Launder and Spalding (1972) [22]. The experimentally determined constant C. is 0.09. Co p k 2 e (2.9) 14T ^ The effects of turbulence are introduced into the momentum equations through the effective viscosity (p eff). The eddy viscosity (pT) is added to the laminar viscosity (p i.) to produce the effective viscosity. I 1 eff = 14 laminar + µ T ^ (2.10) For variables other than the velocity components, a fluctuating component term similar to the Boussinesq approximation is used. -pu - au Iv =rii0 — ay . IT 7 (2.11) The diffusivity term (F 0 ) is determined from the eddy viscosity using the turbulent Prandtl number (a,) for the specific variable. The Prandtl number for diffusion of gaseous species through gaseous host media has been estimated by Launder and Spalding to be of order 1. r4 , _ _ "L eff a0 (2.12) With the equations for k and e defined by Launder and Spalding (1972), the complete generalized governing equation for the TEACH code can be shown to be: 21 a ^ ir as\ + (rp v4, ) =^ at^(4)^(Pu4))^ ar^az variable^r4,^ V a az 0 aaz ( ileff : zu)^aan(riteffaalzr)^ap z ^ eff rr ar)^c, So 1 (mass)^0^ Neff ^i az/^--z." ar‘^arj^-4) U N eff T: a„^ay^,,,^ap r ar^ ar ( r " eff-- 4. " ) k^II eft/ ak^G - pe Neff/ k (C 1 G - C2pe) a t^ C (mass frac.) ^peffia.^ 0 au■2 ^au + G = Per42RT;)^( T:r1 ar ) J^ar ty 7,1 2( 2] The k-e model constants are: CP 0.09 C1 C2 ak K 1.44 1.92 1.0 0.4187 where^a K2 — (C2 - C1 )^/ 2 r ar 2 — eff— 22 2.4 Discretization of Governing Equations The discretization method used by the TEACH code is the control volume formulation. The region of interest is divided into a network of control volumes. Inside each control volume is a nodal point where the values of the dependent variables are to be evaluated. The velocities u and v are computed in control volumes that are staggered in relation to the scalar variables control volumes (C,k,e, and p). The staggered grid prevents the possibility of an unrealistic solution to the pressure field ( Patankar 1980 ). Standard Control Volume (c.v.) • N h- AZ r______ _ ar n _ ,,/A 7 2' E e , Ar 8 1 _,/ __I 8xw axe • r , 7-- S scalar c.v. Figure 2.4.1: standard control volume configuration A typical control volume arrangement is shown in figure 2.4.1. The scalar control volume is shown by the shaded region. The local scalar nodal point is indicated by P, while its adjoining neighbours are indicated by N,S,E,W, referring to northern, southern, eastern, and 23 western locations . The faces of the scalar control volume are indicated by the lower case letters; n,s,e,w. The velocity control volumes have nodal points at locations indicated by the solid arrows in figure 2.4.1. The staggered grid allows the velocity nodal points to be located on the faces of the scalar control volumes, simplifying the calculation of convection through these faces. The integration of the governing equation over time and the control volume is shown in equation (2.13). t+At e a ^vit) 12 — 4)) + 1—a (rp) f^f^ (.13)^az^r ar r dr dz dt t w s t+At e 11 = f f wf (1(r.2) t s .13) 1 — rr — 8 +^rdrdzdt r ark^ar) The first term of the integrated equation can be evaluated as: t+At n e f f at (p4))) dz rdr dt rArAz( (p4)) 1 , - (p4))%) tsw j 2.14 The superscripts denote points in time; 1 indicates the present time, and o indicates the previously calculated point in time. The nodal point value is assumed to prevail over the entire control volume, allowing the integral to be evaluated. The other terms are not as easily integrated over time. Some assumption must be made about the variation of the integrand over time. A common practice is to assume that the result of the integration is equivalent to various contributions of the integrand at the start and at the end of the time step, shown in equation 2.15: ^ 24 f^dt = LEC + ( 1 - f)^A t^ t+At (2.15) where is any dependent variable (u,v,k, e,p,and C).The weighting factor f varies between 0 and 1. If f is chosen to be 0 the discretization becomes fully explicit, and the linear equations can be solved directly from the values at previous time step. If f is chosen to be 1, then the discretization is implicit, and the resulting equations must be solved simultaneously using an appropriate solver. Although the explicit scheme is computationally easier, it becomes unstable if the time step exceeds the limits determined by the Courant condition [3]. For this reason the TEACH code uses the implicit scheme. The remaining integrals can be computed as shown in equation 2.16, again assuming that the integrand is constant over the control volume and omitting the 1 superscripts for present time step. _ f ftf (p u, - r^r a a t+At ^4 en ,^ zdt t w s = {(pu. - P 134 • az ') — e - ( pub - Pe to-aZ) w (2.16) irArAt The integration of the source term can be written in terms of a constant component S c , and a component linearly dependent on the variable, Op, shown in equation (2.17). The discretized governing equations must be linear if they are to be solved easily. Therefore, the discretized source term must be no more complicated than a linear approximation. t+At en f f Sa rdrdzdt =^+ Sptt pirArA t ft w s (2.17) Once the governing equation is integrated, the results must be expressed in terms of the 25 known nodal values. The TEACH code uses the hybrid differencing scheme to express the convective and diffusive fluxes in such terms. Equation 2.16 shows that the fluxes through a scalar control volume face can be written as a convective ( C ) and diffusive ( D ) component. As an example the fluxes at the east face (C e and D e)are shown. Ce = (p e u e. e ) e rAr,^De = (Iltz ) e rAr..^(2.18) The hybrid scheme is a combination of the central and upwind differencing schemes. The central differencing scheme assumes that a linear variation of the dependent variables exists between neighbouring nodes. The convective and diffusive fluxes can be written for the eastern face of a control volume: E ^4 p^4)p - 4 E) • r Ce = (pu) 4 I 2 )rAr,^De = r4 8z. ra . ( ) ) ) (2.19) The central differencing scheme is accurate when diffusion is the important transport mechanism. In flows where convection is dominant and the Peclet number ( Pe = I puSxir I) is greater than 2, the central differencing scheme will become unstable. The upwind differencing scheme has been combined with the central differencing scheme to alleviate the high Peclet number instability of the central differencing scheme. In strongly convective flows the upwind scheme approximates the convective terms with values of 4) from the upwind location. An example is given in equation 2.20. Ce = p e u e4) prAr^for u e > 0 C, = p e u e(b ErAr^for u, < 0 (2.20) 26 The hybrid differencing scheme uses the superior performance of central differencing for low Peclet number flows and upwind differencing for high Peclet number flows. Diffusive terms are always approximated with central differencing. Convective terms use central differencing when IPe 15.2 and upwinding for IPe The complete linearized governing equation for a dependent variable (I) is given in equation 2.21. a ptt p = adO E + a wO w + a NON + a sO s + b + Scz- ArAz a p = a s + a w + a N + a s + a p ° - SprArAz rArAz a p ° - ^ Pp0 At (2.21) b = ap 0 41p A typical coefficient of the hybrid scheme for the east face is given: a E = MAX (^, .D 2 ^ ^2 . (2.22) 2.5 The SIMPLE Algorithm SIMPLE stands for the Semi Implicit Method for Pressure Linked Equations, developed by Patankar and Spalding (1972). In brief, this is an algorithm to solve the momentum equations through an iterative procedure. Pressure and velocity correction equations are developed from the continuity and momentum equations. Complete details of the algorithm 27 are given by Anderson et al [3] and Patankar [26]. As an example of the derivation of the SIMPLE algorithm the axial momentum equation is considered: aeu e = Ea bub + b + (P p -P E )A., ^(2.23) where Ae is the area of the eastern face of a scalar control volume, subscript b is the neighbouring coefficients similar to equation 2.21, and b is the source term. The velocities have subscript e since they are evaluated in the velocity control volume on the east face of a scalar control volume ( see figure 2.4.1 ). To solve the momentum conservation equations the pressure field must be known. If a pressure field is guessed ( guess value denoted by P`) and the momentum equations solved with this incorrect pressure field, then the resulting velocities ( u* and v* ) will not satisfy the continuity equation. The momentum equation based on the guessed pressure field can be written as: a e u* = Ea b u* b + b + (P* p -P* E )A e ^ . (2.24) It is assumed that the velocities and pressures appearing in the discretized axial momentum equation (2.23) are each a sum of a guessed value (0 * ), and a correction (4'), where 4) is either a velocity component or the pressure ( = 4)* + 4' ). If a guessed velocity equation ( 2.24) is subtracted from (2.23) then an equation of velocity and pressure correction terms is obtained (2.25). aeui e = Ea bui b + ( P i p - P I E ) A, .^(2.25) 28 The summation term, la b db , can be dropped from equation (2.25). The path taken to reach a solution is irrelevant so long as a solution is achieved that satisfies momentum and mass conservation. Neglecting this term in effect only alters the route taken in the iterative process, and has no influence on the final solution. Equation (2.25) can then be manipulated to give the velocity correction ue. If this is inserted into the original velocity definition, Ue + ue, then the velocity correction equation (2.26) is obtained. Ue = ^ + cl,(P 1 p - PC)^where d e = a . ^(2.26) A velocity correction equation for v can be obtained by similar methods. The velocity correction equations can be substituted into the discretized continuity equation shown in equation (2.27). (p p -pe p ) rArAz ^ +[(pu).- (pu),,,jrAr +[(pv)^(pv) sjAz = o 2.27 At With some rearrangement of terms the pressure correction equation (2.28) is obtained. a pP i p = a EP I E + a wP I w + a NP 1 N a p = aE + a w + + a sP I s +b aN a E = p ederAr b - (Pp ° - Pp) At rArAz PPu*),,, - ( Pu s ) eirAr + [(pv*) .9 (2.28) - ( pv*),21Az The b term in this relation is the discretized continuity equation of 2.27 with guessed 29 velocity values. The SIMPLE algorithm iteratively guesses and corrects the velocity and pressure fields until convergence is satisfied. Initially, a pressure field is guessed, P* , which is then used to compute the guessed velocity fields with the aid of equations such as 2.24. With the values for guessed velocities, the pressure correction equation 2.28 is solved and the results are added to the guessed pressures to obtain the pressure field ( P=P*+P' ). The velocity field is then computed using equations for axial and radial velocity equivalent to (2.26). Once the updated pressure and velocity fields are computed, relations governing other dependent variables can then be solved. Finally, the newly computed pressure field is used as a guess for the next iteration and the entire procedure is repeated. 30 CHAPTER 3: JET SIMULATION 3.1 General In the previous chapter the governing equations were used in a derivation of fundamental relations use in the TEACH code. This chapter reviews the changes made to the TEACH code and the boundary conditions necessary for in-cylinder engine analysis. The second section of this chapter ( 3.2 ), reviews the basic configuration of the gas injector in the cylinder and how it is modelled. The assumptions upon which the models are based are outlined. The third section ( 3.3 ), shows how the control volume grid is made to contract and expand, simulating the piston motion. Section 3.4 shows how the mass injected into the closed cylinder is accounted for. The final section ( 3.5 ), looks at boundary conditions used in the code. 3.2 CONFIGURATION AND ASSUMPTIONS Four figures will be used in this section to show the location of the jet inside the engine and the manner in which the jet is modelled. Figure 3.2.1 gives a schematic drawing of the injector inside the piston-cylinder assembly of the engine. The injector is centrally located in the cylinder head between the two exhaust valves, and protrudes approximately 1.6 millimetres into the combustion chamber. The jet is an axially symmetric continuous conical sheet emanating from the injector tip inside the cylinder. The moving piston is directly underneath the jet. A detailed view of the gas jet is seen in figure 3.2.2. This figure shows a cut view of the tip portion of the injector. The axial coordinate is on the centerline of the injector with the radial coordinate perpendicular to it. The tip is a poppet valve which moves axially 31 Figure 3.2.1: Conical gas jet in relation to piston and cylinder creating a pathway for the gas. The gap created by the poppet is 0.3 millimetres. The jet is inclined to the radial axis at an angle determined by the injector tip. The dimensions used in the fixed and moving piston models are the same as those in the actual engine combustion chamber. The combustion chamber can be divided into a series of two-dimensional control volumes or cells in the axial and radial plane. Figure 3.2.3 shows the computation grid with respect to the injector, piston and cylinder. The number of control volumes shown does not correspond to the number used in the models. The computation grid is densest near the injector exit since at this point the gradients of all dependent variables are greatest. A higher 32 Figure 3.2.2: Details of the conical gas jet. grid density will better be able to resolve the high gradients at the injector tip. The symmetry axis is at the extreme left edge of the grid at the centerline. A detailed diagram of a portion of the computation grid and its features is given in figure 3.2.4. The computation grid includes 29 cells in the axial direction and 35 in the radial direction, including the boundary and injector cells. Shaded control volumes along the outer edge, are boundary control volumes and do not take part in the analysis. Their purpose is to allow the imposition of boundary conditions which may change from side to side. In the upper left corner of the grid, a group of cells 5 by 9 is segregated by a heavy line. This group represents the injector and has in each cell a value of 0 for all dependent variables. A solid object inside a computational domain is often represented in this fashion. One control volume 33 GRIDDED COMBUSTION CHAMBER INJECTOR COMPUTATION GRID 1 ■ ■■^\ 1.. ■ V/ZI/4111%/-14IZA il 1 V.//■ ',4///:.r/Adw/p///41////4 )1 2.125 in. ^ (N^ SYMMETRY AXIS Figure 3.2.3: Computation grid inside combustion chamber in this region is used to represent the injector outlet. The values of all dependent variables in this control volume are specified; in particular the velocities are set at specific values that will vary with time according to an injection velocity profile. Details of the injection velocity profile will be given in a subsequent section. A simplifying assumption in these models is the axial symmetry. The axial symmetry of the jet allows the use of a two-dimensional grid in the axial-radial plane. The twodimensional grid reduces the computational effort required to analyze the model relative to a three-dimensional grid. Another notable assumption is that of incompressibility. Compressibility effects are present in high speed flows of gasses such as air and natural gas. Though the flow at the injector 34 Computation Grid and Injector tip Region shaded regions are WA r 3 nodes 2 4 nodes i____, 26 nodes boundary cells ./'^I. . A%^ 9 nodes ./A a. Figure 3.2.4: Detailed computation grid exit will be sonic or supersonic, the velocity drop to "low speed" values occurs within a distance of roughly 10% of the cylinder radius ( as described in section 1.2 ). This means that if a suitable equivalence is established for the injector flow the bulk of the fluid may be treated as though it were incompressible. One other significant assumption is that the initial velocities and turbulence quantities are zero. In reality this is not exactly true, at any time inside the cylinder there will be motion, a pressure field and some turbulence. However, the initial fields may be small compared to the fields created by the jet. The effects of the initial fields on the subsequent results may then be insignificant. 35 3.3 Moving Boundary Modification The motion of the piston is simulated in the computation domain through the movement of the lower boundary of the computation grid shown in figure 3.3.1 ( see Watkins 1977) . This face is considered impermeable and is specified to move in a manner consistent with that of a piston inside a cylinder. The computation grid will compress and expand uniformly in the axial direction to simulate the piston motion. This is achieved through the transformation of the axial dimension z, into a nondimensional counterpart , which is specified as: E = Z^ ZH (3.1) The variable zH , is defined as the instantaneous piston location, or the distance from the lower axial boundary to the tip of the injector. The advantage of this transformation is that the values of the new axial variable 4, defined at each node, do not change with piston position. For instance, if each control volume in the compressing region contracted a uniform amount in the axial direction, then the ratio of the new axial node location z to the axial piston location z H will not change. The non-moving injector tip will have z=0 and 4=0, and the piston face will have z=zH and 4=1. The original TEACH code with its non-moving grid formulation can then be used for a moving boundary model if 4 is used instead of z. The moving piston model is divided into compressing and non-compressing grid regions as shown in figure 3.3.1. The two regions are divided by the z = 0 line. The non-compressing grid region is used to prevent the injector tip dimensions from changing as the piston moves. To implement the transformation, the original governing equation must be transformed in 36 radial edge Z=0 ■ Non-Compressing Grid =0 axial Compressing Grid edge boundary moving at piston speed: azH at Figure 3.3.1 Details of compressing and non-compressing grid regions. terms of the new coordinate system. This is achieved by first realizing that the value of any dependent variable (0), must be the same in both systems. 4)(z,r ,t) = 4) 1 ( , r,t) ^ (3.2) The prime represents the variable in the new coordinate system, t is the time, z and r are axial and radial coordinates. The derivatives of these variables must also be equivalent; equation 3.2 is differentiated: (3.3) 37 Observing that the piston location z H is a function of time, the derivative of can be expanded in the following manner: zH = za(t)^= 4zH = (z,t) ( ) (3.4) at dz at,k az^at dt = If this expression is introduced into the right side of 3.3 , then it becomes; a4d r^-`1t) a4dt ar^at az^at^at or (3.5) at' dr +^a dz + adv + atldt ar^at az^at at^at The derivative of the new axial variable can be shown to be, at =^at , ak dzH = _ dzH az^zH ' at^azH dt^ZH ck • (3.6) Introducing equations 3.6 into 3.5 the derivative identities that will, allow the original governing equation to be transformed are: 4 ar = _ dzif 4 1 + aiv^att) = 1 4'^ = ar^az^ at^zH dt at^at at^ • (3.7) The derivative identities shown in 3.7 can be substituted into the original governing equation shown in the previous chapter. The resulting general governing equation of all dependent variables with a moving boundary is: 38 a^ao4))^a 1a (4)p) ^--(pu4)) (rPvC at^zH a ds^at^r ar = a (^–1 a rr — 4) s . zH at zH at^r ar^3r (3.8) 1 The new governing equation contains an additional term compared to its predecessor (equation 2.4). Discretization of this equation can be simplified by defining a new axial velocity: a = u - E^. at (3.9) The new velocity is simply the speed of the fluid relative to the local grid (iI), which is equivalent to the absolute fluid velocity u, less the local speed of the grid ((az H/at)). The new governing equation can be then written in the same form as the original governing equation (2.4): 1 43 18 1a (p4)z,i)^--(pu4))^–r — ar (rPv.) zH at^zH 13E (3.10) r^+^(rr -aA) s = 1 a (A at^r^är^• In a similar fashion the continuity equation can be written in the same format. a^a ar (rpv) = 0 —(po –— 11(pz,i)^ zH at^ r zH (3.11) The new governing equation is discretized after having been multiplied by rz H . The result 39 will be a series of linearized equations similar to the linear equations obtained through the discretization of the non-moving boundary governing equation. Terms such as Az will be replaced by zHA4. The time derivative terms will contain z H at present and past time steps, shown in equations (3.12) and (3.13). For the continuity equation the time derivative term becomes: e n t+At fif W s t - Vz,e) rpArAt . A t at m r— (pz )dtdrdt - (3.12) Here rP is the cylinder radius at the local node p. The density p is the local node value shown in the next section. The superscripts indicated the present point in time (') and the previous point in time (°). For the governing equation the integrated time derivative is: e n t+At A)dtdrelt - (Piz 1+1 Vzm°41?) r ArAt . 111 r—(pz At at z p i (3.13) W S t To construct some of the coefficients of the linear equations resulting from the discretization, the values of the velocity field with respect to the grid (II) will be required. The relative velocity, fi ,should be used for the axial velocity instead of the absolute velocity u. It would then be possible to compute fi directly. To obtain the governing equations for relative velocity the following substitution is made, u= a + dz^ dt (3.14) 40 which is placed in the absolute velocity equation ( from eqn. 2.21 with •:0 as u): A pu = EA nun + A p°u° + S ,^(3.15) to yield, dZH^A A pti +^eTt. -^Vinun dZH + A n &--=') + APIs° +^+ S dt (3.16) p0ZH° - ^r ArA . At P In its final form the linearized axial momentum equation is then: A pt2 = E A nt2 + A p '14° + s, s (EAn. _ Aptp)dzH dt Ap cl itti (3.17) 3.4 DENSITY CHANGE DUE TO PISTON MOTION AND GAS INJECTION The computational domain is a closed system with mass entering through the injector. As mass enters and as the grid compresses and expands, the density and pressure of the gas inside the cylinder will change. The density and pressure fields can be corrected on a global level through an application of the mass conservation relation. A correction term for density and pressure is determined which is then applied to all control volumes inside the domain. The density (p) is defined in terms of a guessed (p*) and correction value (p'). 41 (3.18) The pressure (P) can be corrected by use of the correction term for the density. P =P' + p , ap ap (3.19) The density correction is derived from the summation of the mass conservation equation over each control volume in the computation domain; equation (3.20). (z„pp - zif°pZ) At a^ mmwg. rpArAt - mi =0 (3.20) The present time superscripts (') have been omitted. The term m o is the mass inflow through the injector and is considered to be constant over one time step. Inserting the redefined density (3.18) into (3.20) and rearranging terms the density correction is obtained. P^ La. (zHpp zi;pp) - rpArAt - At vol. zipArA L „„, i ,, .^ (3.21) At This term is used to compute the pressure correction in equation (3.19). The pressure gradient term of equation (3.19) is evaluated assuming isentropic compression and expansion. The pressure gradient can be shown to be: ap = k P .^ aP^P (3.22a) The variable k is the ratio of specific heats. The pressure and density corrections are added to all control volumes during each iteration. These corrections drop to zero as the code 42 reaches convergence. The density is adjusted for the presence of natural gas using the ideal gas law shown in P — P (MWas,g RT ) (3.22b) equation (3.22b). The term MW e„g is the average molecular weight of the gas at the local node ( computed from the gas mass fraction), R is the universal gas constant ( 8.314 L kPa/(K mol)), T is the temperature assuming isentropic compression. 3.5 BOUNDARY CONDITIONS The boundary conditions of this model are the jet inlet velocity and density, and the conditions at the symmetry axis and the solid walls. Figure 3.5.1 shows how the injection speed will initially vary with time. Initially the time step size is 0.1 milliseconds, growing at a rate of 20% per time step. As the poppet valve opens the jet velocity will start at 0 and will increase rapidly until Mach 1 is attained. Mach 1 is the speed of sound in natural gas at the jet orifice. After Mach 1 is reached the jet speed will remain constant. It is assumed that the poppet valve lift is directly proportional to the speed of the jet. This is done to provide a "ramped" starting velocity since impulse staring may prevent convergence of the model. The pressure ratio between the upstream gas tank (P e ) and the engine cylinder (P e) will dictate the exit speed and density at the nozzle. The injector passageway behaves as a converging nozzle, accelerating the gas. If the pressure ratio is above the critical value (above 1.8), the jet will be underexpanded. The pressure at the exit plane P e , will then be higher than 43 Figure 3.5.1: Injection speed as a function of time the cylinder pressure. As the underexpanded jet exits, its diameter increases as it expands to the engine cylinder pressure. Figure 3.5.2 shows the characteristic jet behaviour due to underexpansion. Experiments on underexpanded jets have been conducted by Ewan and Moodie 1986 [11], who show that the underexpanded jet is equivalent to a correctly expanded one emerging from a nozzle of different diameter D om . The equivalent diameter is defined as the diameter of the jet after it has expanded to cylinder pressure and acquired an assumed Mach number of 1. Ewan and Moodie observed that the adjustment region was only a few jet diameters in length so they neglected entrainment during the expansion process, and assumed constant mass flow rate in this rapidly expanding jet. Their experiments showed that with a pressure 44 Figure 3.5.2: Underexpanded jet details. ratio of 5 to 1 the length of the expansion region will only be about 1.3 times the nozzle diameter. With reference to figure 3.5.2, mass conservation implies: PAU. = p egi eq U,4^ (3.22c) Here, subscript e indicates the exit plane and eq indicates the end of the expansion region. Using the ideal gas law and assuming that the temperature at the exit plane is equal to the temperature at the end of the expansion region and that M eg is 1, it can be shown that: p e lf,^P, p eq U 07 P.^A, (3.22d) 45 The equivalent diameter (D eq) is then: (3.23a) Day = The term D 1 , is the injector outlet diameter. The equivalent diameter calculated from eq. 3.23a was used as the initial diameter for the jet calculation. Given that the mathematical flow simulation employed for calculating the jet behaviour implied incompressible flow, it could not exactly represent the effects of Mach number on density and temperature. Hence it was decided to use a starting value of the jet density which was more representative of the downstream conditions than the conditions at the end of the sudden expansion. With this in mind the starting value of the density was taken to be: P Pa MWg R Ta (3.23b) where MWg is the molecular weight of the gas, and P a and T. are the cylinder pressure and temperature. An approximation was required also for the starting value of the jet velocity. Consistent with the idea ( shown to be reasonable by Ewan and Moodie ) that the jet Mach number just downstream of the sudden expansion is nearly sonic, the steady initial value of the jet velocity was taken to be: k R Ta U — ^ eg^MW (3.23c) in which the temperature representative of downstream conditions has been used to compute 46 . the sonic speed. In addition to exit conditions at the nozzle, conditions at the other boundaries must also be specified. There are two boundary conditions of interest on the edge of the computation domain; solid wall and symmetry axis conditions. Figure 3.5.3 shows typical control volumes across a wall and a symmetry boundary. Figure 3.5.3: Boundary control volumes The symmetry axis is the simplest boundary in the model. Since values of all dependent variables are equal on either side, no fluxes take place. Referring to figure 3.5.3, radial velocity V. will be zero. The nodal index (I,J) is respectively the axial and radial position sequentially in the grid. For example site (2,2) is the top left node inside the injector ( referring to figure 3.2.4). To eliminate momentum transfer across the boundary, U. is set 47 equal to Ub. For variables such as gas concentration and turbulence parameters the linearized equation coefficients for example a s at position (1,2), are set equal to zero. The wall boundary condition is more complex than the symmetry axis. The wall is considered impermeable, therefore velocities that penetrate the boundary such as U„ referring to figure 3.5.3, are set equal to zero. Velocities that are parallel to the wall such as V c , will be influenced by the wall. A logarithmic velocity profile for wall boundary layers, described by Launder and Spalding (1972) [22], is used to introduce wall effects. The logarithmic wall velocity (u w) profile for turbulent flow is defined as: uw = --1 1n(Ey . The term u„ is the friction velocity and K ^ (3.25) and E are experimentally determined constants defined as 0.4187 and 9.783 respectively. The friction velocity is a function of the wall shear stress 'T y,: 1 UT = t w)i (3.26) The non dimensional distance from the wall y + , is a function of the friction velocity, the laminar viscosity, and the normal distance from the wall y p . Y^ P Yp u,^ (3.27) Wall effects are introduced into velocities such as V c in figure 3.5.3 through the wall shear stresses. The wall shear stress is computed using the logarithmic law of the wall. Coefficients such as aw, at position (2,J) are set to zero. The wall shear stress is then 48 ^1 ^1 pc 4 k 2 Ku w T w -^it^ ln(Ey+) (3.28) introduced into the source term for the velocities such as V c . The value of c p , is constant and commonly accepted to be 0.09. Values near the wall for turbulent kinetic energy (k) and dissipation (e) are introduced in a manner similar to the velocities. The turbulence parameters near the walls are defined as: 33 c 4 k 2^U P. ' e -^k = 11 2 7 KYp^ I c 2 W (3.29) 49 CHAPTER 4: ROUND FREE JET 4.1 General A computer model of an injection of gas into a combustion chamber is of little value if it is not known whether the phenomena is being simulated reasonably. The validity of a code can be determined through a comparison of its results to experimental measurements of the flow being numerically simulated. However, it is not always feasible or practical to construct an experiment to obtain data to verify a numerical model. Data are available on the subject of the transient and steady-state incompressible round free jets. Even though the round free jet are not exactly the same as the confined jet under study, data on such jets can be used to check the performance of the code. A numerical model of a transient jet inside a large closed chamber is used to compare the performance of the TEACH code to known results ( see appendix A ). Briefly, the jet is formed by two control volumes in the wall of the chamber, giving a jet radius of 2 cm. The wall opposite the jet is 1.2 m from the orifice and the adjacent wall is 0.8 m from the centerline of the jet. Velocities in some areas of the chamber may indicate laminar flow conditions, invalidating the k-e turbulence model in these regions. Profiles are obtained in regions near the orifice of the jet where the flow is expected to be turbulent. Having the jet orifice formed by two control volumes limits the possible velocity profiles at the exit plane to uniform flow. More control volumes would need to be added to simulate a non-uniform velocity profile. However, this is not always possible since very small control volumes may cause numerical stability. 50 It must first be shown (section 4.2) that the walls which enclose the transient jet in the numerical model are far enough from the origin so that the jet will behave, for a limited time, as though it were in a constant pressure environment. This is crucial if steady-state jet theory and results are to be used. The next section (4.3), will show the behaviour of the code as the time step size and grid size is changed. The final section ( 4.4 ), will compare both axial and radial profile data to results produced by the TEACH code. Details of the comparisons and of the grid configuration of these tests are given in appendix A. 4.2 STEADY-STATE JET Before a comparison can be made to steady-state data, it must be shown that the TEACH code can produce such results. In spite of the fact that solid wall boundary conditions were imposed to enclose the jet, the walls were made far enough from the jet so that their influence was not immediate. Figure 4.2.1 shows the velocity on the jet axis as a function of distance from the orifice. A uniform inlet velocity profile of 10 m/s was imposed at time zero. A time step size of 0.02 seconds was used. It is clear that the jet does reach a steady-state; within 0.8 seconds the velocity profile along the jet axis is steady. As expected, the region nearest the jet inlet reaches the steady-state most rapidly. In time the presence of the wall will induce a recirculating flow pattern and invalidate the comparison to free jet data. Radial velocity and concentration profiles taken at a distance of 0.435 m and 0.46 m respectively from the jet orifice are shown at various points in time in figures 4.2.2 and 4.2.3. Both these figures show very much the same behaviour. The profiles for velocity and 51 Figure 4.2.1: Axial velocity along jet axis at various times; 10 m/s initial jet speed concentration reach a steady-state very rapidly after the jet starts. After roughly 0.2 seconds the velocity and concentration fields have stabilized, again establishing the validity of this model for comparison to steady-state data. 4.3 GRID SIZE AND TIME STEP DEPENDENCE If the grid size and time-step size are chosen too large, the results of the analysis may contain considerable error. In a finite-difference discretization it is possible to obtain an estimate of the order of magnitude of the error through truncation terms [3]. Commonly the error would be found to be of the order of Ax or (Ax) 2 , depending on the term, 'x' being any dependent variable. The control volume formulation will yield errors similar to that of the 52 RADIAL CONCENTRATION PROFILES 0.6 0.55 2 6P Time 0.5 • 0.22 sec 0.45 NE 0.18 sec 0.4 0 0.14 sec 0.35 ^ 0.3 A 0.06 sec 0.25 2 0.10 sec 0.2 0.15 0.1 0.05 0 0.' 5 0.'^ 0.^O.^0. 5 0.'^0.105 0.12 RADIAL DISTANCE FROM JET AXIS (M) Figure 4.2.2: Radial concentration profiles at various times; 0.46 m from jet orifice; 10 m/s initial jet speed. finite-difference procedure. Clearly then a reduction of the time-step size grid spacing should reduce the error of the results. It is common practice in numerical analysis to run a model with various sized grids and time steps to determine the sensitivity of the results. Figures 4.3.1 and 4.3.2 show the effects of a reduction of the time step size on the radial profiles of steady-state velocity and concentration. The velocity and concentration profiles appear very similar in nature and exhibit the same behaviour with changes in time step size. It is clear that only the largest time step size, 0.02 seconds, produces any deviation in velocity and concentration profiles. Other time step sizes appear to yield similar results as their curves tend to be equivalent. This suggests that the code tends to be insensitive to time step size as long as the time steps are kept small. 53 RADIAL VELOCITY PROFILES 65.5 5- Time • 0.22 sec 4.543.53.5 a CIE 0.18 sec O 0.14 sec 2.5 2 1.5 1 ^ 0.10 sec A 048 sec 0.5 0-0.5 1 O oxiie 0.b3 0.d46 0.b6 0.075 ate RADIAL DISTANCE FROM JET AXIS 0.105 0.12 (in) Figure 4.2.3: Radial velocity profiles at various times; 0.46 m from jet orifice; 10 m/s initial speed. Figure 4.3.3 shows steady-state radial concentration profiles with three different grid densities. The densest grids , 34 axially by 30 radially and 31 by 27, have relatively similar profiles. This suggests that for these grid densities grid size independence has been attained. The concentration profiles of the least dense grid, 29 by 25, appears to deviate somewhat from the other two. This suggests that for an analysis of this nature the grid density should be higher than 31 axial nodes by 27 radial nodes, to ensure step size independence. 4.4 COMPARISON WITH ANALYTICAL AND EXPERIMENTAL RESULTS The round free jet has been investigated by a number of researchers in the past who have defined the jet behaviour both experimentally and analytically. One of the first to look at the 54 Figure 4.3.1: Radial velocity profiles for various time step sizes; 0.435 m from jet orifice; 10 m/s initial jet speed. round free incompressible jet was Tollmien [2]. From his experimental work he was able to confirm the validity of a simple algebraic expression for the axial velocity profile; U/U 0 =0.96/(ax/R.) where U is the axial velocity, R. is the initial jet radius, a is a constant ( see appendix A ). Based on his analysis it is possible to compare his results to that predicted by the TEACH code. Figure 4.4.1 shows the axial velocity profile of Tollmien's analysis compared to that of the TEACH code. It is clear that there exists good agreement between the numerical results and that predicted by Tollmien. The shapes of the curves appear to match very closely. The difference between the two curves appears uniform over their length, suggesting that it is not a breakdown of the basic theory at fault but rather some assumption about the initial ^ 55 RADIAL CONCENTRATION PROFILES 0.8 0.55 TIME STEP SIZE ^0.5^ 0.4^ • 0.02^sac A 0.013^sec 0 0.01^sac 0.35^ 0 0.0087 sec 0.3^ X 0.004 mm 0.46^ V 0.25 0.2 0.15 0.1 0.05 0 0.'^ a^0.1^0.12 RADIAL DISTANCE (m) Figure 4.3.2: Radial concentration profiles at various time step sizes; 0.46 m from jet orifice; 10 m/s initial jet speed; mass fraction in kg/kg. conditions. To'Innen states that the virtual origin of the jet is 0.0878 m before the x = 0 starting point, for a starting radius of 2 cm. The virtual origin is where the jet apparently starts from as a point source of infinite velocity but with finite momentum. If the virtual origin location is now set to 0.125 m, then the curves tend to collapse to one, shown in figure 4.4.2. The purpose in applying this new virtual origin location is to demonstrate the similarity in shape of the curves in figure 4.4.1. The sharp change in profile expected at the end of the core region of the jet is not well produced in the numerical results. A possible explanation is the tendency for many numerical codes to overestimate diffusion, called false diffusion [26], which would account for the quick 56 1 .2 1 .1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 Figure 4.3.3: Radial concentration profiles at various grid sizes; 0.46 m from jet orifice; 10 m/s initial jet speed; C m is jet centerline concentration. initial drop off of the jet momentum on the axis as it quickly diffuses radially. As well, using only two cells for the jet outlet may lead to lack of resolution in the outlet region. Birch et al [7] have produced similar results for an incompressible jet of methane gas into air. Figure 4.4.3 shows the axial concentration profile for a methane jet erupting into air compared to results provided by the TEACH code for a similar situation. An offset of 0.02 m is applied to the curve, similar to the previous figure. Again there is a reasonably good agreement between the shapes of the two curves. The slow change in the profile at the end of the core region encountered in the previous case of air into air is also present in the case of methane into air, indicating that the problem may be dependent on the numerical 57 AXIAL VELOCITY PROFILE 1 .2 1 .1 1 + TOLLMIEN M9 0 D NUMERICAL 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 • 0.2^0.4^0.6^0.8^1^1.2 AXIAL DISTANCE (m) Figure 4.4.1: Comparison of numerical and analytical axial velocity profiles; Ucl is 10 m/s the initial jet speed; U m is jet centerline velocity. methodology and not specific to the phenomena. 4.5 RADIAL PROFILES AND JET PENETRATION DATA Tollmien has investigated the radial profiles of concentration and velocity fields. He was able to develop tables of non-dimensionalized variables to describe their profiles at any reasonable distance from the jet orifice. His results showed that the radial profiles could be collapsed onto one universal curve [2]. Selecting the velocity profile at 0.435 m away from the jet orifice in the numerical results, a comparison can be made to Tollmien's data. Schlichting [29] also has made an analysis of the round free turbulent jet and developed an analytical model of its steady-state structure. Figure 4.5.1 shows a comparison of radial 58 Figure 4.4.2: Corrected prediction of the axial velocity profile; U 0 is 10 m/s the initial jet speed; U m is the jet centerline velocity. velocities between Tollmien , Schlichting and the numerical analysis. The radial distance has been non-dimensionalized by dividing by the radial distance r c , which is where the velocity has dropped to half its centerline velocity U.. The profiles are reasonably similar in nature and deviate only slightly near the jet edges. This would suggest that there is a slight discrepancy between the numerical and experimental results on the location of the jet boundary. However, in general the overall agreement is reasonably good. Tollmien has also investigated the profiles of concentration in turbulent jets. The results obtained by Tollmien and those for a numerical simulation are compared in figure 4.5.2. Again the results are plotted in a non-dimensional form. The agreement between the 59 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 Figure 4.4.3: Comparison of numerical and experimental axial concentration profiles for a methane jet; 10 m/s initial jet speed; C o is the jet outlet concentration. experimental and numerical analysis in this case is quite good. Very little deviation between the two curves occurs, and in general the profiles are virtually indistinguishable. The concentration profile in a methane jet has been studied by Birch et al. [7]. A numerical model of the round free transient jet of methane into air has been constructed to compare results with Birch's data. The results are shown in figure 4.5.3, where the radial coordinate is non-dimensionalized in this case by division with the axial distance from the jet (0.46 m). The agreement here is not as good as in the previous case. Although the profiles are similar, there is some deviation near the centre of the profiles. It appears as though the numerical code slightly overestimates the diffusion of the methane in the radial direction. 60 RADIAL VELOCITY PROFILES 1.2 1.1 RESULTS 1 • TINIAVEN 0.9 a 0.8 X NUMERICAL 0.7 0 SCHLICTING 0.6 0.5 0.4 0.3 02 0.1 0 -0.1 -0.2 0.25 0.5 0.75 1 1.25 1.5 1.75 2 2.25 2.5 2. DIMENSIONLESS RADIAL DISTANCE r/r C Figure 4.5.1: Comparison of velocity profiles of Tollmien, Schlichting and the numerical analysis; 10 m/s initial jet speed; U m is the jet centerline velocity. This behaviour would be consistent with the high momentum diffusion in the case of the axial jet. The methane and air jets possess different profiles as a result of the density differences between the two. Transient behaviour in free jets has been studied by Kuo and Bracco [21]. Their interest was in the penetration times of a transient round jet erupting into a quiescent medium. They measured the time required for the velocity at a given point on the jet axis to reach a certain percentage of its steady-state velocity. In figure 4.5.4, the non-dimensionalized penetration curve for a jet to reach 99.9% of its steady-state velocity is presented. The time and distance coordinates are non-dimensionalized by the Reynolds number at the orifice and the initial jet 61 CONCENTRATION PROFILES 1.2 1.1 RESULTS 1 • TOLLMIEN'S 0.9 CIE NUMERICAL 0.8 a 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0.25^0.5^0.^1^125^1.5 DIMENSIONLESS RADIAL DISTANCE r/re Figure 4.5.2: Comparison of numerical radial concentration profile with Tollmien's data; 0.46 m from jet orifice; 10 m/s initial jet speed; C m is the jet centerline concentration. diameter, details of which are given in appendix A. Clearly, the curves are very similar in behaviour in their rates of penetration. The overall shape of the plots appears to be in good agreement. It is only when the penetration times of the two sets of data are compared is there some difference. The TEACH code appears to overestimate the penetration time of the free jet compared to the results of Kuo and Bracco. The variation at the far ends of the curves are most likely the effects of the wall. The difference between the two curves is consistent for the most part meaning that on a qualitative basis the TEACH code is justified in such an analysis. In general, from the results presented it is clear that the TEACH code can predict the axial 62 RADIAL CONCENTRATION PROFILE 1.2 1.1 1 • NUMERICAL 0.9 a A( BIRCH 0.8 Li)^0.7 L) 0.6 0.5 0.4 0.3 0.2 0.1 0 a^0.1^0.15^02^0.25^0.3 DIMENSIONLESS RADIAL DISTANCE (r/z) Figure 4.5.3: Comparison of a numerical model and Birch's data of a radial concentration profile of a methane jet; 10 m/s initial jet speed; 0.46 m from jet orifice; C m is the jet centerline concentration. and radial concentration and velocity profiles of round free turbulent jets. In some instances such as the axial profile and some radial profiles of air into air, the agreement between the experimental and analytical results to that of the model is very good. The grid sizing and time -step size necessary to achieve grid and time-step size independence appears to be reasonable. Clearly, the comparisons made in this chapter give a good indication of the validity of using the TEACH code for the applications discussed in this thesis. 63 T star JET PENETRATION 140 130 120 + NUMERICAL 110 KUO & BRACCO 100 90 90 70 80 50 40 30 20 10 0 11 13 15 17 19 21 23 25 27 X star Figure 4.5.4: Comparison of modelled jet penetration times with that of Kuo and Bracco; 10 m/s initial jet speed. 64 CHAPTER 5: FIXED PISTON RESULTS 5.1 GENERAL Before proceeding to the moving piston model, the fixed piston model is studied to identify the variables which have the greatest influence on the jet development. The engine combustion chamber outlined in chapter 3 is modelled with a fixed piston height. Parameters such as the injection angle, the pressure ratio between the fuel tank and the combustion chamber, and the distance between the piston face and the cylinder head ( clearance gap ), are varied. For each of these three parameters there are only two values selected; the injection angle will either be 10 or 30 degrees from the cylinder head, the pressure ratio is either 5:1 or 2:1, and the clearance gap is either 0.3125 inches or two times this value. The gap 0.3125 inches is used since this will be the actual clearance in an engine with a flat piston face, 17 to 1 compression ratio, and a 5 inch piston stroke. A table is shown below which summarizes the various cases. Clearance Gap Jet Angle Pressure Ratio Case 1 0.3125" 10 2: 1 Case 2 0.3125" 30 2: 1 Case 3 0.625" 30 2: 1 Case 4 0.625" 10 2: 1 Case 5 0.3125" 10 5:1 Case 6 0.3125" 30 5:1 Case 7 0.625" 10 5:1 Case 8 0.625" 30 5: 1 65 As an introduction, the development with time of the jet sheet of case 1 is shown in figure 5.1.1 in terms of the mass fraction concentration of the natural gas. Each contour represents a line of equal gas concentration in kilograms of gas per kilogram of air/gas mixture. The jet plume penetrates outward from the injector tip in the top left corner of each plot, along the top wall of the chamber. By 3.59 milliseconds after the start of injection, the jet has encountered the far wall and begins to fill the chamber with natural gas. Figures 5.1.2 and 5.1.3 give concentration and velocity plots for the 30 degree case, showing its development. The thirty degree jet propagates toward the bottom wall and not the top wall as in the previous case. Again the jet encounters the far wall and turns to fill the remaining chamber. The velocity fields show zones of recirculation under the injector, which appear to correspond to zones of high gas concentration. Figures such as the ones just shown are used throughout the next two chapters to demonstrate the jet character. This chapter is divided into six sections. The second section of this chapter (5.2), compares experimental results with numerical results under similar initial and boundary conditions. The third section ( 5.3 ) looks at how a change in jet angle influences jet mixing and penetration. The effects of an increase in pressure ratio are described in section 5.4. A higher pressure ratio injects a greater amount of fluid momentum into the chamber and its effects are shown. Section 5.5 shows the development of a region combustible mixture as the jet penetrates into the chamber. The final section (5.6) investigates the effects of the clearance gap on the jet behaviour. 5.2 COMPARISON TO EXPERIMENT One of the best methods by which the performance of a numerical model can be assessed 66 radial distance (m) time = 1.14 ms 0 o 0 O O O 0 0 O 0 0 O 0 0 01 01 11- 0 0 co cn O O 1111111^111111111 O O O co time --a 2.13 ms o. time EN 3.59 ms 16 0.04 time nu 5.72 ms 0.4 o n4 ^nfut Figure 5.1.1 Concentration field for case 1 showing jet development with time, clearance gap 0.3125 inches, pressure ratio 2:1, jet angle 10 deg., contour interval 0.04 (kg/kg), injector tip in top left corners. 67 radial distance (m) time = 2.13 ms 0 0 O 0 O O O cri^co 0 0 0 0 0 0 0 0 co^ 0 0 co 0 0 0 0 -P• IIIIIIIIITIIIIIIIII O 0 21 0 0 0 co 0 0 time = 3.59 ms time = 5.72 ms Figure 5.1.2 Concentration field for case 8 showing jet development with time, clearance gap 0.625 inches, pressure ratio 5:1, jet angle 30 deg., contour interval 0.08 (kg/kg). ^ a) U) O 0 C-) C/) C) O 01, . 4. 4 4 . a 4 4. 4 ^ 1C, I^4 • 40441 ,4444. 0 4 44 0 4 11 4 MR R a a a a 1, * AMAMI. Pt IR • Ot IP It a a Rt ~MR It 44144444.444.4 i p R a .. 4416 Malan kit 6 6 444 444 4.^a 44 a Qtr g,' R'RMMIRRIR IR IR k 1 466* 64416**^ 464* 46444666^ 444 46064~^ 444 444444461464 V. i * . 6+16446444464***4 44 • at 44 a 64+4 446644444 * * i ... 4441444MN.144 4. 16 It 116 i 14 644 4 4444444***** * * 4. * ^ti Q.) U) 0 O 444 4 .164444444644 # ^ ... 4 + I I P 6µ 4 4444441i1 44444446.4 ^ . 4 4 4} t D IT . • t 4 „ : * . ,, , :, , ii ri, i - -4 --. '. za 4 . . . p 1 T -a .., 7. 7 i I 2• r . - I P 444 4 4444444466 4 • ^ .. 444 4, 4.1444, 4 4 4 41,1. ^ ^ . + 444 4446446664 • 4. 4. + a 4 464 +44444444 ^ 4µ1444464444i 1 ^ 414 4111444444 4 ^ 4 414 4 441444444 4 4 4 4 4 * 4 a 111141111111114 4 4 4 4 11114444 .,..,- -*:%* 14 a A 4 -" _-4 -----_-_,^/m61 --1' ''' )4-- , e --if,. 141^ I. V 4 r e -4--0 ---4 "r4 -- a , ^v les "44**444 4 4^ 4 4-4-4*-64.--- at t 4 +I.?? • a r 4 At I I kb* 4114661460*-46-4 4- 4- e- 4- 4- 4 f ♦ r 6461446444644414.16 ^ “1.1.414,14. ^ ... 4441µ44444444 4t ^ 464 4 466444 ^ ....^• * I I r r P • 4 is ^ 414+ 41414441 414+ 1144444 t 6 4 t 4444 4µ444444M4 tat 44 4 4 4 4 1- c..) +44 4 4444444+ ^ ^ Q.) rf) 44 7 I I ---..., r1 444 4 4444466-44 .6 1 + 4 4 4 * a * ''' E . . t t I 444 444444444 41114^ + 4 4 4 0^ . 4 4 .„,TI 4 1 1 1 4U14": ^ 4414 4144444444 4 4 4 4 4 4^' ' c:- 7:......,,,„ " a 'a , ^,` ), a ), ,t --■ :, : , ..-- :_. —: .,,,,,: 4. - 4_+_,0 . i 1411111114414 4 4 4 * 4 4 4 4 A '''' -' 7 t I ii:ii: 7_ --1 Au^ Lit^4,- C...) 441•1144444•4^* 4 %,,,a1: .3:1,1:6-" - 4 :1) 41,4444kiti.44 . --`,Z * 41:01 4 - ,,,- : -!- : -.,---- 17.7 .... . -.4 ...> 69 time = 1.14 ms time — 2.13 ms Figure 5.2.1 Schlieren image compared to case 8, jet angle 30 deg., pressure ratio 5:1, clearance gap 0.625 inches, concentration contours 0.08 (kg/kg), injector tip in top left corners. 70 Figure 5.2.2 Schlieren image comparison to case 7, jet angle 10 deg., pressure ratio 5:1, clearance gap 0.625 inches, concentration contour interval 0.08 (kg/kg), injector tip in the top left corners. 71 time = 1.14 ms time = 2.13 ms oght Figure 5.2.3 Schlieren image comparison to case 1, jet angle 10 deg., pressure ratio 2:1, clearance gap 0.3125 inches, concentration contour interval 0.08 (kg/kg), injector tip in top left corners. 72 time = 1.14 ms time = 2.13 ms Figure 5.2.4 Schlieren image compared to case 6, jet angle 30 deg., pressure ratio 5:1, clearance gap 0.3125 inches, concentration contour interval 0.08 (kg/kg), injector tip in top left corners. 73 is by comparison to experiments. Work on an experimental model of similar shape and initial conditions as the work done in this chapter was performed concurrently by Patric Ouellette [25]. Some of the results he has obtained are presented. Four cases were selected for the comparison, two 10 degree and two 30 degree injection angle cases. Figures 5.2.1 through 5.2.4 show schlieren photographs of a jet of natural gas corresponding to the cases considered. In each figure there are two schlieren images; the topmost is at 1 millisecond after the start of injection and the bottommost is at 2 milliseconds. The photographs are of the complete chamber while the concentration plots are of half the cylinder with the centerline on the left edge. Figure 5.2.1 shows the experimental result for the 30 degree case, with a larger chamber. The image in the photograph clearly shows the shape and the extent to which the jet has penetrated. The jet appears roughly bell shaped which is closely parallelled by the numerical results shown below. The qualitative agreement between the photograph and the numerical model suggests that the assumptions made in the model were reasonable. Figure 5.2.2 shows the comparison to the 10 degree, large chamber case. The jet is clearly penetrating into the chamber along the top wall which was predicted in the numerical model. The jet penetration appears to be well predicted. The 2 ms images both have penetrations on the order of 4.5 cm. The numerical model predicts a jet held tight against the top wall, while the photograph shows a much looser jet which has penetrated axially into the chamber. It should be considered that the results of the numerical analysis are averaged values, while the photograph gives an image that is one instant in time, which may suggest one reason for the difference. 74 Figure 5.2.3 shows the comparison to the 10 degree, small chamber case. The photographic images are not very well defined in this case, but it is possible to ascertain the edge of the jet. The 1 ms photograph clearly shows the jet clinging to the top wall which is also seen in the numerical model. In the 2 ms photograph it is again apparent that the actual jet has penetrated axially more and radially less than in the numerical model results. This indicates that the 10 degree results may not be as accurate as the 30 degree cases. Figure 5.2.4 shows the comparison to the 30 degree, small chamber case. The image in the photograph does not offer the detail of previous images. However, it is possible to define the edge of the jet. What can be seen is a jet that has reached the bottom wall at 1 ms in the first photograph, unlike the numerical model. The second photograph suggests that the model and the experiment give the same jet penetration. In general, the photographs compare well with the numerical results. In most cases the numerical results were qualitatively similar to the experiment. The resolution was not sufficiently detailed to make a quantitative comparison. However, the results obtained by the experiment indicate that the numerical code simulates the real flow reasonably well. 5.3 INJECTION ANGLE This section investigates the influence of jet angle on the gas jet inside a fixed piston chamber to ascertain if it is of significance. Two comparisons are made; a 10 degree and a 30 degree jet each with similar pressure ratios and chamber sizes. The character of the jet appears radically altered by the initial injection angle. The concentration plots for 2:1 pressure ratio and small clearance gap, shown in figure 5.3.1, demonstrate the magnitude of the change brought about by a 20 degree change in injection 75 angle. The jet in the 30 degree case starts off at 30 degrees but as it penetrates into the chamber, shifts direction until the jet is nearly vertical. The jet then attaches itself to the bottom wall and proceeds to fill the chamber from this position. The velocity plot , shown in figure 5.3.2, confirms the jet path shown in the previous figure. It is clear that a high degree of recirculatory flow is developed just under the injector. A zone of high gas concentration exists under the injector suggesting gas is drawn into the region of recirculation. The 10 degree case behaves in a manner that is in complete contrast to the thirty degree case. The jet now clings to the top wall as it fills the chamber. Figure 5.3.2, showing the velocity field for this case gives no indication of the recirculation zone under the injector that was apparent in the thirty degree case. In general the magnitudes of the velocities in the 10 degree case appear less than in the 30 degree case, suggesting a cause for the difference in the concentration fields. Figure 5.3.3 shows concentration fields for cases with larger clearance gaps and pressure ratios of 5 to 1. It is quite clear that the jet patterns seen in the previous cases are again repeated in these cases. The thirty degree jet still moves towards the bottom wall and the ten degree jet is still fixed to the top wall, in spite of the change of pressure ratio and clearance gap. The velocity fields for these two cases in figure 5.3.4 show the same pattern as seen in cases 1 and 2. Clearly the injection angle has a considerable influence on the jet character regardless of the clearance gap or pressure ratio. It has been established that the jet is sensitive to injection angles, but the range and character of this sensitivity has not been explored. It has yet to be shown if the jet gradually switches between the top and bottom walls or if this change is abrupt. 76 time 2.13 ms jet angle 10 deg. jet angle N- 30 deg. time = 3.59 ms jet angle = 10 deg. jet angle = 30 deg. Figure 5.3.1 Concentration fields for case 1 and 2 showing the jet angle effect, pressure ratio 2:1, clearance gap 0.3125 inches, contour intervals 0.04 (kg/kg), injector tip in top right corners. 77 time am 2.13 ms jet angle 10 deg. VELOCITY SCALE.^100 m/sec jet angle 30 deg. time = 3.59 ms jet angle an 10 deg. VE OCI A E ^100 m/sec . VELOCITY SCALE.^100 m/sec U( T jet angle = 30 deg. Figure 5.3.2 Velocity fields for cases 1 and 2, showing jet angle effects, pressure ratio 2:1, clearance gap 0.3125 inches, injector tip in top left corners. 78 jet angle 10 deg. jet angle = 30 deg. Figure 5.3.3 Concentration fields for cases 7 and 8 showing jet angle effect, pressure ratio 5:1, time 3.59 ms, contour interval 0.04 (kg/kg), injector tip in top left corners. ^ ^ ► TTh s Vi4 • •Trimi ^'---- ^T-4 1*-4- •+ -1-4•44÷ -V•er ^ ++ ::•••• •►■••• • ••••••••■• + + Uq CD ^+A+ CD 44 1,1 CA ■1) 4.) 73 A A A ++ a42 4 a a a ++ CD ++ ^ +++ 4 4 4 ARAN ^ AA ^ ++ 4444 • ^414411.1t t PD A +++++ A A II 4 A A A A4 , 14 • + }44}1444144}4h Cr4 CL C'l 4 1••••••■•■ A A AxAmAN A A %%AAA+ ^ AIMN A A^ C) N ul -t 1 AAAAAAAAAN 4} 4 A 4 AAAAAAAAAAAA. +++++ +44+1.4411+1, 444 . Ci) A A A A A A A A A A A AAAA44, cb o ^4 fa. , ^4 CD 00 0 ^a 4. o .0. t a 0 5' (IQ 5^++ ^ +1.91.1.1.1144 AAR ^AA 491 0 0 n " CD CD P CD (#3 CD q 1=t. 0 • 0, CD 'A A %%AAA% AA 4 S A 9. A A • a A An•%ma aaa444,1 ^++ +++ ^ 444-44449+44 4444 a a 4 r f+ 4 ,4 4++44444+944-44411 - Fr 5 9141■ , f a • • + f • • + f • 414+44144 • A• +++++ 1114,4.111 f A 4 4***AAA +++++ .4444441 80 jet angle 18.75 deg. jet angle = 19.0 deg. jet angle 19.375 deg. Figure 5.3.5 Concentration fields showing jet angle sensitivity, pressure ratio 2:1, time 2.13 ms, contour interval 0.08 (kg/kg), injector in top left corners. 81 Figure 5.3.5 gives a clear indication of the very abrupt nature of the switch. The concentration contours for the 19 degree case looks similar to that of the 30 degree plots. The 18.75 degree field looks like the 10 degree case. The difference in jet patterns between the 19.0 and the 18.75 degree cases is clear evidence of the jet's extreme sensitivity to the jet angle. As a confirmation of these results, an experimental investigation of the 20 degree jet performed by Patric Ouellette [25] revealed that the switch over occurred abruptly near 20 degrees jet angle. In general the path of the gas jets between stationary walls has been shown to be highly dependent on the angle at which they are injected into the chamber. The presence of a wall will prevent the jet from entraining surrounding air as would be the case in a "free" jet. A low pressure zone develops between the wall and the jet, which is shown in the next section. The low pressure zone may bend the jet towards the adjacent wall. The abrupt change in jet behaviour is an indication of the relative intensity of the low pressure regions on either side of the jet. In the next chapter the effect of jet angle with the moving boundary (piston motion) is described. 5.4 PRESSURE RATIO High speed gas jet injection requires a high pressure ratio between the injector and the combustion chamber. It is possible that this pressure ratio will be great enough to make the jet underexpanded. In this case a greater amount of fluid momentum will be injected into the chamber through a higher inlet density than a lower pressure ratio case. Figure 5.4.1 shows the effect of the pressure ratio on concentration fields for two 10 degree jets in the small chamber driven by different pressure ratios. In both cases the gas jet 82 clings to the top wall and induces a clockwise recirculation inside the entire chamber seen in the velocity plot of figure 5.4.2. The jet shapes and penetration distances are much the same, with the higher pressure ratio case having slightly greater penetration . The similarity in concentration and velocity fields together suggest that for the 10 degree cases the pressure ratio is not a strongly influencing factor. Figure 5.4.3 shows the pressure ratio effect on concentration fields for two 30 degree jets in the large chamber driven by different pressure ratios. The higher pressure ratio shows a significantly different concentration and velocity field than the lower pressure ratio case. Here, the jet has a tendency to move out into the chamber before moving towards the bottom wall. The lower pressure ratio jet appears to immediately move underneath the injector and along the symmetry axis. The differences in the two are stressed by the velocity fields in figure 5.4.4 , which shows the dissimilar paths taken initially. It appears that the extra momentum of the high pressure ratio case causes the jet to continue in the path it was initially directed for a longer period of time than the low pressure ratio case. This is evident in the velocity fields where the high momentum jet strikes the bottom wall further radially than the low momentum jet. An investigation into the relative pressure fields produced by these cases may explain the variation in behaviour. Figure 5.4.5 shows the pressure fields near the injector tip 2.13 milliseconds after injection start for the 5:1 and 2:1 pressure ratios respectively. Clearly, the two pressure fields are completely different. As has been mentioned, the jet is a thin sheet. Therefore, it may be sensitive to low pressure regions, which was stated in the last section. The low pressure ratio case has one zone of high pressure in the corner opposite the 83 time= 2.13 ms pressure ratio 5:1 pressure ratio 2:1 time 3.59 ms pressure ratio 5:1 pressure ratio 2:1 0.‘ 0.04 Figure 5.4.1 Concentration fields for cases 1 and 5 showing pressure ratio effects, jet angle 10 deg., clearance gap 0.3125 inches, contour intervals 0.04 (kg/kg), injector tip in the top left corners. 84 time = 2.13 ms V EL CITY,"!§C^100 m/sec 0111KOM pressure ratio 5:1 11111 (i ( pressure ratio 2:1 = I I' 1 I time = 3.59 ms ^ \TEL CITA;^100 m/sec MI I WIN pressure ratio 5:1 pressure ratio 2:1 Figure 5.4.2 Velocity fields of cases 1 and 5 showing pressure ratio effects, jet angle 10 deg., clearance gap 0.3125 inches, injector tip in top left corner. 85 pressure ratio 5:1 pressure ratio 2:1 Figure 5.4.3 Concentration fields for cases 8 and 3 showing pressure ratio effects, jet angle 30 deg., time 2.13 ms, clearance gap 0.625 inches, contour interval 0.08 (kg/kg), injector tip in the top left corners. ^ • ▪ Mm 8 ,48880444.44418 4400S ^ PPPPPPPPP 4 4 4 44,80441/81.88, If V 8 P r, IP S 4 4 4 4 ^1 of4.48 4- +11444.1* 441-4 Off l st. 4 4 i.“ % MINNA a a r 4- Er Or 4 + 4. 4 - 1 4 4 4 w 8 4 4 WW C.) W*0 • ao,66444 a.) ^A4.4. ^4 144 484848^ WWWWW R /8 R R 4 4 4 WW 444~60 .141.444444444 f 4 f 4 4 O 16 4 f C 4 f 444 4444444^ W 4 f* 4 4 W 14,14444444 4 4 4 ^ WWW 4414-11.44444iii*** t I v*. • t 1' T r W 4 L.) 1,444.4,6444.(ii... ^+ cf)^444,/,,a„ii^"14 • saW.e /40C-ICA'. 4.1 ^T W 44+ 4 44444444.4 • 4 W ^ IM+441414++4+464 1 4 4 4 4 4 8 I Cri 4 A r l / 87 injector tip and one zone of low pressure directly underneath it. The values of pressure are given relative to the zero datum. The lowest pressure in this case is approximately -6400 Pa, underneath the injector. Observing the corresponding velocity field of figure 5.4.4, the low pressure zone under the injector is located in a region of intensely recirculating flow. The high pressure region is associated with an area where the velocity field is rapidly changing direction or is stagnating, as in the corner of the chamber. The presence of high pressure zones indicate the existence of stagnation points. The presence of low pressure regions tend to indicate the existence of recirculating flow patterns. The high pressure ratio case has three zones of high pressure (stagnation), one underneath the injector, one in the corner opposite the injector and one located radially one centimetre from the corner high pressure region. Two low pressure regions exist; one directly adjacent to the injector tip radially outward, and one beneath it centred in the middle of the high pressure regions. The very low pressure zone opposite the injector ( with the majority of the contours removed for clarity ), is associated again with a region of recirculatory flow behaviour as seen in the corresponding velocity plot of figure 5.4.4. It appears as though a low pressure region underneath the injector, that may divert the jet, never develops in the high pressure ratio case. It is evident that in the 10 degree jet angle cases the effect of pressure ratio on the jet behaviour is minimal. Only in the case of large injection angle and clearance gap does the pressure ratio have a significant effect. 5.5 COMBUSTIBLE MIXTURE REGION Information on the mass fraction of natural gas present in the chamber at any point in time ^ 88 pressure ratio 5:1 contour interval 250 Pa ^Pp^p^PPP^P^PPS) §^§4§ o § ^2^t2/^ g Ei^2^2 i pressure ratio 2:1 contour interval 1000 Pa Figure 5.4.5 Pressure fields for cases 8 and 3 showing pressure ratio effect, jet angle 30 deg., clearance gap 0.625 inches, injector tip in top left corners. 89 provides an opportunity to investigate the development of a combustible mixture region. In the actual engine combustion will commence soon after a combustible mixture is produced. This will eventually distort the flow field and concentration field to configurations other than what is predicted here. Nevertheless, this section will give an insight into the speed and shape of a combustible region as it begins to develop. For natural gas the combustible range can arbitrarily be defined as the volume lying between the concentration contours of 0.035 to 0.07. These numbers represent the flammability limit at standard temperature and pressure, and correspond to a relative air-fuel ratio of 1.6 and 0.8. The model is run assuming standard temperature and pressure conditions at all times. The results shown have only a qualitative significance since the high temperatures of the actual engine will alter the combustible limits. The higher pressure of the actual engine has little effect on these limits. The 10 degree jet in the large chamber is selected to demonstrate the growth of the combustible region. Figure 5.5.1 shows the region of combustible mixture at various points in time. As is expected, the combustible region moves and expands as the jet does. For most of the jet, the combustible zone is limited to a band of very small thickness. Only when the jet hits the far wall and turns to fill the cavity does the zone attain appreciable thickness. It is clear that the jet entrains air as it expands outward into the cavity. The combustible region remains thin near the jet outlet where considerable entrainment of air is expected. Only when the jet has circulated into a region of relatively slow moving fluid does the combustible zone enlarge. 90 time = 1.14 ms time — 2.13 ms o .07 0.035 time =-- 5.72 ms 0 .02,..............^. .^ . ,.......,I Figure 5.5.1 Combustible region development for case 7, pressure ratio 5:1, jet angle 10 deg., clearance gap 0.625 inches, contours given in mass fraction (kg/kg), injector tip in top left corners. 91 5.6 CLEARANCE GAP An investigation of clearance gap effects will be a useful exercise in gaining a more comprehensive understanding of the jet behaviour. The effect of clearance gap is most apparent when parameters such as the injection angle and pressure ratio are held fixed and the clearance gap is allowed to vary. Figures 5.6.1 and 5.6.2 compare the effects of clearance on the concentration fields for two jets of 30 degrees angle each of a different sized chamber. Clearly, the concentration fields of both cases are similar. The jet, once it has emerged from the injector, attaches itself to the piston face. Both jets appear to have roughly the same penetrations, with the smaller chamber jet having propagated slightly further. In both cases a zone of high gas concentration exists just below the injector tip. The one noticeable difference in behaviour between the two cases is the mixing occurring in the chambers. The smaller chamber is nearly completely occupied by gas of some concentration at 3.59 milliseconds, while this is not the case for the larger chamber. The gas appears to spread upwards, proportionally, at a faster rate in the smaller chamber, referring specifically to figure 5.6.2. If this behaviour were repeated in the comparison of 10 degree jets then it can be inferred that the smaller chamber generally yields slightly higher mixing rates than the larger. Figures 5.6.3 and 5.6.4 show two cases with a jet angle of 10 degrees and a pressure ratio of 2:1. It is clear that the behaviour of the jet in these cases is similar to that of the previous two. The results suggest that the clearance gap has only a slight influence on the overall jet behaviour. In both sets of cases just compared, the jets do not change significantly with clearance gap. However, the proportionally faster spreading seen previously in the smaller 92 clearance gap = 0.3125 inches clearance gap 0.625 inches radial distance (m) Figure 5.6.1 Concentration fields for cases 6 and 8 showing the effect of clearance gap, pressure ratio 5:1, jet angle 30 deg., time 2.13 ms, concentration contours 0.08 (kg/kg), injector tip in the top left corner. 93 clearance gap = 0.3125 inches clearance gap = 0.625 inches Figure 5.6.2 Concentration fields for cases 6 and 8 showing clearance gap effect, pressure ratio 5:1, jet angle 30 deg., time 3.59 ms, concentration contour interval 0.08 (kg/kg), injector tip in the top left corners. 94 chamber is again present in figure 5.6.4. It is clear from the results just shown that the clearance gap has only a small effect on the complete jet character, independent of pressure ratio or jet angle. 95 clearance gap a- 0.3125 inches .......■,.. ,—.., sz pA clearance gap = 0.625 inches Figure 5.6.3 Concentration fields for cases 1 and 4 showing clearance gap effect, pressure ratio 2:1, jet angle 10 deg., time 2.13 ms, concentration contour interval 0.08 (kg/kg), injector is in the top left corners. 96 clearance gap 03125 inches clearance gap 0.625 inches — illiWjRr= ------- IT;:rr;- 0.04. ------ ° J" ------ — ,^ .----o.04 Figure 5.6.4 Concentration fields for cases 1 and 4 showing clearance gap effect, pressure ratio 2:1, jet angle 10 deg., time 3.59 ms, concentration contour interval 0.04 (kg/kg), injector in the top left corners. ^ ^ 97 clearance gap 0.3125 inches VELOCITY SCALE:^100 m/sec •*----7- EFF_FFFFFF_F^F_^f^f 411 1; .e •4,44*. 10.V4‘ss;^\ ,AN clearance gap =1 0.625 inches VELOCITY SCALE. , 100 m/sec I • t i,._^ ^i^f^f •^t^i^i^d•i kr.:-^f_r_f_Efff ; \ [tiff!!! 1 1 1 • I . , ^t 1^.i...).4._ i ^4 ^\ " •fl^tV 1. IP V ,, 4,- 4- +^+^+^4+ + • 0,4 • 4- + +^4-^4. 4- I+ ••^4-^♦ '40‘')), ,^, 4. ^4^4.^44 4 a 4-^4 • *^ \ : , ,,24 \ \; :ie^.,,,^*^4-^4.^4. r 4- 4- 4- ‘). ),^ •^ ,‘ :,., \ \,, \^4 4.. 1^1 \ 4. 4+ 4. Iiirc7!..r•kk4,_*-- , ., F • ^• a *^ A .': y 9 -9 -....„^--.^ -• -4^ • ----• --• --• --■ —•^--. ,,^ , ■^,\^) .. . .^ i',^ ,'s\g, ‘„, 4- 4.^4. • + * 44. 4 4* 4. •^4- F -•^ 4 + * 44- 4 4 • 4 • 1 1 I^1 i^I^a^ 4• 44. 444 4- •^4 0 -0^ -• --4 44 + 4 44 4 * • • • 44+ • 4+ 4 • 4 4 • 4 • V • • • • 4^ -0 4 4 +^ -0 —0^ -4 a l I. 4 444.• • 4- ♦ • • • 44 0. I. • 0, • • • 4 **..144+41-444.4.44.4_,___^,6 .....,4' 4 ---'4 —_. —_, —4 —4 --4 —4 ____. __4 ^—9^-4^—0^—4^ —•^-4 • • 4 9 Figure 5.6.5 Velocity fields for cases 6 and 8 showing clearance gap effect, pressure ratio 5:1, jet angle 30 deg., time 3.59 ms, injector in the top left corner. 98 CHAPTER 6: MOVING PISTON MODEL 6.1 GENERAL The moving piston model investigates the effects of piston motion through a moving boundary. The flat bottom wall of the chamber is made to approach the stationary top wall using the coordinate transformation of the axial spatial variable outlined in chapter 3. A bowl shape could have been imposed on the piston, but was not done in this analysis. The moving boundary characteristically moves in a manner consistent with the piston of a 71 series Detroit Diesel engine. The motion of the piston is obtained through an analysis similar to that of Heywood [15], given the true engine dimensions. Injection of the gas in the engine occurs, at the earliest, 35 degrees crank angle rotation before the top most piston location ( 35 degrees before top dead centre; BTDC ). The injection in the model will be started from this crank position with a velocity profile equivalent to the fixed piston model. In the actual Pressure Ratio Engine Speed Jet Angle Case 1 1200 RPM 30 degrees 2:1 Case 2 1200 RPM 30 degrees 5:1 Case 3 1200 RPM 60 degrees 2:1 Case 4 1200 RPM 60 degrees 5:1 Case 5 600 RPM 60 degrees 2:1 Case 6 600 RPM 30 degrees 2:1 Case 7 1200 RPM 30 degrees 2:1 Case 8 1200 RPM 10 degrees 2:1 Case 9 600 RPM 10 degrees 2:1 Case 10 1200 RPM 10 degrees 5:1 engine the injection starting angle and its duration will depend on engine load and speed. The 99 earliest injection crank angle is selected to give an opportunity to observe all possible piston motion effects. Ten cases are considered, tabulated above. Three parameters are selected to observe their effects on the development of the jet; engine speed , jet angle and pressure. Case 7 is similar to case 1, only the initial conditions of the air and gas in the cylinder are consistent with that of the actual engine at 35 deg. BTDC. This comparison is made, in section 6.2, to demonstrate the effect of variation in gas temperature and pressure. Section 6.3 focuses on the effect of pressure ratio on the jet behaviour, and will compare results to the fixed boundary. The effect of injection angle on the jet behaviour is then investigated in section 6.4. In section 6.5 the engine speed is considered to study its influence on the jet. What is also shown is that the results of the moving piston model are fundamentally different than the fixed piston analysis through a comparison of both models with similar conditions. In the final section ( 6.6 ) the development of a combustible mixture region under the influence of a moving piston is considered. 6.2 EFFECT OF CYLINDER PRESSURE AND TEMPERATURE LEVEL The compression ratio of the engine is approximately 17:1, giving a high temperature and pressure of the host air when the jet commences. In all cases the jet commences at 35 degrees BTDC. At this point the compression ratio is at 6 to 1, giving a host pressure and temperature of 1279 kPa and 590 K. Ideally, the initial temperature and pressure of the moving piston model should be consistent with the engine. However, to make a valid comparison to the fixed piston model the initial conditions must be similar. Therefore, it would be expected that the same standard temperature and pressure used in the fixed piston case will be used in the 100 moving piston model. To justify this approach it is necessary to show that the dynamics of the phenomena do not depend appreciably on the absolute values of the pressures and densities of the air and gas but on their ratios. Figure 6.2.1 shows the concentration fields of a 30 degree jet with both standard temperature and pressure and engine initial conditions at 6.8 and 16.3 degrees BTDC. The results appear almost identical, which suggests that the values of the pressure and density fields are not the controlling factors of the jet dynamics. In fact it appears that the ratio of the densities of the injected gas and the host gas is the important factor. The small increase in the laminar viscosity of the fluids with the increased temperature (proportional to T 314 ) of the actual engine configuration appears to be insignificant. This seems logical when it is realized that for most of the flow field the turbulent viscosity may be much larger than the laminar viscosity. It can then be concluded that for velocity and concentration fields it matters little if the actual engine conditions or that of standard temperature and pressure are used as initial conditions. 6.3 PRESSURE RATIO EFFECTS As in the fixed piston model, the influence of pressure ratio on the jet behaviour is investigated. The higher pressure ratios are expected to impart a higher momentum to the jet, which coupled with the effects of the piston motion may yield results that are significantly different than the fixed boundary model. Previously it was shown that the pressure ratio had only a very limited effect on the jet behaviour. Only in the case of a 30 degree jet with enlarged chamber was there a significant effect. Two values of pressure ratio between the fuel tanks and the combustion chamber are used. 101 C.A. Ins 16.3 deg. BTDC Engine initial conditions STP initial conditions C.A. = 6.8 deg. BTDC Engine initial conditions STP initial conditions Figure 6.2.1 Concentration fields for cases 7 and 1 comparing actual engine and STP initial conditions, pressure ratio 2:1, jet angle 30 deg., engine speed 1200 RPM, contour intervals 0.08 (kg/kg), injector tip in top left corners. 102 The 2 to 1 ratio is considered since it is near this value that the flow through the injector will become choked ( Mach 1 ). The higher pressure ratio will be choked as well, but will also be underexpanded. In the engine, the pressure ratio will change as the cylinder pressure increases with compression. A constant pressure ratio is used in the moving piston model. The grid spacing would need to be adjusted at every time step otherwise to account for the pressure ratio change, based on the effective diameterarguments. Figures 6.3.1 and 6.3.2 show the concentration and velocity field for a 30 degree jet with pressure ratios of 2 and 5 to 1. Clearly the difference between the two cases is considerable. The higher pressure ratio case has a velocity and concentration field that attaches itself to the top wall and moves along this boundary. At no time does this jet attach itself to the lower wall as in the fixed piston model. The lower pressure ratio jet appears to initially move towards the top wall and then as top dead centre is approached, switches to the bottom wall. It is evident from these results that the pressure ratio in conjunction with the effects of piston motion can have a significant effect on the gas propagation in the chamber. In the fixed piston model the 30 degree jet always moved towards the bottom wall with only a minor change resulting from an increase in pressure ratio. There, the higher pressure ratio large chamber case had the 30 degree jet hitting the lower wall further radially than the lower pressure ratio case. In the moving piston model the jets behave quite differently in several ways. First, the jet does not immediately attach itself to the bottom wall. The lower pressure ratio jet eventually attaches itself to the lower wall, but the higher pressure ratio case never does. Secondly, the jets develop in a completely different manner. The high pressure ratio jet is firm against the top wall as it entrains air to fill the far end of the chamber. The low 103 C.A. =. 16.3 deg. BTDC Pres. Ratio am 5:1 Pres. Ratio ma 2:1 C.A.= 3.2 deg. ATDC Pres. Ratio = 5:1 Pres. Ratio in 2:1 Figure 6.3.1 Concentration fields for cases 1 and 2 showing pressure ratio effects, jet angle 30 deg., engine speed 1200 RPM, contour interval 0.08 (kg/kg),injector tip in top left corners. 8 44 U 0 O n II 16 ^ a • • 16 • 11, 16 4 + 4 ar16416 16 It % 4 * 4, 19444414••414:14 4 4 4 st 4 4 4 ‘M41010441644 4 416 4 4 4 1 11140111•14.41616 Wu 4 16 4 4 4 4 II 1 1 III 1 III 11 . tick 8 8 8 16 4 4 4 4481414116 4 8 4 4 161616 4 * la 4 46 4 4 .,k^ , :,,,,..,...., ■i, , .:17 16 4 * tlahr 4 4 16816 11•466 • 66 16 16 4 4 46446- 4 4 4 4 16 16 4 It* 6-k • 4 4 8 4 4 444 4 4 4 k 16 16 4 444 4 8 4 4 8 4 4 4 4 4 8 4 4444 4- 4 • 4 4 • 4 4 4 4 4 4- 4 4 • 4 8 4- 4 4 • • 4. 4- 4 4 4_4_14 46_4_, 4 • 1--I O _J LLJ .4,0440=44•4164164/66 , • 66 • • • 4 4 4414011114144 164444 • 4 4 4 4^+ 41044141114444444444^• 4 • + 4 4H84L4.4M?f4 , 444 4 4 4 • 4 4444414144414444444 66 +^+ 4 444444 4 4 4 4 • 4 4 • •••••••••■••••••• I. • a • a • • 44••••••••••••••• a A 444444111114411440•1444 4 4 46644441111111144*44144 i 61 t M 66 4 44444•11MP/4444*48.^• 4444 * • 4 • * 4044."""" 11 • * • A d 4 4 14 4 4 4 16 4 4 4 t If Id 8 t 4 ark V 6 II * 41- le 16 V rrt t 16 st Ir^11 41 48 8 4 4-^4 4 4- 4-^Iv^4 4r Z4 4 4 t ounItft.4 4.4±± t_t U 1tl CDn o ill1411114844444.4 • 4 4 4 INk.,....„,,,34 3 4 * 4444 4 + 4 lillharterwa *v.* • * 111 I Nkuirve******* 4 a in Ill ^444 4 • • 'If^ 44 4 4 * ‘4 !Ir .44. I I '^ 4 ". 161648 4 4 4 8164 4 4 • 8844 18 4 tan 16446 4 ^ ^ 84 4 ^ 16. 4 44 8 4 4 4 6- 4 4 8 4 44.4 4 4 44 4 11 4 4444 4 ' +4 4 4 li 611444 t• 4 4-.426 •• 105 pressure ratio jet gradually attaches itself to the bottom wall, initially filling the area nearest the injector tip. Finally, the higher pressure ratio adds stability to the jet behaviour. The high pressure ratio jet does not change its path away from the top wall. The lower pressure ratio jet switches movement from the top wall to the bottom wall as the piston approaches. The pressure ratio effects mentioned have only become apparent in the presence of the moving piston. Other cases may not show the same variation in behaviour from pressure ratio effects as is seen in these cases. Figure 6.3.3 shows concentration plots for a 10 degree jet, at 2 and 5 to 1 pressure ratio. In this instance the concentration fields appear to have the same general structure, independent of pressure ratio. The high pressure ratio jet has a slightly greater penetration distance than does the lower pressure ratio case. The similarity in the concentration profiles suggests that the pressure ratio effects seen in the previous case may be influenced by the injection angle. The change in jet behaviour with pressure ratio seen previously is not present in the 10 degree jet. In general it is quite clear that the combination of pressure ratio and moving boundary effects have a profound influence on the jet behaviour, except at 10 degrees injection angle. At 30 degrees angle the higher pressure ratio cases tend to have the jets cling to the top wall at any point in time while the lower pressure cases have the jet near the bottom wall. 6.4 INJECTION ANGLE EFFECTS As has already been shown in the previous section, injection angle is a parameter that can influence the development of the gas jet. In the previous chapter it was shown that the piston motion had a strong influence on the effects of pressure ratio changes. The purpose of this 106 C.A. as 16.3 deg. BTDC Pres. Ratio •N 5:1 0.04 Pres. Ratio •• 2:1 0.04 ^0.04 C.A. 32 deg. ATDC Pres. Ratio - 5:1 Pres. Ratio 2:1 .04 0.04 Figure 6.3.3 Concentration fields for cases 8 and 10 showing pressure ratio effects, jet angle 10 deg., engine speed 1200 RPM, contour interval 0.08 (kg/kg), injector tip in top left corners. 107 section is to determine if a change in the injection angle is also influenced by the moving piston. Figures 6.4.1 and 6.4.2 show concentration and velocity fields for jets of 10, 30 and 60 degrees angle. Similar to the fixed piston model, the moving piston model also exhibits variability in jet behaviour with changes in injection angle. In these cases however the changes due to jet angle are less abrupt than the fixed piston model. The jet behaviour is less well defined as indicated by the 30 degree plot shown in figure 6.4.1. In this case the jet spends a great deal of time between the two boundaries not approaching either wall until well after injection has started . The 60 degree jet firmly adheres to the bottom wall while the 10 degree injection stays on the top wall. It is worth noting that in the fixed piston model, the 30 degree injection always moved immediately towards the bottom wall. In the moving boundary model this is true of the 60 degree case but not necessarily in the 30 degree injection. The effects of the piston motion on the injection angle characteristics of the jet are clear. In contrast to the fixed piston model the immediate attachment of the jet to either the top wall or the bottom wall is not seen in the moving piston case. The 30 degree jet, shown in a previous section, initially moves towards the top wall and switches to the bottom wall well after injection has started. It appears as if the moving piston has some moderating influence on the effects of jet angle. The moving piston effect manifests itself as a net migration of fluid towards the top wall. This obviously causes the jet to be less sensitive to injection angle as compared to the fixed piston case. One final parameter that must be explored in context with the moving piston model is the rate of piston movement on the jet. 108 Injection Angle - 10 deg. 0.04 0.04 ■.......-.... Injection Angle = 30 deg. Injection Angle = 60 deg. Figure 6.4.1 Concentration fields for cases 8,1,and 3 showing injection angle effects, pressure ratio 2:1, engine speed 1200 RPM, crank angle 6.8 deg. BTDC, contour interval 0.08 (kg/kg), injector tip in top left corners. - —- 109 Injection Angle = 10 deg. VEL CIT SCALE (100 m/sec) = rOtarR R 11 It R R IL 11 * ^ PC R R R R IL 4. 4, 4.^# • ^♦ ^ *ff.** • 4. 1. 844014-44.4- • 4. 4- 4. 4- 4- 4- + r + • #^* 4.^• • 4-^4-^• • • •+ • 4 4-^• ♦ •^4-^4. 4-^4-^+ 44- r 4- + 4 e 4. + ( 1 Injection Angle = 30 deg. VELOCIF1 SCALE (100 misec) = Injection Angle = 60 deg. VELOCIF1 SCALE (100 m/ seal = rr 1r r 41, 4," k' r- tr- *r^ rr^Y 4,^* ^ --ft --ft -- ft ft -- L114_4, --t• --4 --P^ 4 Figure 6.4.2 Velocity fields for cases 8,1, and 3 showing jet angle effects, pressure ratio 2:1, engine speed 1200 RPM, crank angle 6.8 deg. BTDC, injector tip in top left corners. 110 6.5 ENGINE SPEED As has been demonstrated in the previous sections, the piston movement has a considerable influence on the jet dynamics in conjunction with jet angle and pressure ratio. This section looks at the effects of the rate of piston motion alone by comparing cases that differ in their engine speeds only. The rate of piston motion will determine the size of some of the source terms in the governing equation for the axial momentum. It is the effect of these terms that will be considered in this section. Figures 6.5.1 and 6.5.2 shows the concentration and velocity fields for a 30 degree jet with a piston moving at 0, 600 and 1200 revolutions per minute (RPM). Figures 6.5.3 and 6.5.4 show similar plots for a 10 degree jet. The figures appear noticeably different in nature. Comparing the moving boundary plots with the fixed boundary plots of both figures, it is clear that the motion of the piston has a substantial effect on jet development. Several points are clear about the effect of the piston speed. First, the behaviour of the 30 degree jet changes with piston speed. The plots for each speed appear different in figure 6.5.1. Secondly, the jet propagation seems inhibited by the motion of the piston. The fixed piston models (0 RPM) have jets that have penetrated into almost all of their chambers compared to the moving piston models. Finally, the 10 degree jet appears sensitive to the moving piston, but not its speed. The difference in the 600 RPM and 1200 RPM plots of the 10 degree jet are insignificant. The fixed piston model shows substantial jet propagation over the two moving piston cases. One possible explanation for the inhibited propagation of the gas jets may be the different amounts of compression between air and natural gas under a similar change in pressure. It 111 0 RPM 4 0.28 a28 600 RPM Fr 's". ' - - ■ —___ ,,, ^....../ - o°' 1200 RPM Figure 6.5.1 Concentration fields for cases 2 (fixed boundary), 6 and 1 showing the effects of engine speed, pressure ratio 2:1, jet angle 30 deg., time 5.3 ms after injection start, contour intervals 0.08 (kg/kg), injector tip in top left corners. ^ 112 0 RPM VELOCITY SCALE: ^100 m/sec 600 RPM VELOCIT r SCALE Ho rn/secl 1H_ 1 I 1 I^ r m • f^ 4 -4^ Ve^I^I^ 4^4^et^ee^f^ a-^•^et^*^•^kJ S^! a^ ^4^•^et^et^f^4-^ el.^•^*^*^• r^it^•^ 0 2. e-^r^ ja^II^ +^to 4^• et f 4 1 +^4,' * • •^• ,,' ---• —4 JO •^• k It * i pr-ttvtorrr‘r •,,, ..... .r.,,4 \,,, y, ,, , ,... ---. --• ..-•^•^A^•^•^et^f^f^ •^••^kg^a 4. -• A^x^•^Or^tx^f^4-^ .^4^k^11 * ^1. • • ••^,T , ..ekelf- ,i— kl 4. ,4 4,^4^•^se^sr^•^4.^ 4-^4^•^r^r 1,•••••••4444.4.4-4_4--4— 4-- 44, 4 --• -• -10^4^4^4^V^4^4^4^4^4 ^4 ...• 4^• :::... ^° ^ r^ •^ '1 ! ..,. ..'° -* —0 tr 1200 RPM v'ELOCI^SCALE (100 misec) = Figure 6.5.2 Velocity fields for cases 2 (fixed boundary), 6 and 1 showing the effects of engine speed, pressure ratio 2:1, jet angle 30 deg., time 5.3 ms after injection start, injector tip in top left corners. 113 0 RPM ^....•■•■•.^m . .16 ^ 0.18 ...^ 0.04... 600 RPM 12^ 0.12 0.04 1200 RPM Figure 6.5.3 Concentration fields for cases 1 (fixed boundary), 9 and 8 showing the effects of engine speed, pressure ratio 2:1, jet angle 10 deg., time 5.3 ms after injection start, contour interval 0.04 (kg/kg), injector tip in top left corners. 114 0 RPM VE OCI^A E• ^ 100 m/sec c. • 71 600 RPM VELl I^SC LE (100 m/sec) = Allir •• • • 'IgiriririCir 01. R It R It R OL t R Or ot Ot it St * 14441. KKKKKKKK 1`4.114.00ra la ft Cleat Pr g F It 0 • IL IL f• =4.4.44t K ^ K^II II F• • 41*^ F ^ 0. • • • • tRIE00••••• • • • • • • • • ((rr^1 1♦( 1 . 4 1 1 I 1200 RPM VEL I CIT , SC LE (100 m/sed = al= ^ ^R .^ i 4'404*^ 4-4-4-4-4-4- 44-4f +4.1.*^ 4-4-4-4-4-4-4-4. 4.44 ^ • • • • • • • :100.0•••••1.4. • • 11. • • • • • • •0-^•^4.^0-^*• • • • • • • • • • • • 4-^• • 4- ••• 4- • • • 44- • Figure 6.5.4 Velocity fields for cases 1 (fixed boundary), 9 and 8 showing the effects of engine speed, pressure ratio 2:1, jet angle 10 deg., time 5.3 ms after injection start, injector tip in top left corners. 115 can be shown that the amount by which a region of gas compresses under a specific change in pressure is proportional to its molecular weight. The relation shown (6.5.1), based on isentropic behaviour, indicates that the ratio of the change of specific volume (u), is proportional to the molecular weight (M) of each species. Since air has a higher molecular weight (28.9 g/mol) than natural gas (17.0 g/mol), 81) CH4 _ 8 V AIR Y AIR ( Y CH4 C 4) (1.4)( MAIR) -^ 1-3^MCH4 (6.5.1) the change in specific volume of the gas will be greater than that for air. The natural gas will therefore be inhibited from expanding outward as the pressure in the chambers increase. This effect is clearly seen in figure 6.5.3 where the moving piston cases have not propagated nearly as far as the stationary piston case. In general it is clear that the piston motion has a significant effect on jet development, especially at higher injection angles. 6.6 COMBUSTIBLE MIXTURE REGION In the fixed piston model the development of a combustible mixture region was computed for the 30 degree jet with a pressure of 2 to 1 and an engine speed of 1200 RPM. The results are shown in figure 6.6.1. Again it must be stated that combustion in an engine environment would distort the images presented in this section. These plots are intended to give an indication of the development of a combustible zone preceding combustion. The combustible limits are taken to be between 0.07 and 0.035 kg/kg mass fraction, which corresponds to a relative air/fuel ratio range of 0.8 to 1.6 respectively. These limits will change with higher temperatures in a real engine. 116 The combustible region development shown in figure 6.6.1 shows several details. First, the combustible zone shows the same inhibited growth as the jet in the previous section. This suggests that the combustible mixture zone exists imbedded in the larger jet concentration field. Secondly, the combustible zone is always a thin region in the moving piston model. The previous chapter showed that in the fixed piston model the combustible region grew considerable after the jet had developed. The compressible effects mentioned in the last section may inhibit the expansion of the combustible zone. 117 C.A. = 23.1 deg. BTDC C.A. = 16.3 deg. BTDC C.A. = 6.8 deg. BTDC 0.07 C.A. - 3.2 deg. ATDC 0•1Ais ..„._ Figure 6.6.1 Combustible region development for case 1, pressure ratio 2:1, jet angle 30 deg., engine speed 1200 RPM, injector tip in top left corners. 118 CHAPTER 7: CONCLUSIONS AND RECOMMENDATIONS The purpose of this study was to investigate the behaviour of a transient injection of natural gas into a simulated combustion chamber. The effects of a number of different parameters including injection angle, engine speed and storage-tank-to-chamber pressure ratio were investigated. A numerical model based on the conservation equations of mass, momentum and species concentration was employed to facilitate this investigation. The TEACH code developed by A.D. Gosman of Imperial College, London was used to model the injection of gas into a fixed piston and moving piston combustion chamber. The numerical code was found to be capable of predicting axial and radial velocity and concentration profiles with acceptable accuracy. The conclusions drawn from the investigation are: 1). The concentration contours calculated from the fixed piston numerical model are in approximate concordance with the photographed shape and penetration of the jet at various times. The 30 degree high pressure ratio jet showed a close correlation in jet shape and penetration to the experimental photograph. 2). The jet angle is the most influential parameter in determining the jet shape. The jet of the fixed piston model is highly sensitive to changes in jet angle. The jet switches attachment from the bottom wall to the top wall near 19 degrees angle. Low pressure zones develop adjacent to the jet that can direct it from the bottom to the top wall. The moving piston model has less sensitivity to jet angle. 3). A higher pressure ratio has significant effects only at high jet angle and in a moving piston chamber. The jet of the fixed piston model showed a small variation with pressure 119 ratio increase in the 30 degree jet and no change with the 10 degree jet. The moving piston model jet showed significant change at 30 degrees jet angle. 4). There is a significant difference in the jet shape and penetration between the moving and fixed piston models. Both the 10 and 30 degree jets showed significant change when piston motion was introduced. The engine speed has a significant effect on the jet shape, only at higher jet angles. The 30 degree jet of the moving piston model showed appreciable change over all engine speeds, while the 10 degree jet showed little variation. 5). The size of the fixed piston chamber had no significant effect on the shape of the jet or its tendency to cling to the top or the bottom wall, for the clearance gaps considered. The recommendations for further research are: 1). A proven compressibility model should be implemented in the code to account for compressibility effects, allowing piston bowl shapes to be studied. 2). Initial air motion in the cylinder created by the scavenging process should be studied through the addition of air inlet ports and exhaust valves. The air motion inside the cylinder was assumed to be non existent at the start of computation. 3). It is recommended that similar calculations be performed on a code with superior accuracy than the. TEACH code, such as KIVA, to validate the results of this study. The present results obtained are only expected to be qualitatively accurate. 120 REFERENCES [1]. Abramovich, S. and Solan, A., "The Initial Development of a Submerged Laminar Round Jet", J. Fluid Mech., Vol. 59, part iv, pp. 791-800, 1973. [2]. Abramovich, G.N., The Theory of Turbulent Jets. MIT press, Cambridge, mass, 1973. [3]. Anderson, D.A., Tannehill, J.C., Pletcher, R.H., Computational Fluid Mechanics and Heat Transfer, Hemisphere Pub. co., New York, 1984. [4]. Bassoli, C. et al., "Two Dimensional Combustion Chamber Analysis of a Direct Injection Diesel". S.A.E. Transactions paper no. 840228, 1984. [5]. Beck, N.J., "Electronic Fuel Injection for Dual Fuel Methane", S.A.E. Transactions paper no. 891652, 1989. [6]. Birch, A.D. et al., "The Structure and Concentration Decay of High Pressure Jets of Natural Gas", Comb. Sci. and Tech. vol 36, pp. 249-261, 1984. [7]. Birch, A.D., et al., "The Turbulent Concentration Field of a Methane Jet", J. Fluid Mech., vol. 88, pp. 431-449, 1978. [8]. Butler, T.D., et al., "Multidimensional Numerical Simulation of Reactive Flow in Internal Combustion Engines", Prog. Energy Comb. Sci., vol. 7, pp. 293-315, 1981. [9]. Canadian Resourcecon Ltd., "Market Assessment: Intensifier Injector for DDEC Controlled Engines" unpublished report, Sept. 1989. [10]. Diwakar, R., "Assessment of the Ability of a Multidimensional Computer Code to Model Combustion in a Homogeneous Charge Engine", S.A.E. Transactions paper no. 840230, 1984. [11]. Ewan, B.C.R., and Moodie, K., "Structure and Velocity Measurements in Underexpanded Jets", Comb. Sci. and Tech., vol. 45, pp.275-288, 1986. [12]. Gaillard, P., "Multidimensional Numerical Study of the Mixing of Unsteady Gaseous Fuel Jets with Air in Free and Confined Situations", S.A.E. Transactions paper no. 840225, 1984. [13]. Gosman, A.D., and Tsui, Y.Y., and Watkins, A.P., "Calculation of the Three Dimensional Air Motion in Model Engines", S.A.E. Transactions paper no. 1989. 121 [14]. Gosman, A.D., Johns R.J.R., Watkins A.P.,"Development of Prediction Methods for in_cylinder Processes in Reciprocating Engines" in Combustion Modelling in Reciprocating Engines,Plenum Press, New York, 1980. [15]. Heywood, J.B., Internal Combustion Engine Fundamentals, McGraw Hill co., New York, 1988. [16]. Henriot, S., and LeCoz, J.F., and Pinchon, P., "Three Dimensional Modelling of Flow and Turbulence in a Four Valve Spark Ignition Engine - Comparison with LDV Measurements", S.A.E. Transactions paper no. 890843, 1989. [17]. Hiroyasu, H., and Arai, M., "Structure of Fuel Sprays in Diesel Engines", S.A.E. Transactions paper no 900475, 1990. [18]. Hutchinson, B.R., and Galpin, P.F., and Raithby, G.D., "Application of Additive Correction Multigrid to the Coupled Fluid Flow Equations", Num. Heat Transfer, pp 133147, 1988. [19]. Kollmann, W., Computational Fluid Dynamics, Hemisphere Pub. co., Washington D.C., 1979. [20]. Kundu, P.K., Fluid Mechanics, Academic Press Inc., San Diego Ca., 1990. [21]. Kuo, T.W., and Bracco, F.V., "On the Scaling of Transient Laminar, and Turbulent Spray Jets", S.A.E. Transactions paper no. 820038, 1982. [22]. Launder, B.E. and Spalding, D.B., Mathematical Models of Turbulence, Academic Press co, London, 1972. [23]. Markatos, N.C., Computer Simulation of Mass Transfer and Combustion in Reciprocating Engines, Hemisphere Pub. co., New York, 1989. [24]. Needham, J.R., et al., "Technology for 1994", S.A.E. Transactions paper no. 891949, 1989. [25]. Ouellette, P., "High Pressure Injection of Natural Gas in Diesel Engines", M.A.Sc. Thesis, Univ. of British Columbia, 1992. [26]. Patankar, S.V., Numerical Heat Transfer and Fluid Flow, Hemisphere Pub. co., New York, 1980. [27]. Reynolds W.C., "Modelling of Fluid Motions in engines - an Introductory Overview" in Combustion Modelling in Reciprocating Engines,Plenum^Press, New York, 1980. 122 [28]. Ricou, F.P., and Spalding, D.B., "Measurements of Entrainment by Axisymetrical Turbulent Jets", J. Fluid Mech., vol. 11, pp. 21-32, 1961. [29]. Schlichting, H., Boundary Layer Theory 7th. ed., McGraw Hill co., New York, 1979. [30]. Tahry, S.H., "k-c Equation for Compressible Reciprocating Engine Flows", J. Energy, vol. 7, pp. 345-353, 1983. [31]. Wakenell, J.F., et al., "High Pressure Late Cycle Direct Injection of Natural Gas in a Rail Medium Speed Diesel Engine", S.A.E. Transactions paper no. 872041, 1987. [32]. Watkins, A.P., "Flow and Heat Transfer in Piston Cylinder Assemblies", Ph.D. Thesis, Univ. of London, London, U.K., 1977. [33]. Weaver, C., "Cost Effectiveness of Alternative Fuels and Conservational Technologies for Reducing Transit Bus Emissions in Santiago Chile", S.A.E. Transactions paper no. 891100, 1989. [34]. Weaver, C., "Natural Gas Vehicles - A Review of the State of the Art", S.A.E. Technical Report no. 892133, 1989. 123 APPENDIX A: GAS JET ANALYSIS A schematic drawing of the computational grid used to compare the numerical results to analytical and experimental results is given below. Two computation nodes at the bottom left corner have their axial velocities fixed to simulate the jet orifice. The jet inlet velocities have a uniform profile, starting at time zero and maintaining a constant speed of 10 m/s. The standard grid used has 29 control volumes in the axial direction and 25 in the radial direction. Three free slip walls and one symmetry axis form the boundaries of the grid. Free slip walls have viscous effects removed by setting the velocity inside the wall equivalent to the velocity adjacent to it. The drawing is not to scale nor does it contain the indicated number of nodes, it is merely given to provide a reference frame. COMPUTATION GRID not to scale jet orifice, 0.02 m radius 25 nodes 0.8 m /t jet orifice 29 nodes 1.2 m ^ The control volume size is set at 1 cm growing at a rate of 10% per control volume from the 124 origin. Tollrnien's and Schlichting's data The following table lists the plotted data obtained from Tollmien's and Schlichting's analysis. The original data tables of Tollmien given in Abramovich [2] have been treated to a regression analysis using the spreadsheet Lotus 123. The resulting values and those determined by Schlichting and those determined numerically are given. Y/Yc Tollmien Numerical Schlichting 0.1 0.983092 1.012758 0.988699 0.2 0.93502 0.966224 0.932514 0.3 0.885116 0.916393 0.876866 0.4 0.833718 0.863756 0.821879 0.5 0.781166 0.808803 0.767678 0.6 0.727798 0.752024 0.714387 0.7 0.673955 0.693911 0.662132 0.8 0.619975 0.634954 0.611037 0.9 0.566197 0.575644 0.561226 1 0.512962 0.516471 0.512826 1.1 0.460607 0.457926 0.465959 1.2 0.409472 0.400499 0.420751 1.3 0.359897 0.344682 0.377328 1.4 0.312221 0.290965 0.335812 1.5 0.266783 0.239838 0.29633 1.6 0.223922 0.191793 0.259006 1.7 0.183977 0.147319 0.223964 1.8 0.147288 0.106907 0.19133 1.9^0.114194 ^ 0.071049 ^ 0.161227 125 2 I 0.085033 I^0.040234^0.133782 H The equations use by Schlichting to obtain the radial velocity profile are listed below. 3 K^1 87c e_x^1 Q^(1 + -7-n 2 ) 4 _ 14Tay - .71^€0 x K=2n f u 2 ydy Experimentally it has been determined that; e ° = 0.0161 Here, y and x are the radial and axial coordinates respectively. Using the data for an axial distance of x = 0.4357 m, the previous data table was produced. Kuo and Bracco iet penetration relations The transient incompressible jet penetration relations which are given, have been used to obtain the following table. The time taken for the jet to reach 70.0% of its steady state speed at a specified stationary point is used in the following scaled parameters. 126 t Uin (R0.053 D t* - x* - ^ D ( RL" 3 Here, t is the time, D the jet starting diameter (0.04 m), U fr, is the jet starting velocity ( 10 mis), and Red is the Reynolds number based on the diameter D. Experimentally, the following relation was found to govern jet penetration. t* = 0.235x` 2 for x* Z 7 numerical^predicted t star x star t star 17.5752 7.0706 11.7486 22.5273 8.2673 16.0617 27.5079 9.6434 21.8538 36.3473 11.2259 29.6151 46.7265 13.0459 39.9956 61.1586 15.1388 53.8577 79.9567 17.5456 72.3443 100.5042 20.3135 96.9697 119.5145 23.4965 129.7403 132.3801 27.1570 173.3135 127 Birch's Radial and Axial Concentration Profiles Birch's axial concentration profile for a steady state methane turbulent jet is given in the following expression. C_ Ca^(z+a) de = d (Rs-) 112 The constant k1 is taken to be 5.05, and the constant 'a' is taken to be 0.02 m, which differs from the value specified by Birch in his experiments. The reason for this may be the inability of the numerical code to adequately model the jet outlet conditions. The densities of the gas ( p 0) and the air (p.) are 0.677 and 1.21 kg/m^3 respectively. Here, d E is the effective diameter, d the initial jet diameter, z the axial coordinate, C the concentration, and Co the initial jet concentration. Data tables obtained are given. axial distance Numerical Predicted 0.026 0.999762 1 0.0381 0.998982 1 0.05141 0.996656 1 0.066051 0.990969 1 0.082156 0.97912 1 0.099872 0.957569 1 0.119359 0.922892 1 0.140795 0.873327 0.939682 0.164374 0.810078 0.819507 0.190312 ......._ .... 0.737445 ___ . 0.718438 /1^./.106 ■•••• •■• 128 0.250227 0.588474 0.559144 0.28475 0.522306 0.495803 0.322725 0.464827 0.440867 0.364497 0.415726 0.39297 0.410447 0.372608 0.351021 0.460992 0.329336 0.314134 0.516591 0.272753 0.281585 0.57775 0.185615 0.252775 0.645025 0.077745 0.227203 0.719028 0.011601 0.204452 0.800431 0.000343 0.184167 Birch's radial concentration profile for a steady state methane jet is given in the following expression. = exp ( -D ( ) 2 ) The constant 'D' was found experimentally to be 73.6 and is used to produce the following data table at an axial distance of z = 0.46 m. Here, r and z are the radial and axial coordinates, C the concentration, and C. is the concentration on the centerline of the jet. r/z numer. C/Cm exper. C/Cm 0.010846 1 0.991379 0.032539 0.93474 0.925034 0.791267 0.814196 0.0564 ^ ^ 0.604872 0.655518 0.082648 129 0.11152 0.474357 0.400377 0.14328 0.280157 0.2207 0.178216 0.109119 0.096559 0.216645 0.029441 0.031605 0.258918 0.005616 0.007198 0.305417 0.000726 0.001043 0.356566 0.000063 0.000086 Tollmien's Axial Velocity Profile The axial velocity profile for a steady state turbulent jet developed by Tollmien and quoted by Abramovich is given as follows; U _ 0.96 U0^ax Ro The constant 'a' is given by Tollmien to be 0.066 for jets with a uniform initial velocity. The axial coordinate 'x' is adjusted to compensate for the undeveloped region of the starting jet. Here, U is the axial velocity and U 0 is the initial axial velocity. The adjustment is specified by Tollmien as: o^ax axax +^where^° - 0. 29 R,^Ro^R, For an initial radius Ro of 0.02 m, this leads to x o = 0.0878 m, confirmed experimentally. For purposes of comparison this value of xo and another at 0.125 m is plotted. The 0.125 m 13 0 value compares best to the numerical tests, and again most like indicates the numerical codes lack of sophistication in jet boundary conditions. Predicted Numerical Axial dist. 1.0000 1.0000 - 0.0050 1.0000 0.9982 0.0150 1.0000 0.9980 0.0265 1.0000 0.9987 0.0397 1.0000 0.9995 0.0549 1.0000 0.9992 0.0724 1.0000 0.9965 0.0925 1.0000 0.9884 0.1157 1.0000 0.9704 0.1423 1.0000 0.9360 0.1729 0.8735 0.8798 0.2080 0.7789 0.8021 0.2485 0.6926 0.7113 0.2950 0.6144 0.6195 0.3485 0.5437 0.5360 0.4100 0.4802 0.4642 0.4808 0.4233 0.4038 0.5622 0.3726 0.3527 0.6558 0.3275 0.3082 0.7634 0.2874 0.2661 0.8871 0.2520 0.2202 1.0294 0.2207 0.1635 1.1931 A schematic drawing of the starting region of the steady state jet showing the initial offset 131 x o is given. Location of Virtual Origin with respect to the jet orifice at x 0. 0 start speed 10 m/s Ro = 0.02 m
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Numerical analysis of high pressure injection of natural...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Numerical analysis of high pressure injection of natural gas into diesel engine combustion chambers Walsh, Paul 1991
pdf
Page Metadata
Item Metadata
Title | Numerical analysis of high pressure injection of natural gas into diesel engine combustion chambers |
Creator |
Walsh, Paul |
Date Issued | 1991 |
Description | The behaviour of a transient turbulent jet of natural gas as it is injected into a simulated combustion chamber of a diesel engine was investigated using numerical techniques. The TEACH code developed by A.D. Gosman of Imperial College, London, was used to investigate the influence of parameters such as injection angle, engine speed, and reservoir tank-to-chamber pressure ratio on the development of the jet. It has been shown that the TEACH code is fully capable of predicting details of jet behaviour such as radial and axial velocity and concentration profiles when compared to known data. The code has been modified to use a compressing and expanding grid to simulate the effects of piston motion. A model of a fixed geometry combustion chamber revealed that the most influential parameter on jet behaviour is the injection angle. The jet had a tendency to adhere to either the top wall of the chamber near the injector tip or to the bottom wall directly opposite depending on the injector angle. The compressing grid simulation showed that the presence of piston motion combined with other parameters such as injection angle and pressure ratio, produced jet characteristics that were dissimilar compared to the fixed boundary model. In general it was shown that the jet was less sensitive to injection angle and strongly influenced by increased pressure ratio as a result of the moving boundary. |
Extent | 5067901 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
File Format | application/pdf |
Language | eng |
Date Available | 2008-09-12 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0080886 |
URI | http://hdl.handle.net/2429/1925 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
Graduation Date | 1992-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
Aggregated Source Repository | DSpace |
Download
- Media
- 831-ubc_1992_spring_walsh_paul.pdf [ 4.83MB ]
- Metadata
- JSON: 831-1.0080886.json
- JSON-LD: 831-1.0080886-ld.json
- RDF/XML (Pretty): 831-1.0080886-rdf.xml
- RDF/JSON: 831-1.0080886-rdf.json
- Turtle: 831-1.0080886-turtle.txt
- N-Triples: 831-1.0080886-rdf-ntriples.txt
- Original Record: 831-1.0080886-source.json
- Full Text
- 831-1.0080886-fulltext.txt
- Citation
- 831-1.0080886.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0080886/manifest