Experimentation and Modeling of Formation in a Small Turbo-Charged Diesel Engine by Matthew Blake Thames A thesis submitted to the Graduate Faculty of Auburn University in partial fulfillment of the requirements for the Degree of Master of Science Auburn, Alabama May 5, 2013 Copyright 2013 by Matthew Blake Thames Approved by: Song-Yul Choe, Professor of Mechanical Engineering David Beale, Professor of Mechanical Engineering Daniel Mackowski, Associate Professor of Mechanical Engineering ii Abstract In this thesis a detailed investigation into the primary formation characteristics and modeling approaches for the formation of in turbo-charged diesel engines is presented. Through this investigation a multi-layered adiabatic core model is developed and implemented into a single-zone internal engine combustion model. The comprehensive engine model includes estimation of heat release rates, heat transfer rates, cylinder pressure, cylinder temperature, and cylinder thermodynamic properties output on a crank angle basis. The model investigation includes an analysis of different reaction mechanisms including the Extended Zeldovich Mechanism, Mechanism, Mechanism, NH/HNO Mechanism, and variations/extensions of these mechanisms. Along with reaction mechanisms, other parameters that are analyzed include: temperature approximations, sensitivity to rate constants, validation of the quasi-steady state approximation, and the effect of multiple injections on simulation results. Common modeling approaches that were tested for their validity include fully mixed models, multi-layered approaches, and lumped modeling approaches. The developed model is a semi-empirical model including: an empirical approximation for the formation of based upon operational characteristics, an exhaust gas recirculation model, and a combustion species equilibrium model. The EGR model is based upon burned gas composition properties which models charge dilution through mixture heat capacity and enthalpy and their respective effects on adiabatic flame temperature. The combustion equilibrium model is developed through the use iii of the minimization of Gibbs Free Energy and is similar to that of NASA?s CEA program. The equilibrium combustion model can produce species concentration of any specified species set with or without EGR, and includes model outputs such as mixture, enthalpy, entropy, frozen specific heat, and equilibrium specific heat. It was found through the development of the model that the addition of an extra 43 reactions had a minimal effect on model accuracy. The rate constant analysis confirmed that the leading reaction in the Extended Zeldovich Mechanism is the most important reaction and serves as the rate limiting reaction for the mechanism. The quasi-steady-state assumption was validated showing only a 1~3% difference in comparison to a kinetically controlled reaction set. Pre-pilot injection scheduling is shown to have a minimal effect on simulation results with the majority of formation originating from pilot and main injection scheduling. The resulting model incorporates the sub-models described above coupled with the use of the Extended Zeldovich Mechanism, an approximated unburned air temperature between the ideal unburned air temperature and the mixture average, full injection scheduling, and a layered adiabatic core approach. The resulting model is shown to predict concentrations, with and without EGR, to within 13% of experimental data. Additional model outputs include relative concentration predictions of carbon dioxide and oxygen based upon complete combustion principles. iv Acknowledgements I am grateful for the love and support of my parents, Grady and Patti Thames. I would not have been able to complete this project without their support and encouragement. I would also like to thank Taylor Wingo for his support in the completion of the test stand and the knowledge he brought to the modeling aspects of internal combustion engines. I would also like to thank the National Center for Asphalt Technology and Dr. Buzz Powell for their support over the duration of my masters program. I would like to extend thanks to Dr. Song-Yul Choe, Dr. David Beale, and Dr. Daniel Mackowski on my committee. I am appreciative of the opportunity to extend my knowledge of internal combustion engines and the guidance that Dr. Choe has given me over the years. Above all I would like to thank God for giving me the ability and opportunity for the advancement of my education. v Table of Contents Abstract ........................................................................................................................................... ii Acknowledgements ........................................................................................................................ iv List of Terms ................................................................................................................................ xiii Abbreviations and Acronyms ...................................................................................................... xiii Equation Symbols ..........................................................................................................................xv List of Figures ................................................................................................................................ ix List of Tables ................................................................................................................................. xi Chapter 1: Introduction ...............................................................................................................1 1.1. Overview and History of Diesel Engine .............................................................................1 1.2. Characterization of Diesel Emissions .................................................................................2 1.3. Oxides of Nitrogen .............................................................................................................4 1.4. Diesel Combustion Process ................................................................................................6 1.5. Oxides of Nitrogen Formation in Diesel Engines ...............................................................8 1.6. Motivation for Research ...................................................................................................11 1.7. Specific Objectives of this Research ................................................................................11 Chapter 2: Literature Review ....................................................................................................13 2.1. Emissions Modeling .......................................................................................................13 2.2. Modeling.................................................................................................................14 2.3. Reaction Mechanisms ............................................................................................14 vi 2.3.1. The Zeldovich Mechanism .............................................................................14 2.3.2. The Prompt NO Mechanism ...........................................................................15 2.3.3. The Mechanism ......................................................................................17 2.3.4. The Mechanism ......................................................................................19 2.3.5. Fuel Bound Nitrogen ......................................................................................20 2.3.6. Other Mechanisms ..................................................................................21 2.4. Modeling Coupled/Decoupled with Combustion ...................................................21 2.5 Detailed Mechanisms and Reduced Mechanisms ...........................................23 2.6 Modeling Formation in Diesel Engines ..................................................................25 Chapter 3: Model Development .................................................................................................30 3.1. Introduction .....................................................................................................................30 3.2. Combustion Model .........................................................................................................30 3.2.1. Governing Equations ..............................................................................................31 3.2.2. Heat Release Model ................................................................................................32 3.2.3. Heat Transfer Model ...............................................................................................32 3.2.4. SOC Model .............................................................................................................33 3.2.5. Supporting Models ..................................................................................................34 3.3. Definition of NO Base Model .........................................................................................35 3.3.1. Equilibrium Composition Model ............................................................................35 3.3.2. Chemical Kinetics Calculation Principle ................................................................49 3.3.3. Adiabatic Layered Core Model ...............................................................................55 3.3.4. EGR Model .............................................................................................................59 3.4. Model Perturbations ........................................................................................................65 vii 3.4.1. Reaction Mechanisms .............................................................................................65 3.4.2. Individual Species Treatment .................................................................................66 3.4.3. Reaction Rate Constants .........................................................................................67 3.4.4. Temperature Approximations .................................................................................69 3.4.5. Multiple Injections ..................................................................................................71 3.5. Model ......................................................................................................................72 3.6. and Model ..........................................................................................................76 3.7. Program Structure ...........................................................................................................77 3.8. Solver for Ordinary Differential Equations ....................................................................80 Chapter 4: Experimental Test Stand ........................................................................................81 4.1. Test Stand and Dynamometer ..........................................................................................81 4.2. Sensors .............................................................................................................................83 4.3. Emission Analyzer ...........................................................................................................84 4.4. Data Monitoring and Recording ......................................................................................86 4.5. EGR Control ....................................................................................................................87 Chapter 5: Experimental Results ..............................................................................................89 5.1. Test Procedure and Test Matrix ........................................................................................89 5.2. Load Sweep Results ..........................................................................................................91 5.3. EGR Sweep Results ..........................................................................................................94 Chapter 6: Parameter Identification for NO Model ................................................................95 6.1. Introduction ......................................................................................................................95 6.2. Temperature Approximation Analysis .............................................................................95 6.3. Rate Constants Analysis ................................................................................................106 viii 6.4. Reaction Mechanism Analysis ......................................................................................107 6.5. Validation of the QSSA Assumption .............................................................................111 6.6. Effect of Multiple Injections on Simulation Results ......................................................113 6.7. Effect of Initial NO Concentration from Residual Gas .................................................116 Chapter 7: Final Model Validation .........................................................................................118 7.1. Introduction ...................................................................................................................118 7.2. Parameters for NO Layered Adiabatic Core Model .....................................................118 7.3. Validation of the NO Layered Adiabatic Core Model ..................................................119 7.4. Validation of the EGR Model .......................................................................................122 7.5. Validation of Model .............................................................................................125 7.6. Analysis of Lumped Model ..................................................................................126 7.7. Validation of and Model ..................................................................................128 7.8. Implementation of Emission Models in Single-Zone Combustion Model ...................130 Chapter 8: Conclusions and Recommendations ....................................................................134 References ....................................................................................................................................137 Appendix A ..................................................................................................................................142 Appendix B ..................................................................................................................................145 ix List of Figures Figure 1-1: Engine Heat Release Rate [9] ......................................................................................7 Figure 2-1: Prompt NO Formation Pathways ...............................................................................17 Figure 2-2: Relative concentrations of NO and in IC engines [9] .........................................19 Figure 2-3: Relative contributions of mechanisms [39] .......................................................26 Figure 3-1: Layered Temperature Profiles ....................................................................................57 Figure 3-2: Layer NO Concentrations ..........................................................................................58 Figure 3-3: Total NO Concentration Profile .................................................................................59 Figure 3-4: Simulated Effects of Charge Dilution on Flame Temperature ..................................63 Figure 3-5: Comparison of Unburned Air Temp. to Average Cylinder Temp .............................70 Figure 3-6: Temperature Approximations ....................................................................................70 Figure 3-7: Heat Release and Injector Current .............................................................................71 Figure 3-8: Lumped Model ...................................................................................................73 Figure 3-9: Data Trend .........................................................................................................75 Figure 4-1: Experimental Test Cell...............................................................................................82 Figure 4-2: Experimental Test Bed Setup .....................................................................................82 Figure 4-3: ECOM J2KN Emissions Analyzer .............................................................................85 Figure 4-4: EGR Control Apparatus .............................................................................................88 Figure 5-1: IMEP vs. NO ..............................................................................................................93 Figure 5-2: IMEP vs. ............................................................................................................93 Figure 6-1: Tested Temperature Profiles ......................................................................................96 Figure 6-2: NO Production with Unburned Air Temp. ................................................................98 Figure 6-3: Approximated Error with Progressive Temperature Profiles ....................................99 x Figure 6-4: Progression of Tested Temperature Profiles ..............................................................99 Figure 6-5: Finite Temperature Analysis ....................................................................................100 Figure 6-6: Progressive Error of Finite Temperature Analysis ..................................................101 Figure 6-7: Approximated Temperature Profile .........................................................................102 Figure 6-8: NO production with Temperature Profile #14 .........................................................103 Figure 6-9: NO Layer Profile Results with Approximated Air Temp. .......................................104 Figure 6-10: NO Layer Profile Results with Unburned Air Temp. ............................................104 Figure 6-11: NO Layer Profile Concentration Results with Approximated Air Temp. .............105 Figure 6-12: NO Layer Profile Concentration Results with Unburned Air Temp. ....................105 Figure 6-13: Kinetically Controlled NO production vs. QSSA Controlled NO production .......112 Figure 6-14: Kinetically Controlled NO production vs. QSSA Controlled NO production .......112 Figure 6-15: NO production with Pilot and Main Injections ......................................................113 Figure 6-16: Heat Release and Injector Current Signal ..............................................................114 Figure 6-17: NO production with Pre-Pilot, Pilot, and Main Injections .....................................115 Figure 6-18: NO Layer Concentrations ......................................................................................115 Figure 7-1: NO Results from Layered Routine ...........................................................................121 Figure 7-2: EGR Results for NO ................................................................................................124 Figure 7-3: Concentration Simulation Results ..................................................................126 Figure 7-4: Concentration Simulation Results ..................................................................129 Figure 7-5: Concentration Simulation Results ....................................................................129 Figure 7-6: NO Results from coupling NO Model and Combustion Model ..............................130 Figure 7-7: Predicted Pressure and Measured Pressure ..............................................................131 Figure 7-8: Predicted Temperature and Measured Temperature ................................................131 Figure 7-9: EGR NO Results from coupling EGR Model and Combustion Model ...................132 Figure 7-10: Results from coupling Model and Combustion Model .........................133 xi List of Tables Table 1-1: Characteristics of NO and .....................................................................................4 Table 1-2: Exposure Limits for NO [12] ........................................................................................5 Table 1-3: Exposure Limits for [12] .......................................................................................5 Table 3-1: Coefficients for lean and rich combustion ..................................................................61 Table 3-2: Various Rate Constants for EZM ................................................................................68 Table 4-1: Experimental Engine Specifications ...........................................................................81 Table 4-2: Mean operating Sensors (Slow Time Base) ................................................................83 Table 4-3: High Time Base Sensors .............................................................................................84 Table 4-4: ECOM J2KN Technical Specifications .......................................................................85 Table 5-1: Test Procedure Outline ................................................................................................90 Table 5-2: Test Matrix ..................................................................................................................90 Table 5-3: Load Sweep Results (NO) ...........................................................................................91 Table 5-4: Load Sweep Results ( ) .........................................................................................92 Table 5-5: NO Load Sweep Results..............................................................................................94 Table 5-6: Load Sweep Results ............................................................................................94 Table 6-1: Unburned Charge Air Temp. NO Results ...................................................................97 Table 6-2: NO results of Temp Profile #14 ................................................................................102 Table 6-3: Scheme 1 NO Results ................................................................................................108 Table 6-4: Scheme 2 NO Results ................................................................................................108 xii Table 6-5: Scheme 4 NO Results ................................................................................................109 Table 7-1: NO Results from Layered Routine ............................................................................120 Table 7-2: NO Results from Fully Mixed Routine .....................................................................120 Table 7-3: Corrected NO Results from Layered Routine ...........................................................121 Table 7-4: EGR Results for NO ..................................................................................................123 Table 7-5: Lumped Model Results .....................................................................................127 xiii List of Terms Abbreviations and Acronyms Oxides of Nitrogen AFBR Apparent Fuel Burning Rate BGE Burned Gas Element CI Compression Ignition CLV Calculated Load Value CN Cetane Number CSP Computational Singular Perturbation DE Diesel Exhaust ... Carried on to the next Line EA Excess Air ECU Engine Control Unit EGR Exhaust Gas Re-Circulation EPA Environmental Protection Agency EVO Exhaust Valve Open EZM Extended Zeldovich Mechanism FARG Fuel-Air Residual Gas Mixture FMBR Fuel Mass Burn Rate G Gibbs Free Energy GPM Gallons Per Minute H Enthalpy HC Hydrocarbon xiv HRR Heat Release Rate IC Internal Combustion ILDM Intrinsic Low-Dimensional Manifold IMEP Indicated Mean Effective Pressure MAF Mass Air Flow MATLAB Numerical Computation Program ODE Ordinary Differential Equation ODE Ordinary Differential Equation ode15s Matlab Built In ODE solver ode23s Matlab Built In ODE solver PAH Poly-Nuclear Aromatic Hydrocarbon PE Partial Equilibrium Approximation PELs Permissible Exposure Limits PM Particulate Matter PPM Parts Per Million PWM Pulse Width Modulation QSSA Quasi-Steady State Approximation S Entropy SEZM Super-Extended Zeldovich Mechanism SI Spark Ignition SOC Start of Combustion SOI Start of Injection T Temperature VGT Variable Geometric Turbocharger ZM Zeldovich Mechanism xv Equation Symbols Ignition Delay Transport Delay ? Mean Piston Speed Activation Energy Clearance Volume Compression Ratio Cylinder Head Surface Area Discharge Coefficient Nozzle Minimum Area Fuel Density Engine Speed Mass Air ? Gibbs Function Number of Moles (species i) nSp Number of Species Mole Fraction Universal Gas Constant Equilibrium Constant ? Specific Heat Equivalence Ratio Number of Moles of (species) nEL Number of Elements Number of Atoms of Element j in species i xvi Lagrangian Multiplier Equilibrium Specific Heat Frozen Specific Heat Reaction Specific Heat Stoichiometric Coefficient of Reactants Stoichiometric Coefficient of Products ? Net Production Rate Progress Variable Forward Rate Constant Enhanced Third Body Coefficients Reverse Rate Constant Enthalpy of Unburned Mixture Enthalpy of Burned Mixture [ ] NO Concentration of Layer j Mass Burned in Layer j Residual Fraction Residual Mass Inducted Mass per Cycle Mass EGR Residual Gas Fraction Molecular Weight Mole Fraction Total Molecular Weight Residual Mole Fraction 1 Chapter 1: Introduction 1.1 Overview and History of Diesel Engine The diesel engine was invented by Rudolph Diesel in 1893. R. Diesel received patent no. 608,845 for the design of a highly efficient compression ignition, internal combustion engine. R. Diesel?s design caused much rivalry among energy conversions systems at this time due to its clear superiority to the steam and coal industry. Diesel engines extract the greatest amount of energy from a specific amount of fuel which makes them superior to other energy conversion systems. Diesels achieve their high performance and excellent fuel economy from compressing air to high pressures then injecting a small amount of fuel. The temperatures created in this region cause the highly atomized fuel droplets to evaporate and mix with surrounding air. This mixture of fuel and air then reaches its auto-ignition temperature and burns releasing the chemical energy stored in the fuel. In terms of piston and cylinder arrangement reciprocating engines there are two types: two-stroke engines and four-stroke engines. Two stroke engines complete a combustion cycle in two strokes, while four-stroke engines complete a combustion cycle in four strokes, hence the name two-stroke and four-stroke. There are two dominant four- stroke diesel engine designs which include indirect-injection and direct-injection engines, or more commonly known as IDI and DI engine designs. In IDI diesel engines fuel is injected into a pre-chamber separated from the main combustion chamber via a narrow passage or nozzle. During the compression process air is forced through this nozzle from the cylinder into the pre- chamber at high velocity. Fuel injected into this highly turbulent flow of air, achieves very high mixing rates. Combustion starts in the pre-chamber and the resulting increase in pressure causes the movement of burning gases, fuel, and air into the main combustion chamber. This rapid mixing then occurs in the main chamber as the burning gases mix with the remaining air and 2 combustion is completed. One advantage of the IDI design is lower noise characteristics but the presence of the passage way between the main combustion chamber and pre-chamber increases pumping losses, thus lowering overall mechanical efficiency. In DI diesel engine design fuel is injected directly into the main combustion chamber. The direct injection of fuel directly into the main combustion chamber places extreme importance on air/fuel mixing characteristics in the role of the overall combustion efficiency. This importance placed on air/fuel mixing characteristics has lead to extensive research and development in these areas (i.e. injection pressure, injector nozzle hole geometry, injection profiles, piston geometry etc.) to optimize engine performance and exhaust emissions. Overall DI diesel engines have 10% to 15% better fuel economy than IDI designs. Advantages of diesel engines in respect to other energy conversion systems include lower friction losses, fuel safety characteristics, durability, reliability, and lower pumping losses [12]. In the early 1900?s diesel engines were used in locomotives, submarines, ships, large trucks, and power generation. The first light duty on road applications were introduced in the 1970?s Wingo [42]. Recent reports from regulatory agencies in the United States, such as the Environmental Protection Agency (EPA), indicate that in spite of very restrictive emission regulations, diesel sales in the U.S. light duty truck segment are growing at the rate of 4% per year [12]. Modern diesel advancements include variable geometric turbochargers (VGT), high pressure common rail injection, and exhaust gas recirculation (EGR), all controlled through the engine control unit (ECU). 1.2 Characterization of Diesel Emissions Diesel engine combustion theoretically produces only , , , and . Although diesel exhaust is primarily composed of these species, the complex combustion process produces a wide variety of pollutants that can have adverse health and environmental effects. The 3 combustion processes responsible for the production of these pollutants include incomplete combustion of fuel, reactions between mixture components under high temperature and pressure, and combustion of non-hydrocarbon components of diesel fuel, such as sulfur compounds or various fuel additives [12]. Due to the potential of adverse health and environmental effects of diesel exhaust, the United States, Europe, and other countries have implemented numerous pollutant restrictions and regulations. The common regulated emission species include diesel particulate matter (PM), Nitrogen Oxides ( , Hydrocarbons (HC), and Carbon monoxide (CO). Due to the differences between the combustion processes of spark ignited engines and diesel engines, the pollutants that result from these combustion processes are also quite different. Generally, CO, HC, and emissions are higher in SI engines while diesel engines have significant higher concentrations of PM emissions which are almost obsolete in SI engine exhaust. The majority of today?s diesel pollutant control technologies primarily focus on reducing and PM. There are also non-regulated pollutants contained in diesel exhaust which usually have concentrations much lower than that of the regulated emission species. These non-regulated emission species include: ? PAH, poly-nuclear aromatic hydrocarbons, which are heavy organic compounds found mostly in the PM but some PAHs are also present in the gas phase. ? SOF, soluble organic fraction, which constitutes part of the PM ? Aldehydes, R-CHO, which are derived from hydrocarbons ? Sulfur dioxide, , from sulfur present in the fuel ? Nitrous oxide, 4 ? Metal oxides: Several engine lubricating oil additives include metalorganic compounds resulting in some metal oxide emissions including phosphorus, zinc, and calcium. 1.3 Oxides of Nitrogen Oxides of Nitrogen ( emissions in diesel engines include both NO and , with concentrations ranging from 50-1000ppm. NO is a colorless and odorless gas that originates from the oxidation of atmospheric nitrogen at high temperature and pressure. is a toxic red-brown gas characterized by an unpleasant odor. has been known to be extremely reactive and exhibit strong oxidation properties. In older technology engines, approximately 95% of nitrogen oxides were composed of NO and only 5% of . The relative percent of in the total concentration in newer, turbocharged diesel engines can be as high as 15% or more [12]. Some general characteristics of these pollutants are shown in Table 1-1. Name Nitric Oxide Nitrogen Dioxide Formula Formula Weight 30.01 46.01 Appearance Colorless Gas Red-Brown Gas Density 1.0367 - Melting Point Boiling Point Table 1-1: Characteristics of NO and measurements are commonly obtained with the Chemiluminescence analyzer. The CLD utilizes the chemical reaction between nitric oxide and ozone. Ozone is produced in the instrument and introduced to the nitric oxide species. The excited nitrogen dioxide molecules spontaneously return to their normal state. This transition is accompanied by the emission of red 5 light with photon emissions. This light is then filtered and measured by a photomultiplier. Intense focus by air pollution boards have led to numerous regulations regarding emissions from internal combustion engines. In present day regulations, revisions have led to a 90% decrease in the concentration of regulated emission species from uncontrolled levels. The following are regulations proposed by the EPA for 2007 through 2010 [15]. Exposure to NO can result in symptoms including irritation of the eyes, nose, throat, unconsciousness, and drowsiness. Exposure to can result in coughing, mucoid frothy sputum, dyspnea, chest pain, pulmonary edema, cyanosis, tachycardia, or eye irritation [12]. A table of exposure limits for , and can be seen in Tables 1-2 and 1-3. Authority Limit Value Occupational Safety and Health Limits U.S. OSHA Permissible Exposure Limit(PEL), 8-h average 25 ppm ACGIH Threshold Limit Value (TLV) 8-h average 25 ppm Table 1-2: Exposure Limits for NO Authority Limit Value Occupational Safety and Health Limits U.S. OSHA PEL, Ceiling Value 5 ppm ACGIH STEL, 15-min average 5ppm ACGIH Thershold Limit Value (TLV), 8-h average 3 ppm Ambient Air Quality Limits EPA (NAAQS) Annual average 0.052 Enviroment Canada Annual average (maximum acceptable) 0.050 Enviroment Canada 24-h average (maximum acceptable) 0.110 Enviroment Canada 1-h average (maximum acceptable) 0.210 Table 1-3: Exposure Limits for 6 Nitrogen oxides are one of the most researched pollutants due to their role in the formation of low level ozone or smog. Smog includes many primary and secondary air pollutants such as nitrogen oxides, sulfur oxides, hydrocarbons, carbon monoxide, and nuclei mode particulates. Photochemical smog also produces ground level ozone as a secondary pollutant from reactions between nitrogen oxides and hydrocarbons. Acid rain is another important secondary pollutant produced by the dissolution of acidic gases in water. These acidic gases include , , and hydrochloric acid. Acid rain is believed to play an important role in acidification of lakes and water ways, death of forests, and damage to building materials. The strong oxidation properties of some diesel exhaust can cause significant reduction in structural integrity and aesthetic quality of buildings and other structures. 1.4 Diesel Combustion Process Before discussing formation in diesel engines an overview of the combustion process will be given first. The diesel engine combustion process is an unsteady, heterogeneous, three-dimensional process which in its very nature is complex. The combustion process can be defined to begin with the start of injection near the end of the compression stroke. The liquid fuel injected through the injector nozzle at high velocity atomizes into small droplets. These droplets then vaporize and mix with the high temperature cylinder air. The temperature of this mixture of fuel and air are above the fuel?s ignition point which in turn causes spontaneous ignition after a few crank angle degrees. The cylinder pressure then increases due to the combustion of the fuel-air mixture. The ensuing compression of the unburned portion of cylinder charge shortens the delay for ignition of the unburned portion which has mixed to within combustible limits. The processes of atomization, vaporization, fuel-air mixing and 7 combustion continue until fundamentally all the fuel has passed through each process. The mixing of air remaining in the cylinder with already burning and burned gases continue through the expansion stroke. To summarize this process the overall combustion of a compression ignition engine can be viewed in relation to the heat release in four stages: ignition delay, premixed or rapid combustion phase, mixing-controlled combustion phase, and late combustion phase. The phases in relation to the heat release can be seen in Figure 2-1. Figure 1-1: Engine Heat Release Rate [9] The ignition delay period is the time between the start of injection and the start of combustion. The ignition delay period is one of the most critical stages affecting the overall combustion process, engine noise, and exhaust emissions. Ignition delay is best identified from the change in slope of the heat-release rate, determined from the cylinder pressure data, which occurs at ignition [9]. The ignition delay period involves both physical and chemical processes. Some of the physical processes also known as the mixture preparation include, spray breakup, 8 atomization, vaporization of atomized liquid fuel, fuel-air mixing, and dispersion of fuel vapor with air to within combustible limits. The chemical processes are the pre-combustion reactions of the fuel, air, and residual gas mixture which lead to auto-ignition [9]. The second phase, premixed or rapid combustion phase, is the combustion of fuel-air mixtures that have mixed to within flammability limits during the ignition delay period. This premixed charge burns at a very high rate and is associated with the characteristic noise of diesel engines. Once the fuel-air mixture within flammability limits in the premixed phase has been consumed, the remaining combustion is controlled by the rate at which the mixture becomes available for burning, also known as the mixing-controlled combustion phase. The rate at which this mixture becomes available for burning is controlled by atomization, vaporization, and the mixing of fuel vapor with available air. After the mixing-controlled combustion phase comes the late combustion phase. This phase is characterized by a lower rate of heat release which continues well into the expansion stroke. Professor John B. Heywood [9] suggested several reasons for this slow burn- rate including: a small fraction of fuel may have not burned, a fraction of the fuel energy is present in soot and fuel-rich combustion products and can still be released, the cylinder charge is non-uniform and mixing during this period promotes more complete combustion and less- dissociated product gases, and the kinetics of the final burnout processes become slower as the temperature of the cylinder gases fall during expansion. 1.5 Oxides of Nitrogen Formation in Diesel Engines It can be concluded from the explanation of the diesel combustion process above that the diesel combustion process is itself a very complex and heterogeneous process. Unlike spark ignition engines whose mixture of fuel and air are predominantly homogenous, the combustion 9 of the diesel engine?s mixture depends on the uneven distribution of fuel-air ratios. There are several different factors that affect emission formation in diesel engines such as the physical and chemical processes during the ignition delay, stratification of combustion temperatures, the residence time of burned elements, and specific engine advancements such as exhaust gas recirculation, variable geometric turbine, injection timing, and exhaust gas treatment systems. Due to the number of factors that affect emission formation, sub-models must be created to replicate the effect of such factors. Unlike S.I engines fuel in D.I diesel engines is injected just before combustion begins. The non-uniformities in the fuel distribution result in non-uniform burned gas temperatures and composition. These non-uniformities are more extensively complicated by the premixed phase following the ignition delay where fuel air mixtures with a composition close to stoichiometric burns due to spontaneous ignition and flame propagation. During the mixing controlled combustion phase the burning mixture is likely to be closer to stoichiometric [9]. Further complexities of diesel combustion occur because of the mixing between already burned gases, air, and lean and rich unburned fuel vapor-air mixture, which changes the composition of any elements that burned at a specific equivalence ratio. Also temperature changes occur due to compression and expansion as cylinder pressure rises and falls. The critical equivalence ratio for formation in high-temperature/pressure burned gases is close to stoichiometric [9]. Due to the compression of the cylinder gases burned early in the combustion process these burned elements increase in temperature which in turn, increase NO formation rate. At peak pressure both the decrease in temperature and mixing with cooler gases cause the chemistry to freeze. Because of the large compression ratios of diesel engines this freezing occurs more rapidly in the diesel engines than spark ignited engines [9]. Professor John 10 B. Heywood [9] performed rapid gas sampling valve in quiescent DI diesel engine and formulated the following conclusions: 1) Local concentrations rise from the residual gas value following SOC, to a peak point where the local burned gas equivalence ratio changes from rich to lean. 2) As local burned gas equivalence ratio becomes leaner due to mixing with excess air, concentrations decrease since formation becomes much slower as dilution occurs. 3) As injection timing is retarded the combustion process is also retarded causing concentrations to be lower due to the lower peak temperatures and also from the retarded formation process. 4) The relevance of temperature to formation can be seen in these observations, thus similarly, at high load, with higher peak pressures and temperatures and larger regions of close-to-stoichiometric burned gases levels increase. 5) The amount of fuel injected decreases proportionally as the overall equivalence ratio is decreased, thus much of the fuel still burns at stoichiometric. This shows that emissions should be roughly proportional to the mass of fuel injected. The formation process and total engine out emissions can be reduced by the addition of dilutents. The effects of diluents such as EGR, nitrogen, and argon primarily have the effect of lowering burned gas temperatures by increasing the heat capacity of the mixture. The addition of oxygen will in turn have the opposite effect of increasing flame temperatures therefore 11 increasing emissions [9]. So it is therefore necessary to include the effects of EGR and other dilutents in the calculation of engine out emissions. 1.6 Motivation for Research Over the years pollutant emissions from automobiles, industrial processes, and power plants have become an increasingly interesting field of study. The control of these pollutants has played major roles in the design of combustion systems. Due to the health and environmental effects summarized in Appendix A, there is a dire need for detailed understanding of the formation processes of primary and secondary pollutants. In the past, diesel emissions predictions were based on purely empirical models that would require many hours of engine testing and operational costs. Although the results from these tests were accurate at the time, many new technological advances such as advanced injection timing control, intake charge cooling, injection pressure control, air management, oil consumption control, variable valve actuation technology, water injection, and engine accessories/subsystem component technology require more detailed formation mechanisms. By using a computer model to simulate the effects of pollutant emissions and emerging technologies as a function of engine operating conditions versus extensive dynamometer testing, new engine design and control methods can be explored at a fraction of the cost and time. 1.7 Specific Objectives of this Research The main objective of this research is to construct an emissions model that can be easily implemented into the existing diesel combustion model. Given the emissions equipment 12 provided and the importance of emission regulations, this thesis will primarily focus on modeling. The primary objective of this research is through a detailed investigation encompassing all aspects affecting both formation and destruction; develop a model that shows good predictability across a wide range of operating conditions while maintaining low empirical input. Some of the techniques used and applied in this study have been used in other studies but not directly parallel with the engine characteristics of this study. The overall model will need to be computationally fast, non-engine specific, and represent the physical and chemical processes of formation accurately. 13 Chapter 2: Literature Review 2.1 Emission Modeling Due to increasingly stringent emission standards much emphasis in recent years has been placed on the predictive capability of emission models. The development of such models capable of predicting relationships between engine design, operating conditions, performance, and emissions has proven to be a difficult goal due to the complexity of the problem, especially if advanced control systems are considered. Other objectives that are hoped to be reached through the study of emission modeling is improved understanding of detailed chemical and physical processes controlling formation and destruction of pollutant species. During the 1960?s the availability of high-speed computers became a reality and the interest in engine cycle simulation grew exponentially. With the emerging emission problem, many researchers grew more and more interested in incorporating predictive emissions capabilities. The first predictive emission models for diesel engines began to appear around the 1970?s. Emissions modeling can be generalized by two approaches phenomenological modeling and detailed multidimensional modeling. Multidimensional modeling is characterized by solving partial differential equations governing fluid, heat, and mass transfer. However the complex kinetic mechanisms that control species emission rates coupled with complex CFD code can increase computational times to unreasonable levels. It is because of this reason that the phenomenological approach has been more intriguing to researchers. Phenomenological models benefit by only focusing on those areas most closely related to the phenomenon of interest. Most researchers in today?s time rely almost exclusively on the phenomenological approach. 14 2.2 Modeling Over the years extensive research has been undertaken in order to develop models that better our understanding of formation in IC engines. Due to the extensive amount of research in formation, many different pathways for formation and destruction have been proposed. In present time these diverse modeling approaches have led to models including from as little as three reactions to over two-hundred reactions. These diverse models are spawn from the complexity of the combustion process which varies between different engines. The complex nature of the diesel combustion process can lead to being formed from incoming oxygen and nitrogen, from nitrogen-bound fuel, from the flame zone, and from recombination reactions between intermediates, radicals, and other various species produced during the combustion process. 2.3 Mechanisms There have been numerous suggested pathways for the formation of in internal combustion engines. The following section provides a review of these different mechanisms including the, thermal mechanism (EZM), nitrous oxide mechanism, nitric dioxide mechanism, HNO mechanism, NNH mechanism and others. 2.3.1 The Zeldovich Mechanism The Zeldovich mechanism otherwise known as the thermal mechanism is by far the most common modeling approach for the formation of NO. This mechanism dominates in high- temperature combustion over a wide range of equivalence ratios. The original Zeldovich mechanism included two chemical reaction equations but has been extended to include a third 15 reaction containing OH radicals thus rendering the name the Extended Zeldovich Mechanism (EZM). The general reaction scheme of this mechanism is illustrated below. The forward and reverse reaction rates have been extensively evaluated in numerous studies over a wide range of temperatures and equivalence ratios. The EZM is coupled with the fuel combustion chemistry through the , , and species but in instances where fuel combustion is complete before NO formation becomes significant the two processes can be assumed uncoupled and , , , and can be assumed to be at equilibrium which greatly simplifies the calculation of NO. Most researchers who rely explicitly on this mechanism use the equilibrium assumption and an associated quasi-steady state approximation (QSSA) for the concentration of N species. The concepts of QSSA and decoupling will be described in detail in a later section. While this mechanism has been proven to provide validated results, some researchers claim that use of the EZM alone can lead to errors in respect to experimental data of up to 50%. 2.3.2 The Prompt NO Mechanism The prompt mechanism otherwise known as the Fenimore mechanism was originally proposed as a possible route for NO formation by C.P. Fenimore [36]. Fenimore noticed that the growth of NO in the post-flame gas is well described by the thermal mechanism but could not describe a faster transient formation of NO in the primary reaction zone. Tests ran by Fenimore suggested since this rapid formation did not occur in hydrogen or carbon monoxide flames which 16 implied that this rapid transient formation of NO must be attributable to the attack of hydrocarbon radicals on nitrogen molecules. The attack of nitrogen by hydrocarbon radicals intimately links the prompt mechanism with fuel combustion chemistry. Ignoring hydrocarbons radicals that initiate the mechanism, the Fenimore mechanism can be written as, where R4 is the rate limiting step in the sequence. Miller and Bowman [35] noted that typical levels of prompt NO can range from a few parts per million by volume to more than 100ppmv. Bowman also suggested that some prompt NO formation occurs due to super-equilibrium concentrations of O and OH, which accelerate the thermal mechanism. There have been numerous studies in order to determine the hydrocarbon radical species responsible for prompt NO which include: Among these hydrocarbon radical species, the principal species driving the prompt mechanism are believed to be The prompt mechanism is primarily attributed to fuel rich combustion due to its dependency on fuel combustion chemistry. Beyond equivalence ratios larger than 1.2, other chain branching routes open up and the chemistry becomes much more complex. Miller and Bowman [35] suggested that beyond several factors combine to cause the NO concentration to decrease relative to HCN including: the conversion of HCN to NO by the prompt mechanism is no longer valid, the recycle of NO to HCN by the prompt mechanism begins to inhibit NO production, and R1 17 actually destroys NO rather than producing it. A reaction path diagram illustrating the major steps in prompt NO formation can be seen below in Figure 2-12, where the bold lines represent the most important reaction paths [35]. Figure 2-1: Prompt NO Formation Pathways 2.3.3 The Mechanism While NO is the predominant oxide of nitrogen produced in the cylinder as defined by ( ), some combustion processes, especially diesel engines, emit significant fractions of their total oxides of nitrogen as . The primary chemical reactions involved in this process during combustion are: 18 Miller and Bowman [35] suggested that in the low-temperature regions of flames, significant concentrations are found that can react with NO formed in the high-temperature regions and transported by diffusion to the low temperature region. The removal reactions are rapid, and in the prescence of high radical concentrations, will be converted back to NO. The principal production reaction is R10 and the rate of this reaction is dependent on the concentration of which in turn is dependent on R13 which is the principal formation reaction in low temperature regions. One must also consider the importance of R10 on the availability of hydrogen, produced and consumed by R14, on the production of . From the chained sequence discussed above it is clear that the production of is heavily dependent on the kinetic parameters of R10-R12 and on the accuracy of species diffusional transport calculations, especially for NO and H [35]. can exist as only a transient species therefore requiring the quenching of formed in the flame. This quenching can be the result of turbulent mixing of hot and cold fluid elements characteristic of diesel combustion. Below is a figure of the comparative concentrations of NO vs. in both a S.I. engine and a diesel engine. 19 Figure 2-2: Relative concentrations of NO and in IC engines [9] 2.3.4 The Mechanism The mechanism becomes of interest at fuel-lean and highly dilute conditions such as in gas-turbine manufacturers and heavily diluted combustion systems. The three steps of the -intermediate mechanism can be written as: Miller and Bowman [35] suggested that studies in laminar premixed flames indicate that is a very short lived species in hot combustion gases and that the principal formation reactions involve NO and various nitrogen-containing radicals. The species although being a short- lived intermediate can still cause an increase in NO production from converting back to NO. Mellor [14] developed a 7 reaction mechanism including the mechanism in his research and suggests that at pressures typical of both the premixed and diffusion portions of the combustion 20 process, the third-body reaction leading to the formation of becomes faster than the leading reaction in the Zeldovich mechanism, thus formation in diesels through the mechanism can thus proceed more efficiently than through the traditional route. Mellor [14] also stressed the importance of the mechanism in high pressures claiming it is pressure sensitive and that decomposition of in the combustion products during the power stroke can also occur by both the reverse EZM and the second order step that produces . Despite this, Heywood [9] included the mechanism in his reaction scheme and compared the results to that of the EZM and concluded that it did little to change the effect of overall production. 2.3.5 Fuel Bound Nitrogen A major source of emissions in the combustion of fossil fuels is the nitrogen chemically bound in the fuel. Typically coal and coal-derived fuels contain 0.5-2.0% nitrogen by weight. In combustion of fuels bound with nitrogen the nitrogen contained in the fuel can be converted to hydrogen cyanide or ammonia. These reactions follow that of the Fenimore mechanism. Miller and Bowman [35] suggests that the extent of conversion of fuel nitrogen to NO is nearly independent of the identity of the model compound, but is strongly dependent on the local combustion environment and on the initial level of nitrogen compound in the fuel-air mixture. Therefore a reaction mechanism for the production of NO from nitrogen bound fuel must contain reactions involving the oxidation of hydrogen cyanide and ammonia. Due to the low nitrogen concentrations in diesel fuel the effect of fuel bound nitrogen in this study will be neglected. 21 2.3.6 Other Mechanisms There are many different possibilities in reaction combinations between intermediates and transient species which will ultimately be converted back to NO. These reaction pathways contain many different combinations of species and routes such as the oxidation of hydrogen cyanide, the oxidation of species, oxidation of NNH, and HNO. A detailed study of these mechanisms can be found in Miller and Bowman [35] a few of which are shown below. 2.4 Modeling Coupled/Decoupled with Combustion Physically speaking formation and combustion are coupled because they occur at the same time. Decoupling the two processes arises due to the complexity of fuel combustion chemistry. Take methane combustion as an example. Methane combustion is by far the most researched area of combustion which leads to much knowledge about the combustion chemical kinetics. Detailed mechanisms for methane combustion usually involve anywhere from 200 to 300 reactions with around 48 species. These fuel combustion kinetics lead to lengthy computational times and are usually limited to simple flame combustion scenarios. For more complex fuels like isooctane this combustion chemistry can reach 800 reactions with over 200 species involved. Thus in much research, especially engine applications, the decoupling assumption helps to simplify the formation problem and shorten computational time 22 dramatically. The assumption of decoupling formation is based on the premise that combustion occurs at high pressure so the flame reaction zone is extremely thin and residence times within this zone are short. This causes the rates of the energy producing reactions in a flame to be sufficiently fast so that burned gases can be assumed to be in thermodynamic equilibrium. Most researchers regard formation as being a slow process compared to the combustion process. With this in mind formation can be treated predominately as a post combustion process by decoupling it from the main combustion process and evaluating the relevant species at equilibrium. Another idea which can be applied to the understanding of decoupling, in terms of post combustion formation, is the temperature stratification developed in the cylinder during combustion. During diesel combustion there exist great temperature gradients. These temperature gradients are caused by the different time histories of gas elements which burn early or late in the combustion process. These time histories can change specifically for engines, due to isentropic compression/expansion of the burned gas elements. Taking note that formation is extremely dependent on temperature one can assume that in most cases any post-flame produced dominates any flame-front produced . Although in some cases the equilibrium assumption is regarded invalid due to rapid freezing of gaseous phase species. This rapid freezing can cause the concentrations of these gaseous species to be several orders of magnitude higher than equilibrium. This is referred to as super-equilibrium. Turns [27] suggested that within flame zones proper and in short-time scale, post-flame processes, the equilibrium assumption is not valid. Super-equilibrium concentrations of O atoms, up to several orders of magnitude greater than equilibrium, greatly increase NO formation rates. Despite that super-equilibrium does occur in some situations, most researchers 23 assume that the processes are decoupled and approximate the relevant species at equilibrium concentrations. 2.5 Detailed Mechanisms and Reduced Mechanisms As research continues in the area of formation, scientists and engineers will continue to add and suggest new chemical pathways for formation. Detailed mechanisms like Miller and Bowman [35] and Kilipinen 97 [38] contain over 200 elementary reactions and include all, and variations of, the mechanisms described in the previous sections. Most of these detailed mechanisms are considered coupled with the fuel combustion kinetics and require detailed fuel burning chemistry in order to determine concentrations of certain species. While the accuracy of these detailed mechanisms are usually better than simple mechanisms, the computational time required for the solution of simulations make them undesirable design tools. In order to have accurate models that have reasonable computational times, researchers must develop reduced mechanisms. Generally the development of reduced mechanisms can be divided into two steps: 1. Establish a skeletal mechanism 2. Develop a reduced mechanism by applying quasi-steady state, partial equilibrium assumptions. To determine a skeletal scheme from a detailed reaction mechanism researchers use sensitivity analysis? to determine the most important formation and destruction steps for NO. The resulting skeletal mechanism is a tradeoff between model complexity, accuracy, and range of applicability. The next step in the development of a reduced mechanism is to apply the quasi- steady state and partial equilibrium assumptions. Both QSSA and PE techniques help the 24 conversion of stiff differential equations to non-stiff differential equations and also reduces the number of unknowns to be solved. These techniques allow the systematic reduction of the skeletal kinetic scheme to an even simpler form known as a semi-global mechanism. Although the above reduction procedure is systematic, it can be algebraically tedious and unwieldy, especially because the algebraic relations for concentrations are frequently nonlinear with imaginary roots and cannot be explicitly expressed. Furthermore, the range of applicability of the reduced mechanism also becomes progressively restricted as more species are assumed to be in steady state and more reactions are assumed to be in partial equilibrium [37]. Other techniques for reducing kinetic mechanisms include the Computational Singular Perturbation (CSP) and Intrinsic Low-Dimensional Manifold techniques (ILDM). The CSP technique is based upon evaluating the time scales of reactions and ordering them from fastest to slowest. Thus as the reaction proceeds the fastest reaction becomes exhausted and its reaction rate approaches zero and it is discarded from the reaction scheme. This procedure is basically another way to apply the partial equilibrium technique on a priori basis. The Intrinsic Low- Dimensional Manifold (ILDM) works on the premise of applying the QSSA and PE techniques automatically over all the possible thermo-chemical states of the system. Taking into account that a reaction mechanism is characterized by fast and slow reactions, it is assumed that the fast reactions can be decoupled from the slow ones, and can be ignored while those with slow time scales are tracked using progress variables. Either which technique is used, the QSSA and PE techniques seem to form the basis of operation. Some examples of reduced mechanisms using the reduction techniques discussed above include, the Goldsworthy [40] mechanism which contains 13 elementary reactions based on the Kilipinen 97 [38] scheme and the Glarborg [41] scheme which is based on the Miller and Bowman [35] detailed mechanism. 25 Quasi-Steady State Approximation The QSSA technique was developed at the start of the last century when the solution of nonlinear reaction systems could not be obtained due to the non-availability of computers. Some researchers suggests that the QSSA technique is an obsolete technique, but as reactions systems become more and more complicated the computational demands on the computer increase greatly. The QSSA technique helps the conversion of stiff differential equations to non-stiff differential equations and also reduces the number of unknowns to be solved. Stiffness arises in the ODE?s due to the varying time constants involved in a reaction system. During a reaction scheme leading from reactants to products, reaction intermediates are produced. These intermediates act as chain carriers between the individual reactions and thus proceed with rapid rates, although the concentrations and the net rate of change of these intermediates are quite low. In other words the rates of formation and consumption of these intermediates can be approximated as equal. 2.6 Modeling Formation in Diesel Engines Most current internal combustion engine models are developed in relation to spark ignition engines. While these models have been shown to give accurate results with IC spark ignition engines, researchers suggests they fail to describe the underlying kinetics of formation in compression ignition engines. This non-linearity between the two engines can be traced to the differences of the combustion processes between the two engines. SI engines are characterized by homogenous combustion processes, whereas diesel combustion is considered a highly heterogeneous process. Due to these differences many researchers have performed 26 detailed chemical kinetic investigations on the mechanisms leading to the formation of in- cylinder in CI engines. Zabetta [39] performed this kinetic investigation for low and medium speed diesel engines using the detailed kinetic model of Kilpinen 97 [38]. Zabetta showed that the proportion of was negligible compared to relative concentrations of and reported that the formation and destruction of NO was mainly dependent on the thermal and mechanism comprising 10 reactions. These 10 reactions included the thermal mechanism, the -intermediate mechanism, and a -extension subset, which included oxidation reactions with HNO and NH. A figure showing the relative contributions of these mechanisms is shown below. Figure 2-14: Relative contributions of mechanisms [39] 27 Zabetta [39] concluded that the improved sub-model produced more accurate predictions than did the detailed Kilipinen 97 mechanism. Zabetta used the PE assumption for the concentrations of radicals O, OH, and H, and used the QSSA approach for the concentration of . Zabetta?s elementary reaction equations can be seen below: Goldsworthy [40] performed several different schemes involving the reactions stated above proposed by Zabetta [39]. Goldsworthy noted that for slow speed diesel engines the pathway added an extra 15% to the predicted by the thermal mechanism at 100% load, reducing to 12% additional at 25% load. Goldsworthy concluded that the relative contribution to NO production from the intermediate reactions R21 and R22 increased as temperature decreased and oxygen content increased, while the contributions from the thermal mechanism and the extension reduced. Prompt NO was considered insignificant in terms of a slow speed diesel engine and fuel bound NO was neglected unless operating on a residual fuel oil. Mellor [14] also suggested the importance of the mechanism in his research and 28 development of a 7 reaction skeletal mechanism. Mellor suggests that at pressures typical of both the premixed and diffusion portions of the combustion process, the third-body reaction leading to the formation of becomes faster than the leading reaction in the Zeldovich mechanism, thus formation in diesels through the mechanism can proceed more efficiently than through the traditional route. In other words the thermal mechanism alone is inadequate to describe NO formation in CI engines. Despite this many researchers rely on the thermal mechanism alone for calculations and correct the predicted values with extensive tuning to engine parameters. Mellor also stated that the decomposition of NO in the combustion products during the power stroke can also occur by both the reverse thermal mechanism and the second order step that produces . The 7 reaction scheme developed by Mellor [14] included the following elementary reactions: Stiesch [24] used a similar model to Mellor [14] that included 10 reactions with variations among the and mechanism. Stiesch claimed that the mechanism provided accurate results and the crucial parameter to accurate predictions was the boundary conditions on the combustion zones. Miller [16] developed a super-extended Zeldovich (SEZM) mechanism for prediction in SI engines. The SEZM includes a set of 67 elementary reactions including 23 chemical species. It was shown that this mechanism alone was capable of predicting concentrations to within 10% of experimental data. The SEZM contains the thermal mechanism, 29 mechanism, mechanism, and oxidation reactions similar to the Zabetta -extension involving HNO, , , and other cyano derivatives. Overall in terms of mechanisms, the thermal mechanism seems to dominate for CI engine formation. Due to the fact that the prompt mechanism is generally characterized by rich combustion and ultra-lean conditions, it is typically neglected for CI engine analysis. The mechanism as shown by the researchers above can be an important source of NO and must be included in the skeletal mechanism. Most of these mechanisms developed by researchers use the techniques discussed in the previous section for mechanism reduction. Rublewski [23] and Ericson [7] suggests while the full kinetic approach improves accuracy in predictions compared to the equilibrium approach, it significantly increases the computational time. Taking into account that is a strong function of local temperature, pressure, the amount of oxygen available, and various temperatures and composition stratifications occurring in the cylinder during the combustion process, another important element to be considered in modeling is the combustion model. Accurately predicting concentrations requires the accurate modeling of the thermodynamic states of burned gases. Many researchers who use single-zone fully mixed models usually incorporate separate models to account for these stratifications for the treatment of emissions. Heywood [23] along with many other researchers use approaches such as layered adiabatic core and zonal models in an attempt to emulate these stratifications and to improve the slope of their models to accurate levels. With new methods of reduction it is imperative that the combustion model includes the effects of control strategies such as EGR. Due to the fact that EGR lowers combustion temperatures, increases intake air temperatures, and dilutes the oxygen concentration it can have a dramatic effect on engine out and must be accounted for. 30 Chapter 3: Model Development 3.1 Introduction This chapter will describe the development of all models and sub-models used in this study. The principles by which the base model is developed will be presented and explained in detail. The perturbations that will be applied and tested for validity parallel with the base model will also be identified. Only a brief description of the Wingo [42] Combustion Model will be presented. The solver used to solve the ODEs will also be explained along with the overall program structure. 3.2 The Combustion Model The combustion model used in this study is a single-zone quasi static model developed by Wingo [42]. The details of this model are given in Wingo [42], only a brief outline of the author?s work is given in the subsequent sections. The following assumptions were made by Wingo [42] in the development of this model: 1. Fuel injected into the cylinder is evaporated immediately so that fuel spray behavior will have no effect on the model. In reality, liquid fuel evaporates and mixes with air to produce a non-uniform fuel/air distribution. 2. The working fluid is a homogeneous ideal gas mixture 3. Crevice flow, including the charge trapped in the volume between piston, rings, and cylinder, is considered to be negligible. 4. Combustion is modeled as a uniformly distributed heat release, which is assumed to be proportional to the AFBR based on the above assumptions. 31 3.2.1 Governing Equations The basis of the single-zone model is a quasi-static application of the first law of thermodynamics. ? (3.1) Where: is the heat-transfer rate through the system boundary is the work done by the system through the boundary displacement, ? is the energy due to mass flow across the system boundary. Equation (3.1) is then modified to only the sensible energy to solve the Apparent Heat Release Rate (AHHR). (3.2) Substituting the differential form of the ideal gas law gives the net heat release rate as a function of differential volume and pressure. ( ) ( ) (3.3) Equations (3.2) and (3.3) are then combined and converted to a crank angle basis. ( ) ( ) (3.4) (3.5) 32 3.2.2 Heat Release Model The heat release model was based on that of the Watson combustion model where the premixed and diffusion-controlled portions of the burned mass rate are controlled by a proportionality factor . (3.6) (3.7) The premixed and diffusion controlled portions are modeled by and respectively, (3.8) (3.9) Where: and are all empirically determined shape parameters and is given by: (3.10) 3.2.3 Heat Transfer Model The heat loss rate term is expressed in the model as, ( ) (3.11) where heat transfer coefficient was given by a correlation developed by Woschni, ( ) (3.12) Where: , , are constants, is the characteristic speed, and L is the characteristic length 33 3.2.4 SOC Model The Start of Injection (SOI) is signaled by an ECU command, which can be measured directly as a current signal and pulse width (PW). For a given fuel injection command, there will be a delay from the start of the command to the needle lift and then a second delay from the needle lift until actual fuel injection, and a third delay from the SOI until the Start of Combustion (SOC). The total delay from the ECU command to SOC is the sum of these delays. Since a needle lift sensor has not been installed on this test stand, the delay from the ECU command to the injector needle lift is estimated from an empirically based equation developed by Gong [45]. The first delay time is related to the magnetic and hydraulic force exerted on the injector plunger. Therefore it is not a function of injection pressure or pulse width. The second delay time is the time between the needle lift and the point of actual fuel injecting into the chamber. This time includes the transport delay for fuel flow and can be assumed to be 0.1 milliseconds. (3.13) [ ] (3.14) Where: is the fuel rail pressure [bar] P is the cylinder pressure [bar] To estimate the third delay from injection to the actual start of combustion an empirical formula found in Heywood [9] which showed good agreement with experimental data in this study is used. [ ] ( ? ) [ ( )( ) ] (3.15) Where p is the cylinder pressure in [bar], T is charge temperature in [Kelvins], ? is the mean piston speed, and is given by: 34 (3.16) 3.2.5 Supporting Models The cylinder volume at any crank position is calculated by, [ ] (3.17) Where: is the clearance volume, is the compression ratio R=2L/S is the ratio of connecting rod length (L) to stroke (S) The combustion chamber surface area at any crank position is given by, (3.18) Where: is the cylinder head surface area is the bore The fuel is estimated using the following equation. ? ? (3.19) Where: is the discharge coefficient is the nozzle minimum area is the fuel density is the pressure drop across the nozzle is the duration of the fuel injection is the engine speed The trapped air mass during the intake stroke is estimated as, (3.20) 35 Where: N is the engine speed in revolutions per minute is the average MAF rate (lb/min) from experimental data The equilibrium model, which is explained in the next section, is used to predict the unburned and burned gas mixture states at different stages of the engine cycle. This model predicts state properties such as specific enthalpy, internal energy, entropy, gas constant, and specific heats. The program uses fitted NASA polynomial curves for the various thermodynamic properties. 3.3 Definition of Base Model In the following sections the treatment and approach of the base model will be identified. The base model consists of the equilibrium calculation model, the adiabatic layered core model, the EGR model, and the chemical kinetics principle. These preceding models will not be modified throughout the model perturbation and sensitivity analysis. With the base model identified, changes in the reaction mechanism, rate constants, treatment of particular radical species, and appropriate temperature approximations will be identified and tested for validity. 3.3.1 Equilibrium Composition Model It has been validated that the equilibrium condition is a good approximation for estimates in internal combustion engines. This equilibrium condition means that the chemical reactions involved produce and remove species at equal rates thus, no net change in species composition results. The following equilibrium method is based on the minimization of Gibbs free energy, G. The Gibbs free energy is defined in terms of thermodynamic properties as, (3.21) The second law can then be expressed as 36 (3.22) which states that the Gibbs function always decreases for a spontaneous, isothermal, isobaric change of a fixed-mass system in the absence of all work effects except boundary work. This relationship allows us to calculate the equilibrium composition of a mixture at a given temperature and pressure. In comparison to Entropy which attains a maximum at equilibrium for fixed-energy and fixed-volume cases, the Gibbs function attains a minimum. Thus, at equilibrium, (3.23) For a mixture of ideal gases, the Gibbs function for the ith species is given by ? ? (3.24) where ? is the Gibbs function of the pure species at the standard-state pressure , which is taken to be 1 atmosphere, and is the partial pressure. When dealing with reacting systems another formulation called the Gibbs function of formation is frequently used. This expression can be written as ? ? ? ? (3.25) For a mixture of ideal gases the Gibbs function change can be expressed as ? ? ? [ ? ] (3.26) where is the number of moles and nSp is the number of species. For fixed temperature and pressure the equilibrium condition becomes (3.27) 37 or ? [ ? ] ? [ ? ] (3.28) The second term of Equation (3.28) is zero since the sum of all partial pressures is zero, i.e. ? , Equation (3.28) reduces to ? [ ? ] (3.29) The Gibbs function of the mixture may also be represented as, (3.30) where is the mixture enthalpy , T is the temperature, and is the mixture entropy. The two quantities of enthalpy and entropy can be represented by, ? ? (3.31) where is the number of moles of species i and ? is the individual species enthalpies which are obtained from, ? ? ? ? . (3.32) and ? ? (3.33) where is the number of moles of species i and ? is the individual species entropies which are obtained from, ? ? ( ) ? ? . (3.34) 38 Substituting Equation (3.32) and Equation (3.34) into Equation (3.31) and Equation (3.33) respectively yields: ?[ ? (? ? ) ] (3.35) ? ? ? ? . (3.36) Substituting Equation (3.35) and Equation (3.36) into Equation (3.30) yields: ? { [ ? ?] ( ) (? ? ) ? ? } (3.37) Values of the Gibbs function are tabulated as a function of temperature at a fixed reference pressure. Therefore for a specified temperature Equation (3.37) becomes: ? [ ? (3.38) Taking the value of to be 1atm Equation (3.38) can be represented as: ? [ ? ] (3.39) Recalling that for an ideal gas , where is the mole fraction defined as the Gibbs Free Energy Function can be defined for each component of the mixture as: ? (3.40) Dividing by yields: (3.41) 39 Equation (3.41) is similar to the form of the NASA thermodynamic database data stored in the form of the polynomial curve fits shown below: ? (3.42) ? (3.43) ? (3.44) Specific NASA data used in this study for was determined from the NASA tabulated values of , where is the equilibrium constant and was determined from the relation: (3.45) This data was fitted using the method of least squares with the 9 coefficient polynomial form shown below: [ ] (3.46) Equation (3.41) will now become the basis of the minimization of Gibbs Free Energy. Considering the reaction equation involving all 19 species that will be considered in this study, (3.47) 40 where is the number of moles of the respective species, is the equivalence ratio, [ ] are the number of atoms of each specific element in the fuel, and is represented by: (3.48) An element balance of this global reaction can be written as: ? (3.49) where is the total molar concentrations of atoms of element j ( j=1to nEL) which is set by the reactant concentrations, nSp is the number of species, nEL is the number of elements, and is the number of atoms of element j in species i. This element balance can be written in terms of the species [C H O N] as: (3.50) (3.51) (3.52) (3.53) 41 In order to minimize the Gibbs function for equilibrium, Equation (3.41) must be minimized subject to the constraints of Equation (3.49). To perform this task the method of Lagrangian multipliers is applied. Let be defined from ( j=1 to nEL) and define F as: ? ( ? ) ? ?( ) (3.54) The solution required for equilibrium is when F is a minimum, i.e. for all i= 1? nSp. ? ? ? ? (3.55) Equation (3.55) reduces to: ? ? (3.56) Equation (3.56) gives nSp equations and Equation (3.49) gives nEL equations. There are nSp unknown values of and nEL unknown Lagrangian multipliers . Since the number of equations and unknowns are equal we can solve the system described above. The systems of equations are regrettably non-linear and must be solved iteratively. There are many different approaches to solve these equations but the method used in this study is described below. Let ? be an estimate of and expand as a Taylor series around this estimate: ( ) ? ?( ( )) ? ( ? ) (3.57) Differentiating Equation (3.56) yields: ? (3.58) 42 ? (3.59) ( ) ? ? ? (3.60) ( ) ? ? ? ? (3.61) Substituting into Equation (3.57) yields: ? ? ? ( ?) ? ? ? ( ?) (3.62) Reducing Equation (3.62) yields: ? ? ? ? ? ? (3.63) where is defined as: ? (3.64) Setting Equation (3.63) equal to zero for the minimum and rearranging yields: ? (? ? ? ? ) (3.65) Where a relationship of can be expressed as: ? ? ? ? (3.66) Substituting Equation (3.65) into Equation (3.49): 43 ? ? ? ? ? ? ? ? ? (3.67) Equations (3.67) and Equation (3.66) yield nEL+1 equations with unknowns and ? ? ? . These equations are used to solve for the unknowns and then substituted into Equation (3.65) to get a better estimate of . Equations (3.67) and (3.66) can be written in the form of a matrix operation: ? (3.68) where: [ ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ] ? [ ? ? ? ] (3.70) 44 [ ? ? ? ? ? ] (3.71) The matrix can be represented as: [ ] (3.72) 45 The equilibrium program follows the following steps in order to determine the concentrations. 1. Provide an initial estimate of the number of moles of each component. 2. For the specified temperature, pressure, equivalence ratio, and species calculate the Gibbs Function for each species using Equation (3.46). 3. Generate the matrices A, B, and ? 4. Solve Equations (3.66) and (3.67) for the unknowns using the matrix form Equation (3.68). 5. Calculate the new estimate for the concentrations from Equation (3.65) 6. Check for Convergence 7. If not converged reset the initial estimates with the corrected concentrations. 8. If converged print result There are several ad hoc adjustments in the equilibrium program to enhance convergence. The concentrations are bound between the highest and lowest possible values depending on element availability. The corrections made by Equation (3.65) are also bound to prevent the corrections from making large steps. The initial conditions are initially set as (0.1/nSp). Once the equilibrium concentrations of the species set is known, the next step is to evaluate the properties of the mixture. For determining the mixture enthalpy, entropy and internal energy the process is fairly straightforward and is obtained from the following relations shown below. For the mixture enthalpy, entropy, and internal energy: ? ? (3.73) ? ? (3.74) (3.75) 46 where is the universal gas constant, is the molecular weight of the mixture, T is the temperature, is the mole fraction of species i, and ? is the species enthalpy determined from Equations (3.42-3.44). As for the specific heat there are two approaches which are quite different due to the fact that the equilibrium specific heat differs from that of the frozen specific heat. Equilibrium specific heat refers to the specific heat of a composition that reaches its equilibrium value instantaneously, while the frozen specific heat is associated with infinitely slow reaction times in which the composition remains frozen. The following method follows that of NASA publication Gordon [43]. The equilibrium specific heat can be expressed as the sum of a frozen contribution and a reaction contribution illustrated below: where the frozen specific heat can be represented as: With the definition of specific heat at constant pressure defined as: we can differentiate Equation(3.73) with respect to temperature to obtain an expression for : (3.76) ? ? (3.77) ( ) (3.78) ? ? ( ) (3.79) 47 After making the appropriate derivatives of with respect to temperature the resulting reduced number of temperature derivative equations (excluding solid phase species possibilities) is, where k=1..nEl. These two index equations give five equations for the five unknown temperature derivatives, The approach to solving these equations is the same that was used in the equilibrium composition calculation. All notation such as , and hold the same relationship as before. Equations (3.80) and (3.81) can be written in matrix form as: ? (3.83) where: ?? ( ) ? ( ) ? ? (3.80) ?? ( ) ? ? (3.81) ( ) ( ) ( ) ( ) ( ) (3.82) 48 [ ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ] ? [ ] (3.85) [ ? ? ? ? ? ? ? ? ? ? ] (3.86) 49 Once the five derivatives are found the equilibrium specific heat can be calculated from the Equation (3.87) shown below. 3.3.2 Chemical Kinetics Calculation Principle Once the state of the species has been specified, i.e., QSSA, PE, FE, or kinetically controlled, we can represent a system of elementary equations or a mechanism in compact notation as: ? ? (3.88) where and are the stoichiometric coefficients on the reactants and products side of the equation, respectively, for the jth species in the ith reaction. The single prime superscript ( represents the forward stoichiometric coefficients while the double prime superscript ( represents the reverse stoichiometric coefficients. The net production rate ? for each species in the multistep mechanism can be represented by the following relations: ? ? (3.89) where KnSp is the kinetically controlled species and: (3.90) ?(? ? )( ) (? ? )( ) ? ? ? ( ? ) (3.87) 50 The rate of progress variable for the ith reaction is given by the difference between the forward and reverse rates as: ?[ ] ?[ ] (3.91) where [ ] is the molar concentration of the jth species present in the ith reaction. For third body reactions the rate of progress variable formulation is of the form: (?( )[ ] )( ?[ ] ?[ ] ) (3.92) where are the enhanced third body coefficients. The forward reaction rate constants used in this study are given in the Arrhenius form of: ? (3.93) where the values of the pre-exponential factor A, the temperature exponent b, and the activation energy are given in Appendix A. The reverse reaction rates are related to the forward reaction rates through the use of the equilibrium constant illustrated by: (3.94) The equilibrium constant can be represented in two forms. The first form is based upon partial pressures. Consider the arbitrary bimolecular reaction: ? ? (3.95) According to the chemical kinetic principles above we can write an expression for the time rate of change of species A as: [ ] [ ][ ] [ ][ ] (3.96) For the equilibrium condition, A+B=C+D, the time rate of change of [A] is zero so Equation (3.96) can be expressed as: 51 [ ][ ] [ ][ ] (3.97) which can be rearranged as: [ ][ ] [ ][ ] (3.98) Using the relationship for partial pressures, [ ] ? (3.99) Equation (3.97) can be written for the equilibrium constant in terms of partial pressures as: ? ? ? ? (3.100) where the superscripts are the corresponding stoichiometric coefficients and . The equilibrium constant can also be represented by use of the Gibbs function: ? (3.101) where [ ] (3.102) Since tabulated values for already exist for the equilibrium program, Equation (3.102) and Equation (3.94) can be used to determine the reverse rate coefficients. Now that the chemical kinetic principles have been clarified along with the determination of equilibrium values and rate coefficients, we can explore the equations generated by Equation (3.88). The Quasi-Steady State Assumption (QSSA) technique was developed at the start of the last century when the solution of nonlinear reaction systems could not be obtained due to the non-availability of computers. Some researchers suggests that the QSSA technique is an obsolete technique but as reactions systems become more and more complex the computational demands on the computer increase greatly. The QSSA technique helps the conversion of stiff differential equations to non-stiff differential equations and also reduces the number of 52 unknowns to be solved. Stiffness arises in the ODE?s due to the varying time constants involved in a reaction system. During a reaction scheme leading from reactants to products, reaction intermediates are produced. These intermediates act as chain carriers between the individual reactions and thus proceed with rapid rates, although the concentrations and the net rate of change of these intermediates are quite low. In other words the rates of formation and consumption of these intermediates can be approximated as equal. Consider the arbitrary reaction below: ? (3.103) ? Consider the first reaction in this pair to be slow and thus rate limiting, while the second reaction is extremely fast. According to the law of mass action stated previously we can represent the net production of A atoms as: [ ] [ ][ ] [ ][ ] (3.104) After a fast transient allowing the buildup of A atoms to some low concentration, the two terms on the right hand side of Equation (3.104) become equal and [ ] approaches zero. With [ ] we can determine the steady state concentration of A atoms: [ ][ ] [ ] [ ] (3.105) where, [ ] [ ][ ] [ ] (3.106) 53 Although the QSSA assumption suggests that [ ] does not change with time, [ ] may change as it rapidly readjusts according to Equation (3.104). To determine the time rate of change of Equation (3.106) requires differentiating Equation (3.106) rather than applying Equation (3.106). In the partial equilibrium (PE) assumption we assume that the combustion mechanism involves both fast and slow reactions, where the fast reactions are rapid in both the forward and reverse routes. These fast reactions are usually chain-propagating or chain-branching steps, while the slow reactions are termolecular recombination reactions [27]. Treating the fast reactions as if they were an equilibrated species eliminates the need to write rate expressions for these species. The PE assumption technique can be illustrated by the following arbitrary mechanism taken from Turns [27]: (3.107) where the reactive intermediates are A, B, and AB, while the stable specie are , and . Equilibrating the forward and reverse reaction rates gives: [ ][ ] [ ][ ] [ ][ ] [ ][ ] (3.108) [ ][ ] [ ][ ] These three equations can be solved simultaneously for the radical species A, B, and AB, in terms of the stable species , and . Turns [27] suggests that the net result of invoking either the partial equilibrium assumption or the QSSA approximation is the same: a radical concentration is determined by an algebraic equation rather than the integration of an ordinary differential equation, but it is important to note that physically the two approximations are quite 54 different in that the PE assumption forces a set of reactions to be essentially equilibrated while the QSSA assumption forces the individual net production rate of one or more species to be zero. The determination of which species is considered for PE treatment is usually based upon experience. To illustrate the above chemical kinetic principles a sample calculation of the rates of reactions of a reaction mechanism will be presented. Consider the reaction mechanism: (3.109) Applying the equilibrium condition and assuming the concentrations of , O, H, OH, and are at their equilibrium values, the one-way rates of reactions can be represented as : [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] (3.110) [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] [ ] 55 Treating , , and as the kinetically controlled species we can develop the reaction rates which can be represented as: [ ] [ ] (3.111) [ ] Where , , and can be represented as: [ ][ ] [ ][ ] (3.112) [ ][ ] Assuming that and qualify as a QSSA species the right side of their rate equations can be approximated as zero and and can be solved in terms of the R?s. These values of and can then be substituted into ?s rate equation which then can be written as: [ ] [ [ ] [ ] ] (3.113) 3.3.3 Adiabatic Layered Core Model The combustion model used in this study outputs a temperature based on the ideal gas law. While this temperature is a valid approximation for the average cylinder temperature it is 56 too low to represent local burned gas temperatures and therefore not a valid temperature to be used in modeling. From the ideal gas law the mean cylinder temperature is given by: (3.114) where is the gas constant, is the cylinder pressure, is the cylinder volume, and is the mixture mass. As discussed in the literature review a valid approximation of representing local and burned gas temperature stratifications is the multi-layered concept. The multi-layered concept is especially suited for diesel combustion due to the great number of temperature gradients that occur in the cylinder. Studies have shown that some of these early burning regions can reach higher temperatures than those of the final burning stages. In the multi-layered concept a new temperature profile is created every time a new burned gas element (BGE) is created from the start of combustion at each crank angle until the burned fraction equals one. This layered concept works on the assumptions that there are no mixing between each of the BGE layers, where as in the actual process the burned zone behavior lies somewhere between mixing and non-mixing. Each BGE has its own temperature history during the compression/expansion processes which is assumed to take place as an isentropic compression/expansion processes. This isentropic compression and expansion can be represented as: ( ) (3.115) Since the enthalpy, pressure, and charge composition is known about the unburned zone, the flame temperature can be calculated by means of an isenthalpic process. Without heat transfer and work this constant enthalpy process can be described as . From this assumption an approximation of the flame temperature is the adiabatic flame temperature calculated from the 57 specified unburned air temperature. The equivalence ratio used in the model was based upon studies done by Heywood [9] and other researchers and was set at 0.85 for all points. The adiabatic temperature is calculated by an iterative process of selecting flame temperatures until the enthalpy of the combustion products equals the enthalpy of the initial gas mixture. The mixture enthalpies are determined by calculating the sum of the enthalpies for each species, i.e., ? . The enthalpies for each species are determined by the polynomial function Equations (3.42-3.44). Since the adiabatic flame temperature is extremely sensitive to residual gas mixture and EGR dilution, the details of how this is handled will be discussed in the next section. Figure 3-1 shown below represents the temperature profiles from the layered adiabatic routine where the black lines indicate each BGE profile and the red line indicates the adiabatic temperature profile calculated from the isenthalpic process. Figure 3-2 shows the corresponding NO layer concentration profiles. Figure 3-1: Layered Temperature Profiles 58 Figure 3-2: Simulated Layer NO Concentrations Figure 3-2 indicates the large amount of NO stratification that can occur in the cylinder due to the various temperature gradients. Figure 3-2 also shows that these layer NO concentrations are sensitive to their equilibrium values indicated by the negative formation rates. To determine the overall NO concentration within the adiabatic layered core, the mass burned during each crank angle is used to mass weight the results according Equation (3.116), [ ] ? [ ] ? (3.116) where and [ ] represents the mass burned and final NO concentration of each layer. Since the burned mass fraction of the first and last layered profiles is small the contribution from these layers to the overall adiabatic core concentration is limited. Once the layers have been 59 mass averaged, the net production of NO during one combustion cycle is represented in Figure 3- 3. Figure 3-3: Simulated Total NO Concentration Profile 3.3.4 EGR Model Exhaust gas recirculation is one of the most common techniques used by engine manufactures to control the relative proportion of engine out emissions. EGR has advantages at low-load conditions where equivalence ratios are high, but re-circulating exhaust gas with high soot concentrations at low or high load conditions can have adverse effects on engine durability and lube oil quality. EGR directly affects two primary variables in formation which include temperature and oxygen availability. EGR both increases the heat capacity of the burned gases, and depletes the oxygen concentration thus lowering combustion temperatures and yields. It is because of the great impact of EGR on yields that EGR 60 must be accounted for in this study. The first step in creating the model is to determine the thermodynamic properties of the fuel-air-residual gas mixture (FARG). The residual fraction can be represented as: ? (3.117) where is the residual mass and is the inducted mass per cycle. A typical residual fraction for diesel engines is usually just a few percent. Once EGR is introduced into the cylinder the residual gas fraction can be represented as: ( ) (3.118) At low temperatures (under 1700K) the overall combustion reaction can be written as: (3.119) Validated approximations for lean and rich combustion are: (3.120) For lean or stoichiometric cases, atom balance equations are sufficient to determine the product composition. For the rich case an equilibrium constant for the reaction below is introduced. (3.121) where the equilibrium constant is: (3.122) In the rich case the variable is given by: 61 ? (3.123) where: (3.124) [ ] (3.125) (3.126) Solutions for the coefficients for both lean and rich cases are given below in Table 3-1. i SPECIES 1 2 3 0.79+ 4 0.21 0 5 0 6 0 Table 3-1: Coefficients for lean and rich combustion The analysis above allows the calculation of the coefficients of the combustion products at low temperatures. The next step is to introduce the burned gas fraction into the combustion equation. Let Equation (3.119) be rewritten in the form: (3.127) where is the reactant coefficient, and is the product coefficient. The mass fractions of the residual mixture can be given in terms of the residual-gas fraction by: (3.128) 62 where where is the molecular weight of the species, is the mole fraction, and is the total molecular weight. Using the same relation as Equation (3.128) the mole fractions can be expressed by: (3.129) where the residual mole fraction is expressed as: [ ( ) ] (3.130) So with the initial concentrations of the products calculated from Table 3-1, the mole fractions of the residual composition can be obtained. The fuel-air residual gas composition procedure above follows the pseudo algorithm shown below. 1. Determine , , , , , from TABLE 2. Determine mole fractions of products: ? and find ? 3. Determine mole fractions and molecular weight of fuel-air: 4. Determine the mole fractions of fuel-air residual gas mixture 5. Calculate the thermodynamic properties of the mixture The thermodynamic properties of the FARG mixture completely define the beginning thermodynamic state of the reactants. An iterative procedure is performed next in order to determine the adiabatic flame temperature. Recall in this study we are evaluating the combustion process at equilibrium therefore in the calculation of the equilibrium product states we do not consider the FARG mixture. The FARG mixture is used in the adiabatic flame temperature calculation only. The enthalpy of the FARG mixture defines the enthalpy to be matched by the adiabatic temperature model i.e. . The effect of EGR/residual gas mixture on flame temperature can be seen Figure 3-4. 63 Figure 3-4: Simulated Effects of Charge Dilution on Flame Temperature Incorporating the influence of residual mixture into the equilibrium composition calculation requires the implementation of Equation (3.129) into the equilibrium model described in section 3.3.2. Recall the original combustion equation, (3.131) Implementation of Equation (3.129) into Equation (3.47) yields for the reactants, 64 (3.132) This implementation effectively changes the element balance previously defined by Equations (3.50-3.53) to, (3.133) (3.134) (3.136) From this point on the equilibrium procedure proceeds as outlined in section 3.3.2 with the exception of the changes in the element balance equations (3.50-3.53). (3.135) 65 3.4 Model Perturbations In the following sections, the parts of the model that will be subjected to variations will be identified. The purpose of this investigation is determine the optimum operational characteristics that will yield the most accurate results given the combustion model and engine type. These variations include reaction mechanisms, individual species treatment, rate constants, and temperature approximations. 3.4.1 Reaction Mechanisms A total of five different schemes will be presented. These mechanisms include the Extended Zeldovich Mechanism, Mechanism, mechanism, NH/HNO mechanism, and variations/extensions of these mechanisms. Prompt NO will be neglected along with fuel-bound NO. The different schemes will contain a total of 19 different chemical species. The chemical species set includes: , , , , , , , , , , , , , , , , , , and . Among this set, the kinetically controlled species are , , , , , and . The SEZM found in Miller [16] will be chosen as the base mechanism in this study. This reaction set includes a total of 67 reactions including all of the major mechanisms recommended by other researchers as an important source of NO. Since there is much research surrounding the thermal or EZM mechanism, the first scheme to be analyzed, Scheme 1, will be the EZM. This scheme will include reactions found in Appendix B. Scheme 2 will include the effects of adding the mechanism. Scheme 2 will include reactions found in Appendix B. Since it has been suggested that concentrations can be significant in diesel engines, Scheme 3 will analyze the effect of adding the mechanism. Scheme 3 will include reactions found in Appendix B. Scheme 4 will include all the 66 mechanisms listed above with the addition of oxidation reactions involving NH and HNO. Scheme 4 will include reactions found in Appendix B. The last scheme, Scheme 5, includes the preceding reactions with more variations and extensions of the NH/HNO mechanism. Scheme 5 will include reactions found in Appendix B. For the analysis of these different mechanisms the temperature profile supplied to the adiabatic temperature calculation will be supplied at the unburned air temperature. The rate constants for these reaction schemes will be those supplied in the Miller [16] mechanism which can be found in Appendix B. 3.4.2 Individual Species Treatment The perturbation that will be tested regarding individual species treatment is the QSSA assumption. Based upon numerous studies performed by other researchers, the species in this study that have been suggested as valid QSSA species are , , , ,and . For all the schemes there are a total of six unknowns. These unknowns include the concentrations of , , , , , and . These can be put in vector form to give: ? [ ] (3.142) where represent the unknown concentrations of , , , , , and . Representing these unknowns with the independent variable time gives: ? ? ? (3.143) The variable ? represents the vector of pre-exponential factors and the known equilibrium concentrations form the coefficients of these equations. In order to represent those species which are assumed to be QSSA species Equation (3.69) becomes: ?? ? ? (3.144) 67 ? ( ? ? ) ? (3.145) where ?? and ? represent the concentrations of the non-QSSA and QSSA species respectively. The basis behind the QSSA assumption for the species , , , , and comes from rough estimates for the time rate of change of , , , , and being several orders of magnitude shorter than that of NO. 3.4.3 Reaction Rate Constants The forward and reverse rate constants for IC modeling have been measured in numerous experimental studies. While the majority of these rate constants are similar, some vary considerably. For the purpose of simplicity and due to the facts that there are not enough variations in rate constants to cover all of the elementary reaction equations from schemes 2-5, and assuming that the thermal mechanism will be the governing reaction set, we will limit the analysis of the rate constants to only the thermal mechanism i.e. Scheme 1. The table below shows various rate constants used in conjunction with the thermal mechanism. These reaction rates have the form, ? ( ) (3.149) 68 Elementary Rxn A b Source REACTION EQUATION 1 1.58E14 0 38031 Mellor [14] 7.6E13 0 38000 Heywood [9] 3.3E12 0.2~0.3 0 Rublewski [23] 1.8E14 0 38370 Zeldovich [44] 3.5E13 0 330 Miller [35] REACTION EQUATION 2 2.65E12 0 3226 Mellor [14] 6.4E9 1 3150 Heywood [9] 6.4E9 1 3160 Rublewski [23] 1.8E10 1 4680 Zeldovich [44] 2.65E12 1 6400 Miller [35] REACTION EQUATION 3 7.33E13 0 564 Mellor [14] 4.1E13 0 0 Heywood [9] 3.8E13 0 0 Rublewski [23] 7.1E13 0 450 Zeldovich [44] 7.33E13 0 1120 Miller [35] Table 3-2: Various Rate Constants for EZM It can be seen from the table above that the majority of these established reaction rates are very similar. Despite the small differences between these rates, each of the above set of rates will be subjected to analysis. Once the optimum set of rates from these five sources is determined a more general sensitivity analysis will be performed on that specific set of reactions. This post-sensitivity analysis of the chosen rates will subject the rate constants to a multiplicative and a reducing factor, first to the set of rate constants, then individually. This will allow an examination of the sensitivity of each individual reaction. 69 3.4.4 Temperature Approximations In the majority of the models that use the adiabatic layered core concept, the temperature profile from which the adiabatic temperature is determined is normally the unburned air temperature. However in diesel engines, the actual unburned air temperature can substantially stray from the modest profile of the ideal unburned air temperature. Unlike SI engines, diesel engine combustion does not necessarily conform to a flame-front approach within the cylinder, which is often the approximation of simple one-zone models, instead diesel combustion as stated in the literature review, burns in a turbulent diffusion flame as fuel and unburned air mixed together. This non-uniform, heterogeneous combustion process leads to the formation of highly radiative soot particles. The radiation from these soot particles is approximately five times the radiation from the gaseous combustion products. It is estimated that 20-35% of the total heat transfer within the diesel combustion chamber is attributed to radiative heat transfer. With this in mind it is apparent that the initial temperature from which the adiabatic flame temperature is calculated from, is worthy of investigation. Figure 3-5 illustrates the differences between the unburned and average cylinder air temperatures. It can also be seen in Figure 3-5 that an approximation with the unburned charge temperature alone would fail to illustrate some of the critical NO formation regions i.e. high temperature and high pressure regimes. From the assumption that the local unburned air temperature will be located somewhere between the mixed and ideal unburned air temperatures, the temperature that is supplied for the adiabatic flame temperature calculation will be varied between these two boundaries and analyzed. This stratification of the unburned charge temperature will look similar to Figure 3-6. This analysis should supply a good approximation as to the appropriate local unburned air temperature. 70 Figure 3-5: Comparison of Unburned Air Temp and Average Cylinder Temp Figure 3-6: Temperature Approximations T D C 10 20 30 40 50 60 70 80 90 600 800 1000 1200 1400 1600 1800 2000 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) U n b u r n e d A i r T e m p e r a t u r e A v e r a g e C y l i n d e r T e m p e r a t u r e T D C 10 20 30 40 50 60 70 80 90 600 800 1000 1200 1400 1600 1800 2000 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) U n b u r n e d A i r T e m p e r a t u r e A v e r a g e C y l i n d e r T e m p e r a t u r e A p p r o x i m a t e d L o c a l U n b u r n e d C h a r g e T e m p 71 3.4.5 Multiple Injections Unfortunately for this study, control over the fuel injection timing, shape, and pulse counts is not available. This lack of control leads to uncertainties in the cycle simulation concerning how to treat the areas between each injection. A figure showing some experimental results from this study of heat release vs. injector current can be seen below. Figure 3-7: Heat Release and Injector Current At high loads (2250 RPM and higher) Figure 3-9 is generally representative of the injection cycle and heat release. At light loads this figure changes to include only two injection pulses. This inconsistency between injection schedules will knowingly introduce errors in the modeling. As seen in the figure, the first injection causes little change in the overall heat release curve in the vicinity surrounding it. This lack of change in this region could be explained by the 72 temperatures and pressures residing near the first injection not supporting auto-ignition. Based upon this assumption, the first injection will first be neglected. To determine the effects of neglecting injection pulses, a sensitivity analysis will be performed, first, including both injections and treating this as one combustion regime, and second, neglecting the second injection and treating the main injection as the only residing injection schedule. 3.5 Model In order to develop a proper model we must consider formation inside the cylinder as well. Chemical equilibrium considerations suggest that at temperatures typical of in cylinder flame temperatures, ratios should be negligibly small. While this assumption is true for spark ignition engines, in diesels can constitute from 10% to 50% of the total oxides of nitrogen emissions. Unlike NO formation, formation is much more complex and is dependent on operating parameters such as IMEP, RPM, and the local concentration of NO. The highest levels of are found at low RPM and low load percentages. The highest levels at low load are consistent with the explanation by Heywood [9] stating that, ?The highest levels of occur at light load when cooler regions which could quench the conversion of back to NO are widespread.? For low RPM the concept would be the same but would be dependent on the mixing rates inside the cylinder. Thus at low RPM?s less mixing occurs therefore large regions of cooler unmixed air reside, which can quench the conversion of back to NO. An additional concept to consider for smaller concentrations at higher loads would be the effect of higher gas temperatures. With the increase of fuel and thus heat release, higher temperatures allow the conversion of and other intermediate species back to NO more easily. It is evident from the 73 explanation above that formation must be expressed in terms of operating parameters, rather than simple temperature dependent reaction mechanisms, as was the case with NO formation. Many researchers tend to either combine both and NO concentrations into one lump sum or use a multiplicative factor to make concentration corrections to determine a total engine out . While this approach may be appropriate for spark ignition engines where concentrations are relatively low, this concept is incorrect when applied to diesel engines due to their large concentrations of . When using the principles in sections 3.3 and 3.4 the resulting solution of [ ] is only valid for concentrations of NO and should not be applied to a total concentration that includes lumped sums of species other than NO. An example of results using this lumped concentration to predict engine out levels is shown in Figure 3-10. Figure 3-8: Lumped Model 70 80 90 100 110 120 130 140 0 200 400 600 800 1000 1200 I M E P ( p s i ) S p e c i e s C o n c e n t r a t i o n ( p p m ) L u m p e d N O x V a l u e P r e d i c t e d L u m p e d N O x V a l u e A c t u a l N O V a l u e A c t u a l N O 2 V a l u e 74 As can be seen in the Figure, the predictions for the total engine out are greatly over- predicted. This error is caused by the rapid drop in the concentration of indicated by the green circle. At light load the relative concentration of to the total lumped value is 41%, but as load is increased the concentration of quickly drops to only contribute a total of about 1% to the total lumped concentration. This rapid decrease in the concentration of should provide explanation to the conclusions of other researchers claiming that the EZM or extensions/variations of this mechanism greatly over-predict concentrations. The concept of a lumped model will be explained in greater detail in a later section. With the analysis presented above it is clear that we must treat the concentrations seperately from that of other species. An attempt at developing a reaction mechanism to predict concentrations was undertaken, but sufficient lack in confidence of the rate coefficients and the lack of a present reaction mechanism that has been validated lead to a more empirical method of predicting concentrations. Taking the concepts presented at the beginning of the section we can represent the concentrations of as a function of the primary controlling variables. [ ] (3.150) By analyzing the Figure 3-11 shown below, we can approximate the general shape of the concentrations curve having a simple cubic shape. This general curve can be represented by, (3.151) 75 Figure 3-9: Data Trend By using a cascading method of least squares we can develop an expression that will predict concentrations for the entire load and RPM range based upon the primary variables. To use this cascading method of least squares we first introduce three different cubic equations based on the primary variables. [ ] (3.153) [ ] (3.152) 70 80 90 100 110 120 130 140 - 4 0 - 2 0 0 20 40 60 80 100 120 140 160 I M E P ( p s i ) N O 2 C o n c e n t r a t i o n ( p p m ) 76 [ ] [ ] [ ] [ ] (3.154) Solving each individual equation using the method of least squares for parameters will yield the appropriate first approximations of tuning coefficients for the complete model equation. The number of each individual primary controlling variable should be proportional to the number of data points in the tuning range i.e. tuning 20 data points will have 20 different values for RPM, IMEP, and [NO]. Once the initial approximations of parameters have been determined, we can perform the same method of least squares for the complete model equation. This equation will have the form, [ ] [ ] [ ] [ ] (3.155) The solution of parameters will provide Equation (3.155) with the proper tuning coefficients to be capable of predicting concentrations over the desired tuning range. 3.6 and Model Diesel engines are superior to other energy conversions systems due to their great combustion characteristics. Combustion efficiency for direct injection diesel engines is usually around 90% which allows the assumption that the combustion process is essentially complete. With this assumption a simple formulation can be made to determine the general percentages of carbon dioxide and oxygen found in the exhaust. The combustion equation for a hydrocarbon fuel can be represented as, 77 ( ) [ ] [ ] [ ] [ ] (3.156) Using the complete combustion assumption correlations for the amounts of , , , and in the exhaust gas can be developed. The relative molar amounts of each species can be represented by, ( ) ( ) (3.157) To find the relative percentages of each species involves simply dividing each species molar amount by the total molar amount. This method is applicable to diesel engines due to their lean operation and great combustion characteristics and thus yields good results. 3.7 Program Structure A model diagram depicting the general flow of the model and its interactions with sub- models is shown below. The outputs from the combustion model are the initial inputs to the main program. The outputs from the combustion model were initially incremented at 0.36 degrees, but this data is shortened to 1 degree once it enters the main program to reduce computational time. Among the first steps in the model include determining the appropriate unburned charge temperature, and whether or not to include heat loss from the flame temperature. Once these parameters have been identified the adiabatic flame temperature of the 78 first layer can be calculated. Once the adiabatic flame temperature is calculated for each layer, the corresponding temperature profile for the remainder of the layer can be calculated according to a reversible process. The next step in the program is to evaluate the equilibrium concentrations and determine the level at each crank angle for each layer from the specified reaction mechanism. The steady state value from each layer is then mass averaged and a final concentration can be obtained. 79 Yes Model Flow Diagram No Yes Yes Yes No No No Combustion Model Outputs Reduce Data Array Length For all Data Entries For [ ] Calculate Unburned Charge Temp For i @ start =SOC Calculate Adiabatic Temperature ; Heat Loss Compute , Woschini Calculate air Properties from Gibbs For [ ] Calculate Layer Temp are now specified For Specify Initial Conditions Compute Equilibrium Concentrations from Gibbs Equilibrium Compute Rate Constants [ ] [ ] ? [ ] Solve, 80 3.8 Solver for Ordinary Differential Equations Numerical computations to solve chemical kinetics of reacting systems for species concentrations often require significant computational time due to the stiffness inherent in chemical kinetics computations. Stiffness arises in the chemical kinetic ODE?s due to the varying time constants involved in a reaction system. In order to solve these ODE?s several different MATLAB built in solvers were used. The first solve built in solver used was the MATLAB built in ode15s. The solver ode15s is a variable order solver based on numerical differentiation formulas. It is a multistep solver that optionally uses the Gear backward differentiation formulas. The other solver ode23s is based on a modified Rosenbrock formula of order2. Because this is a one step solver it is more efficient in some cases than ode15s at crude tolerances. 81 Chapter 4: Experimental Setup 4.1 Test Stand and Dynamometer The test stand uses a 1.1L three cylinder turbocharged diesel engine from Hyundai/Kia of South Korea. The engine specifications can be seen in Table 4-1. Manufacturer KIA/HYUNDAI Type I-3 DOHC Displacement 1120cc Bore x Stroke 75.00mm x 84.50mm Compression Ratio 17.80:1 Fuel System Common Rail Direct Injection Aspiration Turbocharged Max. Output 74.0 bhp at 4000RPM Max Torque 113 lb-ft at 1900 RPM Redline / Idle 4500RPM / 890RPM Table 4-1: Experimental Engine Specifications The engine load is applied by a 13? toroidal dynamometer which is controlled via a servo- operated load valve. The test cell was recently moved from a renovated test cell location into a new laboratory that was designed specifically for engine testing. The new engine test cell was equipped with ventilation systems for both ambient air exchanges and exhaust gas removal. The ventilation system allowed the removal of the old cooling system, which consisted of a shell and tube heat exchanger, and reinstated the stock radiator. The fuel system uses a five gallon racing fuel tank and delivers fuel by gravity feed to the fuel filter and then to the engine. The components of the test stand were for the most part complete but a considerable amount of work was required to adapt the test stand to the new location. 82 Figure 4-1: Experimental Test Cell Figure 4-2: Experimental Test Bed Setup 83 The dynamometer is a 13? single-stator absorber from Land and Sea that operates at 10 PSI dynamic pressure and 30 PSI static pressure for testing purposes. The water flows through a pressure reducing valve from a ?? city main supplied at 90 PSI. A paddle wheel flow sensor is installed between the pressure regulator and the servo load valve to measure water flow going to the dynamometer. A servo motor is used to control the position of the rising stem valve in the load valve to regulate water flow to dynamometer. 4.2 Sensors Two different groups of sensors are used to obtain all the necessary information about the operating condition of the engine. The following sensors are used to measure the mean operating parameters of the engine and are on a relatively slow time base (sampling 2 to 10 Hz). Temperature: Intake (at Filter) Intake (at Manifold) Oil (at crankcase) Coolant (pump outlet) Exhaust (downpipe) Dyno Outlet Omega K type thermocouples National Instruments TC01 Conditioning Modules Intake Manifold Pressure Calibrated Voltage input from OEM Sensor Load Torque Full Bridge Strain Gauge(Land and Sea) National Instruments SG-04 conditioning modules Mass Air Flow Calibrated Voltage input (counter) from OEM sensor Air/Fuel Ratio Innovate Motorsports LC-1 lambda sensor Engine Speed BEI H25D Incremental Encoder Dyno Water Flow Omega FP7001A paddlewheel flow sensor Table 4-2: Mean operating Sensors (Slow Time Base) The fuel cell is placed on a SVI-20B digital scale. This scale has a RS-232 input and a measurement is taken every two seconds to determine the overall fuel consumption. Other 84 sensors that take measurements based on a crank angle basis, and operate at very high sampling rates (50,000 to 150,000 Hz) are shown below. In-cylinder Pressure Kistler 6056A Kistler 5010B Charge Amplifier Injector Current Tektronix A6302 Current Probe Tektronix TM502A Charge Amplifier Fuel Rail Pressure Calibrated Voltage input from OEM Sensor Crank Angle Position BEI H25D Incremental Encoder Table 4-3: High Time Base Sensors 4.3 Emissions Analyzer The emissions analyzer used in this study is the ECOM J2KN. The J2KN is a portable microprocessor controlled, electrochemical sensor based analyzer. The J2KN is suitable for a wide variety of applications such as: Boiler/Burner Flue Gas Analysis, Engine Emissions Analysis, Pollution Source Testing, and Environmental Reporting of Stack Gas Parameters. The J2KN provides measurements for , CO (low range 0-4,000ppm), and/or CO (high range 0- 40,000ppm), NO, , ( ), , and combustibles. This unit will also measure gas and ambient temperatures, stack draft, and perform a standardized smoke test. Calculated parameters of this unit include , combustion efficiency, excess air, and losses. The unit comes with an onboard printer for hard copy printing of vital stack parameters and data-logging is available via on-board memory card. The unit also comes with real time transmission to a PC via wireless interface device. A picture of the device and a table of technical specifications are provided below in Figure 4-3 and Table 4-4. 85 Figure 4-3: ECOM J2KN Emissions Analyzer M E A S U R E M E N T ( K E Y ) R A N G E A C C U R A C Y R E S O L U T IO N S E N S O R L IF E S E N S O R T Y PE O x y g e n ( O ) 0 - 2 1 v o l % 0 . 1 v o l % 0 . 1 v o l % 2 Y e a r s E l e c tr o c h e mi c a l Ca r b o n M o n o x i d e ( C) 0 - 4 , 0 0 0 p p m ? 2 % M e a s u r e d * 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l Ca r b o n M o n o x i d e ( V ) 0 - 4 0 , 0 0 0 p p m ? 2 % M e a s u r e d * 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l N i tr i c O x i d e ( N ) 0 - 4 , 0 0 0 p p m ? 2 % M e a s u r e d * 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l N i tr i c O x i d e ( . N ) 0 - 4 0 0 p p m ? 2 % M e a s u r e d * 0 . 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l N i tr o g e n D i o x i d e ( X ) 0 - 5 0 0 p p m ? 2 % M e a s u r e d * 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l N i tr o g e n D i o x i d e ( . X ) 0 - 5 0 p p m ? 2 % M e a s u r e d * 0 . 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l S u l f u r D i o x i d e ( S ) 0 - 5 , 0 0 0 ? 2 % M e a s u r e d * 1 p p m 3 Y e a r s E l e c tr o c h e mi c a l Co mb u s ti b l e s ( H ) 0 - 6 . 0 0 % v o l . ? 2 % M e a s u r e d * 0 . 0 1 % 5 y e a r s P e e l l i s te r G a s T e mp e r a tu r e 3 2 ? - 1 8 0 0 ? F ? 2 % M e a s u r e d * 1 ? F 1 0 Y e a r s N i Cr N i A mb i e n t T e mp e r a tu r e 0 - 2 5 0 ? F ? 2 % M e a s u r e d * 1 d e g F 1 0 Y e a r s S e mi - Co n d u c tor D r a f t / P r e s s u r e ? 4 0 " H 2 o ? 2 % M e a s u r e d * 0 . 1 % H 2 O 1 0 Y e a r s D M S Ca r b o n D i o x i d e 0 - 2 0 % ? 5 % M e a s u r e d * 0 . 1 % v o l % 4 Y e a r s I n f r a r e d Ca r b o n M o n o x i d e 1 0 0 0 - 6 3 , 0 0 0 p p m ? 5 % M e a s u r e d * 1 0 p p m 4 Y e a r s I n f r a r e d O 2 Co r r e c ti o n 0 - 2 0 % O x y g e n S mo k e S c a l e 0 - 9 Ca r b o n D i o x i d e ( CO2 ) 0 - CO2 ma x o f f u e l Ca l c u l a te d E f f i c i e n c y 0 - 9 9 . 9 % Ca l c u l a te d E x x e s s A i r ( L a mb d a ) 1 - I n f i n i ty Ca l c u l a te d Table 4-4: ECOM J2KN Technical Specifications 86 4.4 Data Monitoring and Recording In combustion related studies parameters like cylinder pressure and other related events are typically referenced with a crankshaft position. The BEI H25D incremental encoder is used to measure the crank angle position and has a resolution of 1000 cycles per revolution. The encoder is mounted directly to the driveshaft with a rigid metal coupling. A National Instruments PCI-6259 Data Acquisition card, along with National Instruments SCB-68 breakout boards, are used to input the measured signals from the high sampling rate sensors into the computer. The slow based time scale sensors use National Instruments SC-2345 boards for data transmission to the computer. Measured signals are scanned at a constant rate from the DAQ card. This scan rate is defined in a user written labVIEW program. The labVIEW programs used in this study were custom designed by Wingo [42]. The PCI-6259 is a 16-bit Multichannel DAQ card capable of 1.25 mega-samples per second [MS/s]. The labVIEW programs designed by Wingo [42] samples the data off of the card at each pulse of the incremental encoder to capture cylinder pressure, injector current, and fuel rail pressure at each resolute crank angle. These LabVIEW programs store the data in arrays and outputs them in Excel format for easy use with Matlab. Separate LabVIEW programs are used for engine fuel and dynamometer load control. A NI PCI-6025E DAQ card and SCB-68 breakout board are used to output the signals to the ECU and servo motor on the load valve. The engine fuel load control sends variable voltage signal to the accelerator position sensor. The dynamometer load control is a pulse with modulation (PWM) signal that determines the position of the servo motor on the load valve. Water flow into the dynamometer is measured by a paddlewheel flow sensor. The frequency of the signal is measured through a counter input and converted to gallons per minute. 87 4.5 EGR Control In this study emission modeling is the primary objective. The effects of EGR on emission levels can be drastic, so it is of great importance that control over this operating parameter is permitted. The original EGR valve is actuated by a PWM signal initiated from an ECU command. The solenoid operated EGR device opens and closes the spring loaded valve depending on the amperage, duty cycle, and frequency of the command sent by the ECU. While this is an effective way of on the fly control over the amount of EGR permitted into the cylinder, it can lead to drastic errors in modeling due to its inconsistency and variable valve opening percentages. In order to control the amount of EGR a simple mechanical device machined from a piece of aluminum is installed in place of the solenoid to allow the user to choose the amount of EGR permitted. This device is shown below in Figure 4-4. 88 Figure 4-4: EGR Control Apparatus By adjusting the bolt on top of the base plate, the user can control the position of the rising spring valve thus allowing EGR control. The relative proportion of EGR% can be calculated by assuming that the total air flow taken at the filter before throttling the EGR valve remains the same. With this assumption the EGR% can be calculated by monitoring the drop in the air flow measured at the filter and comparing this value to the un-throttled value of MAF at the filter. 89 Chapter 5: Experimental Results In this chapter the test procedure along with the test matrix will be presented. Load sweep data is taken at three different load points including 50%, 75%, and 100% load values. This load value is the calculated load value given from the ECU. These load values were tested encompassing the entire effective RPM range (1500-3000 RPM). For EGR sweep testing, specific RPM points were chosen and the EGR sweep was conducted at these RPMs with varying load percentages. The EGR percentages tested were 25%, 50%, and 75% recycled gas. The EGR percentage is calculated as the ratio of intake air flow at the filter with EGR valve closed, to that of the air flow with EGR valve open. 5.1 Test Procedure and Test Matrix In engine testing it is imperative that all tests are to be performed under similar conditions, in the same sequences, and under the same time constraints. It was observed over this testing period that even a small amount of humidity (caused from rain) in the supply air to the engine caused a drastic variation in test results. The test procedure used in this study is presented below. 90 START-UP/WARM-UP 1 Calibrate injector current amplifier, wideband O2 sensor, dynamometer strain gauge, and ECOM J2KN 2 Launch LabVIEW and other engine monitoring programs. 3 Turn on engine exhaust fan, room ventilation, dyno water supply, and additional cooling fans 4 Warm engine to and maintain intake temperature under 5 Shut Down Engine MOTORED TRACE 6 Unplug fuel injector wire for #1 cylinder 7 Start engine and engine monitoring programs. Motor engine at desired RPM 8 Run the high speed LabVIEW program to record 8 runs once engine is at the warm up condition 9 Shut off engine and attach fuel injector wire to #1 cylinder TEST 10 Start engine and achieve desired speed and load. Start engine monitoring program 11 Run high speed LabVIEW program to record 8 runs once engine has reached a steady state condion. 12 With engine at steady state start ECOM analyzer continous log, start wideband O2 sensor log, start main monitoring program, and start fuel consumption program log 13 Continue the previous step for a minimum of 5 minutes 14 Shut down engine 15 Do post calibration test on all equipment to ensure no variations in test equipment has occured over the test period 16 CRUNCH DATA Table 5-1: Test Procedure Outline A test matrix which contains all the points tested is shown below. Those in yellow represent points which were subjected to an EGR sweep. TEST RPM LOAD % EGR % TEST RPM LOAD % EGR % 1 1200 50 12 2250 50 2 1200 75 13 2250 75 3 1500 50 25-75 14 2250 100 4 1500 75 25-75 15 2500 50 25-75 5 1500 100 16 2500 75 25-75 6 1800 50 25-75 17 2500 100 7 1800 75 25-75 18 2750 50 8 1800 100 19 2750 75 9 2100 50 25-75 20 3000 50 10 2100 75 25-75 21 3250 50 11 2100 100 22 3500 50 Table 5-2: Test Matrix 91 5.2 Load Sweep Results A compilation of the load sweep results is shown below. The experimental NO concentrations followed the expected trend of increasing concentrations with increasing IMEP. The table also shows a significant amount of at light loads. NO Concentrations (ppm) Calculated Load Value RPM 50% 75% 100% 1200 263 433 NA 1500 235 387 649 1800 200 367 678 2100 186 339 635 2250 211 453 600 2500 186 358 635 2750 196 394 NA 3000 226 NA NA 3250 271 NA NA 3500 297 NA NA Table 5-3: NO Load Sweep Results (ppm) It can be seen in the table that NO levels generally increase as the calculated load value rises. This increase in NO levels with load can be explained by the higher peak pressures and temperatures that occur as the CLV is increased resulting in larger regions of close to stoichiometric burned gas and hence higher NO levels. This suggests that NO formation is roughly proportional to the injected fuel mass. The relative NO concentrations as a function of engine speed seem to be higher at lower RPMs, decrease to their lowest values around the peak torque RPM (2250 RPM), and then increase again as the engine speed surpasses this peak torque RPM. This could possibly be explained by the amount of NO decomposition that occurs as a function of engine speed. At light loads there exist in the cylinder large pockets of cooler 92 regions that are not affected by the burned gas regions. These cooler regions interact with NO formed from the burned gas regions by freezing the NO chemistry. This freezing of NO chemistry also can occur at higher engine speeds where there are much better mixing rates limiting the amount of NO decomposition that can occur resulting in higher NO concentrations. Concentrations (ppm) Calculated Load Value RPM 50% 75% 100% 1200 154 252 NA 1500 159 272 8 1800 142 135 6 2100 142 109 3 2250 123 76 3 2500 93 96 3 2750 51 20 NA 3000 92 NA NA 3250 56 NA NA 3500 63 NA NA Table 5-4: Load Sweep Results (ppm) concentrations can be seen to be significant in this study, with the highest concentrations occurring at low RPM. These high concentrations at low load conditions are indicative of pockets of cooler air hindering the conversion of back to NO. Figure 5-1 and Figure 5-2 show the relative grouping of these results of NO and as a function of IMEP. 93 Figure 5-1: IMEP vs. NO Figure 5-2: IMEP vs. 200 250 300 350 400 450 100 200 300 400 500 600 700 I M E P ( P S I ) N O p p m I M E P v s . N O 1 5 0 0 R P M 1 8 0 0 R P M 2 1 0 0 R P M 2 2 5 0 R P M 2 5 0 0 R P M 200 250 300 350 400 450 0 50 100 150 200 250 300 I M E P ( P S I ) N O 2 p p m I M E P v s . N O 2 1 5 0 0 R P M 1 8 0 0 R P M 2 1 0 0 R P M 2 2 5 0 R P M 2 5 0 0 R P M 94 5.3 EGR Sweep Results Below is a table showing the results from the EGR sweep tests. NO Concentrations (ppm) EGR % RPM-CLV 0% 10% 20% 30% 1500-50% 235 172 135 90 - - - - - 1800-50% 200 183 137 79 1800-75% 367 290 162 100 2100-50% 186 167 112 67 2100-75% 339 242 162 87 2500-50% 186 146 111 73 2500-75% 358 284 182 111 Table 5-5: NO Load Sweep Results (ppm) (ppm) EGR % RPM-CLV 0% 10% 20% 30% 1500-50% 159 85 67 39 - - - - - 1800-50% 142 85 57 24 1800-75% 135 44 19 3 2100-50% 142 62 46 21 2100-75% 109 23 8 1 2500-50% 93 38 32 16 2500-75% 96 24 7 1 Table 5-6: Load Sweep Results (ppm) The drastic effect of charge dilution can be seen in these tables. Two main processes contribute to the reduction of concentrations through EGR dilution. The first is oxygen displacement. By the addition of already burned gases the initial charge has less oxygen than that of the fresh air, this reduction in oxygen levels leads to lower combustion rates, lower peak pressures, and hence lower temperatures. The second effect explaining this drastic drop in local concentrations is that the diluents act as a heat sink in the combustion process. Diluents with high heat capacities such as and absorb much of the heat produced from the combustion processes causing a drop in combustion temperatures leading to lower formation rates. 95 Chapter 6: Parameter Identification for NO model 6.1 Introduction In this chapter the model perturbations presented in Chapter 3.4 will be analyzed. An in- depth analysis on all perturbations will be presented along with the simulation results for each. The purpose of this chapter is to identify the appropriate parameters that will become the final NO model implemented into the single-zone combustion model. For each perturbation or sensitivity analysis, the parameters that are used as reference or base parameters will be specified. The analysis? to be presented in this chapter include appropriate temperature approximations, rate constants, reaction mechanisms, species treatment, and injection scheduling. 6.2 Temperature Approximation Analysis As has been stated throughout this paper, NO formation and destruction is extremely dependent on local burned gas temperatures. Due to this fact, the determination of an appropriate unburned gas temperature for the calculation of the adiabatic temperature will be analyzed first. Assuming that the thermal or EZM will be the dominant reaction mechanism, this will suffice as the base reaction mechanism. As for the rate constants, those supplied in Heywood [9], will be used as base rate constants. The reaction set will be solved as a kinetically controlled set, and the combustion cycle will be treated as one large cycle beginning at the second injection. As was stated in section 3.4.4 there is a large difference between the average cylinder temperature and that of the unburned charge temperature calculated by from a reversible process. Between these two temperature profiles will serve as a starting point for the analysis. Treating 96 the average mixture temperature as the upper boundary and the unburned temperature as the lower boundary, a set of temperature profiles are created between these two boundaries and analyzed for their validity. For the first analysis a very broad range of temperature profiles were tested which encompassed the entire boundary. Figure 6-1 shown below represents a depiction of the temperature profiles tested. Figure 6-1: Tested Temperature Profiles From the analysis, the unburned charge temperature seemed to suffice for low load conditions but under-predicted NO levels at higher loads. A possible explanation for this phenomenon could be due to the various mixing rates and heat release occurring inside the cylinder. For low fuel load conditions we can expect the greater percentile of the cylinder contents to be at the unburned cylinder charge temperature while at higher fuel loads and thus higher heat release rates, we can expect the local unburned temperature to stray from this unburned charge T D C 20 40 60 80 100 120 140 160 180 200 400 600 800 1000 1200 1400 1600 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) A p p r o x i m a t e T e m p e r a t u r e P r o f i l e U n b u r n e d T e m p e r a t u r e A v e r a g e T e m p e r a t u r e 97 temperature and favor the shape of the cylinder average temperature. A table showing the results of the unburned charge temperature approximation is shown below. The letter A represents the tuning constant used for the rate limiting reaction R1 i.e. ? Load Value 50% 75% 100% A RPM Exp. Est. Exp. Est. Exp. Est. 2.3617E12 1800 200 200 367 283 678 538 2.3617E12 2100 186 196 370 236 632 468 2.3617E12 2250 211 210 453 408 600 436 2.3617E12 2500 186 200 407 257 635 384 Table 6-1: Unburned Charge Air Temp. NO Results (ppm) As can be seen in Table 6-1, with the use of the unburned charge temperature those points at 50% load value are predicted to within 5% while the higher loads are under-predicted by 35% or more. In addition to the prediction errors, the use of the unburned temperature seemed to yield incorrect results with regards to the NO production rate shown in Figure 6-2. 98 Figure 6-2: NO Production with Unburned Air Temp. It can be seen in Figure 6-2 that the NO produced by the second injection is almost equal to the amount of NO produced by the main injection. In the actual combustion process the relative proportions of NO produced by each injection should be proportional to the amount of fuel injected for each injection. To be an applicable model these results will not suffice. As the tested approximate temperature approached the average of the unburned charge temperature and the average cylinder temperature, the collective average error between these twelve points seem to improve. If the approximated temperature surpassed the average of the unburned charge temperature and the cylinder average temperature the collective error seemed to increase. A figure showing the progression of the average error as a function of the approximated temperature is shown below. - 1 0 0 10 20 30 40 50 60 70 80 90 0 5 10 15 20 25 C r a n k A n g l e ( d e g ) N O ( p p m ) 99 Figure 6-3: Approximated Error with Progressive Temperature Profiles The x-label in this figure represents the different temperature profiles being used from the lower boundary ( ) to the upper boundary ( ). The figure below represents the progression of the x-axis. Figure 6-4: Progression of Tested Temperature Profiles 0 2 4 6 8 T a p p 12 14 16 18 20 10 11 12 13 14 15 16 17 18 19 20 A v e r a g e E r r o r o f 1 2 T e s t P o in t s A p p r o x i m a t e d A i r T e m p e r a t u r e P r o f i l e T D C 20 40 60 80 100 120 140 160 180 200 400 600 800 1000 1200 1400 1600 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) A p p r o x i m a t e T e m p e r a t u r e P r o f i l e U n b u r n e d T e m p e r a t u r e A v e r a g e T e m p e r a t u r e 100 As can be seen in Figure 6-3 the collective error falls sharply when approaching which is represented by, (6.1) Based on the results shown in Table 6-1 we can assume that the unburned air temperature is not a valid temperature approximation for all points and that an appropriate unburned air temperature will be located somewhere in the vicinity of . To further this investigation another set of temperature profiles are chosen for analysis but this time, are selected to closely surround . The temperature profiles used in this analysis can be seen below in Figure 6-5 where the thick black line indicates . Figure 6-5: Finite Temperature Analysis T D C 10 20 30 40 50 60 70 80 90 100 400 600 800 1000 1200 1400 1600 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) 101 Upon examining these temperature profiles in Figure 6-5, the best results were obtained for temperature profiles just after . The collective error as a function of the progressive temperature profiles, indicated by the black arrow in Figure 6-5, can be seen in Figure 6-6 shown below. Figure 6-6: Progressive Error of Finite Temperature Analysis In order to stay within the goal of this analysis which is within a 15% error for each individual point, the temperature profile chosen is profile number 14. The simulation results for this temperature profile can be seen in Table 6-2. 0 1 2 3 4 5 6 7 8 T a p p 10 11 12 13 14 15 16 17 18 10 11 12 13 14 15 16 A p p r o x i m a t e d A i r T e m p e r a t u r e P r o f i l e A v e r a g e E r r o r o f 1 2 T e s t P o i n t s 102 These results satisfy the tolerance of 15% error and have a collective average error of 7.84%. This temperature profile was created with the relation shown below, (6.2) This approximated temperature in comparison to the unburned and cylinder average temperature has the profile shown in Figure 6-7 below. Figure 6-7: Approximated Temperature Profile in Comparison with Average and Unburned Air Temp 10 20 30 40 50 60 70 80 90 600 700 800 900 1000 1100 1200 1300 1400 1500 1600 C r a n k A n g l e ( d e g ) T e m p e r a t u r e ( K ) U n b u r n e d A i r T e m p . C y l i n d e r A v e r a g e T e m p A p p r o x i m a t e d A i r T e m p Load Value 50% 75% 100% A RPM Exp. Est. Exp. Est. Exp. Est. 9.0503E11 1800 200 200 367 424 678 715 9.0503E11 2100 186 165 339 298 632 556 9.0503E11 2250 211 228 453 454 600 536 9.0503E11 2500 186 183 358 364 635 538 Table 6-2: NO results of Temp Profile #14 (ppm) 103 As can be seen in Table 6-2 the results predict relatively well at low loads but still under predict NO levels at higher loads. The NO production rate for this temperature profile also displays a more believable result than Figure 6-2 which is used in conjunction with the unburned air temperature. Figure 6-8: NO production with Temperature Profile #14 Figure 6-8 shows a more realistic result than Figure 6-2 with the production of NO relatively proportional to the fuel injected for each injection. Both the layer NO concentration profiles for the chosen approximated air temperature and the unburned air temperature are shown below in Figures 6-9 and 6-10 respectively. - 1 0 0 10 20 30 40 50 60 70 80 90 0 5 10 15 20 25 C r a n k A n g l e ( d e g ) N O ( p p m ) 104 Figure 6-9: NO Layer Concentration Profile Results with Approximated Air Temp. Figure 6-10: NO Layer Concentration Profile Results with Unburned Air Temp. - 1 0 0 10 20 30 40 50 60 70 80 90 0 50 100 150 200 250 300 C r a n k A n g l e ( d e g ) N O ( p p m ) - 1 0 0 10 20 30 40 50 60 70 80 90 0 200 400 600 800 1000 1200 N O ( p p m ) C r a n k A n g l e ( d e g ) 105 The differences between these temperatures can be seen in the figures by the discrepancy of over-lapping NO concentration profiles. These NO concentration profiles can be seen clearer by using a mesh approach shown below. Figure 6-11: NO Layer Profile Concentration Results with Approximated Air Temp. Figure 6-12: NO Layer Profile Concentration Results with Unburned Air Temp. 0 20 40 60 80 020 4060 80 0 50 100 150 200 250 50 100 150 200 250 T D C 50 100 T D C 50 100 0 200 400 600 800 1000 C r a n k A n g l e ( d e g ) N O C o n c e n t r a t io n ( p p m ) 0 100 200 300 400 500 600 700 800 900 1000 106 Using the mesh plots the effect of using the approximated temperature profile can be seen to accent the main injection better than that of the unburned air temperature. From this analysis presented above, temperature profile 14 will be used throughout the remainder of this study. As was presented in Chapter 3, a heat loss analysis was performed surrounding the temperature profiles in Figures 6-4 and 6-5. The results from this study indicated that the relative error of the measurements did not change drastically and that better results were obtained neglecting heat loss. It is for this reason that heat loss from the flame will be neglected for the remainder of this study and that temperature profile 14 will serve as the base temperature in the calculation of the adiabatic flame temperature. 6.3 Rate Constant Analysis Since the thermal mechanism will be the governing reaction set, the thermal mechanism will be used as the base set. Recall the set of reaction rates used in previous studies shown in Table 3-2. From the analysis, those reaction rate sets used by Mellor [14], Heywood [9], and Zeldovich [44] yielded the best results. For the most part the results were very close to each other with the reaction sets presented by Rublewski [44] and Miller [35] over-predicting results a few percent more than the other reaction rate sets. The critical difference between the reaction sets that must be noted is the use of or . For those sets used by Rublewski [23] and Miller [35] the use of is employed, which implies using the right side of Equation (6.3). [ ] [ ] [ ] [ ] (6.3) With the analysis between the different reaction sets now established to be small, the next analysis will look at the multiplicative factors of each set first as a group and then second individually. By multiplying all the reaction rates by a factor of 2 the corresponding NO concentrations increased by a factor of ~1.9. Through the analysis of each individual reaction 107 rate it was determined that reaction was the rate limiting reaction, controlling the majority of NO production. Reactions and have a relatively minor effect on the total NO concentration only affecting the output by 2 or 3%. 6.4 Reaction Mechanism Analysis In this section all of the mechanisms discussed in Chapter 3 will be tested for their validity. These mechanisms include schemes 1-5. The reaction rates used for each reaction can be found in Appendix B, and the corresponding temperature profile used is profile 14 from section 6.1. The results from Scheme 1 have already been presented in 6.2 but to clarify, the results will be reviewed. Scheme 1 included only the thermal mechanism and the results were within the range of 15% for all test points with an average error of 7.84%. Scheme 1 showed good agreement with experimental data at low RPM and low load conditions, but under-predicted NO levels at higher loads. The results for four different RPM values are shown below in Table 6-3. For each specified RPM three different load points are used including 50%, 75%, and 100% calculated load values. 108 Table 6-3: Scheme 1 NO Results (ppm) Scheme 2 is tested next with the addition of the mechanism. The results from this mechanism were very close to that of Scheme 1 with a slight reduction of predicted values. The results of Scheme 2 are shown below in Table 6-4. Table 6-4: Scheme 2 NO Results (ppm) Load Value 50% 75% 100% A RPM Exp. Est. Exp. Est. Exp. Est. 9.0503E11 1800 200 200 367 424 678 715 9.0503E11 2100 186 165 339 298 632 556 9.0503E11 2250 211 228 453 454 600 536 9.0503E11 2500 186 183 358 364 635 538 Load Value 50% 75% 100% - RPM Exp. Est. Exp. Est. Exp. Est. - 1800 200 200 367 418 678 700 - 2100 186 164 339 295 632 548 - 2250 211 226 453 447 600 527 - 2500 186 182 358 360 635 529 109 Compared with that of Scheme 1, Scheme 2 was relatively ineffective at low load conditions, but resulted in a reduction of predicted values for 75% and 100% loads of 1~2%. With the addition of Scheme 3 into the reaction mechanism the predicted values were relatively unaffected and only resulted in reductions of .2-.5%. This relative minor change suggests the minor importance of the mechanism. The next scheme included in the mechanism is Scheme 4 which includes the addition of oxidation reactions involving NH and HNO. The results of the inclusion of this mechanism had similar results to the addition of Scheme 2, with reducing predicted values. The results of this scheme are shown below in Table 6-5. Table 6-5: Scheme 4 Results (ppm) With the addition of Scheme 4 predicted values are again relatively unaffected at low load conditions but under higher load conditions are reduced by ~3.5%. Scheme 5 had similar results to that of Scheme 3 with relatively minor effects and therefore will not be discussed further. Overall the addition of an extra 41 reactions only lead to reducing predicted values by an average of 3%. It cannot be stated conclusively that this reduction in predicted values can be attributed to the extended reaction mechanism due to the Load Value 50% 75% 100% - RPM Exp. Est. Exp. Est. Exp. Est. - 1800 200 200 367 403 678 676 - 2100 186 167 339 288 632 535 - 2250 211 222 453 436 600 512 - 2500 186 182 358 349 635 511 110 small value. When dealing with small numbers such as molar fractions, truncation errors can be introduced through the computer which can lead to small fluctuations in final model results. The reduction seen by the addition of these mechanisms also did not increase the prediction accuracy of the model, but in fact made it worse. Some error could have been introduced in respect to the rate constants used. Most of these rate constants were determined through methane/air combustion experiments held at constant pressures (mostly atmospheric). The lack of knowledge surrounding some of the rate constants in the extended schemes could have played a role in being relatively ineffective. An important observation in this analysis is tuning of the reaction mechanism. For Scheme 1 tuning the model was accomplished by adjusting the parameter A in the Arrhenius expression to the match the NO concentration at the lowest IMEP test point. It was observed that with the addition of more and more reactions the previous rate limiting reaction became less effective and a different tuning approach was adopted. Those reactions that dominated the reaction mechanism included , , , , and . In order to tune the model with larger schemes a single reducing factor is applied to all reaction rates. The dominating behavior displayed by some of the reaction rates could be the result of out of date rate coefficients. Most of the rate constants listed in Appendix B were taken from studies done over thirty years ago from today?s date. In order to test the validity of the rate coefficients used in the extended mechanisms, rate coefficients taken from GRI MECH 3.0 are implemented. These rate coefficients are updated frequently and are the most current representation of rate coefficients for the majority of the reactions included in the extended schemes. The inclusion of the new rate coefficients lead to either minimal change in output or an increase in model error. Therefore those rate constants stated in Heywood [9] along with the EZM prove to yield the best results when coupled with the adiabatic layered core model. 111 6.5 Validation of the QSSA assumption The QSSA assumption as discussed in Chapter 3 is primarily used to help alleviate some of the difficulties that arise when dealing with stiff differential equations. While this assumption works well for fairly short reaction mechanisms, its derivation can get lengthy for complex reaction mechanisms. This section will look at simulation results of mechanisms treated both as a kinetically controlled set and with the QSSA assumption. The first mechanism that will be analyzed is the EZM. For Scheme 1 the QSSA assumption yielded the same result as did the kinetically controlled approach. This result was expected due to the very low concentrations of equilibrium N concentrations. The results are so close that no stipulation can be seen between the two approaches, therefore validating the QSSA assumption for the thermal mechanism. For Scheme 2 the QSSA assumption is very close to the kinetically controlled result with only a 2-5% difference. Figure 6-13 and Figure 6-14 show these differences between the kinetically controlled and QSSA assumption for 1800 RPM 50% CLV and 1800 RPM 100% CLV. As can be seen in the figures, the QSSA assumption tracks well for the majority of the NO trace but lacks some of the formation characteristics of the kinetically controlled NO trace. With the difference between the results of the kinetically controlled trace and the QSSA trace established to be negligible, the QSSA assumption can be considered an appropriate assumption. However, an important note to be considered when using the QSSA assumption that must be considered is while the QSSA assumption works well for small reaction mechanisms the derivation for larger mechanisms can become quite lengthy and have several complex roots resulting in tedious manipulation to achieve desired results. 112 Figure 6-13: Kinetically Controlled NO production vs. QSSA Controlled NO production Figure 6-14: Kinetically Controlled NO production vs. QSSA Controlled NO production - 2 0 0 20 40 60 80 100 0 100 200 300 400 500 600 700 800 C r a n k A n g l e ( d e g ) N O ( p p m ) K i n e t i c a l l y C o n t r o l l e d Q S S A A s s u m p t i o n - 1 0 0 10 20 30 40 50 60 70 80 90 0 50 100 150 200 250 C r a n k A n g l e ( d e g ) N O ( p p m ) K i n e t i c a l l y C o n t r o l l e d Q S S A A s s u m p t i o n 113 6.6 Effect of Multiple Injections on Simulation Results All the analysis above were treated as specified in Chapter 3, meaning, the combustion duration was treated as one single duration starting from the second injection for a three pulse injection and starting from the first injection for a two pulse injection scheduling. As is stated in Chapter 3, the first injection for a three pulse injection schedule is neglected due to the fact that the amount of fuel injected in this pulse was only 2-3% of the entire amount injected over the combustion cycle. In this section different approaches to handling the injection scheduling, and the start of combustion will be analyzed. There are two figures shown below, the first, Figure 6- 15, represents the total NO formation over one combustion cycle and the second figure, Figure 6- 16 represents the injection scheduling. As can be seen in Figure 6-15, there is an initial sharp rise of NO concentration due to the pilot injection. Figure 6-15: NO production with Pilot and Main Injections - 2 0 0 20 40 60 80 100 0 100 200 300 400 500 600 C r a n k A n g l e ( d e g ) T o t a l N O c o n c e n t r a t io n ( p p m ) 114 Figure 6-16: Heat Release and Injector Current Signal As can be seen in Figure 6-16 there are actually three injections while the simulation only includes the second and the main injection. By only including the second and the main injection the combustion cycle simulation begins at the second injection and is carried through the entire combustion cycle until EVO. By including the first injection and beginning the cycle simulation at the first injection we can expect the predicted NO concentrations to increase, due to the higher temperatures reached by some of the early burning elements. As can be seen in Figure 6-17 the inclusion of the first injection can be seen by a small initial rise due to the first injection. Shown in Figure 6-18 are the corresponding first 10 NO layer concentration profiles where the red lines indicate the three pulse schedule and the blue represent the two pulse schedule. The profiles - 2 0 0 - 1 5 0 - 1 0 0 - 5 0 T D C 50 100 150 200 - 1 0 0 0 0 1000 2000 H e a t R e l e a s e 0 - 0 . 1 0 0 . 1 0 . 2 I n j e c t o r C u r r e n t C r a n k A n g l e ( d e g ) 115 created by the three injection simulation completely trump those of the two injection simulation, but are controlled well due to the low amount of mass burned during the first injection. Figure 6-17: NO production with Pre-Pilot, Pilot, and Main Injections Figure 6-18: NO Layer Concentrations (Red Three Pulse Schedule/Blue Two Pulse Schedule) - 4 0 - 2 0 T D C 20 40 60 80 100 0 100 200 300 400 500 600 C r a n k A n g l e ( d e g ) T o t a l N O c o n c e n t r a t io n ( p p m ) T h r e e I n j e c t i o n s T w o I n j e c t i o n s - 4 0 - 2 0 0 20 40 60 80 100 0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 x 1 0 4 C r a n k A n g l e ( d e g ) L a y e r N O c o n c e n t r a t io n s ( p p m ) 116 For this particular load point which is 2500 RPM@100%CLV, the inclusion of the first injection improved the accuracy of the model prediction. The load points that contain the three pulse injection schedules are typically the 100% CLV points. The same trend of increasing model prediction accuracy is observed for the remaining points with the addition of the pre-pilot injection. It is for this reason that for 100% CLV points that the model will incorporate the effects of the pre-pilot injection. 6.7 Effect of Initial NO Concentration from Residual Gas Diesels generally have very small concentrations of residual gas in the cylinder left over from previous combustion cycles. This residual gas concentration is usually around 3~5%. For this study a constant residual concentration of 5% is assumed. A sensitivity analysis was performed to determine the effects of this small NO residual concentration as an initial condition on simulation results. The effect of adding the residual fraction increased the total NO concentration simulation output by only 3~5%. While for low NO concentrations (~200ppm) this small increase is negligible, for larger NO concentrations (~600ppm) this small percentage can lead to a considerable amount of added NO. However, with the initial condition premise added to the NO cycle simulation and tuned to fit the NO concentration at the lowest IMEP, the same results without the initial concentration consideration are achieved. By adding the residual NO concentration a resulting decrease in the tuning constant A will occur which will ultimately lead to the same simulation results that were achieved without considering the residual NO concentration, therefore the residual NO concentration will be omitted. An important observation to consider when adding the residual fraction is the transient nature of effects of residual gas or EGR. When a residual fraction or EGR is considered the cycle simulation must 117 be treated as a iterative process, meaning that the output of the simulation will serve as the new input for the initial concentration of residual NO. The process of re-assigning the output as an input continues until the difference between the output and the residual NO fraction input is negligibly small. An explanation of why this is applicable stems from the fact the measurement that is observed with the emissions equipment is actually a steady state value which is the result of a fast transient response to the initial condition. 118 Chapter 7: Final Model Validation 7.1 Introduction This chapter will begin by clarifying all the parameters that will be used in the final NO model based on the analysis done in the previous chapter. Once the parameters have been identified, a validation of the layered adiabatic core model will be introduced in comparison to a fully mixed adiabatic core model, along with the validation of the EGR model. Until this point in the study the NO model has been used in conjunction with experimental data to obtain the most accurate tuning results. The first 6 sections of this chapter will present the findings of each individual model based on the experimental data. In order to test the accuracy of the coupling of the combustion model and the emission model, the NO model will be implemented into the combustion model and the results will be presented. 7.2 Parameters for the Layered Adiabatic Core NO Model This section will serve to clarify all of the chosen parameters used in the NO model. For the approximation of the unburned temperature used in the calculation of the adiabatic temperature it is found that profile 14 is the most accurate choice across the experimental data range. This approximated temperature is calculated by the expression below. (7.1) Since little or no change is observed between reaction mechanisms, the Extended Zeldovich Mechanism will serve as the reaction mechanism along with the rate constants found in Heywood [9]. For each layer the mass averaging technique discussed in Chapter 3 will be used to determine the total NO concentration for each layer. For tuning the model only the parameter 119 A from the Arrhenius expression is used. This parameter is tuned for the lowest IMEP point within the data set and is not changed over the entire data range. The local equivalence ratio used in determining the adiabatic flame temperature and the equilibrium concentrations is set to . With the availability of computers, solving simultaneous differential kinetic equations is an easy task, thus the QSSA technique will not be considered and the reaction mechanism will be solved as a kinetically controlled species set. For injection scheduling, the simulation includes all injections from pre-pilot to the main injection. The simulation starts at the first SOC point determined from the first injection and continues until EVO. 7.3 Validation of the Layered Adiabatic Core NO Model The validation of the layered adiabatic core model is accomplished in this study by comparing the results to that of a fully mixed adiabatic core model. The fully mixed adiabatic core model is treated with all the base assumptions stated in section 7.2 except for the layered approach. The fully mixed model predicts an adiabatic temperature at each crank angle and solves the governing kinetic differential equations to determine the amount of NO formation for each crank angle. The results from the layered adiabatic core and the fully mixed adiabatic core model are shown in Table 7-1 and Table 7-2 respectively. From the tables shown below, the fully mixed adiabatic core model seems to predict well at low load conditions but seems to lack the predictive capability for higher load conditions. This under-prediction from the fully mixed model is consistent from analysis from other researchers as well. The layered adiabatic core model predicts well at both low and higher load conditions consistent with a 15% error boundary. The single tuning of the parameter A to the lowest IMEP test point is the only tuning parameter used in both the fully mixed and the layered model, limiting the amount of empirical 120 adjustments to only one. There is of course an alternative tuning method to tune both of these NO models by adjusting the parameter A based upon loading conditions. By adjusting this parameter over different loading conditions better results can be obtained, however, this adjustment will limit the predictive capability of the model between different engines and thus cause the corresponding model to be engine specific. Table 7-1: NO Results from Layered Routine (ppm) Table 7-2: NO Results from Fully Mixed Routine (ppm) Load Value 50% 75% 100% A RPM Exp. Est. Err. Exp. Est. Err. Exp. Est. Err. 9.0501E11 1800 200 200 0 367 423 15.2% 678 713 5.1% 9.0501E11 2100 186 167 10.2% 339 298 12% 632 571 9.65% 9.0501E11 2250 211 222 5.2% 453 454 0% 600 546 9% 9.0501E11 2500 186 182 2.15% 358 363 1.4% 635 549 13.5% Load Value 50% 75% 100% A RPM Exp. Est. Err. Exp. Est. Err. Exp. Est. Err. 7.2986E12 1800 200 200 0 367 394 7.35% 678 620 8.56% 7.2986E12 2100 186 160 14% 339 278 18% 632 478 24% 7.2986E12 2250 211 212 0% 453 398 12% 600 465 22% 7.2986E12 2500 186 174 6.45% 358 330 7.8% 635 465 26% 121 The accuracy of the layered model can be improved by considering the accuracy of the emissions equipment. Taking into account that the specified accuracy is +/-2% of the measured value the corrected NO concentrations and prediction errors are shown below in Table 7-3. Table 7-3: Corrected NO Results from Layered Routine (ppm) A graphic depiction of each of these test points starting with 1800 RPM and moving to the right then progressing to the higher RPM values is shown below. Figure 7-1: NO Results from Layered Routine 1 2 3 4 5 6 7 8 9 10 11 12 100 200 300 400 500 600 700 800 T e s t P o i n t N O C o n c e n t r a t io n ( p p m ) E x p e r i m e n t a l P r e d i c t e d Load Value 50% 75% 100% A RPM Exp. Est. Err. Exp. Est. Err. Exp. Est. Err. 9.0501E11 1800 200 200 0% 374 423 13% 691 713 3.1% 9.0501E11 2100 182.3 167 8.38% 332 298 10% 619 571 7.8% 9.0501E11 2250 215 222 3.15% 454 454 0% 588 546 7.14% 9.0501E11 2500 182 182 0% 365 363 .006% 622 549 11.7% 122 7.4 Validation of the EGR Model Upon evaluating the EGR model several aspects of the effects of EGR dilution were noted. The first of these effects is the rise in temperature of the incoming fresh air/EGR mixture. To account for the temperature increase adiabatic mixing is assumed to take place between the re-circulated exhaust gases and the incoming fresh air mixture. The justification of this adiabatic assumption was determined from comparing the adiabatic mixture temperature calculation with measurements taken inside the intake manifold after mixing between fresh air and EGR dilution occurred. The adiabatic approach provides reasonable approximations for the mixture temperature between fresh air and EGR dilution. The second effect of EGR dilution is the rise in the exhaust gas temperature. At low CLV?s the effect of EGR dilution on exhaust gas temperature is negligible but at higher fuel loads the rise in exhaust gas temperature was shown to be noteworthy. For 75% calculated load value it was observed that 10%, 20%, and 30% EGR dilution, caused a corresponding rise in exhaust gas temperature of approximately respectively. The rise in exhaust gas temperature is unexpected since the primary effect of EGR dilution is lowering the adiabatic flame temperature. With the introduction of EGR into the cylinder gaseous species that are highly effective at absorbing radiative heat flux such as are increased, while the oxygen content is displaced. With the displacement of oxygen less efficient combustion takes place which results in an increase of highly radiative soot particles which emit radiative heat fluxes five to ten times that of gaseous species. With an increase of gaseous species such as , more and more of the radiative heat flux created by both the combustion of fuel and soot particles are absorbed by gases in the cylinder causing a corresponding increase in exhaust gas temperature. With this EGR acts as both a heat sink for combustion, effectively lowering the enthalpy of the incoming mixture thus 123 lowering the adiabatic flame temperature, and a radiative heat flux sink, absorbing more of the radiative heat flux from combustion. The rise in exhaust gas temperature can be seen in experimental studies based upon the fact that this rise due to radiative heat transfer is proportional to . Since this rise in EGT temperature is relatively constant a simple empirical correlation was implemented in the EGR model. The last effect that needed special attention is the affect of the residual NO concentration in the cylinder due to the amount of NO residual fraction in the recycled exhaust gas. Better results were obtained by treating the residual NO fraction as an inert species that is non-reacting with other species during the combustion process. With this assumption the residual NO fraction contained in the recycled exhaust gas is added to the overall NO production from the NO/EGR model. The results of the EGR model are shown below in Table 7-4. Table 7-4: EGR Results for NO (ppm) EGR % ~10% ~20% ~30% CLV RPM Exp. Est. Err. Exp. Est. Err. Exp. Est. Err. 50% 1800 186 190 4.3% 135 104 22.9% 73 67 8% 50% 2100 167 167 0% 112 81 27% 90 63 30% 50% 2500 146 172 17.8% 111 83 27% 73 66 9.5% 75% 1800 290 305 5.1% 162 184 13.5% 100 98 2% 75% 2100 242 267 10.33% 162 167 3% 87 79 9% 75% 2500 284 294 3.5% 182 196 7.7% 111 101 9% 124 A graphic depiction of each of these test points starting with 1800 RPM and moving to the right then progressing to the higher RPM values is shown below in Figure 7-2. Figure 7-2: EGR Results for NO It can be seen in Table 7-4 that the EGR model predicts fairly well with only a few results varying considerably from the experimental data. Most of the error in Table 7-4 can be attributed to the approximation of the EGR percentage. As was stated in Chapter 5 the EGR percentage is based upon measuring the air flow at the filter with no EGR and comparing this measurement with that of the air flow at the filter with EGR. The entire EGR percentage calculation is based upon the assumption that the total air flow will remain unchanged from the measurement with 0% EGR. While this technique is valid for the majority of the points some of the experimental data suggests that small changes in the operating conditions occurred due to 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 50 100 150 200 250 300 350 T e s t P o i n t N O C o n c e n t r a t i o n ( p p m ) E x p e r i m e n t a l P r e d i c t e d 125 fluctuations driven by an ECU command. Recall that in this study that no control over the ECU commands is available therefore error can be expected in simulation results due to ECU driven changes in operating characteristics such as injection timing, pulse-width, and pulse-counts. With the exception of the limited availability of ECU control the EGR model developed in this study will serve as an exceptional starting point once total engine control is implemented. An important observation to consider when adding the residual fraction is the transient nature of effects of residual gas or EGR. When a residual fraction or EGR is considered the cycle simulation must be treated as a iterative process, meaning that the output of the simulation will serve as the new input for the initial concentration of residual NO. The process of re-assigning the output as an input continues until the difference between the output and the residual NO fraction input is negligibly small. An explanation of why this is applicable stems from the fact the measurement that is observed with the emissions equipment is actually a steady state value which is the result of a fast transient response to the initial condition. 7.5 Validation of Model After tuning Equation (3.155) the empirical parameters determined for the model are shown below. -1.3E-07 8.4E-04 -1.79 0.00104 -0.3988 45.7736 -5.3E-06 0.00631 -2.0846 The results from executing the cycle simulation are shown below in Figure 7-3. 126 Figure 7-3: Concentration Simulation Results The empirical model shows excellent results with experimental data. The only disadvantage of this approach is the fully empirical tuning constants. An attempt to create a reaction mechanism model for was applied but the results suggested that there are many factors that are important in formation which cannot be represented by a simple mechanism approach. 7.6 Analysis of Lumped model Recalling the discussion in Chapter 3 on lumped modeling, a lumped model is one in which there are no differentiations between the species and NO. The combined concentrations of the species are the objective value of the model. One great disadvantage of 1 2 3 4 5 6 7 8 9 10 11 12 0 20 40 60 80 100 120 140 160 T e s t P o i n t N O 2 C o n c e n t r a t i o n ( p p m ) E x p e r i m e n t a l P r e d i c t e d 127 such a model is that the relative proportions and information regarding the species formation characteristics are insufficient. The same sensitivity and perturbation analysis that was performed for the layered model was also applied to the lumped model. The results indicated that the lumped model is capable of producing results within 15% percent of experimental data. The major difference between the results of the two approaches is that the over-predictions of experimental and predicted data are more common. The temperature profile used for the lumped model is determined from the relationship below. (7.2) Equation (7.2) shows that for the approximated unburned air temperature the lumped model is weighted more towards the ideal unburned air temperature. This temperature profile can explain the conclusions by other researchers of using the ideal unburned air temperature in the calculation of the adiabatic flame temperature. The resulting corrected concentrations for the lumped model are shown below in Table 7-5. Table 7-5: Lumped Model Results (ppm) Load Value (ppm) 50% 75% 100% A RPM Exp. Est. Err. Exp. Est. Err. Exp. Est. Err. 4.143E12 1800 342 342 0% 502 483 3.7% 684 782 14% 4.143E12 2100 328 313 12% 448 380 15% 635 705 11% 4.143E12 2250 309 330 1.3% 529 590 11.5% 603 632 4.19% 4.143E12 2500 279 319 14% 454 434 .2% 638 559 12.3% 128 As can be seen in Table 7-5, the lumped model does predict concentrations fairly well. As compared to the non-lumped model, the lumped model has an increase of average error with respect to experimental data of +3.26%. Again the great disadvantage of the lumped model is the lack of information regarding the relative proportions of and NO. The non-lumped model has proven to be the better choice in regards to experimental error, but if measurements of the relative amounts of and NO were not available the lumped model will suffice to provide good results. 7.7 Validation of and Model Despite the simplicity of the and models an analysis showing their validity will be presented. The results from the and models are shown below in Figure 7-4 and Figure 7-5. To enhance the model some ad hoc adjustments were made to account for a small percentage of combustion inefficiency. The overall results from the model showed very good agreement with experimental data. 129 Figure 7-4: Concentration Simulation Results Figure 7-5: Concentration Simulation Results 1 2 3 4 5 6 7 8 9 10 11 12 5 6 7 8 9 10 11 12 13 14 T e s t P o i n t C O 2 C o n c e n t r a t io n ( M o la r % ) E x p e r i m e n t a l P r e d i c t e d 1 2 3 4 5 6 7 8 9 10 11 12 0 2 4 6 8 10 12 14 T e s t P o i n t O 2 C o n c e n t r a t i o n ( M o l a r % ) E x p e r i m e n t a l P r e d i c t e d 130 7.8 Implementation of Emission Models in Single-Zone CIDI combustion Model Up until this point in the study all of the above analysis has been done in conjunction and tuned with experimental data. All of the inputs to each emission model have been determined from solving the governing thermodynamic equations stated at the beginning of Chapter 3. This section will discuss the results from coupling the emission models with the combustion model briefly described at the beginning of Chapter 3. Six points from this analysis are tested in conjunction with the combustion model. The first of the models to be implemented is the NO model. The results for the NO concentrations with the coupled NO model can be seen in Figure 7-6. The good agreement with experimental data proves that the combustion model accurately predicts the average cylinder temperature and pressure. Figure 7-6: NO Results from coupling NO Model and Combustion Model 1 2 3 4 5 6 50 100 150 200 250 300 350 400 450 N O T e s t P o i n t N O C o n c e n t r a t i o n ( p p m ) E x p e r i m e n t a l P r e d i c t e d 131 The temperature and pressure approximations for point 2 are shown below in Figure 7-7 and Figure 7-8. The good match between predicted and experimental data shows the effectiveness of the combustion model and the relationship between good prediction results with the NO model as well. Figure 7-7: Predicted Pressure and Measured Pressure Figure 7-8: Predicted Temperature and Measured Temperature - 2 0 0 - 1 5 0 - 1 0 0 - 5 0 0 50 100 150 200 0 1 2 3 4 5 6 x 1 0 6 C r a n k A n g l e ( d e g ) P r e s s u r e ( P a ) P R E D I C T E D D A T A - 2 0 0 - 1 5 0 - 1 0 0 - 5 0 0 50 100 150 200 200 400 600 800 1000 1200 1400 1600 1800 C r a n k A n g l e ( d e g ) T e m p e ra tu re ( K ) P R E D I C T E D D A T A 132 The next model to be tested coupled with the combustion model is the EGR model. Good predictions of NO concentrations are a good indication that the EGR model will be accurate as well. A sample of nine EGR test points is shown below in Figure 7-9. Figure 7-9: EGR NO Results from coupling EGR Model and Combustion Model The results of coupling the combustion model and the EGR model seemed to improve the predictions from testing results used in conjunction with experimental data. This improvement in predictions can be attributed to slight variations in the predicted temperature profiles. Due to the fact that the equivalence ratio is treated as an input to the combustion model the predictions of the and models will show the same results depicted in section 7.6. For this reason an analysis with the combustion model will not be presented. The next model to be coupled with the combustion model and tested for validity is the model. Once again with 1 2 3 4 5 6 7 8 9 50 100 150 200 250 300 E G R T e s t P o i n t N O C o n c e n t r a t i o n ( p p m ) E x p e r i m e n t a l P r e d i c t e d 133 great confidence in the predicted NO concentrations the expected trend of the model should display the same results. The results from the six sample points are shown below in Figure 7-10. Figure 7-10: Results from coupling Model and Combustion Model The model predicts fairly well with some over and under predictions with small error. The error in these predictions are directly related to the predicted NO concentrations. Because the model is an empirical one small variations in equation parameters can have a larger effect on the simulation results. 1 2 3 4 5 6 40 60 80 100 120 140 160 N O 2 T e s t P o i n t N O 2 C o n c e n t r a t i o n ( p p m ) E x p e r i m e n t a l P r e d i c t e d 134 Chapter 8: Conclusions and Recommendations Conclusions The primary objective of this thesis was to develop an emissions model through an extensive and thorough investigation. The emission model not only had to be accurate with low empirical input but also easily incorporated into the current combustion model. Due to the limitations on current emissions equipment the number of emission species modeled was limited to include emissions. Throughout the model development every aspect involved in the modeling process was tested and analyzed. The main emission model development also includes several sub-models to replicate the effects of changing operating characteristics and modeling approaches. It was shown that the layered adiabatic core approach is valid for diesel combustion only if coupled with the appropriate approximation of the local unburned air temperature. The adiabatic layered approach coupled with the simple Extended Zeldovich Mechanism alone is shown to provide concentrations predictions within 13% of experimental data. The developed layered model requires only the adjustment of the parameter A in the Arrhenius expression for the entire effective RPM range. The value for this tuning constant used in this study was set at 9.0501E11. Some supplementary conclusions developed during the perturbation analysis are presented below. Fully Mixed Adiabatic Core models can present large error in predictions due to their lack of representing local burned gas elements accurately. Fully mixed models also require constant adjustment of the rate limiting tuning constant A. Adding additional reactions to the EZM, including the majority of well known contributors, affected the model output only by about 1.5%. The unburned gas temperature calculated from a reversible process is ineffective at representing local unburned gas temperatures. 135 A gas temperature between the average mixture temperature and the unburned gas temperature represents the local unburned gas temperatures more efficiently. The QSSA assumption is a valid approximation for the EZM, with differences from the kinetically controlled approach of only +/-5ppm. Small variations between rate constants do not affect model output by a considerable amount. The most accurate rate constants reside close to those presented in Heywood [9]. All injections must be accounted for in modeling. Small burning elements produced from pre-pilot or pilot injections can lead to a significant amount of concentrations due to the fact they are compressed to higher temperatures. The effect of adding an initial NO concentration to the cycle simulation affected the output of the simulation by a small percent. The effect of adding an initial concentration is cumbersome due to the fact that this will only lower the rate limiting tuning constant, A, and thus the same results will be achieved. It is critical to develop separate models for both NO and , due to the differences in the formation characteristics of these species. Lumped models can perform to show errors within 15% using the methods described in this thesis, but the more accurate model is one with separate treatment of each individual species. Lumped models also lose the characteristics of individual species formation. EGR increases exhaust gas temperatures significantly at higher load values. With the introduction of EGR into the cylinder gaseous species that are highly effective at absorbing radiative heat flux such as are increased, while the oxygen content is displaced. With the displacement of oxygen less efficient combustion takes place which results in an increase of highly radiative soot particles which emit radiative heat fluxes five to ten times that of gaseous species. With an increase of gaseous species such as , more and more of the radiative heat flux 136 created by both the combustion of fuel and soot particles are absorbed by gases in the cylinder causing a corresponding increase in exhaust gas temperature. With this EGR acts as both a heat sink for combustion, effectively lowering the enthalpy of the incoming mixture thus lowering the adiabatic flame temperature, and a radiative heat flux sink, absorbing more of the radiative heat flux from combustion. Recommendations One of the most critical aspects in engine testing is the control of operating parameters. If inadequate control of operational parameters is present, a great amount of error can be introduced when constructing models. In this study the primary concern when constructing the engine models was the limited control of some operational parameters. The critical parameters that were in question during this study were the control over the amount of fuel burned and the injection scheduling. Although relative close approximations of burned mass and SOC points can be produced using the pressure derived heat release and empirical correlations, there will still be some error introduced due to the critical importance of these parameters. The current test stand and the corresponding engine models can be improved if this test stand were equipped with such components as an injector needle lift sensor and an injector driver. With these components the effects of these operating parameters can be analyzed more accurately, resulting in more accurately designed engine models. It is of the utmost importance that these critical operating parameters be controllable. 137 References 1. Andersson, M., Johansson, B., Hultqvist, A., and N?hre, C., "A Real Time NOx Model for Conventional and Partially Premixed Diesel Combustion," SAE Technical Paper 2006-01-0195, 2006. 2. Bagal, N., Rutland, C., Foster, D., Narayanaswamy, K., and He, Y., "CO Emission Model for an Integrated Diesel Engine, Emissions, and Exhaust Aftertreatment System Level Model," SAE Technical Paper 2009-01-1511, 2009. 3. Baratta, M., Catania, A., Spessa, E., and Vassallo, A., "Development and Assessment of a Multizone Combustion Simulation Code for SI Engines Based on a Novel Fractal Model," SAE Technical Paper 2006-01-0048, 2006. 4. Costa, M., Merola, S., and Vaglieco, B., "Modelling of soot formation in diesel engines exploiting measurements of soot volume fraction and diameter," SAE Technical Paper 2001-24-0011, 2001. 5. Dodge, L., Kubesh, J., Naegeli, D., and Campbell, P., "Modeling NOx Emissions from Lean-Burn Natural Gas Engines," SAE Technical Paper 981389, 1998. 6. Duffy, K., and Mellor, A., "Further Developments on a Characteristic Time Model for NOx Emissions from Diesel Engines," SAE Technical Paper 982460, 1998. 7. Ericson, C., Westerberg, B., Andersson, and M.Egnell, R., ?Modelling Diesel Engine Combustion and NOx Formation for Model Based Control and Simulation of Engine and Exhaust Aftertreatment Systems? SAE Technical paper 2006-01-0687, 2006. 138 8. Ferguson, C.R., Internal Combustion Enignes, 1986, New York, NY: John Wiley and Sons, Inc. 9. Heywood, J.B., Internal Combustion Enigne Fundamentals, 1988, New York, NY: McGraw-Hill, Inc. 10. Jay, S., B?ard, P., and Pires da Cruz, A., "Modeling Coupled Processes of CO and Soot Formation and Oxidation for Conventional and HCCI Diesel Combustion," SAE Technical Paper 2007-01-0162, 2007. 11. Kim, H., and Sung, N., "Combustion and Emission Modeling for a Direct Injection Diesel Engine," SAE Technical Paper 2004-01-0104, 2004. 12. Majewski, W.A., and Khair, M.K., Diesel Emissions and Their Control, 2006, Warrendale, PA: SAE International. 13. Mattavi, J.N., and Amann, C.A., Combustion Modeling in Reciprocating Engines, 1980, New York, NY: Plenum Press. 14. Mellor, A., Mello, J., Duffy, K., Easley, W., and Faulkner, J., "Skeletal Mechanism for NOx Chemistry in Diesel Engines," SAE Technical Paper 981450, 1998. 15. Merrion, D., "Heavy Duty Diesel Emission Regulations - Past, Present, and Future," SAE Technical Paper 2003-01-0040, 2003. 16. Miller, R., Davis, G., Lavoie, G., Newman, C., and Gardner, T., "A Super-Extended Zel'dovich Mechanism for Nox Modeling and Engine Calibration," SAE Technical Paper 980781, 1998. 17. Miller, R., Russ, S., Weaver, C., Kaiser, E., Newman, C., Davis, G., and Lavoie, G., "Comparison of Analytically and Experimentally Obtained Residual Fractions and NOX Emissions in Spark-Ignited Engines," SAE Technical Paper 982562, 1998. 139 18. Olikara, Cherian and Borman, Gary L., ? A Computer Program for Calculating Properties of Equilibrium Combustion Products with Some Applications to I.C. Engines,? SAE Technical Paper 750468, 1975. 19. Papagiannakis, R., Hountalas, D., and Kotsiopoulos, P., "Experimental and Theoretical Analysis of the Combustion and Pollutants Formation Mechanisms in Dual Fuel DI Diesel Engines," SAE Technical Paper 2005-01-1726, 2005. 20. Payri, F., Arr?gle, J., L?pez, J., and Mochol?, E., "Diesel NOx Modeling with a Reduction Mechanism for the Initial NOx Coming from EGR or Re-entrained Burned Gases," SAE Technical Paper 2008-01-1188, 2008. 21. Perini, F., Mattarelli, E., and Paltrinieri, F., "Development and Validation of Predictive Emissions Schemes for Quasi-Dimensional Combustion Models," SAE Technical Paper 2010-01-0148, 2010. 22. Ramos, J.I., Intenal Combustion Engine Modeling, 1989, New York, NY: Hemisphere Publishing Corporation. 23. Rublewski, M., and Heywood, J., "Modeling NO Formation in Spark Ignition Engines with a Layered Adiabatic Core and Combustion Inefficiency Routine," SAE Technical Paper 2001-01-1011, 2001. 24. Stiesch, G. and Merker, G., "A Phenomenological Model for Accurate and Time Efficient Prediction of Heat Release and Exhaust Emissions in Direct-Injection Diesel Engines," SAE Technical Paper 1999-01-1535, 1999. 25. Tan, P., Deng, k., Lu, J., Lou, D., and Wan, G., "A Composition-Based Model for Particulate Matter Emission of Direct Injection Diesel Engines," SAE Technical Paper 2005-01-3463, 2005. 140 26. Tan, P. and Lu, J., "Particulate Matter Emission Modeling of Diesel Engines," SAE Technical Paper 2003-01-1904, 2003. 27. Turns, S.R., An Introduction to Combustion, 2000, New York, NY: Mcgraw-Hill Companies, Inc. 28. Val?rio, M., Raggi, K., and Sodr?, J., "Model for Kinetic Formation of CO Emissions in Internal Combustion Engines," SAE Technical Paper 2003-01-3138, 2003. 29. Westlund, A., and ?ngstr?m, H., "Fast Physical Prediction of NO and Soot in Diesel Engines," SAE Technical Paper 2009-01-1121, 2009. 30. Morgan, W.K.C., Reger, R.B., and Tucker, D.M., (1997) ?Health Effects of Diesel Emissions.? Annals of Occupational Hygiene 41, pp.643-658. 31. Lloyd, A.C., and Cackette, T.A. 2001 Critical Review-Diesel Enignes: Enviromental Impact and Control. EM 2001, 8 (6), pp.34-41. 32. Kagawa, J., ?Health effects of diesel exhaust emissions-a mixture of air pollutants of worldwide concern,? Toxicology, 181-182(2002), pp.349-353. 33. McClellan, R.O., ?Health effects of exposures to diesel exhaust particles,? Annu Rev Pharmacol Toxicol 27: pp.279-300(1987). 34. Buttsworth, David R. (2002) Spark ignition internal combustion engine modeling using Matlab. Technical Report. University of Southern Queensland, Faculty of Engineering and Surveying, Toowoomba, Australia. 35. Miller, J. A., and Bowman, C. T., ?Mechanism and Modeling of Nitrogen Chemistry in Combustion,? Progress in Energy and Combustion Science, 15: pp.287-338 (1989). 141 36. Fenimore, C. P., ?Formation of Nitric Oxide in Premixed Hydrocarbon Flames,? Thirteenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, PA, pp373-380, 1970. 37. Law, C.K., Combustion Physics, 2006, New York, NY: Cambridge University Press. 38. Coda Zabetta, E., Kilpinen,P., Hupa, M., Stahl, K., Leppalahti, J., Cannon, M., Nieminen, J., ?Kinetic Modeling Study on the Potential of Staged Combustion in Gas Turbines for the Reduction of Nitrogen Oxide Emissions from Biomass IGCC Plants,? Energy and Fuels 2000, 14, pp.751-761. 39. Zabetta, E. C., Kilpinen, P., ?Improved NOx Submodel for In-Cylinder CFD Simulation of Low-and Medium-Speed Compression Ignition Engines,? Energy and Fuels, 2001, 15, pp.1425-1433. 40. Goldsworthy, L., ?Reduced Kinetics Schemes for Oxides of Nitrogen Emissions from a Slow-Speed Marine Diesel Engine,? Energy and Fuels, 2002, 17, pp.450-456. 41. Glarborg, P., ?A reduced mechanism for nitrogen chemistry in methane combustion,? Twenty-Fourth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, PA, pp.889-898, 1992. 42. Wingo, T. J., Construction of a Diesel Engine Test Stand and a Crank Angle Based Heat Release Model. Ms Thesis., Auburn University. 43. Gordon, S., McBride J.B., ?Computer Program for Calculation of Complex Chemical Equilibrium Compositions and Applications?, NASA Reference Publication 1311,1994. 44. Zeldovich, Y.B. (1946), Acta Physiochem, URSS, 21, 577. 45. Kim, Chung-Gong. A Crank Angle Resolved CIDI Engine Combustion Model with Arbitrary Fuel Injection for Control Purpose. PhD diss., Ohio State University. 142 Appendix A BIOLOGICAL IMPACT OF DIESEL EMISSION COMPONENTS Emission Component Atmospheric Reaction Products Biological Impact Gas Phase Carbon Monoxide --- Highly toxic to humans; blocks oxygen uptake Nitrogen Oxides Nitric acid, ozone Nitrogen dioxide is a respiratory track irritant and major ozone precursor. Nitric acid contributes to acid rain. Sulfur Dioxide Sulfuric Acid Respiratory tract irritation. Contributor to acid rain. Carbon Dioxide --- Major Contributor to global warming Saturated Hydrocarbons Aldehydes, alkyl nitrates, ketones Respiratory tract irritation. Reaction products are ozone precursors ( in the presence of Nox) Unsaturated Hydrocarbons Aldehydes, ketones Respiratory tract irritation. Some alkenes are mutagenic and carcinogenic. Reaction products are ozone precursors (in the presence of Nox). Formaldehyde Carbon Monoxide, Hydroperoxyl radicals Formaldehyde is a probable human carcinogen and an ozone precursor (in the presence of Nox). Higher Aldehydes Peroxyacyl nitrates Respiratory tract and eye irritation; causes plant damage. Monocyclic Aromatic Compounds Hydroxylated and Hydroxylated- nitro derivatives Benzene is toxic and carcinogenic in humans. mutagenic in bacteria. PAHs nitro-PAHs Some of these PAHs and nitro- PAHs are known as mutagens and carcinogens 143 Nitro-PAHs Quinones and hydroxylated-nitro derivatives Some reaction products are mutagenic in bacteria Particulate Phase Elemental Carbon --- Nuclei adsorb organic compounds; size permits transport deep into the lungs. Inorganic Sulfates --- Respiratory tract irritation. Aliphatic Hydrocarbons (C14- C35) Little information; possibly aldehydes, ketones, and alkyl nitrates. Unkown PAHs (4 rings and more) (e.g., pyrene, benzo(a)pyrene) Nitro-PAHs(4 rings or more); nitro-PAH lactones Larger PAHs are major contributors of carcinogens in combustion engines. Many nitro-PAHs are potent mutagens and carcinogens Nitro-PAHs(3 rings and more) (e.g., nitropyrenes) Hydroxylated-nitro derivatives Many nitro-PAHs are potent mutagens and carcinogens. Some reaction products are mutagenic in bacteria(Ames assay). Health Effects of DE [12] 144 Appendix B Tabulation of published reaction for the H-O-N reaction system used in reaction mechanisms. Rate coefficients of the form . Units are moles, cubic centimeters, seconds, Kelvin. Reaction A b 9.03E11 0 38031 6.4E9 1 3160 4.6E13 0 0 6.9E13 0 13400 1E14 0 14100 7.58E13 0 7600 6.9E23 -2.5 32170 1E13 0 10000 1E14 0 25000 ? 6E14 -1.5 4985 ? 1E13 0 300 ? 3.467E14 0 740 ? 1.9E12 0 13500 ? 5E12 0 0 ? 2E12 0 -240 ? 2.3E14 0 18000 ? 2.730E13 0 53800 145 ? 2.7E14 -1.0 63140 ? 3.7E14 -0.5 40280 ? 6E14 -1.5 8000 ? 5.9E10 0.5 21550 ? 1.29E14 0 2165 ? 1.32E15 0 11230 ? 3.59E15 0 18430 ? 1.6E12 0.56 755 ? 5E11 0.5 0 ? 6.44E11 0 1460 ? 4.38E12 0 6546 ? 1.9E16 0 25179 ? 2.1E12 0.5 0.0 ? 1.4E13 0 1510 ? 5E11 0.5 0.0 ? 2.24E13 0 10600 ? 1.995E12 0 3000 ? 1E11 0.5 2000 ? 1E13 0 1000 ? 1E13 0 1000 ? 1.9E11 0 1000 ? 3.16E11 0 1000 ? 6.14E12 0 17100 146 ? 1E13 0 1000 ? 2E10 0.5 2230 ? 6.310E11 0.5 0 ? 3.16E21 -2.0 42000