Stochastic Simulation of Soot Formation Evolution in Counterflow Diffusion Flames

Soot generally refers to carbonaceous particles formed during incomplete combustion of hydrocarbon fuels. A typical simulation of soot formation and evolution contains two parts: gas chemical kinetics, which models the chemical reaction from hydrocarbon fuels to soot precursors, that is, polycyclic aromatic hydrocarbons or PAHs, and soot dynamics, which models the soot formation from PAHs and evolution due to gas-soot and soot-soot interactions. In this study, two detailed gas kinetic mechanisms (ABF and KM2) have been compared during the simulation (using the solver Chemkin II) of ethylene combustion in counterflow diffusion flames. Subsequently, the operator splitting Monte Carlo method is used to simulate the soot dynamics. Both the simulated data from the twomechanisms for gas and soot particles are compared with experimental data available in the literature. It is found that both mechanisms predict similar profiles for the gas temperature and velocity, agreeing well with measurements. However, KM2 mechanism provides much closer prediction compared to measurements for soot gas precursors. Furthermore, KM2 also shows much better predictions for soot number density and volume fraction than ABF..e effect of nozzle exit velocity on soot dynamics has also been investigated. Higher nozzle exit velocity renders shorter residence time for soot particles, which reduces the soot number density and volume fraction accordingly.


Introduction
Around 80%-85% world energy comes from combustion of fossil fuel [1]. e formation and evolution of soot (i.e., carbon particles resulting from incomplete combustion of hydrocarbons) is an important and constantly studied field in combustion due to its practical significance in the production of technical carbon (such as filler in rubber, component of printing paints), as well as in the combustion efficiency and human health [2,3]. Understanding the mechanism of soot formation is a long-standing challenge in combustion research. Quantitative knowledge of soot formation has been largely derived from three types of work [4]: measurement of soot volume fraction, number density, and particle size distributions (PSDs); development of detailed chemical mechanisms for the formation of polycyclic aromatic hydrocarbons; and development of soot population dynamics models to describe the evolution of the particle ensemble.
Diffusion flames, where fuel and air are introduced separately into the combustion chamber and form an ignitable mixture by diffusion, are widely used in practical combustion systems for safety reasons. Counterflow diffusion flames are frequently used in experimental and theoretical research because they represent essentially a one-dimensional structure of diffusion flame, which provides valuable data with respect to optimizing combustion processes and experimental data for validation of flame modelling.
Recently, the authors [13] developed an efficient operator splitting Monte Carlo method for the simulation of general aerosol dynamics.
e Monte Carlo method has been coupled with the chemical kinetics solver Chemkin II to simulate soot formation and evolution in a counter ow di usion ame [14]. is work is to further explore the simulation framework of coupling the kinetics solver with the stochastic method and to provide more detailed simulation results on the soot dynamics in counter ow di usion ames, so as to investigate the e ects of di erent popular chemical kinetic mechanisms and nozzle exit velocity on the soot formation and evolution.

Methodology
e simulation of soot formation and evolution in di usion ames are accomplished in two steps, that is, gaseous chemical kinetics simulation to determine the concentration of gaseous soot precursors (i.e., polycyclic aromatic hydrocarbons, PAHs) and stochastic simulation for soot particle dynamics.

Gaseous Chemical Kinetics.
Gaseous chemical kinetics are handled by the open source software Chemkin II [15]. In the current setting of counter ow di usion ame (Figure 1), fuel gas (C 2 H 4 ) and oxidizer eject from two opposing nozzles, respectively. When fuel and oxidizer mix in the middle, a stable sheet of ame forms. Along the nozzle axial direction, the model equations (derived from mass, momentum, and energy conservations) can be reduced to one dimension, and gas velocity, temperature, and reactant product concentration can be obtained through solving the 1D model equations with the Newton iteration method [16]. In Chemkin II, the chemical and physical properties of materials are stored in a library le, which can be adapted according to a user's need. A user provides the kinetic mechanism so as to determine the reaction route.
Here, two kinetic mechanisms (i.e., ABF and KM2) are used to investigate the e ects of the di erent mechanisms. e ABF mechanism [17] contains 101 species and 543 reactions, leading to the maximum product pyrene (C 16 H 10 ). e KM2 mechanism [18,19] contains 202 species and 1351 reactions, leading to the maximum product coronene (C 24 H 12 ).

Soot Dynamics.
e particle size distribution, n, of soot particle satis es the following general dynamic equation [20]: e right most three terms refer to the nucleation, growth, and coagulation, respectively. e di usion term can be neglected owing to large Schmidt number for soot particles [20]. e velocity on the left hand side is the gas velocity corrected by the thermophoretic e ect for particles [21]. e convection term can be implicitly solved with the introduction of the Lagrangian time t * [22,23].
en, the general dynamic equation is simpli ed as (the star notion in the Lagrangian time t * has been omitted) is transformation converts the Eulerian point of view in (1) to the Lagrangian point of view in (3).
Nucleation is the process of a large number of gas molecules forming a stable nucleus. e nucleation term is modelled as In the KM2 mechanism, there are 8 PAHs products. Any two PAHs molecules may nucleate to form a nascent soot particle, which results in 36 di erent nucleation processes. e above nucleation model is derived from a simple particle-particle collision model. A collision coe cient is assigned to every PAH, which is [0.006, 0.01, 0.01, 0.01, 0.011, 0.011, 0.014, and 0.02] for the 8 PAHs (ordered from the smallest to the biggest, see the legend in Figure 2). e pair collision coe cient is chosen as the smaller one between two collision PAHs.
Coagulation is the process that two particles collide and coalesce into a bigger particle. e coagulation term is modelled by the well-known Smoluchowski's equation Surface growth includes surface chemical reactions and physical condensation. e surface chemical reaction includes oxidation and the hydrogen-abstraction/acetyleneaddition (HACA) mechanism [24]. Oxidation is modelled by the reaction of soot particles with OH and O 2 molecules. e HACA process can be described as Condensation modelling is similar to nucleation (4), except that one of the two PAHs concentrations is replaced by the soot particle concentration.

Stochastic Simulation of Soot Dynamics.
e stochastic simulation adopts the operator splitting Monte Carlo method developed recently by the authors [13]. is method is numerically highly e cient and quite exible in accommodating various models for the general aerosol dynamics. Details can be found in the original paper.

Veri cation and Validation.
In principle, this simulation work used two methods: Chemkin for gas chemical kinetics and stochastic simulation tool for soot dynamics. Chemkin is the de facto standard for the simulation of chemical kinetics. A converged numerical solution is generally believed to be true for the model equations, and it requires no veri cation. Since the simulation results greatly depend on the kinetic mechanism used, most of the validation work aims to make accurate prediction as possible for quantities available in experimental measurements. e stochastic method used here for soot dynamics was developed by the authors, which has been thoroughly veri ed against classical testing cases [13], and has been successfully used to simulate soot dynamics in combustion ames [14,25] and general aerosol in turbulent ows [26,27].
In order to compare simulation results from the two kinetic mechanisms, that is, ABF and KM2, a nonpremixed ame of ethylene-oxygen in the canonical counter ow con guration is investigated. On the fuel side, the gas is a mixture of 75% C 2 H 4 and 25% Ar. On the oxidizer side, it is a mixture of 22% O 2 and 78% Ar. e nozzle gas exit velocities are 13.16 and 16.12 cm/s on the fuel side and oxidizer side, respectively. e corresponding strain rates are 9.4 and 11.51 #/s. ese settings are the same as those in the experiment [28]. Figure 3 shows the temperature pro les obtained from the two mechanisms. e measurement data from [28] is also included in the plot. Both mechanisms predict almost the same temperature pro le, very close to the measurement data. e peak position of the simulated temperature pro le di ers from the measurement only a little bit. Both mechanisms predict the temperature pro le quite satisfactorily.
Soot particles are formed from precursor PAH molecules. In the ABF mechanism, it contains four PAHs, denoted as A1 to A4, corresponding to benzene C 6 H 6 , naphthalene C 10 H 8 , phenanthrene C 14 H 10 , and pyrene C 16 H 10 , respectively. In the KM2 mechanism, much larger PAH molecules are included, up to coronene (C 24 H 12 ). e precursors mostly determine the soot particle number density through the nucleation process, that is, conversion of PAH gas to particles. Figure 4 compares the pro les of mole fractions of PAHs from A1 to A4 from simulations with the two mechanisms along with the measurement [29]. e KM2 mechanism gives better prediction than the ABF mechanism for all PAHs, when compared with the experiment data. e ABF mechanism is found to underpredict A4 considerably. It is worth noting that the comparison with experiment data for KM2 modelling is not completely satisfactory either.
ere are large uncertainties in combustion kinetic models [30], which may come from the extrapolation of knowledge of smaller species reaction, missing reaction pathways, uncertainty in the reaction  coe cient, and others. Figure 2 shows the mole fraction of larger PAHs, which is not available in the ABF mechanism. Coronene concentration is much higher than that of A4.
is fact has a very large impact on the soot volume faction. It is shown in Figure 5 that the ABF mechanism underpredicts the soot volume fraction considerably, while the result from KM2 compares rather well with the measurement data. Neglect of the larger PAHs in the ABF mechanism is the underlying reason for such underprediction. Overall, the KM2 mechanism gives much better predictions than the ABF mechanism, although with higher numerical cost due to much higher number of species in the simulation. In Figure 5, it is observed, both in measurement data and simulation, that the volume fraction drops from the peak much faster on the fuel side than on the oxidizer side. at is because nucleation and condensation renders the volume fraction to increase quickly towards the peak. e decrease of volume fraction o the peak on the oxidizer side is mostly due to gas transport (convection and di usion), which is a process much slower than the soot dynamics.

E ects of Nozzle Exit Velocity.
In the last subsection, it is shown that the KM2 mechanism gives much better predictions than the ABF mechanism. So from now on, only the simulation results from the KM2 mechanism are presented. Figure 6 shows the pro les of soot volume fraction under various nozzle exit velocity, from 20 to 40 cm/s (the same on  both sides). It is clear that higher exit velocity renders much lower soot volume fraction, due to the reduction of residence time of soot particle, which causes soot to have less time to form and grow. On the other hand, higher nozzle exit velocity means higher strain rate, and it is found [31] that increasing the strain rate could reduce concentrations of PAHs, which are the precursors to soot formation. When the nozzle exit velocity doubles from 20 to 40 cm/s, the magnitude of the corresponding volume fraction decreases sharply, by a factor nearly 10. e combined e ect of less newly nucleated particles and less growth time makes the nozzle exit velocity a crucial factor in determining the soot volume fraction. Figure 7 shows the average soot size, which is calculated thorough the third and zeroth order of moments (the moment is de ned with respect to the soot diameter) of the soot particle size distribution. When the nozzle exit velocity is smaller, soot particles are larger on average. is comes from the same fact of shorter residence time for higher exit velocity. When the nozzle exit velocity increases from 20 to 40 cm/s, the peak average soot diameter decreases from 5.1 to 2.5 nm, around a factor of 2, which corresponds to a factor of 8 on the decrease of volume fraction. From Figure 6, it can be noted that the magnitude of volume fraction for the case of 20 cm/s is nearly 10 times higher than that of 40 cm/s; hence, it can be concluded that the nozzle exit velocity has a strong e ect on the particle growth process, while it has only a relatively mild e ect on the nucleation process.
Soot particle size distribution is the most informative statistical description on the collective behavior of soot particles. Figure 8 shows the size distribution at various height locations from the fuel nozzle. Since the ame is stable, the distribution at various heights re ects the evolution of the distribution along the axial direction. From the temperature pro le in Figure 3, it is known that the ame sheet is at the height H 0.82 cm. Figure 8 shows the distribution at three locations, H 0.775, 0.75, and 0.675 cm. Further away from the ame sheet, nascent soot particles grow bigger due to surface growth and coagulation. A noticeable trough in the distribution is formed away from the ame sheet, when particles of intermediate size are greatly consumed due to coagulation while nascent nano-size particles are not able to grow quickly enough to make up the depletion. Such distribution has been observed very often in experiments.

Conclusions
Soot formation and evolution in counter ow di usion ames have been investigated by coupling the chemical kinetics Distance from fuel nozzle (cm) Soot volume fraction ×10 -6 ABF mechanism KM2 mechanism B. C. Choi solver Chemkin II with the operator splitting Monte Carlo method, which is an e cient stochastic method for simulating aerosol dynamics.
In this study, two detailed gas kinetic mechanisms (ABF [17] and KM2 [18,19]) have been compared during the simulation (using the solver Chemkin II) of ethylene combustion in counter ow di usion ames. Subsequently, the operator splitting Monte Carlo method is used to simulate the soot dynamics. Both the simulated data from the two mechanisms for gas and soot particles are compared with the experimental data available in the literature. It is found that both mechanisms predict similar pro les for the gas temperature and velocity, agreeing well with the measurements.
However, the KM2 mechanism provides much closer prediction for soot gas precursors, that is, PAHs, when compared with the measurements. In the ABF mechanism, the maximum PAH molecule is pyrene, while in the KM2 mechanism, the maximum PAH molecule is coronene. On one hand, coronene is found to have much higher propensity to nucleate than pyrene. On the other hand, coronene concentration observed in KM2 is found to be much higher than the pyrene concentration (which is the PAH of highest nucleation propensity in ABF). When compared with measurements on the soot number density and volume fraction, KM2 also shows much better predictions than ABF. e e ect of nozzle exit velocity on soot dynamics has also been investigated. Higher nozzle exit velocity renders shorter residence time for soot particles to form and grow, which reduces the soot number density and volume fraction accordingly. However, the effect of nozzle velocity on the number density is far more mild than that on the volume fraction, which means nucleation is slightly affected by the nozzle velocity, while surface growth is greatly affected by that. e residence time is a crucial factor to determine the soot particle size. Nomenclature n: Soot particle distribution function u → : Gas flow velocity D: Diffusion coefficient P vdw : Van der Waals constant of PAH k B : Boltzmann constant m 0 : Mass of the gas mixture N AVO : Avogadro constant d PAH : Diameter of PAH T: Temperature x: Volume fraction of PAH β(v, u): Collision kernel function, i.e., the rate of particle size v, coagulates with u particle size k s : Reactant gas concentration C g : Fraction of active surface α: Diffusion coefficient x s : Active site number density in surface reaction Δ: Variation of soot mass in elementary reaction S i : Surface area of soot particle.

Conflicts of Interest
e authors declare that they have no conflicts of interest.