- Email: [email protected]
air mixtures: two dimensional detailed modeling and laser based diagnostics
Spark Ignited Hydrogen/Air Mixtures: Two Dimensional Detailed Modeling and Laser Based Diagnostics M. THIELE and J. WARNATZ
Interdisziplina ¨res Zentrum fu ¨r Wissenschaftliches Rechnen, Universita ¨t Heidelberg, Im Neuenheimer Feld 368, D-69120 Heidelberg, Germany
A. DREIZLER*
Energie- und Kraftwerkstechnik, Tu. Darmstadt, Petersenstrae 30, 64287 Darmstadt, Germany
S. LINDENMAIER, R. SCHIEßL and U. MAAS
Institut fu ¨r Technische Verbrennung, Universita ¨t Stuttgart, Pfaffenwaldring 12, D-70569 Stuttgart, Germany
and A. GRANT and P. EWART
Clarendon Laboratory, Oxford University, Parks Road, Oxford OX1 3PU, UK This study reports a detailed 2-D model which describes the spark ignition of an initially quiescent hydrogen/air mixture. The model includes the compressible Navier-Stokes equations, detailed chemistry and molecular transport in the gas phase as well as heat conduction to the electrodes. The spark is modeled for the phases subsequent to breakdown using the Maxwell equations for quasi-stationary conditions for the electric field. Initial and boundary conditions necessary for the simulations are chosen in accordance with experimental values. Heat conduction to the electrodes and different electrode shapes are investigated. The influence of these parameters on the shape of the initial flame kernel is discussed and compared qualitatively to experimental results. Spark ignition experiments are performed using a highly reproducible ignition system. Shapes of the early flame kernels are monitored by 2-D laser-induced fluorescence (PLIF) imaging of OH radicals produced during the ignition and the combustion process. The investigations are performed for different equivalence ratios. In addition, for a central position within the flame kernel, temperatures are measured at different times after ignition using vibrational coherent anti-Stokes Raman spectroscopy (CARS) of nitrogen. © 2002 by The Combustion Institute
INTRODUCTION Spark ignition processes are of enormous practical importance. In technical combustion devices, for example, engines or gas turbines [1, 2], spark ignition initiates combustion and influences the temporal development of this process in a critical way. Safety concepts in atomic power plants rely on controlled burning of hydrogen resigned from the reactor [3]. Therefore a deeper understanding of spark ignition is required. In future, a detailed theoretical understanding may improve current concepts of spark plugs in connection with improved electrical ignition devices. To further optimize such models a comparison with experimental data is still necessary. *Corresponding author. E-mail: [email protected] 0010-2180/02/$–see front matter PII 0010-2180(01)00333-9
Spark ignition constitutes a very complex interplay between plasma kinetics, chemical kinetics, molecular transport processes and fluid dynamics. A complete mathematical simulation of a spark ignition is a difficult task because of the enormous numerical problems, due to the stiffness and high dimensionality of the problem (each chemical species introduces an additional conservation equation). Furthermore, experimental investigations of spark kernels and their transition to flame kernels are rendered difficult because of very short process times, extremely high core temperatures and large gradients in the refractive index. In general, only flame propagation subsequent to spark ignition can be studied by laser diagnostic methods. Therefore, internal structures of the plasma core remain mostly unknown. In the literature various simplified models of spark ignition have been presented. The followCOMBUSTION AND FLAME 128:74 – 87 (2002) © 2002 by The Combustion Institute Published by Elsevier Science Inc.
SPARK IGNITED HYDROGEN/AIR MIXTURES ing examples give a brief survey of the field, although this is far from being a comprehensive list. In a 2-D setting the effect of the electrical spark on the flow field has been studied for an air atmosphere [4] and for one step propane chemistry [5]. Heat conduction to the electrodes was studied by Pischinger et al. [6]. A detailed reaction mechanism for a methane-air mixture has been employed by Sher et al. [7, 8] employing rectangular electrodes. Detailed chemistry in conjunction with detailed molecular transport has been used in a one-dimensional geometry, neglecting heat conduction to the electrodes [9]. None of the above studies combine all aspects in one model. In this study an improved model is discussed which considers heat conduction from the gas phase to the electrodes, detailed chemistry and molecular transport as well as the coupling of the gas dynamics to the properties of the electrical discharge through heating by the electrical current. The computational model is reduced to two dimensions by assuming rotational symmetry along the axis of the electrodes. Different scenarios are considered for a lean hydrogenoxygen mixture. In general, the model is used to study the transition from the initial phase (governed by the shock wave formation) to the flame propagation. The influence of the electrode as well as the influence of heat conduction from the hot gas to the electrode are investigated. A qualitative comparison is made with experimental data. For this comparison initial and boundary conditions necessary for the numerical simulation are taken in accordance with the experimental settings. Optical diagnostic methods have become a standard tool for the investigation of chemically reacting flows. Traditional techniques like Schlieren [10, 11] photography and flame emission [12] filming are line-of-sight integrated. Therefore, two-dimensional laser-based techniques were developed to investigate the internal structures of combustion processes [13]. Laser-induced fluorescence (LIF) counts among the most popular techniques because of its experimental simplicity. The spatial distribution of minority species can be detected using commercially available equipment [14, 15]. The extraction of quantitative concentration measurements from LIF signals is often rendered
75
difficult because of collisionally induced energy transfer [16, 17]. This problem can be overcome by application of various strategies, as shown for example in [15, 16, 18]. However, in many cases, qualitative information is sufficient. Two-dimensional LIF of OH radicals is an ideal method for the comparison of early sparkignited flame kernels with corresponding numerical simulations. For this reason, the OH distribution in a plane parallel to the spark has been investigated. Coherent anti-Stokes Raman scattering (CARS) is a well-established laser-based spectroscopic technique for the measurement of species concentration and temperature. CARS has proved its value in a large variety of combustion devices, from furnaces [19] to internal combustion engines [20]. The CARS signal generation process and its dependence on physical properties, such as temperature and pressure, is complex and is discussed in detail in [16]. For time-varying combustion processes such as spark ignition, single shot multiplex CARS [21] measurements are necessary. In this study, temperatures in the center of the electrode gap have been determined 5 and 10 ms after the breakdown using vibrational nitrogen CARS. For all experimental studies a well characterized and highly reproducible ignition system presented in a previous study [22] has been used. In a parametric study, gas mixture, spark parameters, and electrode geometry have all been varied.
MODELING PRINCIPLES Equations The model for the simulation of the early development of a flame kernel initiated by an electrical spark is based on the conservation equations for reacting flows, including detailed chemistry. To make the problem treatable the following assumptions are made: ● ● ●
The plasma channel formed after the breakdown is cylindrical. Local thermal equilibrium exists The influence of the magnetic field is negligible.
76
M. THIELE ET AL. 1 1 ⭸T ⫹ vgradT ⫹ div j q ⫹ pញ gradv ⭸t cv cv ⫹
1 cv
冘 冋 u 共M ⫺ div j 兲册 ⫽ 1c 共q ⫹ q 兲 ns i
i
i
i
i
v
s
r
(1) a source term qs is introduced to include the properties of the electric spark. Combining Ohms Law for the electric current flow (jel) with the electrical power input leads to j el ⫽ 兩E兩 2.
Fig. 1. Illustration of the computational area. The conducting plasma channel is assumed to be cylindrical. For details concerning initial and boundary conditions see text.
● ●
The plasma is quasi-neutral. The plasma is transparent to its own radiation
In modeling spark ignition, these assumptions are commonly made [4, 7, 8, 23]. The first assumption allows the geometry to be reduced to two dimensions because of rotational symmetry around the plasma channel as shown in Fig. 1. The second assumption is based on the fact that internal relaxation processes in molecules are faster than the time scales of the investigated processes. The last two assumptions are needed to derive the equation for the electric potential [4, 24]. Temperature- and pressure-dependent transport coefficients for the mixture can be calculated from binary data using mixing rules. The gas phase equations consist of the compressible Navier-Stokes equations (conservation of momentum and total mass) together with energy conservation. For each of the chemical species in the reaction mechanism (9 species and 38 elementary reactions [25], includes dissociation), a mass conservation equation is added to the system, which is closed by the ideal gas law. All equations are solved in their primitive variable formulation as presented in [9]. In the energy equation [26]
(2)
Assuming a cylindrical channel, the radial field component can be neglected because the current I flows along the z-direction only [4, 23]. The source term can then be calculated from the current strength by integrating over the plasma. As soon as the cylindrical shape gets distorted this assumption is no longer valid. From the conservation of charge, Ohms Law and Maxwell’s equations of electrostatics, a time-dependent differential equation for the ជ ⫽⫺grad⌽ results. scalar electric potential E Quasistationarity leads to div(⫺grad⌽) ⫽ 0
(3)
because of the small relaxation times compared with the time-scales of the reacting flow. This equation is added to the system and via ⌽ one ជ and the electric source qs. can compute E The electrical conductivity is derived from temperature and pressure by interpolation of the data of Yos [27] for ionized air. Notice, that the influence of the ions on the chemical mechanism with regard to flame propagation is negligible and is not accounted for in the present approach. Initial and Boundary Conditions Initial and boundary conditions have to be specified for the solution. In this work, we restrict ourselves to simulations for initially quiescent mixtures. The model is, however, able to treat initial velocity fields as well. The properties of the channel at the end of the breakdown are used as initial temperature profile and are taken from [28]. Species profiles have been calculated according to the spark characteris-
SPARK IGNITED HYDROGEN/AIR MIXTURES tics, that is, the initial temperature and pressure profile is used for a homogenous computation [29]. After some computation time, the breakdown time, species profiles according to these temperatures and pressures result. The axis of symmetry (r ⫽ 0), as well as the equipotential line between the electrodes (z ⫽ Z), are treated as symmetry lines. The other two boundaries are treated as walls (compare Fig. 1). For symmetry lines, the normal velocity component is set to zero as are all other gas dynamics variables. For wall boundaries no slip conditions for the velocity and vanishing normal gradients of all other gas dynamics variables are specified. For the equation for the electric potential vanishing normal gradients are set on all boundaries except for the equipotential line, where ⌽ ⫽ 0 is specified. At the electrode vanishing gradients for the species mass fractions are assumed. Adiabatic conditions as well as heat conduction have been implemented for the energy equation via:
⭸T ⫺ T E兲 ⫽ 0, ⫹ h共T ⭸nជ
再
冉冑
⫽ ⫺ a 2 ⌬共r,z兲 2
1 共⌬:⌬兲 2
冊冎
.
(6)
Herein ⌬ denotes the stress tensor, a is a damping factor of the artificial viscosity, “:” denotes the two-fold contraction and ⌬(r,z) accounts for the local mesh size in both spatial directions r and z. Because of the stiffness introduced by the chemical kinetics, we used an implicit solution method. This signifies in this context, that all equations except the potential equation are solved implicitly. The electrical source term is computed from the conductivity and the electrical field of the preceding time step. The potential equation is solved after each time step using a Newtonian scheme. The resulting ordinary differential equations are solved by the extrapolation solver LIMEX [33], including step size and order control. The regular discretization scheme leads to block-structured Jacobi- and iteration matrices which can be solved by block-ILU decomposition, reducing computation time and storage requirements [29, 34].
(4)
where nជ denotes the outwards directed normal, h the heat exchange coefficient and the heat conductivity. Setting h ⫽ 0, ⫽ 1 ones obtain adiabatic conditions while using both, and h, mixed conditions for heat conduction result. The current input is treated by setting ⫺ I共t兲 ⭸⌽ ⫽ 2 ⭸z r elec
77
(5)
where relec denotes the radius of the electrodes. Numerical Solution For the numerical solution the method of lines is employed. Spatial discretization is performed by finite differences [26]. Artificial viscosity is added to resolve the pressure wave emitted from the hot channel. An artificial viscosity coefficient proportional to the flux of momentum is added to the viscosity coefficient following [30] and applying the equations of Ostwaldde-Waele [31]. It is set to [32]
EXPERIMENTAL SETUP Ignition Vessel A spark ignition system has been designed at the Institut fu ¨r Technische Verbrennung for the generation of highly reproducible sparks. This allows comparison between experimental results and the predictions of the model for various gas mixtures and flow conditions. The reproducibility of the system can be proved by measurements of flame propagation. For this purpose, planar laser-induced fluorescence measurements of the spatial distribution of OH radicals have been performed in a previous study [22]. By imaging the temporal evolution of several individual ignition events using a novel high speed LIF image acquisition system [35] variations of the spark characteristics with respect to the flame initiation have been investigated. The cycle-to-cycle variations in flame propagation for the employed ignition system were found to be less than 2%. Details concerning the test of the reproducibility as well as the
78
M. THIELE ET AL.
Fig. 2. Sketch of the ignition constant volume vessel. The device allows full optical access. Electrodes are made from 1-mm thick tungsten wire. In its standard configuration the tips of the electrodes are sharpened (conical shape) as shown in the insert. For PLIF applications a plane approximately 0.55 mm apart from the electrodes is used. CARS measurements are performed at the center of the electrode gap.
design of the ignition system can be found in [22]. The ignition system is attached to a constant volume ignition chamber (cylindrical geometry, volume 2.8 l, diameter 160 mm) with full optical access (Fig. 2). The ignition is initiated in the geometric center of the vessel. The diameter of the tungsten electrodes was 1 mm, and the spark gap was also set to 1 mm. The electrode geometry generally used is shown in Fig. 2. Typical current and voltage traces of the spark used during this parameter study are presented in Fig. 3. The spark consists of a short breakdown phase (150 ns duration) and a long user-selectable arc phase. In the present study this was set to 55 and 100 s, respectively. The electrical
energy of this spark was adjusted to 8.7 mJ well above minimum ignition energy [36]. Because of imperfect conversion thermal energy input to the gas phase is estimated to be 4.5 mJ based on calorimetric measurements conducted for similar spark conditions [22]. The temporal evolution of the pressure inside the vessel was monitored with a piezo-electric pressure transducer. Following a standardized purging procedure, the cell was filled with a lean mixture of hydrogen/air (5, 15, or 20% hydrogen in air). The partial pressures of the gases were monitored with a calibrated baratron (MKS 122) and could be adjusted to an accuracy of approximately 0.5 mbar. A pressure of 1 bar before ignition was used for all the results presented below. The time between ignition events was sufficiently long to ensure the cylinder walls had equilibrated to 305 K. OH Planar Laser-Induced Fluorescence (PLIF) Measurements and Data Reduction
Fig. 3. Current and voltage traces of a typical spark event. The peaks at the leading edges of the traces are because of the breakdown which can not be temporally resolved by the probes.
OH radicals produced during the combustion of the hydrogen/air mixture were probed using the Q1(8) line in the (1– 0)-band of the A2⌺⫹-X2⌸ transition. The excitation wavelength of around 283 nm was generated by a dye laser (Rhodamine 590) which was pumped with a XeCl excimer laser (Lambda Physik). After type I frequency doubling in a BBO crystal pulse energies of 0.2 mJ were obtained. The laser
SPARK IGNITED HYDROGEN/AIR MIXTURES beam was formed into a sheet (height 4 mm, thickness 50 m as measured by use of the “razor blade” method) intersecting the flame kernel in a vertical plane displaced approximately 0.55 mm from the electrodes as shown in Fig. 2. Fluorescence light at around 310 nm was detected perpendicular to the laser propagation direction using appropriate dielectric filters and an intensified CCD camera (Princeton Instruments, 1024 ⫻ 256 pixel, 16 bit). The ignition system, the laser and the camera system were triggered by two pulse/delay generators (Stanford Research Systems DG-535). The laser pulse and the exposure of the camera system were delayed relative to the beginning of the spark. The following data reduction procedure was used to extract the flame radius from the OH PLIF images: 1.
The 16 bit gray-scale images were converted to binary images using a simple thresholding technique. This is possible because of a very high signal-to-noise ratio of the original intensity images. The center of the flame and the approximate position of the flame front were extracted from these binary images. 2. The original gray-scale images were filtered by a 2-D Gaussian-filter to reduce noise. This linear filter preserves the position of steep gradients in the image. 3. A one-dimensional intensity profile was produced along a horizontal line through the center of the flame (established in Step 1) through the flame front. 4. The exact position of the flame front is located within a predefined area in the vicinity of the rough flame front obtained from the binary image (Step 1). The flame front position is defined as the absolute maximum of the respective gradients which lie within this predefined area. 5. The horizontal flame diameter was taken as the distance between the opposite flame front positions. Ensemble averaged horizontal flame diameters and standard deviations were derived from a number of single-shot events performed for identical boundary conditions.
79
Vibrational CARS Measurements A single mode, frequency-doubled Nd:YAG laser (Continuum PowerLite 8000) is used to provide the pump beams for the CARS process; part of the output is split off to pump the broadband Stokes laser (Mode-x ML-2). This laser is based on the modeless laser [37], which has been shown to provide the most precise single-shot CARS measurements when used in conjunction with a single-mode pump laser [38]. The central wavelength of the Stokes laser is tuned to around 607 nm; the spectral width is sufficient to simultaneously access the cold band and the first hot band of the N2 vibrational CARS spectrum, allowing multiplex measurements to be made. A delay line is used to ensure temporal overlap of the pump and Stokes beams. Beam splitter BS1 (compare Fig. 4) produces two pump beams disposed vertically above each other. The energy of the pump beams is adjusted by rotation of beam splitter BS2 which transmitts a variable fraction of the incident energy depending on the angle of incidence. The Stokes beam energy is set by selecting a given area of the beam using aperture A. The glass block GB in the Stokes beam is rotated to adjust the spatial overlap of the beams and hence to optimize the CARS signal. The Stokes and pump beams display very similar collimation and beam size: no significant signal increase is attained by adjusting the collimation and diameter of the Stokes beam. The HeliumNeon laser is used to provide a mock signal beam for alignment purposes. The beams are crossed by a lens of focal length 250 mm into the forward folded boxcars geometry [39]; the resulting interaction region is approximately 2 mm long and 100 m wide. Care has been taken to ensure that the interaction region is situated in the center of the spark gap. The ignition cell is filled with a mixture of 20% hydrogen and 80% air and then ignited. CARS measurements are timed relative to the breakdown by the use of programmable pulse generators (Stanford Research Systems DG535). The generated CARS signal is collimated, filtered by suitable interference filters and then focused onto the entrance slit of a 1 m Czerny-
80
M. THIELE ET AL.
Fig. 4. Experimental set-up of the CARS spectrometer.
Turner type spectrometer (Hilger Monospek 1000). The spectra are recorded on a CCD camera (Princeton Instruments 512TKB). Camera operation and data acquisition are controlled by a PC. Because of the high reproducibility of the ignition cell, the single shot CARS spectra were averaged in software. Background levels are subtracted and the CARS spectra are then normalized using the non-resonant spectrum generated in propane. Non-resonant spectra are recorded in situ by filling the ignition cell with pure propane. This normalization takes account of the Stokes laser spectral profile and spatial non-uniformities in the camera CCD chip. For CARS data analysis the Sandia Laboratory code CARSFIT [40] is used to find a best-fit temperature to the averaged spectra. RESULTS AND DISCUSSION Experimental Observations As standard conditions the parameters of the arc phase are set to a duration of 100 s and a current of 1.5 A corresponding to an electrical energy of 8.7 mJ as calculated from current and voltage traces. Figure 5 shows qualitatively the
evolution of the spark kernel via relative OH concentration fields in a plane shifted 0.55 mm with respect to the electrodes. Correction of the laser intensity profile has not been performed. For an ignitable mixture of 20% H2 in air (left column of Fig. 5) the flame develops from an initially oval kernel. Note, that for process times larger than 250 s in vertical direction only the central part of the kernel is monitored because of the limited height of the laser light sheet. However, for a mixture of 5% H2 in air (below the ignition limit) a substantial amount of OH is produced as well. The sequence on the right column of Fig. 5 shows a “dying” plasma kernel where no transition to a self sustaining flame propagation is observed because of the lean fuel/air mixture. For this case remarkable OH concentrations can be monitored up to 1 ms after the end of the spark. By using these raw data, OH profiles perpendicular to the axis of the electrodes are extracted using the image processing procedure described above. Figure 6 shows qualitatively radial OH concentrations profiles from 50 to 1,300 s. The influence of high core temperatures of approximately 5,000 K [41] causes high peak OH concentrations (50 s after break-
SPARK IGNITED HYDROGEN/AIR MIXTURES
81
Fig. 5. Relative OH concentration profiles recorded by PLIF. These sequences show the temporal evolution of a flame kernel for an ignitable mixture of 20% H2 in air and a mixture of 5% H2 in air (below the ignition limit). Notice, that the beam profile has not been corrected for inhomogeneous intensity distribution. For process times later than 250 s the flame kernel exceeds in vertical direction the height of the laser light sheet.
down) which decrease with time when the spark is switched off. 80 s after breakdown the generation and evolution of a flame front can already be observed propagating radially out-
wards. At this time the flame kernel starts to separate from the plasma kernel. Subsequently, the flame front propagates into the fresh gases while in the vicinity of the elec-
82
Fig. 6. Radial OH profiles obtained in a direction perpendicular to the pointing of the electrodes (radial direction). Intensities are given in relative units; spatial variation of local quenching conditions was not taken into account. Experimental conditions were: 20% hydrogen in air, sharpened (conical) electrode tips, 100 s long arc phase and 1.5 A current during arc phase. High OH concentration levels are observed during the spark which level out for later process times. Heat conduction to the electrodes is manifested by decreasing OH PLIF signal intensities in the vicinity of the electrodes.
trodes, heat transfer to the electrodes and heat conduction in the gas phase cools down the burnt gas as manifested by decreasing OH concentrations and therefore slows down flame propagation in vertical direction parallel to the electrodes. Flame front positions were extracted using the image processing strategy outlined above. In this approach the position of the flame front is defined via the steepest gradient in the rising edge of the flame propagating into the unburned gases. Figure 7 shows radii of flame kernels from 50 to 300 s. Error bars denote standard deviations. Mixtures were varied from 20% ( ⫽ 1.68), 15% ( ⫽ 2.38) to 5% ( ⫽ 7.98) H2 in air. In addition to the standard spark conditions (duration 100 s, current 1.5 A), a short arc phase (55 s) with low current (0.5 A) as well as different electrode geometry (stump electrodes, standard spark parameters) were investigated. The observed propagation can be divided into two phases: in a first phase up to approximately 80 s, the propagation is dominated by expansion caused by the high temperatures in the kernel. In a transition regime the propagation velocity slows down and meets the second phase which covers times larger than 100
M. THIELE ET AL.
Fig. 7. Radii of the flame kernel at different times after breakdown. In a parametric study the gas mixture, the spark duration and intensity as well as the electrode geometry is changed. For the purpose of comparison, the radius of the flame kernel obtained from numerical simulations is included to this figure. Notice, that the experimental data stem from a observation plane shifted 0.55 mm from the electrodes to reduce scattered light, while the numerical results are taken from the central plane through the electrodes. Clearly, the observed kernel size depends on the plane of observation at early stages of the ignition process (compare experimental data points at t ⫽ 50 s). For later times, however, this is of minor influence.
s. This phase is governed by the speed of the flame propagation. Consequently, flame propagation of the mixture with ⫽ 1.68 shows the highest velocity of flame propagation of 7.3 msec⫺1. The parametric variation for stump electrodes (see below) as well as shorter and less intense sparks clearly shows that these parameters are of minor influence on the radial flame propagation perpendicular to the electrodes. The propagation velocity parallel to the electrodes and the cross-section of the flame kernel is strongly influenced by different electrode geometries, as shown in Fig. 8. This figure compares flame kernels 190 s after breakdown. For electrodes with sharpened tips (conical shape of electrode tips, compare insert of Fig. 2) an almost oval flame kernel is observed with small indentations in the vicinity of the electrodes. For this case the center of the flame kernel always appeared at the same location indicating that the spark is anchored by the tips of the electrodes. Using stump electrodes (cylindrical shape of electrode tips, electrode diameter 1 mm) for these early process times a toroidal shape is observed. The asymmetry ob-
SPARK IGNITED HYDROGEN/AIR MIXTURES
83
Fig. 8. Influence of the electrode geometry on the crosssection of the flame kernel at 190 s after breakdown. For sharpened electrodes an oval shape is observed (top). In case of stump electrodes (bottom) a toroidal shape can be identified.
served for the stump electrodes occurred frequently indicating varying locations of the spark. As a consequence of an off-center position of the spark a part of the flame has to propagate through the electrode gap and is more severely disturbed by wall effects, resulting in the toroidal shape. To investigate the influence of heat conduction to the electrodes experimentally, temperature measurements inside the flame kernel are carried out using CARS. These measurements are performed in the center of the electrode gap. Because of comparable size of the interaction volume (2 mm long, 100 m in diameter) and the flame kernel at early process times, temperature measurements were restricted to late stages of ignition (t ⬎ 5 ms) to ensure almost homogeneous temperature distributions within the measurement volume. Spectra recorded at too early times (less than 2 ms) clearly display superimposed contributions from gases at widely different temperatures. Such spectra are not readily amenable to analysis. We note that so far only provisional research has been conducted on the problem of pronounced temperature variations throughout the interaction region [42– 44]. In addition, at early times of the ignition process the signal-to-noise ratio is too low to give reliable temperature results. This problem can be attributed to misalignment of the signal by beam steering: at these early delays, the flame kernel is comparable in size to the length of the interaction region and the beams therefore experience strong refractive index gradients. Therefore, only CARS data
Fig. 9. Top graph shows a CARS spectrum averaged for 10 shots at a delay of 5 ms with best-fit temperature of 1,353 K. Bottom graph is a 10 shot average at a delay of 10 ms with best-fit temperature of 1,248 K.
measured at times larger than 5 ms are evaluated for temperature determination. Exemplarily Figure 9 shows experimental data (dashed line) in comparison to the best fit result. For 5 and 10 ms after breakdown the respective temperatures are 1353 K and 1248 K, which indicates the temperature decrease because of heat losses to the elctrodes. A global cooling rate has been estimated using the CARS data and the temperature information from the 2-D simulation 300 s after breakdown. At 300 s in the vicinity of the electrode gap a homogeneous temperature distribution of approximately 2,000 K is observed. At these times effects of the spark on the temperature distribution are already equalized. Using this temperature and the measured CARS temperatures at 5 and 10 ms after breakdown an exponential decay curve of the form T(t) ⫽ T0 ⫹ A exp(-t/) is calculated from the data. From this fit, values of T0 ⬇ 1200 K, ⬇ 2.5 ms and A ⬇ 900 K can be estimated. The error limits in T0, A and are because of
84 uncertainties in the time t⬘ when the adiabatic flame temperature is reached (a range of t⬘ between 10 s and 1 ms was taken into account), and because of the estimated experimental errors in the determination of the temperatures at 5 ms and 10 ms by CARS (⫾40 K). For the decay constant , the error is approximately 30%. Even though these values can not be compared to numerical simulations currently the CARS results are valuable for future extensions of the calculations. Numerical Simulation Numerical simulation has been performed for comparison to experimental data. Because of the computational costs these numerical studies are limited to process times of 300 s. In accordance with experimental parameters, a spark duration of 100 s and an arc phase current of 1.5 A has been employed. The geometry of the electrodes has been used in close agreement to that of the experimental investigation (sharpened tips, for simulations concerning stump tips see [45]). Figure 10 shows the evolution of the spark and flame kernel from 1 to 300 s via OH concentration distributions. Notice that within this representation the rotationary symmetry is used to construct the 3-D appearance of the kernel from 2-D calculations. To allow a comparison with the measurements a cut through the 3-D flame kernel has been performed at a plane 0.5 mm apart from the electrodes. In agreement with experimental observations the kernel shows an almost oval geometry at early times (up to 100 s). The propagation velocity perpendicular to the electrodes is larger than parallel to the electrodes resulting in an oval shape of the flame kernel. This is in qualitative agreement to the experimental findings presented in Fig. 5. However, comparing the ratio of the horizontal to the vertical extension, within the simulation the oval flame kernel is much more flat. This behavior can be explained by a larger energy content during the breakdown phase compared to experimental conditions (note that the higher energy content during the breakdown phase is an initial condition for the simulations and is taken from the literature [28]). The higher the energy of the break-
M. THIELE ET AL. down the larger is the extension of the initial plasma kernel. Because of the influence of the electrodes this initial plasma kernel is not spherical but shows larger extensions in a direction perpendicular to the electrodes leading to an oval shape too flat compared to the experimental data. For times beyond 150 s a toroidal shape develops from an initially observed oval shape. This has not been observed that clearly in the experimental investigation. To investigate the influence of heat conduction to the electrodes in more detail, simulations with heat conduction switched off are compared to those where heat conduction is taken into account. A substantial difference to the sequence presented in Fig. 10 (heat conduction switched on) is observed for times after 130 s. This indicates that, at least for the conditions investigated in this work, heat losses to the electrodes are of no significant importance during the initial phase and that the shape observed up to 130 s is determined by other parameters such as energy content of the spark or electrode geometry. Radial profiles of absolute OH concentrations are shown in Fig. 11. During the spark extremely high OH concentrations are observed which level out when the spark is switched off. This is in accordance to experimental findings shown in Fig. 6 where radial OH profiles are shown 50 and 80 s after breakdown. At these early times in the gap between the electrodes OH radical concentrations are a factor of 2 to 3 above the maximum OH concentrations observed 200 s after breakdown. From Fig. 11 it can be deduced, that the flame front separates from the plasma kernel at approximately 50 to 70 s after breakdown. This “birth” of the flame kernel compares well with experimental findings shown in Fig. 6 where a flame front is separating from the core between 50 and 80 s as obvious from the corresponding radial profiles. Flame growth in radial direction in the symmetry plane between the electrodes is included in Fig. 7 for comparison to experimental data. Similar to the image post-processing procedure for experimental data, the flame front is defined via the steepest spatial gradient of the OH concentration outside the spark-activated volume. Qualitatively an identical behavior is ob-
SPARK IGNITED HYDROGEN/AIR MIXTURES
85
Fig. 10. Results from numerical simulation showing the evolution of the spark kernel from 1 to 300 s via OH concentration fields. OH concentrations are given in mole fractions. Notice, that a plane 0.5 mm apart from the electrodes has been constructed by the assumption of radial symmetry. Axis scalings are in millimeters.
86
M. THIELE ET AL.
Fig. 11. Calculated radial profiles of OH mass fractions from 1 to 200 s. High OH concentrations are apparent during the spark and level out when the spark is switched off. The generation and development of a flame front leaving the plasma kernel can be observed between 50 to 70 s after breakdown.
served as in the experiments: in a first phase for times up to 60 s the propagation is fast and is driven by plasma expansion. In a subsequent transition region the propagation velocity slows down and approaches an expansion velocity in a second phase which is—in accordance to the experimental data— governed by the laminar flame propagation speed. In the calculations, however, propagation starts from kernels which are larger by approximately 400 m. This observation is because of a higher energy content during the breakdown. Notice, that temperatures at the end of the breakdown phase are taken as initial condition from [28] and have not been varied for closer accordance with experimental results. A better description of the initial conditions would improve the results. CONCLUSIONS Within this study the early stages of ignition and subsequent flame propagation have been investigated for lean hydrogen/air mixtures. Flame
kernel shapes have been investigated experimentally using LIF of OH radicals as well as numerically using a 2-D model. In general a satisfying agreement between experimental and numerical results has been found. A parametric variation of spark characteristics, equivalence ratio of the gas mixture and electrode geometry has been carried out. It was shown that the length of the arc phase as well as the geometrical shape of the electrodes does not influence the temporal development of the radial extension of the kernel significantly. In vertical direction parallel to the electrodes both in the experiment and the simulation flame propagation is slower because of pronounced heat losses to the electrode material. At the early stages of the ignition process investigated in this study the effect of heat conduction to the electrodes is more pronounced for stump electrode geometries. Therefore the overall kernel cross section in case of stump electrodes appears as a torus while for sharpened electrodes the kernel is close to an oval shape
SPARK IGNITED HYDROGEN/AIR MIXTURES with small indentations in the vicinity of the electrodes. Results from experiments and simulations suggest the birth of a self sustaining flame propagation for process times between 50 to 70 s after breakdown. M. Thiele is grateful for financial support of the Max Buchner Stiftung. Financing of this project by the DAAD is gratefully acknowledged. We are grateful for support by the Forschungszentrum Karlsruhe and the Physikalisch-Technische Bundesanstalt. REFERENCES 1.
2. 3. 4. 5.
6. 7. 8. 9. 10.
11. 12. 13.
14.
15. 16.
17.
Maly, R. R., Twenty-fifth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1994, p. 111. Dale, J. D., Checkel, M. D., and Smy, P. R., Prog. Energy Combust. Sci. 23:379 (1997). Breitung, W., FZK Nachrichten 29:347 (1997). Akram, M., AIAA J. 34:1835 (1996). Kono, M., Kumagai, S., and Sakai, T., Sixteenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1976, p. 757. Pischinger, S., and Heywood, J. B., SAE Technical Paper No. 900021. Kravchik, T., and Sher, E., Combust. Flame 99: 635 (1994). Kravchik, T., Sher, E., and Heywood, J. B., Combust. Sci. Tech. 108:1 (1995). Maas, U., and Warnatz, J., Combust. Flame 74:53 (1988). Boston, P. M., Bradley, D., Lung, F. K.-K., Vince, I. M., and Weinberg, F. J., Twentieth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1984, p. 141. Lim, M. T., Anderson, R. W., and Arpaci, V. S., Combust. Flame 69:303 (1987). Akindele, O. O., Bradley, D., Mak, P. W., and McMahon, M., Combust. Flame 47:129 (1982). Hanson, R. K., Twenty-first Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1986, p. 1677. Wolfrum, J., Twenty-seventh Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1998, p. 1. Kohse-Ho ¨inghaus, K., Prog. Energy Combust. Sci. 20: 203 (1994). Eckbreth, A. C., Laser Diagnostics for Combustion Temperature and Species, 2nd Edition, Abacus Press, Tunbridge Wells, G. B., 1987. Dreizler, A., Tadday, R., Monkhouse, P., and Wolfrum, J., Appl. Phys. B 57:85 (1993).
18. 19. 20. 21. 22.
23. 24.
25. 26. 27.
28. 29. 30. 31.
32. 33. 34. 35. 36.
37. 38. 39. 40. 41.
42. 43. 44. 45.
87 Stepowski, D., and Cotterau, M. J., Appl. Opt. 18:354 (1979). Ferrario, A., et al., Paper WD2, CLEO Meeting, Baltimore, 1983. Stenhouse, I. A., et al., Appl. Opt. 18:3819 (1979). Snelling, D. R., et al., Appl. Opt., 26:99 (1987). Dreizler, A., Lindenmaier, S., Maas, U., Hult, J., Alde´n, M., and Kaminski, C. F., Appl. Phys. B. 70:287 (2000). Sher, E., Ben-Yai’sch, J., and Kravchik, T., Combust. Flame 89:186 (1992). Scha¨fer, M., Schmidt, R., and Ko ¨hler, J., Twenty-sixth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1996, p. 2701. Warnatz, J., Maas, U., and Dibble, R. W., Combustion, 2nd edition, Springer Verlag, Berlin, 1999. Maas, U., and Warnatz, J., Impact Compu. Sci. Eng. 1:394 (1989). Yos, J. M., Transport Properties of Nitrogen, Hydrogen and Air to 30000 K, Technical Memorandum RADTM-63–7, AVCO Corporation, Wilmington (MA), March., 1963. Scha¨fer, M., Der Zu ¨ndfunke, Dissertation, Universita¨t Stuttgart, 1997. Maas, U., Dissertation, Heidelberg University, 1988. Noh, W. F., J. Comp. Phys. 72:78 (1987). Bird, R. B., Stewart, W. E., and Lightfoot, E. N., Transport Phenomena, John Wiley and Sons, New York, 1960. Thiele, M., Dissertation, University of Heidelberg, 1999. Deufflhard, P., Nowak, U., Tech. Rep. SFB 123, University of Heidelberg (1985). Hackbusch, W., Iterative Lo ¨sung groer schwachbesetzter Gleichungssysteme, Teubner Verlag, Stuttgart, 1992. Kaminski, C., Hult, J., and Alde´n, M., Appl. Phys. B 68:757 (1999). Lewis, B., and von Elbe, G., Combustion, Flames and Explosions of Gases, 3rd edition, Academic Press, Orlando., 1987. Ewart, P., Opt. Commun. 55:124 (1985). Snelling, D. R., et al., Appl. Opt. 33:8295 (1994). Eckbreth, A. C., Appl. Phys. Lett. 32:421 (1978). Farrow, R. L., Private communication. Maly, R. R., Vogel, M., Seventeenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1978, p. 821. Boquillon, J. P., et al., Opt. Lett. 13:722 (1988). Thumann, A., et al., Appl. Opt. 34:3313 (1995). Parameswaran, T., et al., Appl. Opt. 35:5461 (1996). Thiele, M., Warnatz, J., and Maas, U., SAE Technical Paper Series No. 01-1178, 1999.
Received 16 October 2000; revised 16 May 2001; accepted 18 August 2001
Interdisziplina ¨res Zentrum fu ¨r Wissenschaftliches Rechnen, Universita ¨t Heidelberg, Im Neuenheimer Feld 368, D-69120 Heidelberg, Germany
A. DREIZLER*
Energie- und Kraftwerkstechnik, Tu. Darmstadt, Petersenstrae 30, 64287 Darmstadt, Germany
S. LINDENMAIER, R. SCHIEßL and U. MAAS
Institut fu ¨r Technische Verbrennung, Universita ¨t Stuttgart, Pfaffenwaldring 12, D-70569 Stuttgart, Germany
and A. GRANT and P. EWART
Clarendon Laboratory, Oxford University, Parks Road, Oxford OX1 3PU, UK This study reports a detailed 2-D model which describes the spark ignition of an initially quiescent hydrogen/air mixture. The model includes the compressible Navier-Stokes equations, detailed chemistry and molecular transport in the gas phase as well as heat conduction to the electrodes. The spark is modeled for the phases subsequent to breakdown using the Maxwell equations for quasi-stationary conditions for the electric field. Initial and boundary conditions necessary for the simulations are chosen in accordance with experimental values. Heat conduction to the electrodes and different electrode shapes are investigated. The influence of these parameters on the shape of the initial flame kernel is discussed and compared qualitatively to experimental results. Spark ignition experiments are performed using a highly reproducible ignition system. Shapes of the early flame kernels are monitored by 2-D laser-induced fluorescence (PLIF) imaging of OH radicals produced during the ignition and the combustion process. The investigations are performed for different equivalence ratios. In addition, for a central position within the flame kernel, temperatures are measured at different times after ignition using vibrational coherent anti-Stokes Raman spectroscopy (CARS) of nitrogen. © 2002 by The Combustion Institute
INTRODUCTION Spark ignition processes are of enormous practical importance. In technical combustion devices, for example, engines or gas turbines [1, 2], spark ignition initiates combustion and influences the temporal development of this process in a critical way. Safety concepts in atomic power plants rely on controlled burning of hydrogen resigned from the reactor [3]. Therefore a deeper understanding of spark ignition is required. In future, a detailed theoretical understanding may improve current concepts of spark plugs in connection with improved electrical ignition devices. To further optimize such models a comparison with experimental data is still necessary. *Corresponding author. E-mail: [email protected] 0010-2180/02/$–see front matter PII 0010-2180(01)00333-9
Spark ignition constitutes a very complex interplay between plasma kinetics, chemical kinetics, molecular transport processes and fluid dynamics. A complete mathematical simulation of a spark ignition is a difficult task because of the enormous numerical problems, due to the stiffness and high dimensionality of the problem (each chemical species introduces an additional conservation equation). Furthermore, experimental investigations of spark kernels and their transition to flame kernels are rendered difficult because of very short process times, extremely high core temperatures and large gradients in the refractive index. In general, only flame propagation subsequent to spark ignition can be studied by laser diagnostic methods. Therefore, internal structures of the plasma core remain mostly unknown. In the literature various simplified models of spark ignition have been presented. The followCOMBUSTION AND FLAME 128:74 – 87 (2002) © 2002 by The Combustion Institute Published by Elsevier Science Inc.
SPARK IGNITED HYDROGEN/AIR MIXTURES ing examples give a brief survey of the field, although this is far from being a comprehensive list. In a 2-D setting the effect of the electrical spark on the flow field has been studied for an air atmosphere [4] and for one step propane chemistry [5]. Heat conduction to the electrodes was studied by Pischinger et al. [6]. A detailed reaction mechanism for a methane-air mixture has been employed by Sher et al. [7, 8] employing rectangular electrodes. Detailed chemistry in conjunction with detailed molecular transport has been used in a one-dimensional geometry, neglecting heat conduction to the electrodes [9]. None of the above studies combine all aspects in one model. In this study an improved model is discussed which considers heat conduction from the gas phase to the electrodes, detailed chemistry and molecular transport as well as the coupling of the gas dynamics to the properties of the electrical discharge through heating by the electrical current. The computational model is reduced to two dimensions by assuming rotational symmetry along the axis of the electrodes. Different scenarios are considered for a lean hydrogenoxygen mixture. In general, the model is used to study the transition from the initial phase (governed by the shock wave formation) to the flame propagation. The influence of the electrode as well as the influence of heat conduction from the hot gas to the electrode are investigated. A qualitative comparison is made with experimental data. For this comparison initial and boundary conditions necessary for the numerical simulation are taken in accordance with the experimental settings. Optical diagnostic methods have become a standard tool for the investigation of chemically reacting flows. Traditional techniques like Schlieren [10, 11] photography and flame emission [12] filming are line-of-sight integrated. Therefore, two-dimensional laser-based techniques were developed to investigate the internal structures of combustion processes [13]. Laser-induced fluorescence (LIF) counts among the most popular techniques because of its experimental simplicity. The spatial distribution of minority species can be detected using commercially available equipment [14, 15]. The extraction of quantitative concentration measurements from LIF signals is often rendered
75
difficult because of collisionally induced energy transfer [16, 17]. This problem can be overcome by application of various strategies, as shown for example in [15, 16, 18]. However, in many cases, qualitative information is sufficient. Two-dimensional LIF of OH radicals is an ideal method for the comparison of early sparkignited flame kernels with corresponding numerical simulations. For this reason, the OH distribution in a plane parallel to the spark has been investigated. Coherent anti-Stokes Raman scattering (CARS) is a well-established laser-based spectroscopic technique for the measurement of species concentration and temperature. CARS has proved its value in a large variety of combustion devices, from furnaces [19] to internal combustion engines [20]. The CARS signal generation process and its dependence on physical properties, such as temperature and pressure, is complex and is discussed in detail in [16]. For time-varying combustion processes such as spark ignition, single shot multiplex CARS [21] measurements are necessary. In this study, temperatures in the center of the electrode gap have been determined 5 and 10 ms after the breakdown using vibrational nitrogen CARS. For all experimental studies a well characterized and highly reproducible ignition system presented in a previous study [22] has been used. In a parametric study, gas mixture, spark parameters, and electrode geometry have all been varied.
MODELING PRINCIPLES Equations The model for the simulation of the early development of a flame kernel initiated by an electrical spark is based on the conservation equations for reacting flows, including detailed chemistry. To make the problem treatable the following assumptions are made: ● ● ●
The plasma channel formed after the breakdown is cylindrical. Local thermal equilibrium exists The influence of the magnetic field is negligible.
76
M. THIELE ET AL. 1 1 ⭸T ⫹ vgradT ⫹ div j q ⫹ pញ gradv ⭸t cv cv ⫹
1 cv
冘 冋 u 共M ⫺ div j 兲册 ⫽ 1c 共q ⫹ q 兲 ns i
i
i
i
i
v
s
r
(1) a source term qs is introduced to include the properties of the electric spark. Combining Ohms Law for the electric current flow (jel) with the electrical power input leads to j el ⫽ 兩E兩 2.
Fig. 1. Illustration of the computational area. The conducting plasma channel is assumed to be cylindrical. For details concerning initial and boundary conditions see text.
● ●
The plasma is quasi-neutral. The plasma is transparent to its own radiation
In modeling spark ignition, these assumptions are commonly made [4, 7, 8, 23]. The first assumption allows the geometry to be reduced to two dimensions because of rotational symmetry around the plasma channel as shown in Fig. 1. The second assumption is based on the fact that internal relaxation processes in molecules are faster than the time scales of the investigated processes. The last two assumptions are needed to derive the equation for the electric potential [4, 24]. Temperature- and pressure-dependent transport coefficients for the mixture can be calculated from binary data using mixing rules. The gas phase equations consist of the compressible Navier-Stokes equations (conservation of momentum and total mass) together with energy conservation. For each of the chemical species in the reaction mechanism (9 species and 38 elementary reactions [25], includes dissociation), a mass conservation equation is added to the system, which is closed by the ideal gas law. All equations are solved in their primitive variable formulation as presented in [9]. In the energy equation [26]
(2)
Assuming a cylindrical channel, the radial field component can be neglected because the current I flows along the z-direction only [4, 23]. The source term can then be calculated from the current strength by integrating over the plasma. As soon as the cylindrical shape gets distorted this assumption is no longer valid. From the conservation of charge, Ohms Law and Maxwell’s equations of electrostatics, a time-dependent differential equation for the ជ ⫽⫺grad⌽ results. scalar electric potential E Quasistationarity leads to div(⫺grad⌽) ⫽ 0
(3)
because of the small relaxation times compared with the time-scales of the reacting flow. This equation is added to the system and via ⌽ one ជ and the electric source qs. can compute E The electrical conductivity is derived from temperature and pressure by interpolation of the data of Yos [27] for ionized air. Notice, that the influence of the ions on the chemical mechanism with regard to flame propagation is negligible and is not accounted for in the present approach. Initial and Boundary Conditions Initial and boundary conditions have to be specified for the solution. In this work, we restrict ourselves to simulations for initially quiescent mixtures. The model is, however, able to treat initial velocity fields as well. The properties of the channel at the end of the breakdown are used as initial temperature profile and are taken from [28]. Species profiles have been calculated according to the spark characteris-
SPARK IGNITED HYDROGEN/AIR MIXTURES tics, that is, the initial temperature and pressure profile is used for a homogenous computation [29]. After some computation time, the breakdown time, species profiles according to these temperatures and pressures result. The axis of symmetry (r ⫽ 0), as well as the equipotential line between the electrodes (z ⫽ Z), are treated as symmetry lines. The other two boundaries are treated as walls (compare Fig. 1). For symmetry lines, the normal velocity component is set to zero as are all other gas dynamics variables. For wall boundaries no slip conditions for the velocity and vanishing normal gradients of all other gas dynamics variables are specified. For the equation for the electric potential vanishing normal gradients are set on all boundaries except for the equipotential line, where ⌽ ⫽ 0 is specified. At the electrode vanishing gradients for the species mass fractions are assumed. Adiabatic conditions as well as heat conduction have been implemented for the energy equation via:
⭸T ⫺ T E兲 ⫽ 0, ⫹ h共T ⭸nជ
再
冉冑
⫽ ⫺ a 2 ⌬共r,z兲 2
1 共⌬:⌬兲 2
冊冎
.
(6)
Herein ⌬ denotes the stress tensor, a is a damping factor of the artificial viscosity, “:” denotes the two-fold contraction and ⌬(r,z) accounts for the local mesh size in both spatial directions r and z. Because of the stiffness introduced by the chemical kinetics, we used an implicit solution method. This signifies in this context, that all equations except the potential equation are solved implicitly. The electrical source term is computed from the conductivity and the electrical field of the preceding time step. The potential equation is solved after each time step using a Newtonian scheme. The resulting ordinary differential equations are solved by the extrapolation solver LIMEX [33], including step size and order control. The regular discretization scheme leads to block-structured Jacobi- and iteration matrices which can be solved by block-ILU decomposition, reducing computation time and storage requirements [29, 34].
(4)
where nជ denotes the outwards directed normal, h the heat exchange coefficient and the heat conductivity. Setting h ⫽ 0, ⫽ 1 ones obtain adiabatic conditions while using both, and h, mixed conditions for heat conduction result. The current input is treated by setting ⫺ I共t兲 ⭸⌽ ⫽ 2 ⭸z r elec
77
(5)
where relec denotes the radius of the electrodes. Numerical Solution For the numerical solution the method of lines is employed. Spatial discretization is performed by finite differences [26]. Artificial viscosity is added to resolve the pressure wave emitted from the hot channel. An artificial viscosity coefficient proportional to the flux of momentum is added to the viscosity coefficient following [30] and applying the equations of Ostwaldde-Waele [31]. It is set to [32]
EXPERIMENTAL SETUP Ignition Vessel A spark ignition system has been designed at the Institut fu ¨r Technische Verbrennung for the generation of highly reproducible sparks. This allows comparison between experimental results and the predictions of the model for various gas mixtures and flow conditions. The reproducibility of the system can be proved by measurements of flame propagation. For this purpose, planar laser-induced fluorescence measurements of the spatial distribution of OH radicals have been performed in a previous study [22]. By imaging the temporal evolution of several individual ignition events using a novel high speed LIF image acquisition system [35] variations of the spark characteristics with respect to the flame initiation have been investigated. The cycle-to-cycle variations in flame propagation for the employed ignition system were found to be less than 2%. Details concerning the test of the reproducibility as well as the
78
M. THIELE ET AL.
Fig. 2. Sketch of the ignition constant volume vessel. The device allows full optical access. Electrodes are made from 1-mm thick tungsten wire. In its standard configuration the tips of the electrodes are sharpened (conical shape) as shown in the insert. For PLIF applications a plane approximately 0.55 mm apart from the electrodes is used. CARS measurements are performed at the center of the electrode gap.
design of the ignition system can be found in [22]. The ignition system is attached to a constant volume ignition chamber (cylindrical geometry, volume 2.8 l, diameter 160 mm) with full optical access (Fig. 2). The ignition is initiated in the geometric center of the vessel. The diameter of the tungsten electrodes was 1 mm, and the spark gap was also set to 1 mm. The electrode geometry generally used is shown in Fig. 2. Typical current and voltage traces of the spark used during this parameter study are presented in Fig. 3. The spark consists of a short breakdown phase (150 ns duration) and a long user-selectable arc phase. In the present study this was set to 55 and 100 s, respectively. The electrical
energy of this spark was adjusted to 8.7 mJ well above minimum ignition energy [36]. Because of imperfect conversion thermal energy input to the gas phase is estimated to be 4.5 mJ based on calorimetric measurements conducted for similar spark conditions [22]. The temporal evolution of the pressure inside the vessel was monitored with a piezo-electric pressure transducer. Following a standardized purging procedure, the cell was filled with a lean mixture of hydrogen/air (5, 15, or 20% hydrogen in air). The partial pressures of the gases were monitored with a calibrated baratron (MKS 122) and could be adjusted to an accuracy of approximately 0.5 mbar. A pressure of 1 bar before ignition was used for all the results presented below. The time between ignition events was sufficiently long to ensure the cylinder walls had equilibrated to 305 K. OH Planar Laser-Induced Fluorescence (PLIF) Measurements and Data Reduction
Fig. 3. Current and voltage traces of a typical spark event. The peaks at the leading edges of the traces are because of the breakdown which can not be temporally resolved by the probes.
OH radicals produced during the combustion of the hydrogen/air mixture were probed using the Q1(8) line in the (1– 0)-band of the A2⌺⫹-X2⌸ transition. The excitation wavelength of around 283 nm was generated by a dye laser (Rhodamine 590) which was pumped with a XeCl excimer laser (Lambda Physik). After type I frequency doubling in a BBO crystal pulse energies of 0.2 mJ were obtained. The laser
SPARK IGNITED HYDROGEN/AIR MIXTURES beam was formed into a sheet (height 4 mm, thickness 50 m as measured by use of the “razor blade” method) intersecting the flame kernel in a vertical plane displaced approximately 0.55 mm from the electrodes as shown in Fig. 2. Fluorescence light at around 310 nm was detected perpendicular to the laser propagation direction using appropriate dielectric filters and an intensified CCD camera (Princeton Instruments, 1024 ⫻ 256 pixel, 16 bit). The ignition system, the laser and the camera system were triggered by two pulse/delay generators (Stanford Research Systems DG-535). The laser pulse and the exposure of the camera system were delayed relative to the beginning of the spark. The following data reduction procedure was used to extract the flame radius from the OH PLIF images: 1.
The 16 bit gray-scale images were converted to binary images using a simple thresholding technique. This is possible because of a very high signal-to-noise ratio of the original intensity images. The center of the flame and the approximate position of the flame front were extracted from these binary images. 2. The original gray-scale images were filtered by a 2-D Gaussian-filter to reduce noise. This linear filter preserves the position of steep gradients in the image. 3. A one-dimensional intensity profile was produced along a horizontal line through the center of the flame (established in Step 1) through the flame front. 4. The exact position of the flame front is located within a predefined area in the vicinity of the rough flame front obtained from the binary image (Step 1). The flame front position is defined as the absolute maximum of the respective gradients which lie within this predefined area. 5. The horizontal flame diameter was taken as the distance between the opposite flame front positions. Ensemble averaged horizontal flame diameters and standard deviations were derived from a number of single-shot events performed for identical boundary conditions.
79
Vibrational CARS Measurements A single mode, frequency-doubled Nd:YAG laser (Continuum PowerLite 8000) is used to provide the pump beams for the CARS process; part of the output is split off to pump the broadband Stokes laser (Mode-x ML-2). This laser is based on the modeless laser [37], which has been shown to provide the most precise single-shot CARS measurements when used in conjunction with a single-mode pump laser [38]. The central wavelength of the Stokes laser is tuned to around 607 nm; the spectral width is sufficient to simultaneously access the cold band and the first hot band of the N2 vibrational CARS spectrum, allowing multiplex measurements to be made. A delay line is used to ensure temporal overlap of the pump and Stokes beams. Beam splitter BS1 (compare Fig. 4) produces two pump beams disposed vertically above each other. The energy of the pump beams is adjusted by rotation of beam splitter BS2 which transmitts a variable fraction of the incident energy depending on the angle of incidence. The Stokes beam energy is set by selecting a given area of the beam using aperture A. The glass block GB in the Stokes beam is rotated to adjust the spatial overlap of the beams and hence to optimize the CARS signal. The Stokes and pump beams display very similar collimation and beam size: no significant signal increase is attained by adjusting the collimation and diameter of the Stokes beam. The HeliumNeon laser is used to provide a mock signal beam for alignment purposes. The beams are crossed by a lens of focal length 250 mm into the forward folded boxcars geometry [39]; the resulting interaction region is approximately 2 mm long and 100 m wide. Care has been taken to ensure that the interaction region is situated in the center of the spark gap. The ignition cell is filled with a mixture of 20% hydrogen and 80% air and then ignited. CARS measurements are timed relative to the breakdown by the use of programmable pulse generators (Stanford Research Systems DG535). The generated CARS signal is collimated, filtered by suitable interference filters and then focused onto the entrance slit of a 1 m Czerny-
80
M. THIELE ET AL.
Fig. 4. Experimental set-up of the CARS spectrometer.
Turner type spectrometer (Hilger Monospek 1000). The spectra are recorded on a CCD camera (Princeton Instruments 512TKB). Camera operation and data acquisition are controlled by a PC. Because of the high reproducibility of the ignition cell, the single shot CARS spectra were averaged in software. Background levels are subtracted and the CARS spectra are then normalized using the non-resonant spectrum generated in propane. Non-resonant spectra are recorded in situ by filling the ignition cell with pure propane. This normalization takes account of the Stokes laser spectral profile and spatial non-uniformities in the camera CCD chip. For CARS data analysis the Sandia Laboratory code CARSFIT [40] is used to find a best-fit temperature to the averaged spectra. RESULTS AND DISCUSSION Experimental Observations As standard conditions the parameters of the arc phase are set to a duration of 100 s and a current of 1.5 A corresponding to an electrical energy of 8.7 mJ as calculated from current and voltage traces. Figure 5 shows qualitatively the
evolution of the spark kernel via relative OH concentration fields in a plane shifted 0.55 mm with respect to the electrodes. Correction of the laser intensity profile has not been performed. For an ignitable mixture of 20% H2 in air (left column of Fig. 5) the flame develops from an initially oval kernel. Note, that for process times larger than 250 s in vertical direction only the central part of the kernel is monitored because of the limited height of the laser light sheet. However, for a mixture of 5% H2 in air (below the ignition limit) a substantial amount of OH is produced as well. The sequence on the right column of Fig. 5 shows a “dying” plasma kernel where no transition to a self sustaining flame propagation is observed because of the lean fuel/air mixture. For this case remarkable OH concentrations can be monitored up to 1 ms after the end of the spark. By using these raw data, OH profiles perpendicular to the axis of the electrodes are extracted using the image processing procedure described above. Figure 6 shows qualitatively radial OH concentrations profiles from 50 to 1,300 s. The influence of high core temperatures of approximately 5,000 K [41] causes high peak OH concentrations (50 s after break-
SPARK IGNITED HYDROGEN/AIR MIXTURES
81
Fig. 5. Relative OH concentration profiles recorded by PLIF. These sequences show the temporal evolution of a flame kernel for an ignitable mixture of 20% H2 in air and a mixture of 5% H2 in air (below the ignition limit). Notice, that the beam profile has not been corrected for inhomogeneous intensity distribution. For process times later than 250 s the flame kernel exceeds in vertical direction the height of the laser light sheet.
down) which decrease with time when the spark is switched off. 80 s after breakdown the generation and evolution of a flame front can already be observed propagating radially out-
wards. At this time the flame kernel starts to separate from the plasma kernel. Subsequently, the flame front propagates into the fresh gases while in the vicinity of the elec-
82
Fig. 6. Radial OH profiles obtained in a direction perpendicular to the pointing of the electrodes (radial direction). Intensities are given in relative units; spatial variation of local quenching conditions was not taken into account. Experimental conditions were: 20% hydrogen in air, sharpened (conical) electrode tips, 100 s long arc phase and 1.5 A current during arc phase. High OH concentration levels are observed during the spark which level out for later process times. Heat conduction to the electrodes is manifested by decreasing OH PLIF signal intensities in the vicinity of the electrodes.
trodes, heat transfer to the electrodes and heat conduction in the gas phase cools down the burnt gas as manifested by decreasing OH concentrations and therefore slows down flame propagation in vertical direction parallel to the electrodes. Flame front positions were extracted using the image processing strategy outlined above. In this approach the position of the flame front is defined via the steepest gradient in the rising edge of the flame propagating into the unburned gases. Figure 7 shows radii of flame kernels from 50 to 300 s. Error bars denote standard deviations. Mixtures were varied from 20% ( ⫽ 1.68), 15% ( ⫽ 2.38) to 5% ( ⫽ 7.98) H2 in air. In addition to the standard spark conditions (duration 100 s, current 1.5 A), a short arc phase (55 s) with low current (0.5 A) as well as different electrode geometry (stump electrodes, standard spark parameters) were investigated. The observed propagation can be divided into two phases: in a first phase up to approximately 80 s, the propagation is dominated by expansion caused by the high temperatures in the kernel. In a transition regime the propagation velocity slows down and meets the second phase which covers times larger than 100
M. THIELE ET AL.
Fig. 7. Radii of the flame kernel at different times after breakdown. In a parametric study the gas mixture, the spark duration and intensity as well as the electrode geometry is changed. For the purpose of comparison, the radius of the flame kernel obtained from numerical simulations is included to this figure. Notice, that the experimental data stem from a observation plane shifted 0.55 mm from the electrodes to reduce scattered light, while the numerical results are taken from the central plane through the electrodes. Clearly, the observed kernel size depends on the plane of observation at early stages of the ignition process (compare experimental data points at t ⫽ 50 s). For later times, however, this is of minor influence.
s. This phase is governed by the speed of the flame propagation. Consequently, flame propagation of the mixture with ⫽ 1.68 shows the highest velocity of flame propagation of 7.3 msec⫺1. The parametric variation for stump electrodes (see below) as well as shorter and less intense sparks clearly shows that these parameters are of minor influence on the radial flame propagation perpendicular to the electrodes. The propagation velocity parallel to the electrodes and the cross-section of the flame kernel is strongly influenced by different electrode geometries, as shown in Fig. 8. This figure compares flame kernels 190 s after breakdown. For electrodes with sharpened tips (conical shape of electrode tips, compare insert of Fig. 2) an almost oval flame kernel is observed with small indentations in the vicinity of the electrodes. For this case the center of the flame kernel always appeared at the same location indicating that the spark is anchored by the tips of the electrodes. Using stump electrodes (cylindrical shape of electrode tips, electrode diameter 1 mm) for these early process times a toroidal shape is observed. The asymmetry ob-
SPARK IGNITED HYDROGEN/AIR MIXTURES
83
Fig. 8. Influence of the electrode geometry on the crosssection of the flame kernel at 190 s after breakdown. For sharpened electrodes an oval shape is observed (top). In case of stump electrodes (bottom) a toroidal shape can be identified.
served for the stump electrodes occurred frequently indicating varying locations of the spark. As a consequence of an off-center position of the spark a part of the flame has to propagate through the electrode gap and is more severely disturbed by wall effects, resulting in the toroidal shape. To investigate the influence of heat conduction to the electrodes experimentally, temperature measurements inside the flame kernel are carried out using CARS. These measurements are performed in the center of the electrode gap. Because of comparable size of the interaction volume (2 mm long, 100 m in diameter) and the flame kernel at early process times, temperature measurements were restricted to late stages of ignition (t ⬎ 5 ms) to ensure almost homogeneous temperature distributions within the measurement volume. Spectra recorded at too early times (less than 2 ms) clearly display superimposed contributions from gases at widely different temperatures. Such spectra are not readily amenable to analysis. We note that so far only provisional research has been conducted on the problem of pronounced temperature variations throughout the interaction region [42– 44]. In addition, at early times of the ignition process the signal-to-noise ratio is too low to give reliable temperature results. This problem can be attributed to misalignment of the signal by beam steering: at these early delays, the flame kernel is comparable in size to the length of the interaction region and the beams therefore experience strong refractive index gradients. Therefore, only CARS data
Fig. 9. Top graph shows a CARS spectrum averaged for 10 shots at a delay of 5 ms with best-fit temperature of 1,353 K. Bottom graph is a 10 shot average at a delay of 10 ms with best-fit temperature of 1,248 K.
measured at times larger than 5 ms are evaluated for temperature determination. Exemplarily Figure 9 shows experimental data (dashed line) in comparison to the best fit result. For 5 and 10 ms after breakdown the respective temperatures are 1353 K and 1248 K, which indicates the temperature decrease because of heat losses to the elctrodes. A global cooling rate has been estimated using the CARS data and the temperature information from the 2-D simulation 300 s after breakdown. At 300 s in the vicinity of the electrode gap a homogeneous temperature distribution of approximately 2,000 K is observed. At these times effects of the spark on the temperature distribution are already equalized. Using this temperature and the measured CARS temperatures at 5 and 10 ms after breakdown an exponential decay curve of the form T(t) ⫽ T0 ⫹ A exp(-t/) is calculated from the data. From this fit, values of T0 ⬇ 1200 K, ⬇ 2.5 ms and A ⬇ 900 K can be estimated. The error limits in T0, A and are because of
84 uncertainties in the time t⬘ when the adiabatic flame temperature is reached (a range of t⬘ between 10 s and 1 ms was taken into account), and because of the estimated experimental errors in the determination of the temperatures at 5 ms and 10 ms by CARS (⫾40 K). For the decay constant , the error is approximately 30%. Even though these values can not be compared to numerical simulations currently the CARS results are valuable for future extensions of the calculations. Numerical Simulation Numerical simulation has been performed for comparison to experimental data. Because of the computational costs these numerical studies are limited to process times of 300 s. In accordance with experimental parameters, a spark duration of 100 s and an arc phase current of 1.5 A has been employed. The geometry of the electrodes has been used in close agreement to that of the experimental investigation (sharpened tips, for simulations concerning stump tips see [45]). Figure 10 shows the evolution of the spark and flame kernel from 1 to 300 s via OH concentration distributions. Notice that within this representation the rotationary symmetry is used to construct the 3-D appearance of the kernel from 2-D calculations. To allow a comparison with the measurements a cut through the 3-D flame kernel has been performed at a plane 0.5 mm apart from the electrodes. In agreement with experimental observations the kernel shows an almost oval geometry at early times (up to 100 s). The propagation velocity perpendicular to the electrodes is larger than parallel to the electrodes resulting in an oval shape of the flame kernel. This is in qualitative agreement to the experimental findings presented in Fig. 5. However, comparing the ratio of the horizontal to the vertical extension, within the simulation the oval flame kernel is much more flat. This behavior can be explained by a larger energy content during the breakdown phase compared to experimental conditions (note that the higher energy content during the breakdown phase is an initial condition for the simulations and is taken from the literature [28]). The higher the energy of the break-
M. THIELE ET AL. down the larger is the extension of the initial plasma kernel. Because of the influence of the electrodes this initial plasma kernel is not spherical but shows larger extensions in a direction perpendicular to the electrodes leading to an oval shape too flat compared to the experimental data. For times beyond 150 s a toroidal shape develops from an initially observed oval shape. This has not been observed that clearly in the experimental investigation. To investigate the influence of heat conduction to the electrodes in more detail, simulations with heat conduction switched off are compared to those where heat conduction is taken into account. A substantial difference to the sequence presented in Fig. 10 (heat conduction switched on) is observed for times after 130 s. This indicates that, at least for the conditions investigated in this work, heat losses to the electrodes are of no significant importance during the initial phase and that the shape observed up to 130 s is determined by other parameters such as energy content of the spark or electrode geometry. Radial profiles of absolute OH concentrations are shown in Fig. 11. During the spark extremely high OH concentrations are observed which level out when the spark is switched off. This is in accordance to experimental findings shown in Fig. 6 where radial OH profiles are shown 50 and 80 s after breakdown. At these early times in the gap between the electrodes OH radical concentrations are a factor of 2 to 3 above the maximum OH concentrations observed 200 s after breakdown. From Fig. 11 it can be deduced, that the flame front separates from the plasma kernel at approximately 50 to 70 s after breakdown. This “birth” of the flame kernel compares well with experimental findings shown in Fig. 6 where a flame front is separating from the core between 50 and 80 s as obvious from the corresponding radial profiles. Flame growth in radial direction in the symmetry plane between the electrodes is included in Fig. 7 for comparison to experimental data. Similar to the image post-processing procedure for experimental data, the flame front is defined via the steepest spatial gradient of the OH concentration outside the spark-activated volume. Qualitatively an identical behavior is ob-
SPARK IGNITED HYDROGEN/AIR MIXTURES
85
Fig. 10. Results from numerical simulation showing the evolution of the spark kernel from 1 to 300 s via OH concentration fields. OH concentrations are given in mole fractions. Notice, that a plane 0.5 mm apart from the electrodes has been constructed by the assumption of radial symmetry. Axis scalings are in millimeters.
86
M. THIELE ET AL.
Fig. 11. Calculated radial profiles of OH mass fractions from 1 to 200 s. High OH concentrations are apparent during the spark and level out when the spark is switched off. The generation and development of a flame front leaving the plasma kernel can be observed between 50 to 70 s after breakdown.
served as in the experiments: in a first phase for times up to 60 s the propagation is fast and is driven by plasma expansion. In a subsequent transition region the propagation velocity slows down and approaches an expansion velocity in a second phase which is—in accordance to the experimental data— governed by the laminar flame propagation speed. In the calculations, however, propagation starts from kernels which are larger by approximately 400 m. This observation is because of a higher energy content during the breakdown. Notice, that temperatures at the end of the breakdown phase are taken as initial condition from [28] and have not been varied for closer accordance with experimental results. A better description of the initial conditions would improve the results. CONCLUSIONS Within this study the early stages of ignition and subsequent flame propagation have been investigated for lean hydrogen/air mixtures. Flame
kernel shapes have been investigated experimentally using LIF of OH radicals as well as numerically using a 2-D model. In general a satisfying agreement between experimental and numerical results has been found. A parametric variation of spark characteristics, equivalence ratio of the gas mixture and electrode geometry has been carried out. It was shown that the length of the arc phase as well as the geometrical shape of the electrodes does not influence the temporal development of the radial extension of the kernel significantly. In vertical direction parallel to the electrodes both in the experiment and the simulation flame propagation is slower because of pronounced heat losses to the electrode material. At the early stages of the ignition process investigated in this study the effect of heat conduction to the electrodes is more pronounced for stump electrode geometries. Therefore the overall kernel cross section in case of stump electrodes appears as a torus while for sharpened electrodes the kernel is close to an oval shape
SPARK IGNITED HYDROGEN/AIR MIXTURES with small indentations in the vicinity of the electrodes. Results from experiments and simulations suggest the birth of a self sustaining flame propagation for process times between 50 to 70 s after breakdown. M. Thiele is grateful for financial support of the Max Buchner Stiftung. Financing of this project by the DAAD is gratefully acknowledged. We are grateful for support by the Forschungszentrum Karlsruhe and the Physikalisch-Technische Bundesanstalt. REFERENCES 1.
2. 3. 4. 5.
6. 7. 8. 9. 10.
11. 12. 13.
14.
15. 16.
17.
Maly, R. R., Twenty-fifth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1994, p. 111. Dale, J. D., Checkel, M. D., and Smy, P. R., Prog. Energy Combust. Sci. 23:379 (1997). Breitung, W., FZK Nachrichten 29:347 (1997). Akram, M., AIAA J. 34:1835 (1996). Kono, M., Kumagai, S., and Sakai, T., Sixteenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1976, p. 757. Pischinger, S., and Heywood, J. B., SAE Technical Paper No. 900021. Kravchik, T., and Sher, E., Combust. Flame 99: 635 (1994). Kravchik, T., Sher, E., and Heywood, J. B., Combust. Sci. Tech. 108:1 (1995). Maas, U., and Warnatz, J., Combust. Flame 74:53 (1988). Boston, P. M., Bradley, D., Lung, F. K.-K., Vince, I. M., and Weinberg, F. J., Twentieth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1984, p. 141. Lim, M. T., Anderson, R. W., and Arpaci, V. S., Combust. Flame 69:303 (1987). Akindele, O. O., Bradley, D., Mak, P. W., and McMahon, M., Combust. Flame 47:129 (1982). Hanson, R. K., Twenty-first Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1986, p. 1677. Wolfrum, J., Twenty-seventh Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1998, p. 1. Kohse-Ho ¨inghaus, K., Prog. Energy Combust. Sci. 20: 203 (1994). Eckbreth, A. C., Laser Diagnostics for Combustion Temperature and Species, 2nd Edition, Abacus Press, Tunbridge Wells, G. B., 1987. Dreizler, A., Tadday, R., Monkhouse, P., and Wolfrum, J., Appl. Phys. B 57:85 (1993).
18. 19. 20. 21. 22.
23. 24.
25. 26. 27.
28. 29. 30. 31.
32. 33. 34. 35. 36.
37. 38. 39. 40. 41.
42. 43. 44. 45.
87 Stepowski, D., and Cotterau, M. J., Appl. Opt. 18:354 (1979). Ferrario, A., et al., Paper WD2, CLEO Meeting, Baltimore, 1983. Stenhouse, I. A., et al., Appl. Opt. 18:3819 (1979). Snelling, D. R., et al., Appl. Opt., 26:99 (1987). Dreizler, A., Lindenmaier, S., Maas, U., Hult, J., Alde´n, M., and Kaminski, C. F., Appl. Phys. B. 70:287 (2000). Sher, E., Ben-Yai’sch, J., and Kravchik, T., Combust. Flame 89:186 (1992). Scha¨fer, M., Schmidt, R., and Ko ¨hler, J., Twenty-sixth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1996, p. 2701. Warnatz, J., Maas, U., and Dibble, R. W., Combustion, 2nd edition, Springer Verlag, Berlin, 1999. Maas, U., and Warnatz, J., Impact Compu. Sci. Eng. 1:394 (1989). Yos, J. M., Transport Properties of Nitrogen, Hydrogen and Air to 30000 K, Technical Memorandum RADTM-63–7, AVCO Corporation, Wilmington (MA), March., 1963. Scha¨fer, M., Der Zu ¨ndfunke, Dissertation, Universita¨t Stuttgart, 1997. Maas, U., Dissertation, Heidelberg University, 1988. Noh, W. F., J. Comp. Phys. 72:78 (1987). Bird, R. B., Stewart, W. E., and Lightfoot, E. N., Transport Phenomena, John Wiley and Sons, New York, 1960. Thiele, M., Dissertation, University of Heidelberg, 1999. Deufflhard, P., Nowak, U., Tech. Rep. SFB 123, University of Heidelberg (1985). Hackbusch, W., Iterative Lo ¨sung groer schwachbesetzter Gleichungssysteme, Teubner Verlag, Stuttgart, 1992. Kaminski, C., Hult, J., and Alde´n, M., Appl. Phys. B 68:757 (1999). Lewis, B., and von Elbe, G., Combustion, Flames and Explosions of Gases, 3rd edition, Academic Press, Orlando., 1987. Ewart, P., Opt. Commun. 55:124 (1985). Snelling, D. R., et al., Appl. Opt. 33:8295 (1994). Eckbreth, A. C., Appl. Phys. Lett. 32:421 (1978). Farrow, R. L., Private communication. Maly, R. R., Vogel, M., Seventeenth Symposium (International) on Combustion, The Combustion Institute, Pittsburgh, 1978, p. 821. Boquillon, J. P., et al., Opt. Lett. 13:722 (1988). Thumann, A., et al., Appl. Opt. 34:3313 (1995). Parameswaran, T., et al., Appl. Opt. 35:5461 (1996). Thiele, M., Warnatz, J., and Maas, U., SAE Technical Paper Series No. 01-1178, 1999.
Received 16 October 2000; revised 16 May 2001; accepted 18 August 2001